1 /* -*- mode: C -*- */
2 /*
3 IGraph library.
4 Copyright (C) 2009-2012 Gabor Csardi <csardi.gabor@gmail.com>
5 334 Harvard st, Cambridge MA, 02139 USA
6
7 This program is free software; you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation; either version 2 of the License, or
10 (at your option) any later version.
11
12 This program is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with this program; if not, write to the Free Software
19 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
20 02110-1301 USA
21
22 */
23
24 #define NCOMPLEX /* to make it compile with MSVC on Windows */
25
26 #include <cs/cs.h>
27 #include <igraph.h>
28
check_solution(const igraph_sparsemat_t * A,const igraph_vector_t * x,const igraph_vector_t * b)29 igraph_bool_t check_solution(const igraph_sparsemat_t *A,
30 const igraph_vector_t *x,
31 const igraph_vector_t *b) {
32
33 long int dim = igraph_vector_size(x);
34 igraph_vector_t res;
35 int j, p;
36 igraph_real_t min, max;
37
38 igraph_vector_copy(&res, b);
39
40 for (j = 0; j < dim; j++) {
41 for (p = A->cs->p[j]; p < A->cs->p[j + 1]; p++) {
42 long int from = A->cs->i[p];
43 igraph_real_t value = A->cs->x[p];
44 VECTOR(res)[from] -= VECTOR(*x)[j] * value;
45 }
46 }
47
48 igraph_vector_minmax(&res, &min, &max);
49 igraph_vector_destroy(&res);
50
51 return fabs(min) < 1e-15 && fabs(max) < 1e-15;
52 }
53
main()54 int main() {
55
56 igraph_sparsemat_t A, B, C;
57 igraph_vector_t b, x;
58 long int i;
59
60 /* lsolve */
61
62 #define DIM 10
63 #define EDGES (DIM*DIM/6)
64 igraph_sparsemat_init(&A, DIM, DIM, EDGES + DIM);
65 for (i = 0; i < DIM; i++) {
66 igraph_sparsemat_entry(&A, i, i, RNG_INTEGER(1, 3));
67 }
68 for (i = 0; i < EDGES; i++) {
69 long int r = RNG_INTEGER(0, DIM - 1);
70 long int c = RNG_INTEGER(0, r);
71 igraph_real_t value = RNG_INTEGER(1, 5);
72 igraph_sparsemat_entry(&A, r, c, value);
73 }
74 igraph_sparsemat_compress(&A, &B);
75 igraph_sparsemat_destroy(&A);
76 igraph_sparsemat_dupl(&B);
77
78 igraph_vector_init(&b, DIM);
79 for (i = 0; i < DIM; i++) {
80 VECTOR(b)[i] = RNG_INTEGER(1, 10);
81 }
82
83 igraph_vector_init(&x, DIM);
84 igraph_sparsemat_lsolve(&B, &b, &x);
85
86 if (! check_solution(&B, &x, &b)) {
87 return 1;
88 }
89
90 igraph_vector_destroy(&b);
91 igraph_vector_destroy(&x);
92 igraph_sparsemat_destroy(&B);
93
94 #undef DIM
95 #undef EDGES
96
97 /* ltsolve */
98
99 #define DIM 10
100 #define EDGES (DIM*DIM/6)
101 igraph_sparsemat_init(&A, DIM, DIM, EDGES + DIM);
102 for (i = 0; i < DIM; i++) {
103 igraph_sparsemat_entry(&A, i, i, RNG_INTEGER(1, 3));
104 }
105 for (i = 0; i < EDGES; i++) {
106 long int r = RNG_INTEGER(0, DIM - 1);
107 long int c = RNG_INTEGER(0, r);
108 igraph_real_t value = RNG_INTEGER(1, 5);
109 igraph_sparsemat_entry(&A, r, c, value);
110 }
111 igraph_sparsemat_compress(&A, &B);
112 igraph_sparsemat_destroy(&A);
113 igraph_sparsemat_dupl(&B);
114
115 igraph_vector_init(&b, DIM);
116 for (i = 0; i < DIM; i++) {
117 VECTOR(b)[i] = RNG_INTEGER(1, 10);
118 }
119
120 igraph_vector_init(&x, DIM);
121 igraph_sparsemat_ltsolve(&B, &b, &x);
122
123 igraph_sparsemat_transpose(&B, &A, /*values=*/ 1);
124 if (! check_solution(&A, &x, &b)) {
125 return 2;
126 }
127
128 igraph_vector_destroy(&b);
129 igraph_vector_destroy(&x);
130 igraph_sparsemat_destroy(&B);
131 igraph_sparsemat_destroy(&A);
132
133 #undef DIM
134 #undef EDGES
135
136 /* usolve */
137
138 #define DIM 10
139 #define EDGES (DIM*DIM/6)
140 igraph_sparsemat_init(&A, DIM, DIM, EDGES + DIM);
141 for (i = 0; i < DIM; i++) {
142 igraph_sparsemat_entry(&A, i, i, RNG_INTEGER(1, 3));
143 }
144 for (i = 0; i < EDGES; i++) {
145 long int r = RNG_INTEGER(0, DIM - 1);
146 long int c = RNG_INTEGER(0, r);
147 igraph_real_t value = RNG_INTEGER(1, 5);
148 igraph_sparsemat_entry(&A, r, c, value);
149 }
150 igraph_sparsemat_compress(&A, &B);
151 igraph_sparsemat_destroy(&A);
152 igraph_sparsemat_dupl(&B);
153 igraph_sparsemat_transpose(&B, &A, /*values=*/ 1);
154
155 igraph_vector_init(&b, DIM);
156 for (i = 0; i < DIM; i++) {
157 VECTOR(b)[i] = RNG_INTEGER(1, 10);
158 }
159
160 igraph_vector_init(&x, DIM);
161 igraph_sparsemat_usolve(&A, &b, &x);
162
163 if (! check_solution(&A, &x, &b)) {
164 return 3;
165 }
166
167 igraph_vector_destroy(&b);
168 igraph_vector_destroy(&x);
169 igraph_sparsemat_destroy(&B);
170 igraph_sparsemat_destroy(&A);
171
172 #undef DIM
173 #undef EDGES
174
175 /* utsolve */
176
177 #define DIM 10
178 #define EDGES (DIM*DIM/6)
179 igraph_sparsemat_init(&A, DIM, DIM, EDGES + DIM);
180 for (i = 0; i < DIM; i++) {
181 igraph_sparsemat_entry(&A, i, i, RNG_INTEGER(1, 3));
182 }
183 for (i = 0; i < EDGES; i++) {
184 long int r = RNG_INTEGER(0, DIM - 1);
185 long int c = RNG_INTEGER(0, r);
186 igraph_real_t value = RNG_INTEGER(1, 5);
187 igraph_sparsemat_entry(&A, r, c, value);
188 }
189 igraph_sparsemat_compress(&A, &B);
190 igraph_sparsemat_destroy(&A);
191 igraph_sparsemat_dupl(&B);
192 igraph_sparsemat_transpose(&B, &A, /*values=*/ 1);
193 igraph_sparsemat_destroy(&B);
194
195 igraph_vector_init(&b, DIM);
196 for (i = 0; i < DIM; i++) {
197 VECTOR(b)[i] = RNG_INTEGER(1, 10);
198 }
199
200 igraph_vector_init(&x, DIM);
201 igraph_sparsemat_utsolve(&A, &b, &x);
202
203 igraph_sparsemat_transpose(&A, &B, /*values=*/ 1);
204 if (! check_solution(&B, &x, &b)) {
205 return 4;
206 }
207
208 igraph_vector_destroy(&b);
209 igraph_vector_destroy(&x);
210 igraph_sparsemat_destroy(&B);
211 igraph_sparsemat_destroy(&A);
212
213 #undef DIM
214 #undef EDGES
215
216 /* cholsol */
217 /* We need a positive definite matrix, so we create a full-rank
218 matrix first and then calculate A'A, which will be positive
219 definite. */
220
221 #define DIM 10
222 #define EDGES (DIM*DIM/6)
223 igraph_sparsemat_init(&A, DIM, DIM, EDGES + DIM);
224 for (i = 0; i < DIM; i++) {
225 igraph_sparsemat_entry(&A, i, i, RNG_INTEGER(1, 3));
226 }
227 for (i = 0; i < EDGES; i++) {
228 long int from = RNG_INTEGER(0, DIM - 1);
229 long int to = RNG_INTEGER(0, DIM - 1);
230 igraph_real_t value = RNG_INTEGER(1, 5);
231 igraph_sparsemat_entry(&A, from, to, value);
232 }
233 igraph_sparsemat_compress(&A, &B);
234 igraph_sparsemat_destroy(&A);
235 igraph_sparsemat_dupl(&B);
236 igraph_sparsemat_transpose(&B, &A, /*values=*/ 1);
237 igraph_sparsemat_multiply(&A, &B, &C);
238 igraph_sparsemat_destroy(&A);
239 igraph_sparsemat_destroy(&B);
240
241 igraph_vector_init(&b, DIM);
242 for (i = 0; i < DIM; i++) {
243 VECTOR(b)[i] = RNG_INTEGER(1, 10);
244 }
245
246 igraph_vector_init(&x, DIM);
247 igraph_sparsemat_cholsol(&C, &b, &x, /*order=*/ 0);
248
249 if (! check_solution(&C, &x, &b)) {
250 return 5;
251 }
252
253 igraph_vector_destroy(&b);
254 igraph_vector_destroy(&x);
255 igraph_sparsemat_destroy(&C);
256
257 #undef DIM
258 #undef EDGES
259
260 /* lusol */
261
262 #define DIM 10
263 #define EDGES (DIM*DIM/4)
264 igraph_sparsemat_init(&A, DIM, DIM, EDGES + DIM);
265 for (i = 0; i < DIM; i++) {
266 igraph_sparsemat_entry(&A, i, i, RNG_INTEGER(1, 3));
267 }
268 for (i = 0; i < EDGES; i++) {
269 long int from = RNG_INTEGER(0, DIM - 1);
270 long int to = RNG_INTEGER(0, DIM - 1);
271 igraph_real_t value = RNG_INTEGER(1, 5);
272 igraph_sparsemat_entry(&A, from, to, value);
273 }
274 igraph_sparsemat_compress(&A, &B);
275 igraph_sparsemat_destroy(&A);
276 igraph_sparsemat_dupl(&B);
277
278 igraph_vector_init(&b, DIM);
279 for (i = 0; i < DIM; i++) {
280 VECTOR(b)[i] = RNG_INTEGER(1, 10);
281 }
282
283 igraph_vector_init(&x, DIM);
284 igraph_sparsemat_lusol(&B, &b, &x, /*order=*/ 0, /*tol=*/ 1e-10);
285
286 if (! check_solution(&B, &x, &b)) {
287 return 6;
288 }
289
290 igraph_vector_destroy(&b);
291 igraph_vector_destroy(&x);
292 igraph_sparsemat_destroy(&B);
293
294 #undef DIM
295 #undef EDGES
296
297 return 0;
298 }
299