1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
|
#ifndef LUDECOMP_H
#define LUDECOMP_H
/*
* Copyright 2000 Graeme W. Gill
* All rights reserved.
*
* This material is licenced under the GNU AFFERO GENERAL PUBLIC LICENSE Version 3 :-
* see the License.txt file for licencing details.
*/
#ifdef __cplusplus
extern "C" {
#endif
/* NOTE:- lu decomp rearanges the rows of the matrix */
/* by swapping pointers rather than exchanging data, */
/* so the matrix must be addressed by the **pointer */
/* if it is re-used after an ludecomp!!! */
/* Solve the simultaneous linear equations A.X = B */
/* Return 1 if the matrix is singular, 0 if OK */
int
solve_se(
double **a, /* A[][] input matrix, returns LU decimposition of A */
double *b, /* B[] input array, returns solution X[] */
int n /* Dimensionality */
);
/* Solve the simultaneous linear equations A.X = B, with polishing */
/* Return 1 if the matrix is singular, 0 if OK */
int
polished_solve_se(
double **a, /* A[][] input matrix, returns LU decimposition of A */
double *b, /* B[] input array, returns solution X[] */
int n /* Dimensionality */
);
/* Decompose the square matrix A[][] into lower and upper triangles */
/* NOTE that rows get swaped by swapping matrix pointers! */
/* Return 1 if the matrix is singular. */
int
lu_decomp(
double **a, /* A input array, output upper and lower triangles. */
int n, /* Dimensionality */
int *pivx, /* Return pivoting row permutations record */
double *rip /* Row interchange parity, +/- 1.0, used for determinant */
);
/* Solve a set of simultaneous equations from the */
/* LU decomposition, by back substitution. */
void
lu_backsub(
double **a, /* A[][] LU decomposed matrix */
int n, /* Dimensionality */
int *pivx, /* Pivoting row permutations record */
double *b /* Input B[] vecor, return X[] */
);
/* Polish a solution for equations */
void
lu_polish(
double **a, /* Original A[][] matrix */
double **lua, /* LU decomposition of A[][] */
int n, /* Dimensionality */
double *b, /* B[] vector of equation */
double *x, /* X[] solution to be polished */
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(
double **a, /* A[][] input matrix, returns inversion of A */
int n /* Dimensionality */
);
/* 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 */
int n /* Dimensionality */
);
/* Pseudo-Invert matrix A using lu decomposition */
/* Return 1 if the matrix is singular, 0 if OK */
int
lu_psinvert(
double **out, /* Output[0..N-1][0..M-1] */
double **in, /* Input[0..M-1][0..N-1] input matrix */
int m, /* In Rows */
int n /* In Columns */
);
/* - - - - - - - - - - - - - - - - - - - */
/* Use Cholesky decomposition on a symetric positive-definite matrix. */
/* Only the upper triangle of the matrix A is accessed. */
/* L returns the decomposition */
/* Return nz if A is not positive-definite */
int llt_decomp(double **L, double **A, int n);
/* Solve a set of simultaneous equations A.x = b from the */
/* LLt decomposition, by back substitution. */
void llt_backsub(
double **L, /* A[][] LLt decomposition in lower triangle */
int n, /* Dimensionality */
double *b, /* Input B[] */
double *x /* Return X[] (may be same as B[]) */
);
#ifdef __cplusplus
}
#endif
#endif /* LUDECOMP_H */
|