ROBOTC.net forums
http://robotc.net/forums/

recursive calculation of matrix functions (eg. determinants)
http://robotc.net/forums/viewtopic.php?f=1&t=7547
Page 1 of 1

Author:  Ford Prefect [ Tue Dec 31, 2013 4:47 am ]
Post subject:  recursive calculation of matrix functions (eg. determinants)

hey,
I'm curious how to port the following recursive C code into RobotC (exemplarily: the determinant of square matrices):

Code:
// matrix determinant

double MatrixDet(int N, double A[N][N]) // or float instead
{
    int i,j,i_count,j_count, count=0;
    double Asub[N-1][N-1], det=0;

    if(N==1) return A[0][0];
    if(N==2) return (A[0][0]*A[1][1] - A[0][1]*A[1][0]);

    for(count=0; count<N; count++)
    {
        i_count=0;
        for(i=1; i<N; i++)
        {
            j_count=0;
            for(j=0; j<N; j++)
            {
                if(j == count) continue;
                Asub[i_count][j_count] = A[i][j];
                j_count++;
            }
            i_count++;
        }
        det += pow(-1, count) * A[0][count] * MatrixDet(N-1,Asub);
    }
    return det;
}

int main() {
 
double A[6][6]; // or float instead
 
 A[0][0]=1;
 A[0][1]=2;
 A[0][2]=3;
 A[0][3]=4;
 A[0][4]=5;
 A[0][5]=6;

 A[1][0]=10;
 A[1][1]=11;
 A[1][2]=12;
 A[1][3]=13;
 A[1][4]=14;
 A[1][5]=15;     

 A[2][0]=19;
 A[2][1]=20;
 A[2][2]=21;
 A[2][3]=22;
 A[2][4]=23;
 A[2][5]=24;

 A[3][0]=32;
 A[3][1]=33;
 A[3][2]=34;
 A[3][3]=35;
 A[3][4]=36;
 A[3][5]=37;

 A[4][0]=39;
 A[4][1]=40;
 A[4][2]=41;
 A[4][3]=42;
 A[4][4]=43;
 A[4][5]=44;

 A[5][0]=51;
 A[5][1]=52;
 A[5][2]=53;
 A[5][3]=54;
 A[5][4]=55;
 A[5][5]=56;

 int n=6;
 double det=MatrixDet(n, A);

 printf("%9.3f\n", det);
 return 1;
}

Author:  mightor [ Tue Dec 31, 2013 5:43 am ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

You need to find a way to get rid of the runtime sized array and then it would be fine.

= Xander

Author:  Ford Prefect [ Tue Dec 31, 2013 5:49 am ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

yes, that's the point. IIRC: you already tried it, too, didn't you?
But you can't "get rid of it", it's simply the way the methematical algorithm works.

ps,
the recursive computation of submatrices is not only an issue for determinants but also for matrix adjugates or inverse matrices !

Author:  mightor [ Tue Dec 31, 2013 9:01 am ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

Well, I didn't try very hard :) Mostly because I don't know enough about the algorithm or the actual operation that this is. I think there may be a way to do this, I'll see what I can come up with.

= Xander

Author:  mightor [ Tue Dec 31, 2013 9:27 am ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

Does this do what it's supposed to?

Code:
#define MAX_N 6
float MatrixDet(int N, float *A) // or float instead
{
   writeDebugStreamLine("MatrixDet: N=%d, A: %p", N, A);
   int i,j,i_count,j_count, count=0;
   float Asub[MAX_N * MAX_N];
   float det = 0;

   if(N==1) return A[0];
   if(N==2) return (A[0*N+0]*A[1*N+1] - A[0*N+1]*A[1*N+0]);

   for(count=0; count<N; count++)
   {
      i_count=0;
      for(i=1; i<N; i++)
      {
         j_count=0;
         for(j=0; j<N; j++)
         {
            if(j == count) continue;
            Asub[i*(N-1)+j] = A[i*N+j];
            j_count++;
         }
         i_count++;
      }
      det += pow(-1, count) * A[0*N+count] * MatrixDet(N-1, &Asub[0]);
   }
   return det;
}

= Xander

Author:  Ford Prefect [ Tue Dec 31, 2013 10:44 am ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

Do I see it correctly that you have to linearize the matrix?

what result do you get with the 6x6 matrix above?

---

and then as a test this one:

Code:
  A[0][0]=-1;
  A[0][1]= 2;
  A[0][2]=-3;
  A[0][3]= 0;
  A[0][4]= 0;
  A[0][5]= 0;

  A[1][0]= 2;
  A[1][1]= 1;
  A[1][2]= 0;
  A[1][3]= 0;
  A[1][4]= 0;
  A[1][5]= 0;

  A[2][0]= 4;
  A[2][1]=-2;
  A[2][2]= 0;
  A[2][3]= 0;
  A[2][4]=23;
  A[2][5]=24;

  A[3][0]= 0;
  A[3][1]= 0;
  A[3][2]= 0;
  A[3][3]= 2;
  A[3][4]= 1;
  A[3][5]=-1;

  A[4][0]= 0;
  A[4][1]= 0;
  A[4][2]=10;
  A[4][3]= 2;
  A[4][4]= 1;
  A[4][5]= 0;

  A[5][0]= 0;
  A[5][1]= 0;
  A[5][2]= 0;
  A[5][3]= 5;
  A[5][4]= 2;
  A[5][5]=-3;

Author:  Ford Prefect [ Tue Dec 31, 2013 11:14 am ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

for test reasons:
this is the program for the EV3:
Code:
// hw MatrixLib


#include <stdio.h>
#include <math.h>
#include <fcntl.h>
#include <string.h>
#include <sys/mman.h>
#include <sys/ioctl.h>

#include "ev3_button.h"
#include "ev3_lcd.h"
#include "ev3_constants.h"
#include "ev3_command.h"
#include "ev3_timer.h"
#include "ev3_lcd.h"
#include "ev3_sound.h"
#include "ev3_output.h"


double MatrixDet(int N, double A[N][N]) // or float instead
{
    int i,j,i_count,j_count, count=0;
    double Asub[N-1][N-1], det=0;

    if(N==1) return A[0][0];
    if(N==2) return (A[0][0]*A[1][1] - A[0][1]*A[1][0]);

    for(count=0; count<N; count++)
    {
        i_count=0;
        for(i=1; i<N; i++)
        {
            j_count=0;
            for(j=0; j<N; j++)
            {
                if(j == count) continue;
                Asub[i_count][j_count] = A[i][j];
                j_count++;
            }
            i_count++;
        }
        det += pow(-1, count) * A[0][count] * MatrixDet(N-1,Asub);
    }
    return det;
}

int main() {


  ButtonLedInit();
  LcdInit();
  LcdClean();


  double A[6][6]; // or float instead

  A[0][0]=-1;
  A[0][1]= 2;
  A[0][2]=-3;
  A[0][3]= 0;
  A[0][4]= 0;
  A[0][5]= 0;

  A[1][0]= 2;
  A[1][1]= 1;
  A[1][2]= 0;
  A[1][3]= 0;
  A[1][4]= 0;
  A[1][5]= 0;

  A[2][0]= 4;
  A[2][1]=-2;
  A[2][2]= 0;
  A[2][3]= 0;
  A[2][4]=23;
  A[2][5]=24;

  A[3][0]= 0;
  A[3][1]= 0;
  A[3][2]= 0;
  A[3][3]= 2;
  A[3][4]= 1;
  A[3][5]=-1;

  A[4][0]= 0;
  A[4][1]= 0;
  A[4][2]=10;
  A[4][3]= 2;
  A[4][4]= 1;
  A[4][5]= 0;

  A[5][0]= 0;
  A[5][1]= 0;
  A[5][2]= 0;
  A[5][3]= 5;
  A[5][4]= 2;
  A[5][5]=-3;

  int n=6;
  double det=MatrixDet(n, A);
  char buf[100];

  sprintf(buf, "%9.3f", det); LcdText(1, 0,20, buf);

  while(ButtonWaitForAnyPress(100) != BUTTON_ID_LEFT);
 

  LcdExit();
  ButtonLedExit();
  return 1;
}


this matrix has a det of 74, the det of the matrix of the TOP is 0.
(very exact results!)

Author:  mightor [ Tue Dec 31, 2013 1:19 pm ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

This program seems to produce the required results:

Code:
// matrix determinant

#define MAX_N 6

// taken from
float Determinant(int n, float *a)
{
  int i,j,j1,j2;
  float det = 0;
  float m[MAX_N * MAX_N];

  if (n < 1)
  {
    /* Error */
  }
  else if (n == 1)
  { /* Shouldn't get used */
    det = a[0*n+0];
  } else if (n == 2) {
    det = a[0*n+0] * a[1*n+1] - a[1*n+0] * a[0*n+1];
  } else {
    det = 0;
    for (j1=0;j1<n;j1++)
    {
      for (i=1;i<n;i++)
      {
        j2 = 0;
        for (j=0;j<n;j++)
        {
          if (j == j1)
            continue;
          m[(i-1)*(n-1)+j2] = a[i*n+j];
          j2++;
        }
      }
      det += pow(-1.0,1.0+j1+1.0) * a[0*n+j1] * Determinant(n-1, m);
    }
  }
  return(det);
}

task main() {
  float A[6][6];

  A[0][0]=1;
  A[0][1]=2;
  A[0][2]=3;
  A[0][3]=4;
  A[0][4]=5;
  A[0][5]=6;

  A[1][0]=10;
  A[1][1]=11;
  A[1][2]=12;
  A[1][3]=13;
  A[1][4]=14;
  A[1][5]=15;

  A[2][0]=19;
  A[2][1]=20;
  A[2][2]=21;
  A[2][3]=22;
  A[2][4]=23;
  A[2][5]=24;

  A[3][0]=32;
  A[3][1]=33;
  A[3][2]=34;
  A[3][3]=35;
  A[3][4]=36;
  A[3][5]=37;

  A[4][0]=39;
  A[4][1]=40;
  A[4][2]=41;
  A[4][3]=42;
  A[4][4]=43;
  A[4][5]=44;

  A[5][0]=51;
  A[5][1]=52;
  A[5][2]=53;
  A[5][3]=54;
  A[5][4]=55;
  A[5][5]=56;

  float det1 = Determinant(6, A);


  A[0][0]=-1;
  A[0][1]= 2;
  A[0][2]=-3;
  A[0][3]= 0;
  A[0][4]= 0;
  A[0][5]= 0;

  A[1][0]= 2;
  A[1][1]= 1;
  A[1][2]= 0;
  A[1][3]= 0;
  A[1][4]= 0;
  A[1][5]= 0;

  A[2][0]= 4;
  A[2][1]=-2;
  A[2][2]= 0;
  A[2][3]= 0;
  A[2][4]=23;
  A[2][5]=24;

  A[3][0]= 0;
  A[3][1]= 0;
  A[3][2]= 0;
  A[3][3]= 2;
  A[3][4]= 1;
  A[3][5]=-1;

  A[4][0]= 0;
  A[4][1]= 0;
  A[4][2]=10;
  A[4][3]= 2;
  A[4][4]= 1;
  A[4][5]= 0;

  A[5][0]= 0;
  A[5][1]= 0;
  A[5][2]= 0;
  A[5][3]= 5;
  A[5][4]= 2;
  A[5][5]=-3;

  float det2 = Determinant(6, A);

  writeDebugStreamLine("%9.3f\n", det1);
  writeDebugStreamLine("%9.3f\n", det2);
}

Author:  Ford Prefect [ Tue Dec 31, 2013 2:28 pm ]
Post subject:  Re: recursive calculation of matrix functions (eg. determina

yes, great, thanks, and that's actually almost exactly my code except this original ** thing (and your linearized adaption)! :)

happy new year everybody

Page 1 of 1 All times are UTC - 5 hours [ DST ]
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group
http://www.phpbb.com/