summaryrefslogtreecommitdiff
path: root/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
diff options
context:
space:
mode:
Diffstat (limited to 'lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c')
-rw-r--r--lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c27
1 files changed, 22 insertions, 5 deletions
diff --git a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
index de3cc98..842b5ad 100644
--- a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
+++ b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
@@ -820,9 +820,8 @@ int lu_l_C(KCMAT(a), DVEC(ipiv), CMAT(r)) {
820 820
821 821
822//////////////////// LU substitution ///////////////////////// 822//////////////////// LU substitution /////////////////////////
823char charN = 'N';
824 823
825int luS_l_R(DMAT(a), KDVEC(ipiv), KDMAT(b), DMAT(x)) { 824int luS_l_R(KDMAT(a), KDVEC(ipiv), KDMAT(b), DMAT(x)) {
826 integer m = ar; 825 integer m = ar;
827 integer n = ac; 826 integer n = ac;
828 integer mrhs = br; 827 integer mrhs = br;
@@ -836,10 +835,28 @@ int luS_l_R(DMAT(a), KDVEC(ipiv), KDMAT(b), DMAT(x)) {
836 } 835 }
837 integer res; 836 integer res;
838 memcpy(xp,bp,mrhs*nrhs*sizeof(double)); 837 memcpy(xp,bp,mrhs*nrhs*sizeof(double));
839 integer ldb = mrhs; 838 dgetrs_ ("N",&n,&nrhs,(/*no const (!?)*/ double*)ap,&m,auxipiv,xp,&mrhs,&res);
840 integer lda = n;
841 dgetrs_ (&charN,&n,&nrhs,ap,&lda,auxipiv,xp,&ldb,&res);
842 CHECK(res,res); 839 CHECK(res,res);
843 free(auxipiv); 840 free(auxipiv);
844 OK 841 OK
845} 842}
843
844int luS_l_C(KCMAT(a), KDVEC(ipiv), KCMAT(b), CMAT(x)) {
845 integer m = ar;
846 integer n = ac;
847 integer mrhs = br;
848 integer nrhs = bc;
849
850 REQUIRES(m==n && m==mrhs && m==ipivn,BAD_SIZE);
851 integer* auxipiv = (integer*)malloc(n*sizeof(integer));
852 int k;
853 for (k=0; k<n; k++) {
854 auxipiv[k] = (integer)ipivp[k];
855 }
856 integer res;
857 memcpy(xp,bp,mrhs*nrhs*sizeof(doublecomplex));
858 zgetrs_ ("N",&n,&nrhs,(doublecomplex*)ap,&m,auxipiv,(doublecomplex*)xp,&mrhs,&res);
859 CHECK(res,res);
860 free(auxipiv);
861 OK
862}