[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
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
> arbitrary answers.
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/