1 /*********************************************************/
2 /* TAUCS                                                 */
3 /* Author: Sivan Toledo                                  */
4 /*********************************************************/
5 
6 #include <taucs_config_tests.h>
7 #include <taucs_config_build.h>
8 
9 /*********************************************************/
10 /* Cilk-related stuff                                    */
11 /*********************************************************/
12 
13 #ifdef TAUCS_CILK
14 #undef TAUCS_C99_COMPLEX /* cilk2c can't process complex.h */
15 #endif
16 
17 #ifdef TAUCS_CORE_CILK
18 #ifdef TAUCS_CILK
19 /* We are compiling a Cilk source with a Cilk compiler */
20 
21 
22 #include <cilk.h>
23 #include <cilk-lib.h>
24 
25 #define taucs_cilk   cilk
26 #define taucs_spawn  spawn
27 #define taucs_sync   sync
28 #define taucs_inlet  inlet
29 #define taucs_Self   Self
30 #define taucs_Cilk_active_size Cilk_active_size
31 
32 #else
33 /* We are compiling a Cilk source, but with a C compiler */
34 #define cilk
35 #define spawn
36 #define sync
37 #define inlet
38 #define Self 0
39 #define Cilk_active_size 1
40 
41 #define taucs_cilk
42 #define taucs_spawn
43 #define taucs_sync
44 #define taucs_inlet
45 #define taucs_Self 0
46 #define taucs_Cilk_active_size 1
47 #endif
48 #else /* not CORE_CILK */
49 #define taucs_cilk
50 #define taucs_spawn
51 #define taucs_sync
52 #define taucs_inlet
53 #define taucs_Self 0
54 #define taucs_Cilk_active_size 1
55 #endif
56 
57 /*********************************************************/
58 /* other stuff                                           */
59 /*********************************************************/
60 
61 #ifdef TAUCS_CONFIG_DREAL
62 #define TAUCS_DOUBLE_IN_BUILD
63 #endif
64 #ifdef TAUCS_CONFIG_SREAL
65 #define TAUCS_SINGLE_IN_BUILD
66 #endif
67 #ifdef TAUCS_CONFIG_DCOMPLEX
68 #define TAUCS_DCOMPLEX_IN_BUILD
69 #endif
70 #ifdef TAUCS_CONFIG_SCOMPLEX
71 #define TAUCS_SCOMPLEX_IN_BUILD
72 #endif
73 
74 #if   defined(TAUCS_BLAS_UNDERSCORE)
75 #define taucs_blas_name(x) (x##_)
76 #elif defined(TAUCS_BLAS_NOUNDERSCORE)
77 #define taucs_blas_name(x) (x)
78 #else
79 #error "taucs_blas_[no]underscore_test: linking with the BLAS failed both attempts"
80 #endif
81 
82 #ifdef OSTYPE_win32
83 typedef unsigned long ssize_t;
84 typedef int mode_t;
85 typedef int perm_t;
86 #define random    rand
87 #define srandom   srand
88 #endif
89 
90 #define TAUCS_SUCCESS                       0
91 #define TAUCS_ERROR                        -1
92 #define TAUCS_ERROR_NOMEM                  -2
93 #define TAUCS_ERROR_BADARGS                -3
94 #define TAUCS_ERROR_INDEFINITE             -4
95 #define TAUCS_ERROR_MAXDEPTH               -5
96 
97 #define TAUCS_INT       1024
98 #define TAUCS_DOUBLE    2048
99 #define TAUCS_SINGLE    4096
100 #define TAUCS_DCOMPLEX  8192
101 #define TAUCS_SCOMPLEX 16384
102 
103 #define TAUCS_LOWER      1
104 #define TAUCS_UPPER      2
105 #define TAUCS_TRIANGULAR 4
106 #define TAUCS_SYMMETRIC  8
107 #define TAUCS_HERMITIAN  16
108 #define TAUCS_PATTERN    32
109 
110 #define TAUCS_METHOD_LLT  1
111 #define TAUCS_METHOD_LDLT 2
112 #define TAUCS_METHOD_PLU  3
113 
114 #define TAUCS_VARIANT_SNMF 1
115 #define TAUCS_VARIANT_SNLL 2
116 
117 typedef double    taucs_double;
118 typedef float     taucs_single;
119 
120 /* The macro TAUCS_C99_COMPLEX is defined in */
121 /* build/OSTYPE/taucs_config_tests if the    */
122 /* test program progs/taucs_c99_complex_test */
123 /* compiles, links, and runs.                */
124 
125 /*#if defined(__GNUC__) && !defined(TAUCS_CONFIG_GENERIC_COMPLEX)*/
126 #ifdef TAUCS_C99_COMPLEX
127 /*
128 typedef __complex__ double taucs_dcomplex;
129 typedef __complex__ float  taucs_scomplex;
130 */
131 
132 #include <complex.h>
133 
134 #undef I
135 #ifdef _Imaginary_I
136 #define TAUCS_IMAGINARY_I _Imaginary_I
137 #else
138 #define TAUCS_IMAGINARY_I _Complex_I
139 #endif
140 
141 typedef _Complex double taucs_dcomplex;
142 typedef _Complex float  taucs_scomplex;
143 
144 #define taucs_complex_create(r,i)  ((r)+TAUCS_IMAGINARY_I*(i))
145 #define taucs_ccomplex_create(r,i) ((r)+TAUCS_IMAGINARY_I*(i))
146 #define taucs_zcomplex_create(r,i) ((r)+TAUCS_IMAGINARY_I*(i))
147 
148 #define taucs_add(x,y) ((x)+(y))
149 #define taucs_sub(x,y) ((x)-(y))
150 #define taucs_mul(x,y) ((x)*(y))
151 #define taucs_div(x,y) ((x)/(y))
152 #define taucs_neg(x)   (-(x))
153 
154 #define taucs_dadd(x,y) ((x)+(y))
155 #define taucs_dsub(x,y) ((x)-(y))
156 #define taucs_dmul(x,y) ((x)*(y))
157 #define taucs_ddiv(x,y) ((x)/(y))
158 #define taucs_dneg(x)   (-(x))
159 #define taucs_dconj(x)  (x)
160 #define taucs_dimag(x)    0.0
161 #define taucs_dreal(x)    (x)
162 #define taucs_dminusone -1.0
163 #define taucs_done      1.0
164 #define taucs_dzero     0.0
165 #define taucs_dabs(x)   (fabs(x))
166 #define taucs_dsqrt(x)  (sqrt(x))
167 
168 #define taucs_sadd(x,y) ((x)+(y))
169 #define taucs_ssub(x,y) ((x)-(y))
170 #define taucs_smul(x,y) ((x)*(y))
171 #define taucs_sdiv(x,y) ((x)/(y))
172 #define taucs_sneg(x)   (-(x))
173 #define taucs_sconj(x)  (x)
174 #define taucs_simag(x)    0.0f
175 #define taucs_sreal(x)    (x)
176 #define taucs_sminusone -1.0f
177 #define taucs_sone      1.0f
178 #define taucs_szero     0.0f
179 #define taucs_sabs(x)   ((taucs_single) fabs(x))
180 #define taucs_ssqrt(x)  ((taucs_single) sqrt(x))
181 
182 #define taucs_zadd(x,y) ((x)+(y))
183 #define taucs_zsub(x,y) ((x)-(y))
184 #define taucs_zmul(x,y) ((x)*(y))
185 #define taucs_zdiv(x,y) ((x)/(y))
186 #define taucs_zneg(x)   (-(x))
187 #define taucs_zconj(x)  (conj(x))
188 #define taucs_zimag(x)    (cimag(x))
189 #define taucs_zreal(x)    (creal(x))
190 #define taucs_zminusone -1.0+0.0*TAUCS_IMAGINARY_I
191 #define taucs_zone      1.0+0.0*TAUCS_IMAGINARY_I
192 #define taucs_zzero     0.0+0.0*TAUCS_IMAGINARY_I
193 #define taucs_zabs(x)   (cabs(x))
194 #define taucs_zsqrt(x)  (csqrt(x))
195 
196 #define taucs_cadd(x,y) ((x)+(y))
197 #define taucs_csub(x,y) ((x)-(y))
198 #define taucs_cmul(x,y) ((x)*(y))
199 #define taucs_cdiv(x,y) ((x)/(y))
200 #define taucs_cneg(x)   (-(x))
201 #define taucs_cconj(x)  (conjf(x))
202 #define taucs_cimag(x)    (cimagf(x))
203 #define taucs_creal(x)    (crealf(x))
204 #define taucs_cminusone -1.0f+0.0f*TAUCS_IMAGINARY_I
205 #define taucs_cone      1.0f+0.0f*TAUCS_IMAGINARY_I
206 #define taucs_czero     0.0f+0.0f*TAUCS_IMAGINARY_I
207 #define taucs_cabs(x)   (cabsf(x))
208 #define taucs_csqrt(x)  (csqrt(x))
209 
210 #if defined(TAUCS_CORE_DOUBLE)
211 
212 #define taucs_conj(x)  (x)
213 #define taucs_im(x)    0.0
214 #define taucs_re(x)    (x)
215 #define taucs_minusone -1.0
216 #define taucs_one      1.0
217 #define taucs_zero     0.0
218 #define taucs_abs(x)   (fabs(x))
219 #define taucs_sqrt(x)  (sqrt(x))
220 
221 #elif defined(TAUCS_CORE_GENERAL)
222 #define taucs_im(x)    0.0
223 #define taucs_re(x)    (x)
224 #define taucs_minusone -1.0
225 #define taucs_one      1.0
226 #define taucs_zero     0.0
227 /*
228 #define taucs_conj(x)  (x)
229 #define taucs_abs(x)   (fabs(x))
230 #define taucs_sqrt(x)  (sqrt(x))
231 */
232 #elif defined(TAUCS_CORE_SINGLE)
233 
234 #define taucs_conj(x)  (x)
235 #define taucs_im(x)    0.0f
236 #define taucs_re(x)    (x)
237 #define taucs_minusone -1.0f
238 #define taucs_one      1.0f
239 #define taucs_zero     0.0f
240 #define taucs_abs(x)   ((taucs_single) fabs(x))
241 #define taucs_sqrt(x)  ((taucs_single) sqrt(x))
242 
243 #elif defined(TAUCS_CORE_DCOMPLEX)
244 /*
245 #define taucs_conj(x)  (~(x))
246 #define taucs_im(x)    (__imag__ (x))
247 #define taucs_re(x)    (__real__ (x))
248 #define taucs_minusone -1.0+0.0i
249 #define taucs_one      1.0+0.0i
250 #define taucs_zero     0.0+0.0i
251 #define taucs_abs(x)   taucs_zabs_fn(x)
252 #define taucs_sqrt(x)  taucs_zsqrt_fn(x)
253 */
254 
255 #define taucs_conj(x)  (conj(x))
256 #define taucs_im(x)    (cimag(x))
257 #define taucs_re(x)    (creal(x))
258 #define taucs_minusone -1.0+0.0*TAUCS_IMAGINARY_I
259 #define taucs_one      1.0+0.0*TAUCS_IMAGINARY_I
260 #define taucs_zero     0.0+0.0*TAUCS_IMAGINARY_I
261 #define taucs_abs(x)   (cabs(x))
262 #define taucs_sqrt(x)  (csqrt(x))
263 
264 #elif defined(TAUCS_CORE_SCOMPLEX)
265 /*
266 #define taucs_conj(x)  (~(x))
267 #define taucs_im(x)    (__imag__ (x))
268 #define taucs_re(x)    (__real__ (x))
269 #define taucs_minusone -1.0f+0.0fi
270 #define taucs_one      1.0f+0.0fi
271 #define taucs_zero     0.0f+0.0fi
272 #define taucs_abs(x)   taucs_cabs_fn(x)
273 #define taucs_sqrt(x)  taucs_csqrt_fn(x)
274 */
275 #define taucs_conj(x)  (conjf(x))
276 #define taucs_im(x)    (cimagf(x))
277 #define taucs_re(x)    (crealf(x))
278 #define taucs_minusone -1.0f+0.0f*TAUCS_IMAGINARY_I
279 #define taucs_one      1.0f+0.0f*TAUCS_IMAGINARY_I
280 #define taucs_zero     0.0f+0.0f*TAUCS_IMAGINARY_I
281 #define taucs_abs(x)   (cabsf(x))
282 #define taucs_sqrt(x)  (csqrtf(x))
283 
284 #endif
285 
286 #else /* C99 */
287 
288 typedef struct {double r,i;} taucs_dcomplex;
289 typedef struct {float  r,i;} taucs_scomplex;
290 
291 #define taucs_zcomplex_create(r,i) taucs_zcomplex_create_fn(r,i)
292 #define taucs_ccomplex_create(r,i) taucs_ccomplex_create_fn(r,i)
293 
294 #define taucs_dadd(x,y) ((x)+(y))
295 #define taucs_dsub(x,y) ((x)-(y))
296 #define taucs_dmul(x,y) ((x)*(y))
297 #define taucs_ddiv(x,y) ((x)/(y))
298 #define taucs_dneg(x)   (-(x))
299 #define taucs_dconj(x)  (x)
300 #define taucs_dabs(x)   (fabs(x))
301 #define taucs_dsqrt(x)  (sqrt(x))
302 #define taucs_dimag(x)   0.0
303 #define taucs_dreal(x)   (x)
304 #define taucs_dminusone -1.0
305 #define taucs_done     1.0
306 #define taucs_dzero    0.0
307 
308 #define taucs_sadd(x,y) ((x)+(y))
309 #define taucs_ssub(x,y) ((x)-(y))
310 #define taucs_smul(x,y) ((x)*(y))
311 #define taucs_sdiv(x,y) ((x)/(y))
312 #define taucs_sneg(x)   (-(x))
313 #define taucs_sconj(x)  (x)
314 #define taucs_sabs(x)   ((taucs_single) fabs(x))
315 #define taucs_ssqrt(x)  ((taucs_single) sqrt(x))
316 #define taucs_sim(x)   0.0f
317 #define taucs_sre(x)   (x)
318 #define taucs_sminusone -1.0f
319 #define taucs_sone     1.0f
320 #define taucs_szero    0.0f
321 
322 #define taucs_zadd(x,y) taucs_zadd_fn(x,y)
323 #define taucs_zsub(x,y) taucs_zsub_fn(x,y)
324 #define taucs_zmul(x,y) taucs_zmul_fn(x,y)
325 #define taucs_zdiv(x,y) taucs_zdiv_fn(x,y)
326 #define taucs_zneg(x)   taucs_zneg_fn(x)
327 #define taucs_zconj(x)  taucs_zconj_fn(x)
328 #define taucs_zabs(x)   taucs_zabs_fn(x)
329 #define taucs_zsqrt(x)  taucs_zsqrt_fn(x)
330 #define taucs_zimag(x)    ((x).i)
331 #define taucs_zreal(x)    ((x).r)
332 #define taucs_zminusone taucs_zminusone_const
333 #define taucs_zone      taucs_zone_const
334 #define taucs_zzero     taucs_zzero_const
335 
336 #define taucs_cadd(x,y) taucs_cadd_fn(x,y)
337 #define taucs_csub(x,y) taucs_csub_fn(x,y)
338 #define taucs_cmul(x,y) taucs_cmul_fn(x,y)
339 #define taucs_cdiv(x,y) taucs_cdiv_fn(x,y)
340 #define taucs_cneg(x)   taucs_cneg_fn(x)
341 #define taucs_cconj(x)  taucs_cconj_fn(x)
342 #define taucs_cabs(x)   taucs_cabs_fn(x)
343 #define taucs_csqrt(x)  taucs_csqrt_fn(x)
344 #define taucs_cimag(x)    ((x).i)
345 #define taucs_creal(x)    ((x).r)
346 #define taucs_cminusone taucs_cminusone_const
347 #define taucs_cone      taucs_cone_const
348 #define taucs_czero     taucs_czero_const
349 
350 #if defined(TAUCS_CORE_DOUBLE)
351 
352 #define taucs_add(x,y) ((x)+(y))
353 #define taucs_sub(x,y) ((x)-(y))
354 #define taucs_mul(x,y) ((x)*(y))
355 #define taucs_div(x,y) ((x)/(y))
356 #define taucs_neg(x)   (-(x))
357 #define taucs_conj(x)  (x)
358 #define taucs_abs(x)   (fabs(x))
359 #define taucs_sqrt(x)  (sqrt(x))
360 
361 #define taucs_im(x)   0.0
362 #define taucs_re(x)   (x)
363 #define taucs_minusone -1.0
364 #define taucs_one     1.0
365 #define taucs_zero    0.0
366 
367 #elif defined(TAUCS_CORE_GENERAL)
368 /*
369 #define taucs_add(x,y) ((x)+(y))
370 #define taucs_sub(x,y) ((x)-(y))
371 #define taucs_mul(x,y) ((x)*(y))
372 #define taucs_div(x,y) ((x)/(y))
373 #define taucs_neg(x)   (-(x))
374 #define taucs_conj(x)  (x)
375 #define taucs_abs(x)   (fabs(x))
376 #define taucs_sqrt(x)  (sqrt(x))
377 */
378 #define taucs_im(x)   0.0
379 #define taucs_re(x)   (x)
380 #define taucs_minusone -1.0
381 #define taucs_one     1.0
382 #define taucs_zero    0.0
383 
384 #elif defined(TAUCS_CORE_SINGLE)
385 
386 #define taucs_add(x,y) ((x)+(y))
387 #define taucs_sub(x,y) ((x)-(y))
388 #define taucs_mul(x,y) ((x)*(y))
389 #define taucs_div(x,y) ((x)/(y))
390 #define taucs_neg(x)   (-(x))
391 #define taucs_conj(x)  (x)
392 #define taucs_abs(x)   ((taucs_single) fabs(x))
393 #define taucs_sqrt(x)  ((taucs_single) sqrt(x))
394 
395 #define taucs_im(x)   0.0f
396 #define taucs_re(x)   (x)
397 #define taucs_minusone -1.0f
398 #define taucs_one     1.0f
399 #define taucs_zero    0.0f
400 
401 #elif defined(TAUCS_CORE_DCOMPLEX)
402 
403 #define taucs_complex_create(r,i) taucs_zcomplex_create_fn(r,i)
404 
405 #define taucs_add(x,y) taucs_zadd_fn(x,y)
406 #define taucs_sub(x,y) taucs_zsub_fn(x,y)
407 #define taucs_mul(x,y) taucs_zmul_fn(x,y)
408 #define taucs_div(x,y) taucs_zdiv_fn(x,y)
409 #define taucs_neg(x)   taucs_zneg_fn(x)
410 #define taucs_conj(x)  taucs_zconj_fn(x)
411 #define taucs_abs(x)   taucs_zabs_fn(x)
412 #define taucs_sqrt(x)  taucs_zsqrt_fn(x)
413 
414 #define taucs_im(x)    ((x).i)
415 #define taucs_re(x)    ((x).r)
416 #define taucs_minusone taucs_zminusone_const
417 #define taucs_one      taucs_zone_const
418 #define taucs_zero     taucs_zzero_const
419 
420 #elif defined(TAUCS_CORE_SCOMPLEX)
421 
422 #define taucs_complex_create(r,i) taucs_ccomplex_create_fn(r,i)
423 
424 #define taucs_add(x,y) taucs_cadd_fn(x,y)
425 #define taucs_sub(x,y) taucs_csub_fn(x,y)
426 #define taucs_mul(x,y) taucs_cmul_fn(x,y)
427 #define taucs_div(x,y) taucs_cdiv_fn(x,y)
428 #define taucs_neg(x)   taucs_cneg_fn(x)
429 #define taucs_conj(x)  taucs_cconj_fn(x)
430 #define taucs_abs(x)   taucs_cabs_fn(x)
431 #define taucs_sqrt(x)  taucs_csqrt_fn(x)
432 
433 #define taucs_im(x)    ((x).i)
434 #define taucs_re(x)    ((x).r)
435 #define taucs_minusone taucs_cminusone_const
436 #define taucs_one      taucs_cone_const
437 #define taucs_zero     taucs_czero_const
438 
439 #endif /* SCOMPLEX */
440 
441 #endif
442 
443 extern taucs_double taucs_dzero_const    ;
444 extern taucs_double taucs_done_const     ;
445 extern taucs_double taucs_dminusone_const;
446 
447 extern taucs_single taucs_szero_const    ;
448 extern taucs_single taucs_sone_const     ;
449 extern taucs_single taucs_sminusone_const;
450 
451 extern taucs_dcomplex taucs_zzero_const    ;
452 extern taucs_dcomplex taucs_zone_const     ;
453 extern taucs_dcomplex taucs_zminusone_const;
454 
455 extern taucs_scomplex taucs_czero_const    ;
456 extern taucs_scomplex taucs_cone_const     ;
457 extern taucs_scomplex taucs_cminusone_const;
458 
459 #define taucs_isnan(x) (isnan((double)(taucs_re(x))) || isnan((double)(taucs_im(x))))
460 #define taucs_isinf(x) (isinf((double)(taucs_re(x))) || isinf((double)(taucs_im(x))))
461 
462 #ifdef TAUCS_CORE_SINGLE
463 #define taucs_zero_const     taucs_szero_const
464 #define taucs_one_const      taucs_sone_const
465 #define taucs_minusone_const taucs_sminusone_const
466 
467 #define taucs_zero_real_const     taucs_szero_const
468 #define taucs_one_real_const      taucs_sone_const
469 #define taucs_minusone_real_const taucs_sminusone_const
470 
471 #define taucs_gemm  taucs_blas_name(sgemm)
472 #define taucs_potrf taucs_blas_name(spotrf)
473 #define taucs_herk  taucs_blas_name(ssyrk)
474 #define taucs_trsm  taucs_blas_name(strsm)
475 #endif
476 
477 #ifdef TAUCS_CORE_DOUBLE
478 #define taucs_zero_const     taucs_dzero_const
479 #define taucs_one_const      taucs_done_const
480 #define taucs_minusone_const taucs_dminusone_const
481 
482 #define taucs_zero_real_const     taucs_dzero_const
483 #define taucs_one_real_const      taucs_done_const
484 #define taucs_minusone_real_const taucs_dminusone_const
485 
486 #define taucs_gemm  taucs_blas_name(dgemm)
487 #define taucs_potrf taucs_blas_name(dpotrf)
488 #define taucs_herk  taucs_blas_name(dsyrk)
489 #define taucs_trsm  taucs_blas_name(dtrsm)
490 #endif
491 
492 /*
493 #ifdef TAUCS_CORE_GENERAL
494 #define taucs_zero_const     taucs_dzero_const
495 #define taucs_one_const      taucs_done_const
496 #define taucs_minusone_const taucs_dminusone_const
497 #define taucs_zero_real_const     taucs_dzero_const
498 #define taucs_one_real_const      taucs_done_const
499 #define taucs_minusone_real_const taucs_dminusone_const
500 #endif
501 */
502 
503 #ifdef TAUCS_CORE_SCOMPLEX
504 #define taucs_zero_const     taucs_czero_const
505 #define taucs_one_const      taucs_cone_const
506 #define taucs_minusone_const taucs_cminusone_const
507 
508 #define taucs_zero_real_const     taucs_szero_const
509 #define taucs_one_real_const      taucs_sone_const
510 #define taucs_minusone_real_const taucs_sminusone_const
511 
512 #define taucs_gemm  taucs_blas_name(cgemm)
513 #define taucs_potrf taucs_blas_name(cpotrf)
514 #define taucs_herk  taucs_blas_name(cherk)
515 #define taucs_trsm  taucs_blas_name(ctrsm)
516 #endif
517 
518 #ifdef TAUCS_CORE_DCOMPLEX
519 #define taucs_zero_const     taucs_zzero_const
520 #define taucs_one_const      taucs_zone_const
521 #define taucs_minusone_const taucs_zminusone_const
522 
523 #define taucs_zero_real_const     taucs_dzero_const
524 #define taucs_one_real_const      taucs_done_const
525 #define taucs_minusone_real_const taucs_dminusone_const
526 
527 #define taucs_gemm  taucs_blas_name(zgemm)
528 #define taucs_potrf taucs_blas_name(zpotrf)
529 #define taucs_herk  taucs_blas_name(zherk)
530 #define taucs_trsm  taucs_blas_name(ztrsm)
531 #endif
532 
533 /*********************************************************/
534 /*                                                       */
535 /*********************************************************/
536 
537 typedef struct
538 {
539   int     n;    /* columns                      */
540   int     m;    /* rows; don't use if symmetric   */
541   int     flags;
542   int*    colptr; /* pointers to where columns begin in rowind and values. */
543                   /* 0-based. Length is (n+1). */
544   int*    rowind; /* row indices */
545 
546   union {
547     void*           v;
548     taucs_double*   d;
549     taucs_single*   s;
550     taucs_dcomplex* z;
551     taucs_scomplex* c;
552   } values;
553 
554 } taucs_ccs_matrix;
555 
556 typedef struct {
557   int   type;
558   int   nmatrices;
559   void* type_specific;
560 
561   /* the following may change! do not rely on them. */
562   double nreads, nwrites, bytes_read, bytes_written, read_time, write_time;
563 } taucs_io_handle;
564 
565 /* generate all the prototypes */
566 
567 #define taucs_datatype taucs_double
568 #define taucs_real_datatype taucs_double
569 #define taucs_dtl(X) taucs_d##X
570 #include "taucs_private.h"
571 #undef taucs_real_datatype
572 #undef taucs_datatype
573 #undef taucs_dtl
574 
575 #define taucs_datatype taucs_single
576 #define taucs_real_datatype taucs_single
577 #define taucs_dtl(X) taucs_s##X
578 #include "taucs_private.h"
579 #undef taucs_real_datatype
580 #undef taucs_datatype
581 #undef taucs_dtl
582 
583 #define taucs_datatype taucs_dcomplex
584 #define taucs_real_datatype taucs_double
585 #define taucs_dtl(X) taucs_z##X
586 #include "taucs_private.h"
587 #undef taucs_real_datatype
588 #undef taucs_datatype
589 #undef taucs_dtl
590 
591 #define taucs_datatype taucs_scomplex
592 #define taucs_real_datatype taucs_single
593 #define taucs_dtl(X) taucs_c##X
594 #include "taucs_private.h"
595 #undef taucs_real_datatype
596 #undef taucs_datatype
597 #undef taucs_dtl
598 
599 /*********************************************************/
600 /*                                                       */
601 /*********************************************************/
602 
603 /* now define the data type for the file that we compile now */
604 
605 #ifdef TAUCS_CORE_DOUBLE
606 #define TAUCS_CORE
607 #define TAUCS_CORE_REAL
608 #define TAUCS_CORE_DATATYPE TAUCS_DOUBLE
609 typedef taucs_double taucs_datatype;
610 #define taucs_dtl(X) taucs_d##X
611 #define taucs_values values.d
612 #define taucs_iszero(x) ((x) == 0.0)
613 typedef double taucs_real_datatype; /* omer: this is the datatype of the real and imaginary part of the datatype*/
614 #endif
615 
616 #ifdef TAUCS_CORE_GENERAL
617 #define TAUCS_CORE
618 #define TAUCS_CORE_DATATYPE TAUCS_DOUBLE
619 typedef taucs_double taucs_datatype;
620 typedef double taucs_real_datatype;
621 /*
622 #define TAUCS_CORE_REAL
623 #define TAUCS_CORE_DATATYPE TAUCS_DOUBLE
624 #define taucs_values values.d
625 #define taucs_dtl(X) taucs_g##X
626 #define taucs_iszero(x) ((x) == 0.0)
627 */
628 #endif
629 
630 #ifdef  TAUCS_CORE_SINGLE
631 #define TAUCS_CORE
632 #define TAUCS_CORE_REAL
633 #define TAUCS_CORE_DATATYPE TAUCS_SINGLE
634 typedef taucs_single taucs_datatype;
635 #define taucs_dtl(X) taucs_s##X
636 #define taucs_values values.s
637 #define taucs_iszero(x) ((x) == 0.0f)
638 typedef float taucs_real_datatype; /* omer: this is the datatype of the real and imaginary part of the datatype*/
639 #endif
640 
641 #ifdef  TAUCS_CORE_DCOMPLEX
642 #define TAUCS_CORE
643 #define TAUCS_CORE_COMPLEX
644 #define TAUCS_CORE_DATATYPE TAUCS_DCOMPLEX
645 typedef taucs_dcomplex taucs_datatype;
646 #define taucs_dtl(X) taucs_z##X
647 #define taucs_values values.z
648 #define taucs_iszero(x) (taucs_re(x) == 0.0 && taucs_im(x) == 0.0)
649 typedef double taucs_real_datatype; /* omer: this is the datatype of the real and imaginary part of the datatype*/
650 #endif
651 
652 #ifdef  TAUCS_CORE_SCOMPLEX
653 #define TAUCS_CORE
654 #define TAUCS_CORE_COMPLEX
655 #define TAUCS_CORE_DATATYPE TAUCS_SCOMPLEX
656 typedef taucs_scomplex taucs_datatype;
657 #define taucs_dtl(X) taucs_c##X
658 #define taucs_values values.c
659 #define taucs_iszero(x) (taucs_re(x) == 0.0f && taucs_im(x) == 0.0f)
660 typedef float taucs_real_datatype; /* omer: this is the datatype of the real and imaginary part of the datatype*/
661 #endif
662 
663 #ifndef TAUCS_CORE_DATATYPE
664 typedef taucs_double taucs_datatype;
665 typedef double taucs_real_datatype;
666 #endif
667 
668 /*********************************************************/
669 /*                                                       */
670 /*********************************************************/
671 
672 double taucs_get_nan(void);
673 
674 /*
675    routines for testing memory allocation.
676    Mostly useful for testing programs
677    that hunt for memory leaks.
678 */
679 
680 double taucs_allocation_amount(void);
681 int    taucs_allocation_count(void);
682 int    taucs_allocation_attempts(void);
683 void   taucs_allocation_assert_clean(void);
684 void   taucs_allocation_mark_clean(void);
685 void   taucs_allocation_induce_failure(int i);
686 
687 /*
688    these are meant to allow allocation
689    and more importantly, deallocation,
690    within the testing programs.
691 */
692 
693 #include <stdlib.h>
694 
695 void* taucs_malloc (size_t size)              ;
696 void* taucs_calloc (size_t nmemb, size_t size);
697 void* taucs_realloc(void* ptr, size_t size)   ;
698 void  taucs_free   (void* ptr)                ;
699 
700 #if defined(TAUCS_CORE)
701 
702 #if defined(TAUCS_MEMORY_TEST_yes)
703 
704 #include <stdlib.h>
705 
706 void* taucs_internal_calloc(size_t nmemb, size_t size,char* file, int line);
707 void* taucs_internal_malloc(size_t size,              char* file, int line);
708 void* taucs_internal_realloc(void *ptr, size_t size,   char* file, int line);
709 void  taucs_internal_free(void *ptr,                   char* file, int line);
710 
711 /*
712  #define realloc(x,y) taucs_internal_realloc(x,y,__FILE__,__LINE__)
713  #define malloc(x)    taucs_internal_malloc(x,__FILE__,__LINE__)
714  #define calloc(x,y)  taucs_internal_calloc(x,y,__FILE__,__LINE__)
715  #define free(x)      taucs_internal_free(x,__FILE__,__LINE__)
716 */
717 
718 #define taucs_realloc(x,y) taucs_internal_realloc(x,y,__FILE__,__LINE__)
719 #define taucs_malloc(x)    taucs_internal_malloc(x,__FILE__,__LINE__)
720 #define taucs_calloc(x,y)  taucs_internal_calloc(x,y,__FILE__,__LINE__)
721 #define taucs_free(x)      taucs_internal_free(x,__FILE__,__LINE__)
722 
723 #define realloc(x,y) taucs_must_not_call_realloc_directly(x,y)
724 #define malloc(x)    taucs_must_not_call_malloc_directly(x)
725 #define calloc(x,y)  taucs_must_not_call_calloc_directly(x,y)
726 #define free(x)      taucs_must_not_call_free_directly(x)
727 
728 #else /* TAUCS_CORE, but not memory testing */
729 
730 void* taucs_calloc_stub(size_t nmemb, size_t size);
731 void* taucs_malloc_stub(size_t size);
732 void* taucs_realloc_stub(void *ptr, size_t size);
733 void  taucs_free_stub(void *ptr);
734 
735 #define realloc(x,y) taucs_must_not_call_realloc_directly(x,y)
736 #define malloc(x)    taucs_must_not_call_malloc_directly(x)
737 #define calloc(x,y)  taucs_must_not_call_calloc_directly(x,y)
738 #define free(x)      taucs_must_not_call_free_directly(x)
739 
740 #define taucs_realloc(x,y) taucs_realloc_stub(x,y)
741 #define taucs_malloc(x)    taucs_malloc_stub(x)
742 #define taucs_calloc(x,y)  taucs_calloc_stub(x,y)
743 #define taucs_free(x)      taucs_free_stub(x)
744 
745 #endif
746 #endif
747 
748 /*********************************************************/
749 /*                                                       */
750 /*********************************************************/
751 
752 #ifndef max
753 #define max(x,y) ( ((x) > (y)) ? (x) : (y) )
754 #endif
755 
756 #ifndef min
757 #define min(x,y) ( ((x) < (y)) ? (x) : (y) )
758 #endif
759 
760 /*********************************************************/
761 /*                                                       */
762 /*********************************************************/
763 
764 /* externs */
765 extern int ireadhb_(char*, char*, int*, int*, int*);
766 extern int creadhb_(char*, int*, int*, int*, int*, int*, taucs_scomplex*);
767 extern int dreadhb_(char*, int*, int*, int*, int*, int*, taucs_double*);
768 extern int sreadhb_(char*, int*, int*, int*, int*, int*, taucs_single*);
769 extern int zreadhb_(char*, int*, int*, int*, int*, int*, taucs_dcomplex*);
770 
771 extern int amdexa_(int*, int*, int*, int*, int*, int*, int*, int*, int*,
772 			int*, int*, int*, int*, int*, int*);
773 extern int amdtru_(int*, int*, int*, int*, int*, int*, int*, int*, int*,
774 			int*, int*, int*, int*, int*, int*);
775 extern int amdbar_(int*, int*, int*, int*, int*, int*, int*, int*, int*,
776 			int*, int*, int*, int*, int*, int*);
777 extern int genmmd_(int*, int*, int*, int*, int*, int*, int*, int*, int*,
778 			int*, int*, int*);
779 
780 /*********************************************************/
781 /*                                                       */
782 /*********************************************************/
783 
784 #if (defined(OSTYPE_irix) || defined(OSTYPE_solaris))
785 
786 #include <math.h>
787 #include <ieeefp.h>
788 #define isinf(x) (!finite((x)) && !isnan((x)))
789 
790 #elif defined(OSTYPE_win32)
791 
792 #include <float.h>
793 #define isnan(x)  (_isnan(x))
794 #define isinf(x)  (!(_finite(x)) && !(_isnan(x)))
795 #define finite(x) (_finite(x))
796 
797 #endif
798 
799 /* If these are mactors (e.g., gcc -std=c99), do not declare   */
800 /* otherwise, declare them, since they are not always declared */
801 /* in math.h (e.g., gcc -std=c89 -pedantic); these are for     */
802 /* gcc 3.3.1                                                   */
803 
804 #ifndef isnan
805 extern int isnan(double);
806 #endif
807 #ifndef finite
808 extern int finite(double);
809 #endif
810 #ifndef isinf
811 extern int isinf(double);
812 #endif
813 
814 extern int taucs_potrf(char*, int*, taucs_datatype*, int*, int*);
815 extern int taucs_trsm(char *, char *, char *, char *,
816 			int*, int*, taucs_datatype*, taucs_datatype*, int *,
817 			taucs_datatype*, int *);
818 extern int taucs_gemm(char *, char *, int*, int*, int *,
819 			taucs_datatype*, taucs_datatype*, int *, taucs_datatype*, int *,
820 			taucs_datatype*, taucs_datatype*, int*);
821 extern int taucs_herk(char *, char *,
822 		      int *, int *,
823 		      taucs_real_datatype*,
824 		      taucs_datatype*, int *,
825 		      taucs_real_datatype*,
826 		      taucs_datatype*, int *);
827 
828 taucs_double taucs_blas_name(dnrm2)(int*, taucs_double*, int*);
829 taucs_single taucs_blas_name(snrm2)(int*, taucs_single*, int*);
830 taucs_double taucs_blas_name(dznrm2)(int*, taucs_dcomplex*, int*);
831 taucs_single taucs_blas_name(scnrm2)(int*, taucs_scomplex*, int*);
832 
833 /*********************************************************/
834 /*                                                       */
835 /*********************************************************/
836 
837 
838