1 /******************************************************************************
2 * Copyright 1998-2019 Lawrence Livermore National Security, LLC and other
3 * HYPRE Project Developers. See the top-level COPYRIGHT file for details.
4 *
5 * SPDX-License-Identifier: (Apache-2.0 OR MIT)
6 ******************************************************************************/
7
8 #include <math.h>
9 #include <stdlib.h>
10 #include <stdio.h>
11
12 #include "fortran_matrix.h"
13 #include "_hypre_utilities.h"
14
15 utilities_FortranMatrix*
utilities_FortranMatrixCreate(void)16 utilities_FortranMatrixCreate(void) {
17
18 utilities_FortranMatrix* mtx;
19
20 mtx = hypre_TAlloc(utilities_FortranMatrix, 1, HYPRE_MEMORY_HOST);
21 hypre_assert( mtx != NULL );
22
23 mtx->globalHeight = 0;
24 mtx->height = 0;
25 mtx->width = 0;
26 mtx->value = NULL;
27 mtx->ownsValues = 0;
28
29 return mtx;
30 }
31
32 void
utilities_FortranMatrixAllocateData(HYPRE_BigInt h,HYPRE_BigInt w,utilities_FortranMatrix * mtx)33 utilities_FortranMatrixAllocateData( HYPRE_BigInt h, HYPRE_BigInt w,
34 utilities_FortranMatrix* mtx ) {
35
36 hypre_assert( h > 0 && w > 0 );
37 hypre_assert( mtx != NULL );
38
39 if ( mtx->value != NULL && mtx->ownsValues )
40 hypre_TFree( mtx->value ,HYPRE_MEMORY_HOST);
41
42 mtx->value = hypre_CTAlloc(HYPRE_Real, h*w, HYPRE_MEMORY_HOST);
43 hypre_assert ( mtx->value != NULL );
44
45 mtx->globalHeight = h;
46 mtx->height = h;
47 mtx->width = w;
48 mtx->ownsValues = 1;
49 }
50
51
52 void
utilities_FortranMatrixWrap(HYPRE_Real * v,HYPRE_BigInt gh,HYPRE_BigInt h,HYPRE_BigInt w,utilities_FortranMatrix * mtx)53 utilities_FortranMatrixWrap( HYPRE_Real* v, HYPRE_BigInt gh, HYPRE_BigInt h, HYPRE_BigInt w,
54 utilities_FortranMatrix* mtx ) {
55
56 hypre_assert( h > 0 && w > 0 );
57 hypre_assert( mtx != NULL );
58
59 if ( mtx->value != NULL && mtx->ownsValues )
60 hypre_TFree( mtx->value ,HYPRE_MEMORY_HOST);
61
62 mtx->value = v;
63 hypre_assert ( mtx->value != NULL );
64
65 mtx->globalHeight = gh;
66 mtx->height = h;
67 mtx->width = w;
68 mtx->ownsValues = 0;
69 }
70
71
72 void
utilities_FortranMatrixDestroy(utilities_FortranMatrix * mtx)73 utilities_FortranMatrixDestroy( utilities_FortranMatrix* mtx ) {
74
75 if ( mtx == NULL )
76 return;
77
78 if ( mtx->ownsValues && mtx->value != NULL )
79 hypre_TFree(mtx->value,HYPRE_MEMORY_HOST);
80
81 hypre_TFree(mtx,HYPRE_MEMORY_HOST);
82 }
83
84 HYPRE_BigInt
utilities_FortranMatrixGlobalHeight(utilities_FortranMatrix * mtx)85 utilities_FortranMatrixGlobalHeight( utilities_FortranMatrix* mtx ) {
86
87 hypre_assert( mtx != NULL );
88
89 return mtx->globalHeight;
90 }
91
92 HYPRE_BigInt
utilities_FortranMatrixHeight(utilities_FortranMatrix * mtx)93 utilities_FortranMatrixHeight( utilities_FortranMatrix* mtx ) {
94
95 hypre_assert( mtx != NULL );
96
97 return mtx->height;
98 }
99
100 HYPRE_BigInt
utilities_FortranMatrixWidth(utilities_FortranMatrix * mtx)101 utilities_FortranMatrixWidth( utilities_FortranMatrix* mtx ) {
102
103 hypre_assert( mtx != NULL );
104
105 return mtx->width;
106 }
107
108 HYPRE_Real*
utilities_FortranMatrixValues(utilities_FortranMatrix * mtx)109 utilities_FortranMatrixValues( utilities_FortranMatrix* mtx ) {
110
111 hypre_assert( mtx != NULL );
112
113 return mtx->value;
114 }
115
116 void
utilities_FortranMatrixClear(utilities_FortranMatrix * mtx)117 utilities_FortranMatrixClear( utilities_FortranMatrix* mtx ) {
118
119 HYPRE_BigInt i, j, h, w, jump;
120 HYPRE_Real* p;
121
122 hypre_assert( mtx != NULL );
123
124 h = mtx->height;
125 w = mtx->width;
126
127 jump = mtx->globalHeight - h;
128
129 for ( j = 0, p = mtx->value; j < w; j++ ) {
130 for ( i = 0; i < h; i++, p++ )
131 *p = 0.0;
132 p += jump;
133 }
134 }
135
136 void
utilities_FortranMatrixClearL(utilities_FortranMatrix * mtx)137 utilities_FortranMatrixClearL( utilities_FortranMatrix* mtx ) {
138
139 HYPRE_BigInt i, j, k, h, w, jump;
140 HYPRE_Real* p;
141
142 hypre_assert( mtx != NULL );
143
144 h = mtx->height;
145 w = mtx->width;
146
147 if ( w > h )
148 w = h;
149
150 jump = mtx->globalHeight - h;
151
152 for ( j = 0, p = mtx->value; j < w - 1; j++ ) {
153 k = j + 1;
154 p += k;
155 for ( i = k; i < h; i++, p++ )
156 *p = 0.0;
157 p += jump;
158 }
159 }
160
161
162 void
utilities_FortranMatrixSetToIdentity(utilities_FortranMatrix * mtx)163 utilities_FortranMatrixSetToIdentity( utilities_FortranMatrix* mtx ) {
164
165 HYPRE_BigInt j, h, w, jump;
166 HYPRE_Real* p;
167
168 hypre_assert( mtx != NULL );
169
170 utilities_FortranMatrixClear( mtx );
171
172 h = mtx->height;
173 w = mtx->width;
174
175 jump = mtx->globalHeight;
176
177 for ( j = 0, p = mtx->value; j < w && j < h; j++, p += jump )
178 *p++ = 1.0;
179
180 }
181
182 void
utilities_FortranMatrixTransposeSquare(utilities_FortranMatrix * mtx)183 utilities_FortranMatrixTransposeSquare( utilities_FortranMatrix* mtx ) {
184
185 HYPRE_BigInt i, j, g, h, w, jump;
186 HYPRE_Real* p;
187 HYPRE_Real* q;
188 HYPRE_Real tmp;
189
190 hypre_assert( mtx != NULL );
191
192 g = mtx->globalHeight;
193 h = mtx->height;
194 w = mtx->width;
195
196 hypre_assert( h == w );
197
198 jump = mtx->globalHeight - h;
199
200 for ( j = 0, p = mtx->value; j < w; j++ ) {
201 q = p;
202 p++;
203 q += g;
204 for ( i = j + 1; i < h; i++, p++, q += g ) {
205 tmp = *p;
206 *p = *q;
207 *q = tmp;
208 }
209 p += ++jump;
210 }
211 }
212
213 void
utilities_FortranMatrixSymmetrize(utilities_FortranMatrix * mtx)214 utilities_FortranMatrixSymmetrize( utilities_FortranMatrix* mtx ) {
215
216 HYPRE_BigInt i, j, g, h, w, jump;
217 HYPRE_Real* p;
218 HYPRE_Real* q;
219
220 hypre_assert( mtx != NULL );
221
222 g = mtx->globalHeight;
223 h = mtx->height;
224 w = mtx->width;
225
226 hypre_assert( h == w );
227
228 jump = mtx->globalHeight - h;
229
230 for ( j = 0, p = mtx->value; j < w; j++ ) {
231 q = p;
232 p++;
233 q += g;
234 for ( i = j + 1; i < h; i++, p++, q += g )
235 *p = *q = (*p + *q)*0.5;
236 p += ++jump;
237 }
238 }
239
240 void
utilities_FortranMatrixCopy(utilities_FortranMatrix * src,HYPRE_Int t,utilities_FortranMatrix * dest)241 utilities_FortranMatrixCopy( utilities_FortranMatrix* src, HYPRE_Int t,
242 utilities_FortranMatrix* dest ) {
243
244 HYPRE_BigInt i, j, h, w;
245 HYPRE_BigInt jp, jq, jr;
246 HYPRE_Real* p;
247 HYPRE_Real* q;
248 HYPRE_Real* r;
249
250 hypre_assert( src != NULL && dest != NULL );
251
252 h = dest->height;
253 w = dest->width;
254
255 jp = dest->globalHeight - h;
256
257 if ( t == 0 ) {
258 hypre_assert( src->height == h && src->width == w );
259 jq = 1;
260 jr = src->globalHeight;
261 }
262 else {
263 hypre_assert( src->height == w && src->width == h );
264 jr = 1;
265 jq = src->globalHeight;
266 }
267
268 for ( j = 0, p = dest->value, r = src->value; j < w; j++, p += jp, r += jr )
269 for ( i = 0, q = r; i < h; i++, p++, q += jq )
270 *p = *q;
271 }
272
273 void
utilities_FortranMatrixIndexCopy(HYPRE_Int * index,utilities_FortranMatrix * src,HYPRE_Int t,utilities_FortranMatrix * dest)274 utilities_FortranMatrixIndexCopy( HYPRE_Int* index,
275 utilities_FortranMatrix* src, HYPRE_Int t,
276 utilities_FortranMatrix* dest ) {
277
278 HYPRE_BigInt i, j, h, w;
279 HYPRE_BigInt jp, jq, jr;
280 HYPRE_Real* p;
281 HYPRE_Real* q;
282 HYPRE_Real* r;
283
284 hypre_assert( src != NULL && dest != NULL );
285
286 h = dest->height;
287 w = dest->width;
288
289 jp = dest->globalHeight - h;
290
291 if ( t == 0 ) {
292 hypre_assert( src->height == h && src->width == w );
293 jq = 1;
294 jr = src->globalHeight;
295 }
296 else {
297 hypre_assert( src->height == w && src->width == h );
298 jr = 1;
299 jq = src->globalHeight;
300 }
301
302 for ( j = 0, p = dest->value; j < w; j++, p += jp ) {
303 r = src->value + (index[j]-1)*jr;
304 for ( i = 0, q = r; i < h; i++, p++, q += jq )
305 *p = *q;
306 }
307 }
308
309 void
utilities_FortranMatrixSetDiagonal(utilities_FortranMatrix * mtx,utilities_FortranMatrix * vec)310 utilities_FortranMatrixSetDiagonal( utilities_FortranMatrix* mtx,
311 utilities_FortranMatrix* vec ) {
312
313 HYPRE_BigInt j, h, w, jump;
314 HYPRE_Real* p;
315 HYPRE_Real* q;
316
317 hypre_assert( mtx != NULL && vec != NULL );
318
319 h = mtx->height;
320 w = mtx->width;
321
322 hypre_assert( vec->height >= h );
323
324 jump = mtx->globalHeight + 1;
325
326 for ( j = 0, p = mtx->value, q = vec->value; j < w && j < h;
327 j++, p += jump, q++ )
328 *p = *q;
329
330 }
331
332 void
utilities_FortranMatrixGetDiagonal(utilities_FortranMatrix * mtx,utilities_FortranMatrix * vec)333 utilities_FortranMatrixGetDiagonal( utilities_FortranMatrix* mtx,
334 utilities_FortranMatrix* vec ) {
335
336 HYPRE_BigInt j, h, w, jump;
337 HYPRE_Real* p;
338 HYPRE_Real* q;
339
340 hypre_assert( mtx != NULL && vec != NULL );
341
342 h = mtx->height;
343 w = mtx->width;
344
345 hypre_assert( vec->height >= h );
346
347 jump = mtx->globalHeight + 1;
348
349 for ( j = 0, p = mtx->value, q = vec->value; j < w && j < h;
350 j++, p += jump, q++ )
351 *q = *p;
352
353 }
354
355 void
utilities_FortranMatrixAdd(HYPRE_Real a,utilities_FortranMatrix * mtxA,utilities_FortranMatrix * mtxB,utilities_FortranMatrix * mtxC)356 utilities_FortranMatrixAdd( HYPRE_Real a,
357 utilities_FortranMatrix* mtxA,
358 utilities_FortranMatrix* mtxB,
359 utilities_FortranMatrix* mtxC ) {
360
361 HYPRE_BigInt i, j, h, w, jA, jB, jC;
362 HYPRE_Real *pA;
363 HYPRE_Real *pB;
364 HYPRE_Real *pC;
365
366 hypre_assert( mtxA != NULL && mtxB != NULL && mtxC != NULL );
367
368 h = mtxA->height;
369 w = mtxA->width;
370
371 hypre_assert( mtxB->height == h && mtxB->width == w );
372 hypre_assert( mtxC->height == h && mtxC->width == w );
373
374 jA = mtxA->globalHeight - h;
375 jB = mtxB->globalHeight - h;
376 jC = mtxC->globalHeight - h;
377
378 pA = mtxA->value;
379 pB = mtxB->value;
380 pC = mtxC->value;
381
382 if ( a == 0.0 ) {
383 for ( j = 0; j < w; j++ ) {
384 for ( i = 0; i < h; i++, pA++, pB++, pC++ )
385 *pC = *pB;
386 pA += jA;
387 pB += jB;
388 pC += jC;
389 }
390 }
391 else if ( a == 1.0 ) {
392 for ( j = 0; j < w; j++ ) {
393 for ( i = 0; i < h; i++, pA++, pB++, pC++ )
394 *pC = *pA + *pB;
395 pA += jA;
396 pB += jB;
397 pC += jC;
398 }
399 }
400 else if ( a == -1.0 ) {
401 for ( j = 0; j < w; j++ ) {
402 for ( i = 0; i < h; i++, pA++, pB++, pC++ )
403 *pC = *pB - *pA;
404 pA += jA;
405 pB += jB;
406 pC += jC;
407 }
408 }
409 else {
410 for ( j = 0; j < w; j++ ) {
411 for ( i = 0; i < h; i++, pA++, pB++, pC++ )
412 *pC = *pA * a + *pB;
413 pA += jA;
414 pB += jB;
415 pC += jC;
416 }
417 }
418 }
419
420 void
utilities_FortranMatrixDMultiply(utilities_FortranMatrix * vec,utilities_FortranMatrix * mtx)421 utilities_FortranMatrixDMultiply( utilities_FortranMatrix* vec,
422 utilities_FortranMatrix* mtx ) {
423
424 HYPRE_BigInt i, j, h, w, jump;
425 HYPRE_Real* p;
426 HYPRE_Real* q;
427
428 hypre_assert( mtx != NULL && vec != NULL );
429
430 h = mtx->height;
431 w = mtx->width;
432
433 hypre_assert( vec->height == h );
434
435 jump = mtx->globalHeight - h;
436
437 for ( j = 0, p = mtx->value; j < w; j++ ) {
438 for ( i = 0, q = vec->value; i < h; i++, p++, q++ )
439 *p = *p * (*q);
440 p += jump;
441 }
442
443 }
444
445 void
utilities_FortranMatrixMultiplyD(utilities_FortranMatrix * mtx,utilities_FortranMatrix * vec)446 utilities_FortranMatrixMultiplyD( utilities_FortranMatrix* mtx,
447 utilities_FortranMatrix* vec ) {
448
449 HYPRE_BigInt i, j, h, w, jump;
450 HYPRE_Real* p;
451 HYPRE_Real* q;
452
453 hypre_assert( mtx != NULL && vec != NULL );
454
455 h = mtx->height;
456 w = mtx->width;
457
458 hypre_assert( vec->height == w );
459
460 jump = mtx->globalHeight - h;
461
462 for ( j = 0, q = vec->value, p = mtx->value; j < w; j++, q++ ) {
463 for ( i = 0; i < h; i++, p++)
464 *p = *p * (*q);
465 p += jump;
466 }
467
468 }
469
470 void
utilities_FortranMatrixMultiply(utilities_FortranMatrix * mtxA,HYPRE_Int tA,utilities_FortranMatrix * mtxB,HYPRE_Int tB,utilities_FortranMatrix * mtxC)471 utilities_FortranMatrixMultiply( utilities_FortranMatrix* mtxA, HYPRE_Int tA,
472 utilities_FortranMatrix* mtxB, HYPRE_Int tB,
473 utilities_FortranMatrix* mtxC ) {
474 HYPRE_BigInt h, w;
475 HYPRE_BigInt i, j, k, l;
476 HYPRE_BigInt iA, kA;
477 HYPRE_BigInt kB, jB;
478 HYPRE_BigInt iC, jC;
479
480 HYPRE_Real* pAi0;
481 HYPRE_Real* pAik;
482 HYPRE_Real* pB0j;
483 HYPRE_Real* pBkj;
484 HYPRE_Real* pC0j;
485 HYPRE_Real* pCij;
486
487 HYPRE_Real s;
488
489 hypre_assert( mtxA != NULL && mtxB != NULL && mtxC != NULL );
490
491 h = mtxC->height;
492 w = mtxC->width;
493 iC = 1;
494 jC = mtxC->globalHeight;
495
496 if ( tA == 0 ) {
497 hypre_assert( mtxA->height == h );
498 l = mtxA->width;
499 iA = 1;
500 kA = mtxA->globalHeight;
501 }
502 else {
503 l = mtxA->height;
504 hypre_assert( mtxA->width == h );
505 kA = 1;
506 iA = mtxA->globalHeight;
507 }
508
509 if ( tB == 0 ) {
510 hypre_assert( mtxB->height == l );
511 hypre_assert( mtxB->width == w );
512 kB = 1;
513 jB = mtxB->globalHeight;
514 }
515 else {
516 hypre_assert( mtxB->width == l );
517 hypre_assert( mtxB->height == w );
518 jB = 1;
519 kB = mtxB->globalHeight;
520 }
521
522 for ( j = 0, pB0j = mtxB->value, pC0j = mtxC->value; j < w;
523 j++, pB0j += jB, pC0j += jC )
524 for ( i = 0, pCij = pC0j, pAi0 = mtxA->value; i < h;
525 i++, pCij += iC, pAi0 += iA ) {
526 s = 0.0;
527 for ( k = 0, pAik = pAi0, pBkj = pB0j; k < l;
528 k++, pAik += kA, pBkj += kB )
529 s += *pAik * (*pBkj);
530 *pCij = s;
531 }
532 }
533
534 HYPRE_Real
utilities_FortranMatrixFNorm(utilities_FortranMatrix * mtx)535 utilities_FortranMatrixFNorm( utilities_FortranMatrix* mtx ) {
536
537 HYPRE_BigInt i, j, h, w, jump;
538 HYPRE_Real* p;
539
540 HYPRE_Real norm;
541
542 hypre_assert( mtx != NULL );
543
544 h = mtx->height;
545 w = mtx->width;
546
547 jump = mtx->globalHeight - h;
548
549 norm = 0.0;
550
551 for ( j = 0, p = mtx->value; j < w; j++ ) {
552 for ( i = 0; i < h; i++, p++ )
553 norm += (*p) * (*p);
554 p += jump;
555 }
556
557 norm = sqrt(norm);
558 return norm;
559 }
560
561 HYPRE_Real
utilities_FortranMatrixValue(utilities_FortranMatrix * mtx,HYPRE_BigInt i,HYPRE_BigInt j)562 utilities_FortranMatrixValue( utilities_FortranMatrix* mtx,
563 HYPRE_BigInt i, HYPRE_BigInt j ) {
564
565 HYPRE_BigInt k;
566
567 hypre_assert( mtx != NULL );
568
569 hypre_assert( 1 <= i && i <= mtx->height );
570 hypre_assert( 1 <= j && j <= mtx->width );
571
572 k = i - 1 + (j - 1)*mtx->globalHeight;
573 return mtx->value[k];
574 }
575
576 HYPRE_Real*
utilities_FortranMatrixValuePtr(utilities_FortranMatrix * mtx,HYPRE_BigInt i,HYPRE_BigInt j)577 utilities_FortranMatrixValuePtr( utilities_FortranMatrix* mtx,
578 HYPRE_BigInt i, HYPRE_BigInt j ) {
579
580 HYPRE_BigInt k;
581
582 hypre_assert( mtx != NULL );
583
584 hypre_assert( 1 <= i && i <= mtx->height );
585 hypre_assert( 1 <= j && j <= mtx->width );
586
587 k = i - 1 + (j - 1)*mtx->globalHeight;
588 return mtx->value + k;
589 }
590
591 HYPRE_Real
utilities_FortranMatrixMaxValue(utilities_FortranMatrix * mtx)592 utilities_FortranMatrixMaxValue( utilities_FortranMatrix* mtx ) {
593
594 HYPRE_BigInt i, j, jump;
595 HYPRE_BigInt h, w;
596 HYPRE_Real* p;
597 HYPRE_Real maxVal;
598
599 hypre_assert( mtx != NULL );
600
601 h = mtx->height;
602 w = mtx->width;
603
604 jump = mtx->globalHeight - h;
605
606 maxVal = mtx->value[0];
607
608 for ( j = 0, p = mtx->value; j < w; j++ ) {
609 for ( i = 0; i < h; i++, p++ )
610 if ( *p > maxVal )
611 maxVal = *p;
612 p += jump;
613 }
614
615 return maxVal;
616 }
617
618 void
utilities_FortranMatrixSelectBlock(utilities_FortranMatrix * mtx,HYPRE_BigInt iFrom,HYPRE_BigInt iTo,HYPRE_BigInt jFrom,HYPRE_BigInt jTo,utilities_FortranMatrix * block)619 utilities_FortranMatrixSelectBlock( utilities_FortranMatrix* mtx,
620 HYPRE_BigInt iFrom, HYPRE_BigInt iTo,
621 HYPRE_BigInt jFrom, HYPRE_BigInt jTo,
622 utilities_FortranMatrix* block ) {
623
624 if ( block->value != NULL && block->ownsValues )
625 hypre_TFree( block->value ,HYPRE_MEMORY_HOST);
626
627 block->globalHeight = mtx->globalHeight;
628 if ( iTo < iFrom || jTo < jFrom ) {
629 block->height = 0;
630 block->width = 0;
631 block->value = NULL;
632 return;
633 }
634 block->height = iTo - iFrom + 1;
635 block->width = jTo - jFrom + 1;
636 block->value = mtx->value + iFrom - 1 + (jFrom - 1)*mtx->globalHeight;
637 block->ownsValues = 0;
638 }
639
640 void
utilities_FortranMatrixUpperInv(utilities_FortranMatrix * u)641 utilities_FortranMatrixUpperInv( utilities_FortranMatrix* u ) {
642
643 HYPRE_BigInt i, j, k;
644 HYPRE_BigInt n, jc, jd;
645 HYPRE_Real v;
646 HYPRE_Real* diag; /* diag(i) = u(i,i)_original */
647 HYPRE_Real* pin; /* &u(i-1,n) */
648 HYPRE_Real* pii; /* &u(i,i) */
649 HYPRE_Real* pij; /* &u(i,j) */
650 HYPRE_Real* pik; /* &u(i,k) */
651 HYPRE_Real* pkj; /* &u(k,j) */
652 HYPRE_Real* pd; /* &diag(i) */
653
654 n = u->height;
655 hypre_assert( u->width == n );
656
657 diag = hypre_CTAlloc(HYPRE_Real, n, HYPRE_MEMORY_HOST);
658 hypre_assert( diag != NULL );
659
660 jc = u->globalHeight;
661 jd = jc + 1;
662
663 pii = u->value;
664 pd = diag;
665 for ( i = 0; i < n; i++, pii += jd, pd++ ) {
666 v = *pd = *pii;
667 *pii = 1.0/v;
668 }
669
670 pii -= jd;
671 pin = pii - 1;
672 pii -= jd;
673 pd -= 2;
674 for ( i = n - 1; i > 0; i--, pii -= jd, pin--, pd-- ) {
675 pij = pin;
676 for ( j = n; j > i; j--, pij -= jc ) {
677 v = 0;
678 pik = pii + jc;
679 pkj = pij + 1;
680 for ( k = i + 1; k <= j; k++, pik += jc, pkj++ ) {
681 v -= (*pik) * (*pkj);
682 }
683 *pij = v/(*pd);
684 }
685 }
686
687 hypre_TFree( diag ,HYPRE_MEMORY_HOST);
688
689 }
690
691 HYPRE_Int
utilities_FortranMatrixPrint(utilities_FortranMatrix * mtx,const char * fileName)692 utilities_FortranMatrixPrint( utilities_FortranMatrix* mtx, const char *fileName) {
693
694 HYPRE_BigInt i, j, h, w, jump;
695 HYPRE_Real* p;
696 FILE* fp;
697
698 hypre_assert( mtx != NULL );
699
700 if ( !(fp = fopen(fileName,"w")) )
701 return 1;
702
703 h = mtx->height;
704 w = mtx->width;
705
706 hypre_fprintf(fp,"%ld\n",h);
707 hypre_fprintf(fp,"%ld\n",w);
708
709 jump = mtx->globalHeight - h;
710
711 for ( j = 0, p = mtx->value; j < w; j++ ) {
712 for ( i = 0; i < h; i++, p++ )
713 hypre_fprintf(fp,"%.14e\n",*p);
714 p += jump;
715 }
716
717 fclose(fp);
718 return 0;
719 }
720
721