lemon-project-template-glpk
comparison deps/glpk/src/glplux.h @ 11:4fc6ad2fb8a6
Test GLPK in src/main.cc
author | Alpar Juttner <alpar@cs.elte.hu> |
---|---|
date | Sun, 06 Nov 2011 21:43:29 +0100 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:e472f9972deb |
---|---|
1 /* glplux.h (LU-factorization, bignum arithmetic) */ | |
2 | |
3 /*********************************************************************** | |
4 * This code is part of GLPK (GNU Linear Programming Kit). | |
5 * | |
6 * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, | |
7 * 2009, 2010, 2011 Andrew Makhorin, Department for Applied Informatics, | |
8 * Moscow Aviation Institute, Moscow, Russia. All rights reserved. | |
9 * E-mail: <mao@gnu.org>. | |
10 * | |
11 * GLPK is free software: you can redistribute it and/or modify it | |
12 * under the terms of the GNU General Public License as published by | |
13 * the Free Software Foundation, either version 3 of the License, or | |
14 * (at your option) any later version. | |
15 * | |
16 * GLPK is distributed in the hope that it will be useful, but WITHOUT | |
17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY | |
18 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public | |
19 * License for more details. | |
20 * | |
21 * You should have received a copy of the GNU General Public License | |
22 * along with GLPK. If not, see <http://www.gnu.org/licenses/>. | |
23 ***********************************************************************/ | |
24 | |
25 #ifndef GLPLUX_H | |
26 #define GLPLUX_H | |
27 | |
28 #include "glpdmp.h" | |
29 #include "glpgmp.h" | |
30 | |
31 /*---------------------------------------------------------------------- | |
32 // The structure LUX defines LU-factorization of a square matrix A, | |
33 // which is the following quartet: | |
34 // | |
35 // [A] = (F, V, P, Q), (1) | |
36 // | |
37 // where F and V are such matrices that | |
38 // | |
39 // A = F * V, (2) | |
40 // | |
41 // and P and Q are such permutation matrices that the matrix | |
42 // | |
43 // L = P * F * inv(P) (3) | |
44 // | |
45 // is lower triangular with unity diagonal, and the matrix | |
46 // | |
47 // U = P * V * Q (4) | |
48 // | |
49 // is upper triangular. All the matrices have the order n. | |
50 // | |
51 // The matrices F and V are stored in row/column-wise sparse format as | |
52 // row and column linked lists of non-zero elements. Unity elements on | |
53 // the main diagonal of the matrix F are not stored. Pivot elements of | |
54 // the matrix V (that correspond to diagonal elements of the matrix U) | |
55 // are also missing from the row and column lists and stored separately | |
56 // in an ordinary array. | |
57 // | |
58 // The permutation matrices P and Q are stored as ordinary arrays using | |
59 // both row- and column-like formats. | |
60 // | |
61 // The matrices L and U being completely defined by the matrices F, V, | |
62 // P, and Q are not stored explicitly. | |
63 // | |
64 // It is easy to show that the factorization (1)-(3) is some version of | |
65 // LU-factorization. Indeed, from (3) and (4) it follows that: | |
66 // | |
67 // F = inv(P) * L * P, | |
68 // | |
69 // V = inv(P) * U * inv(Q), | |
70 // | |
71 // and substitution into (2) gives: | |
72 // | |
73 // A = F * V = inv(P) * L * U * inv(Q). | |
74 // | |
75 // For more details see the program documentation. */ | |
76 | |
77 typedef struct LUX LUX; | |
78 typedef struct LUXELM LUXELM; | |
79 typedef struct LUXWKA LUXWKA; | |
80 | |
81 struct LUX | |
82 { /* LU-factorization of a square matrix */ | |
83 int n; | |
84 /* the order of matrices A, F, V, P, Q */ | |
85 DMP *pool; | |
86 /* memory pool for elements of matrices F and V */ | |
87 LUXELM **F_row; /* LUXELM *F_row[1+n]; */ | |
88 /* F_row[0] is not used; | |
89 F_row[i], 1 <= i <= n, is a pointer to the list of elements in | |
90 i-th row of matrix F (diagonal elements are not stored) */ | |
91 LUXELM **F_col; /* LUXELM *F_col[1+n]; */ | |
92 /* F_col[0] is not used; | |
93 F_col[j], 1 <= j <= n, is a pointer to the list of elements in | |
94 j-th column of matrix F (diagonal elements are not stored) */ | |
95 mpq_t *V_piv; /* mpq_t V_piv[1+n]; */ | |
96 /* V_piv[0] is not used; | |
97 V_piv[p], 1 <= p <= n, is a pivot element v[p,q] corresponding | |
98 to a diagonal element u[k,k] of matrix U = P*V*Q (used on k-th | |
99 elimination step, k = 1, 2, ..., n) */ | |
100 LUXELM **V_row; /* LUXELM *V_row[1+n]; */ | |
101 /* V_row[0] is not used; | |
102 V_row[i], 1 <= i <= n, is a pointer to the list of elements in | |
103 i-th row of matrix V (except pivot elements) */ | |
104 LUXELM **V_col; /* LUXELM *V_col[1+n]; */ | |
105 /* V_col[0] is not used; | |
106 V_col[j], 1 <= j <= n, is a pointer to the list of elements in | |
107 j-th column of matrix V (except pivot elements) */ | |
108 int *P_row; /* int P_row[1+n]; */ | |
109 /* P_row[0] is not used; | |
110 P_row[i] = j means that p[i,j] = 1, where p[i,j] is an element | |
111 of permutation matrix P */ | |
112 int *P_col; /* int P_col[1+n]; */ | |
113 /* P_col[0] is not used; | |
114 P_col[j] = i means that p[i,j] = 1, where p[i,j] is an element | |
115 of permutation matrix P */ | |
116 /* if i-th row or column of matrix F is i'-th row or column of | |
117 matrix L = P*F*inv(P), or if i-th row of matrix V is i'-th row | |
118 of matrix U = P*V*Q, then P_row[i'] = i and P_col[i] = i' */ | |
119 int *Q_row; /* int Q_row[1+n]; */ | |
120 /* Q_row[0] is not used; | |
121 Q_row[i] = j means that q[i,j] = 1, where q[i,j] is an element | |
122 of permutation matrix Q */ | |
123 int *Q_col; /* int Q_col[1+n]; */ | |
124 /* Q_col[0] is not used; | |
125 Q_col[j] = i means that q[i,j] = 1, where q[i,j] is an element | |
126 of permutation matrix Q */ | |
127 /* if j-th column of matrix V is j'-th column of matrix U = P*V*Q, | |
128 then Q_row[j] = j' and Q_col[j'] = j */ | |
129 int rank; | |
130 /* the (exact) rank of matrices A and V */ | |
131 }; | |
132 | |
133 struct LUXELM | |
134 { /* element of matrix F or V */ | |
135 int i; | |
136 /* row index, 1 <= i <= m */ | |
137 int j; | |
138 /* column index, 1 <= j <= n */ | |
139 mpq_t val; | |
140 /* numeric (non-zero) element value */ | |
141 LUXELM *r_prev; | |
142 /* pointer to previous element in the same row */ | |
143 LUXELM *r_next; | |
144 /* pointer to next element in the same row */ | |
145 LUXELM *c_prev; | |
146 /* pointer to previous element in the same column */ | |
147 LUXELM *c_next; | |
148 /* pointer to next element in the same column */ | |
149 }; | |
150 | |
151 struct LUXWKA | |
152 { /* working area (used only during factorization) */ | |
153 /* in order to efficiently implement Markowitz strategy and Duff | |
154 search technique there are two families {R[0], R[1], ..., R[n]} | |
155 and {C[0], C[1], ..., C[n]}; member R[k] is a set of active | |
156 rows of matrix V having k non-zeros, and member C[k] is a set | |
157 of active columns of matrix V having k non-zeros (in the active | |
158 submatrix); each set R[k] and C[k] is implemented as a separate | |
159 doubly linked list */ | |
160 int *R_len; /* int R_len[1+n]; */ | |
161 /* R_len[0] is not used; | |
162 R_len[i], 1 <= i <= n, is the number of non-zero elements in | |
163 i-th row of matrix V (that is the length of i-th row) */ | |
164 int *R_head; /* int R_head[1+n]; */ | |
165 /* R_head[k], 0 <= k <= n, is the number of a first row, which is | |
166 active and whose length is k */ | |
167 int *R_prev; /* int R_prev[1+n]; */ | |
168 /* R_prev[0] is not used; | |
169 R_prev[i], 1 <= i <= n, is the number of a previous row, which | |
170 is active and has the same length as i-th row */ | |
171 int *R_next; /* int R_next[1+n]; */ | |
172 /* R_prev[0] is not used; | |
173 R_prev[i], 1 <= i <= n, is the number of a next row, which is | |
174 active and has the same length as i-th row */ | |
175 int *C_len; /* int C_len[1+n]; */ | |
176 /* C_len[0] is not used; | |
177 C_len[j], 1 <= j <= n, is the number of non-zero elements in | |
178 j-th column of the active submatrix of matrix V (that is the | |
179 length of j-th column in the active submatrix) */ | |
180 int *C_head; /* int C_head[1+n]; */ | |
181 /* C_head[k], 0 <= k <= n, is the number of a first column, which | |
182 is active and whose length is k */ | |
183 int *C_prev; /* int C_prev[1+n]; */ | |
184 /* C_prev[0] is not used; | |
185 C_prev[j], 1 <= j <= n, is the number of a previous column, | |
186 which is active and has the same length as j-th column */ | |
187 int *C_next; /* int C_next[1+n]; */ | |
188 /* C_next[0] is not used; | |
189 C_next[j], 1 <= j <= n, is the number of a next column, which | |
190 is active and has the same length as j-th column */ | |
191 }; | |
192 | |
193 #define lux_create _glp_lux_create | |
194 #define lux_decomp _glp_lux_decomp | |
195 #define lux_f_solve _glp_lux_f_solve | |
196 #define lux_v_solve _glp_lux_v_solve | |
197 #define lux_solve _glp_lux_solve | |
198 #define lux_delete _glp_lux_delete | |
199 | |
200 LUX *lux_create(int n); | |
201 /* create LU-factorization */ | |
202 | |
203 int lux_decomp(LUX *lux, int (*col)(void *info, int j, int ind[], | |
204 mpq_t val[]), void *info); | |
205 /* compute LU-factorization */ | |
206 | |
207 void lux_f_solve(LUX *lux, int tr, mpq_t x[]); | |
208 /* solve system F*x = b or F'*x = b */ | |
209 | |
210 void lux_v_solve(LUX *lux, int tr, mpq_t x[]); | |
211 /* solve system V*x = b or V'*x = b */ | |
212 | |
213 void lux_solve(LUX *lux, int tr, mpq_t x[]); | |
214 /* solve system A*x = b or A'*x = b */ | |
215 | |
216 void lux_delete(LUX *lux); | |
217 /* delete LU-factorization */ | |
218 | |
219 #endif | |
220 | |
221 /* eof */ |