Finding the determinant of a matrix

**arnaudk** · Aug 15 '08, 09:52 AM

As noted in a recent thread, finding the determinant by a minors expansion is very slow (O(n!) for an nxn matrix), there are much better methods like the finding the LU decomposition which is O(n^3). Also, there are a number of libraries which can do all this for you, such as GSL.

**curiously enough** · Aug 15 '08, 10:18 AM

Originally posted by arnaudk

As noted in a recent thread, finding the determinant by a minors expansion is very slow (O(n!) for an nxn matrix), there are much better methods like the finding the LU decomposition which is O(n^3). Also, there are a number of libraries which can do all this for you, such as GSL.

I know there are other methods , but I want to use this one, so please read the program and try it and maybe you'll know what's wrong with it.

**arnaudk** · Aug 15 '08, 10:50 AM

What kind of problems are you having with your code? Does it compile? Does it give you the wrong answer?

**curiously enough** · Aug 15 '08, 11:03 AM

Originally posted by arnaudk

What kind of problems are you having with your code? Does it compile? Does it give you the wrong answer?

it compiles,but It gives the wrong answer for square matrices of order greater than or equal to 4

**arnaudk** · Aug 15 '08, 12:09 PM

Since there is probably an arithmetic error in your code, why not get the function DET to print out the matrix that it is given in a nice tabular form, as well as the cofactor so that you can check that everything is working as expected and you are getting the correct cofactor and associated minor?

Some general suggestions:

The return type of main is int, not void.
Use more intuitive variable names, if i represents a column, call it col.
Use brief comments to explain what's going on if you ever want anyone else to lay eyes on your code.
Scrunching up code like if(i==0&&DIM==d im)for(j=0;j<10 0;j++)M[j]=100; makes it hard to read.
Use const wherever you don't expect a variable to be changed, this is to prevent accidentally changing it and to enable some compiler optimizations: float DET(const float A[100][100], const int dim)
Using static variables in a recursively called function is arcane. Rather than using static variables, why not pass the variables to each new function call, in this case, pass det and DIM as arguments to the function DET(), the first call to DET will be called with arguments det=0 and DIM=dim.
You don't need a special case for dim==2; you can do a cofactor expansion on a 2x2 matrix as well, the minors in this case are just rank 1 matrices (i.e. numbers).

**newb16** · Aug 15 '08, 12:45 PM

Why is det int while A is float?

**curiously enough** · Aug 15 '08, 02:09 PM

Originally posted by arnaudk

Since there is probably an arithmetic error in your code, why not get the function DET to print out the matrix that it is given in a nice tabular form, as well as the cofactor so that you can check that everything is working as expected and you are getting the correct cofactor and associated minor?

Some general suggestions:

The return type of main is int, not void.
Use more intuitive variable names, if i represents a column, call it col.
Use brief comments to explain what's going on if you ever want anyone else to lay eyes on your code.
Scrunching up code like if(i==0&&DIM==d im)for(j=0;j<10 0;j++)M[j]=100; makes it hard to read.
Use const wherever you don't expect a variable to be changed, this is to prevent accidentally changing it and to enable some compiler optimizations: float DET(const float A[100][100], const int dim)
Using static variables in a recursively called function is arcane. Rather than using static variables, why not pass the variables to each new function call, in this case, pass det and DIM as arguments to the function DET(), the first call to DET will be called with arguments det=0 and DIM=dim.
You don't need a special case for dim==2; you can do a cofactor expansion on a 2x2 matrix as well, the minors in this case are just rank 1 matrices (i.e. numbers).

Thanks for your first advice, I printed all the minors like you said but still can't figure out what's wrong, and sorry if the randomly chosen variable names cuased you a headache, so if there is something unclear to you I will clear it up. Here's the program:
[code=c]
#include<stdio. h>
#include<conio. h>
#include<math.h >
int M[100];
float DET(float A[][100],int dim)
{
int i=100,j,k,l,m;
static float det=0;
static int DIM=dim;
if(dim>2)
for(i=0;i<DIM;i ++)
{
if(i==0&&DIM==d im)for(j=0;j<10 0;j++)M[j]=100;
for(j=DIM-1;j>=2;j--)
if(M[j]==i)continue;
/*here printing starts*/
for(k=0;k<dim;k ++)
{
printf("\n");
for(j=0;j<DIM;j ++)
{
m=0;
for(l=DIM-1;l>=dim;l--)
if(j!=M[l])m++;
if(m==DIM-dim)printf("%0. 0f ",A[k][j]);
}
}
printf("\n\n");
/*here printing ends*/
M[dim-1]=i;
det+=pow(-1,i+dim-1)*A[dim-1][i]*DET(A,dim-1);
}
if(i==DIM-1)M[dim-1]=100;
if(dim==2)
{
/*here printing starts*/
for(k=0;k<dim;k ++)
{
printf("\n");
for(j=0;j<DIM;j ++)
{
m=0;
for(l=DIM-1;l>=dim;l--)
if(j!=M[l])m++;
if(m==DIM-dim)printf("%0. 0f ",A[k][j]);
}
}
printf("\n\n");
/*here printing ends*/
for(i=0;i<DIM;i ++)
{
k=0;
for(j=DIM-1;j>=2;j--)
if(i!=M[j])k++;
if(k==DIM-2)break;
}
for(j=0;j<DIM;j ++)
{
k=0;
for(l=DIM-1;l>=2;l--)
if(j!=M[l])k++;
if(k==DIM-2&&j!=i)break ;
}
return(A[0][i]*A[1][j]-A[0][j]*A[1][i]);
}
if(dim==1)retur n(A[0][0]);
return(det);
}

void main()
{
float A[100][100];
int n,k,l;
puts("Input size of square matrix:");
scanf("%d",&n);
puts("Input square matrix:");
for(k=0;k<n;k++ )
{
printf("\n");
for(l=0;l<n;l++ )
{
printf("Input a%d%d: ",k+1,l+1);
scanf("%f",&A[k][l]);
}
}
printf("\n\nThe determenant of the matrix is %0.2f",DET(A,n) );
getch();
}
[/code]

**arnaudk** · Aug 15 '08, 04:20 PM

So... were the minors displayed correct or not? For example, run your program on this invertible matrix and check the minors and cofactors are correct and have the right sign

Code:

.
   (1 2 3 4)
det(1 2 1 1) = -12
   (1 1 2 1)
   (5 6 7 8)

and [code]...[/code] tags are good. Please use them.

**JosAH** · Aug 15 '08, 04:53 PM

Originally posted by arnaudk

So... were the minors displayed were correct or not? For example, run your program on this invertible matrix and check the minors and cofactors are correct and have the right sign

Code:

.
   (1 2 3 4)
det(1 2 1 1) = -12
   (1 1 2 1)
   (5 6 7 8)

and [code]...[/code] tags are good. Please use them.

Yup correct (using my LUP function) ;-)

kind regards,

Jos

**arnaudk** · Aug 15 '08, 05:09 PM

Originally posted by JosAH

Yup correct (using my LUP function) ;-)

Phew! I usually prefer to work in special units where 1 = -1 = 2 = 1/2. In those units I never have to hunt down minus signs or factors of 2 and my answers are always correct. :-)

**JosAH** · Aug 15 '08, 05:33 PM

Originally posted by arnaudk

Phew! I usually prefer to work in special units where 1 = -1 = 2 = 1/2. In those units I never have to hunt down minus signs or factors of 2 and my answers are always correct. :-)

Well, my LUP implementation isn't the most efficient one when it comes to pivot
selection but then again it has to deal with symbolics as well and it does a fine
job when it comes to that e.g.

Code:

| a b |
| c d | = a*d-b*c

It can do any larger dimensions as well ;-)

kind regards,

Jos

**curiously enough** · Aug 17 '08, 06:24 AM

Originally posted by arnaudk

So... were the minors displayed correct or not? For example, run your program on this invertible matrix and check the minors and cofactors are correct and have the right sign

Code:

.
   (1 2 3 4)
det(1 2 1 1) = -12
   (1 1 2 1)
   (5 6 7 8)

and [code]...[/code] tags are good. Please use them.

OK, now I've perfected the program, and all the minors are displayed correctly, but It still gives the wrong answer, and I'm still trying to find out what's wrong.
here's the program, test it on C++(V5.02) if you have it:

Code:

#include<stdio.h>
#include<conio.h>
#include<math.h>
int M[100];
float DET(float A[][100],int dim)
{
int i,j,k,l,m;
static float det=0;
static int DIM=dim;
if(dim>2)
	for(i=0;i<DIM;i++)
	{
	if(i==0&&DIM==dim)for(j=0;j<100;j++)M[j]=100;
    m=0;
	for(j=DIM-1;j>=2;j--)
	if(M[j]==i)m++;
    if(m!=0)continue;
   	/*here printing starts*/
   	for(k=0;k<dim;k++)
   	{
	printf("\n");
   	for(j=0;j<DIM;j++)
		{
		m=0;
		for(l=DIM-1;l>=2;l--)
		if(j!=M[l])m++;
		if(m==DIM-2)printf("%0.0f   ",A[k][j]);
		}
    }
   	printf("\n\n");
    /*here printing ends*/
	M[dim-1]=i;
	det+=pow(-1,i+dim-1)*A[dim-1][i]*DET(A,dim-1);
    M[dim-1]=100;
	}
if(dim==2)
 {
    /*here printing starts*/
   	for(k=0;k<dim;k++)
   	{
	printf("\n");
   	    for(j=0;j<DIM;j++)
		{
		m=0;
		for(l=DIM-1;l>=2;l--)
		if(j!=M[l])m++;
		if(m==DIM-2)printf("%0.0f   ",A[k][j]);
		}
    }
   	printf("\n\n");
    /*here printing ends*/
	for(i=0;i<DIM;i++)
	{
	k=0;
  	for(j=DIM-1;j>=2;j--)
	if(i!=M[j])k++;
	if(k==DIM-2)break;
	}
	for(j=0;j<DIM;j++)
	{
	k=0;
	for(l=DIM-1;l>=2;l--)
	if(j!=M[l])k++;
	if(k==DIM-2&&j!=i)break;
	}
	return(A[0][i]*A[1][j]-A[0][j]*A[1][i]);
 }
if(dim==1)return(A[0][0]);
return(det);
}


void main()
{
float A[100][100];
int n,k,l;
puts("Input size of square matrix:");
scanf("%d",&n);
puts("Input square matrix:");
for(k=0;k<n;k++)
{
printf("\n");
for(l=0;l<n;l++)
{
printf("Input a%d%d: ",k+1,l+1);
scanf("%f",&A[k][l]);
}
}
printf("\n\nThe determenant of the matrix is %0.2f",DET(A,n));
getch();
}

I just have one question:
When I write
det+=pow(-1,i+dim-1)*A[dim-1][i]*DET(A,dim-1); where the return value of the DET(A,dim-1) function is det, is the value added to det just (sign)*A[dim-1][i] multiplied by the determinant of its minor?
P.S: I tried your matrix, and got -8550.

**JosAH** · Aug 17 '08, 10:48 AM

Originally posted by curiously enough

When I write
det+=pow(-1,i+dim-1)*A[dim-1][i]*DET(A,dim-1); where the return value of the DET(A,dim-1) function is det, is the value added to det just (sign)*A[dim-1][i] multiplied by the determinant of its minor?

When you write X+= Y then the entire expression Y is evaluated and its value
is added to X for whatever the values of X and Y here the new value of X is
Y larger than it was before.

kind regards,

Jos

ps. you definitely should use other variable names and use a proper indentation
and spaces to make your code more readable.

**curiously enough** · Aug 17 '08, 11:14 AM

Originally posted by JosAH

When you write X+= Y then the entire expression Y is evaluated and its value
is added to X for whatever the values of X and Y here the new value of X is
Y larger than it was before.

kind regards,

Jos

ps. you definitely should use other variable names and use a proper indentation
and spaces to make your code more readable.

You didn't understand what I meant.
When I add to the value of the det static variable a return value which is that very same static variable(which is happening in the function as long as dim-1 is not equal to 2), will it add itself twice? Or will it just add the return value associated with the dim-1 part of the recursivity to the initial value of det in the part of the recursivity where I made the statement det+=pow(-1,i+dim-1)*A[dim-1][i]*DET(A,dim-1); ?

Finding the determinant of a matrix

Finding the determinant of a matrix

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment

Comment