help-gsl
[Top][All Lists]

## Re: [Help-gsl] Determinant of a matrix:

 From: Joakim Hove Subject: Re: [Help-gsl] Determinant of a matrix: Date: Fri, 11 Feb 2005 09:55:01 +0100 User-agent: Gnus/5.1006 (Gnus v5.10.6) Emacs/21.2 (gnu/linux)

```Nomesh Bolia <address@hidden> writes:

> Hi,
>
> I am using the following program to compute the determinant of a matrix A:
>
> double get_det(gsl_matrix * A) {
>   int sign=0; double det=0.0; int row_sq = mat_ptr->size1;
>   gsl_permutation * p = gsl_permutation_calloc(row_sq);
>   gsl_matrix * tmp_ptr = gsl_matrix_calloc(row_sq, row_sq);
>   int * signum = &sign;
>   gsl_matrix_memcpy(tmp_ptr, mat_ptr);
>   gsl_linalg_LU_decomp(tmp_ptr, p, signum);
>   det = gsl_linalg_LU_det(tmp_ptr, *signum);
>   gsl_permutation_free(p);
>   gsl_matrix_free(tmp_ptr);
>   return det;
> }
>
> Basically, the functions:
>
> gsl_linalg_LU_decomp(tmp_ptr, p, signum); and
> det = gsl_linalg_LU_det(tmp_ptr, *signum);
>
> have been used to compute the determinant. However, I get completely

Well,

that seems reasonable. You are determining the LU decomposition and det of
the newly allocated matrix pointed to by tmp_ptr, which is seemingly
not related to A - are the arguments of the gsl_matrix_memcpy() call
correct?

Why don't you just use A - I assume that is the matrix you are
interested in. Alternative code (not tested):

double det_get(gsl_matrix *A, int inPlace) {

/*
inPlace = 1 => A is replaced with the LU decomposed copy.
inPlace = 0 => A is retained, and a copy is used for LU.
*/

double det;
int signum;
gsl_permutation *p = gsl_perumtation_alloc(A->size1);
gsl_matrix *tmpA;

if (inPlace)
tmpA = A;
else {
gsl_matrix *tmpA = gsl_matrix_alloc(A->size1, A->size2);
gsl_matrix_memcpy(tmpA , A);
}

gsl_linalg_LU_decomp(tmpA , p , &signum);
det = gsl_linalg_LU_det(tmpA , signum);
gsl_permutation_free(p);
if (! inPlace)
gsl_matrix_free(tmpA);

return det;
}

HTH - Joakim

--
Joakim Hove
hove AT ift uib no
Tlf: +47 (55 5)8 27 90
Fax: +47 (55 5)8 94 40
http://www.ift.uib.no/~hove/

```