summaryrefslogtreecommitdiff
path: root/lib/Numeric/LinearAlgebra/LAPACK
diff options
context:
space:
mode:
authorReiner Pope <unknown>2008-10-24 08:59:02 +0000
committerReiner Pope <unknown>2008-10-24 08:59:02 +0000
commitd8639b28ec9e83b54b45c987508d270d5469451c (patch)
treeed660ba16b4dd21bd9608fbbde7cc4dcba93f31c /lib/Numeric/LinearAlgebra/LAPACK
parentce8fed3a3558468b128a03cc4c96aa6c11357b4d (diff)
added dgetrs support for solving linear systems from LU factorization
Diffstat (limited to 'lib/Numeric/LinearAlgebra/LAPACK')
-rw-r--r--lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c29
-rw-r--r--lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h2
2 files changed, 31 insertions, 0 deletions
diff --git a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
index 310f6ee..de3cc98 100644
--- a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
+++ b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
@@ -20,6 +20,9 @@
20#define OK return 0; 20#define OK return 0;
21#endif 21#endif
22 22
23#define TRACEMAT(M) {int q; printf(" %d x %d: ",M##r,M##c); \
24 for(q=0;q<M##r*M##c;q++) printf("%.1f ",M##p[q]); printf("\n");}
25
23#define CHECK(RES,CODE) MACRO(if(RES) return CODE;) 26#define CHECK(RES,CODE) MACRO(if(RES) return CODE;)
24 27
25#define BAD_SIZE 2000 28#define BAD_SIZE 2000
@@ -814,3 +817,29 @@ int lu_l_C(KCMAT(a), DVEC(ipiv), CMAT(r)) {
814 free(auxipiv); 817 free(auxipiv);
815 OK 818 OK
816} 819}
820
821
822//////////////////// LU substitution /////////////////////////
823char charN = 'N';
824
825int luS_l_R(DMAT(a), KDVEC(ipiv), KDMAT(b), DMAT(x)) {
826 integer m = ar;
827 integer n = ac;
828 integer mrhs = br;
829 integer nrhs = bc;
830
831 REQUIRES(m==n && m==mrhs && m==ipivn,BAD_SIZE);
832 integer* auxipiv = (integer*)malloc(n*sizeof(integer));
833 int k;
834 for (k=0; k<n; k++) {
835 auxipiv[k] = (integer)ipivp[k];
836 }
837 integer res;
838 memcpy(xp,bp,mrhs*nrhs*sizeof(double));
839 integer ldb = mrhs;
840 integer lda = n;
841 dgetrs_ (&charN,&n,&nrhs,ap,&lda,auxipiv,xp,&ldb,&res);
842 CHECK(res,res);
843 free(auxipiv);
844 OK
845}
diff --git a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h
index 79e52be..e98e9dc 100644
--- a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h
+++ b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h
@@ -84,3 +84,5 @@ int schur_l_C(KCMAT(a), CMAT(u), CMAT(s));
84 84
85int lu_l_R(KDMAT(a), DVEC(ipiv), DMAT(r)); 85int lu_l_R(KDMAT(a), DVEC(ipiv), DMAT(r));
86int lu_l_C(KCMAT(a), DVEC(ipiv), CMAT(r)); 86int lu_l_C(KCMAT(a), DVEC(ipiv), CMAT(r));
87
88int luS_l_R(DMAT(a), KDVEC(ipiv), KDMAT(b), DMAT(x));