summaryrefslogtreecommitdiff
path: root/numlib/ludecomp.c
diff options
context:
space:
mode:
authorJörg Frings-Fürst <debian@jff-webhosting.net>2017-12-03 20:38:41 +0100
committerJörg Frings-Fürst <debian@jff-webhosting.net>2017-12-03 20:38:41 +0100
commitba627dd9ecb578e9852c7b9cce67ec63199d1acf (patch)
tree27c4258311ca8c8ed7ff67a8a0bc7280e8fcae79 /numlib/ludecomp.c
parent69aec3b712232e93600ecd741269fed1f90b412a (diff)
parent3abb40d43649adb3807180692d8579c405524675 (diff)
Merge branch 'release/2.0.0+repack-1'2.0.0+repack-1
Diffstat (limited to 'numlib/ludecomp.c')
-rwxr-xr-x[-rw-r--r--]numlib/ludecomp.c16
1 files changed, 12 insertions, 4 deletions
diff --git a/numlib/ludecomp.c b/numlib/ludecomp.c
index 72e014a..b756bfe 100644..100755
--- a/numlib/ludecomp.c
+++ b/numlib/ludecomp.c
@@ -162,6 +162,8 @@ int n /* Dimensionality */
}
/* Decompose the square matrix A[][] into lower and upper triangles */
+/* NOTE that it returns transposed inverse by normal convention. */
+/* Use sym_matrix_trans() to fix this. */
/* Return 1 if the matrix is singular. */
int
lu_decomp(
@@ -339,6 +341,8 @@ int *pivx /* Pivoting row permutations record */
/* Invert a matrix A using lu decomposition */
+/* NOTE that it returns transposed inverse by normal convention. */
+/* Use sym_matrix_trans() to fix this, or use matrix_trans_mult() */
/* Return 1 if the matrix is singular, 0 if OK */
int
lu_invert(
@@ -385,7 +389,10 @@ int n /* Dimensionality */
return 0;
}
-#ifdef NEVER /* It's not clear that this is correct */
+/* Invert a matrix A using lu decomposition, and polish it. */
+/* NOTE that it returns transposed inverse by normal convention. */
+/* Use sym_matrix_trans() to fix this, or use matrix_trans_mult() */
+/* Return 1 if the matrix is singular, 0 if OK */
int
lu_polished_invert(
double **a, /* A[][] input matrix, returns inversion of A */
@@ -413,8 +420,8 @@ int n /* Dimensionality */
return i;
}
- for (k = 0; k < 10; k++) {
- matrix_mult(t1, n, n, aa, n, n, a, n, n);
+ for (k = 0; k < 20; k++) {
+ matrix_trans_mult(t1, n, n, aa, n, n, a, n, n);
for (i = 0; i < n; i++) {
for (j = 0; j < n; j++) {
t2[i][j] = a[i][j];
@@ -432,9 +439,10 @@ int n /* Dimensionality */
free_dmatrix(t2, 0, n-1, 0, n-1);
return 0;
}
-#endif
/* Pseudo-Invert matrix A using lu decomposition */
+/* NOTE that it returns transposed inverse by normal convention. */
+/* Use matrix_trans() to fix this, or use matrix_trans_mult(). */
/* Return 1 if the matrix is singular, 0 if OK */
int
lu_psinvert(