/[escript]/trunk-mpi-branch/paso/src/SystemMatrix_saveHB.c
ViewVC logotype

Contents of /trunk-mpi-branch/paso/src/SystemMatrix_saveHB.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1136 - (show annotations)
Wed May 9 04:29:15 2007 UTC (12 years, 3 months ago) by ksteube
File MIME type: text/plain
File size: 7563 byte(s)
MPI branch merged with trunk revisions 1126 to 967

1 /* $Id$ */
2
3 /*
4 ********************************************************************************
5 * Copyright 2006 by ACcESS MNRF *
6 * *
7 * http://www.access.edu.au *
8 * Primary Business: Queensland, Australia *
9 * Licensed under the Open Software License version 3.0 *
10 * http://www.opensource.org/licenses/osl-3.0.php *
11 ********************************************************************************
12 */
13
14 /**************************************************************/
15
16 /* Paso: SystemMatrix is saved to Harwell-Boeing format */
17
18 /**************************************************************/
19
20 /* Copyright: ACcESS Australia 2005 */
21 /* Author: imran@esscc.uq.edu.au */
22
23 /**************************************************************/
24
25 #include "Paso.h"
26 #include "SystemMatrix.h"
27
28 /* TODO: Refactor the stuff in here into hbio, just like mmio! */
29
30 static dim_t M, N, nz;
31
32 static int calc_digits( int );
33 static void fmt_str( int, int, int*, int*, int*, char*, char* );
34 static void print_data( FILE*, int, int, int, char*, void*, int, int );
35 static void generate_HB( FILE*, dim_t*, dim_t*, double* );
36
37 /* function to get number of digits in an integer */
38 int calc_digits( int var )
39 {
40 int digits = 1;
41 while( (var/=10) )
42 digits++;
43
44 return digits;
45 }
46
47 /* function to generate the format string.
48 *
49 * use maxlen to determine no. of entries per line
50 * use nvalues to determine no. of lines
51 */
52 void fmt_str( int nvalues, int integer, int *width, int *nlines, int *nperline, char *pfmt, char *fmt )
53 {
54 int per_line;
55 int maxlen = *width;
56
57 if( integer && maxlen < 10 )
58 maxlen = 10;
59 else
60 maxlen = 13;
61
62 per_line = 80 / maxlen;
63 *nlines = nvalues / per_line;
64 if( nvalues % per_line )
65 (*nlines)++;
66
67 *nperline = per_line;
68 if( integer )
69 {
70 sprintf( pfmt, "(%dI%d)", per_line, maxlen );
71 sprintf( fmt, "%%%dd", maxlen );
72 }
73 else
74 {
75 sprintf( pfmt, "(1P%dE%d.6)", per_line, maxlen );
76 sprintf( fmt, "%%%d.6E", maxlen );
77 }
78 *width = maxlen;
79 }
80
81 /* function to print the actual data in the right format */
82 void print_data( FILE *fp, int n_perline, int width, int nval, char *fmt, void *ptr, int integer, int adjust )
83 {
84 double *data = ptr;
85 int entries_done = 0;
86 int padding, i, j;
87 char pad_fmt[10];
88
89 padding = 80 - n_perline*width;
90 sprintf( pad_fmt, "%%%dc", padding );
91
92 if( adjust != 1 )
93 adjust = 0;
94
95 if( integer )
96 {
97 dim_t *data = ptr;
98 for(i=0; i<nval; i++ )
99 {
100 fprintf( fp, fmt, data[i]+adjust );
101 entries_done++;
102 if( entries_done == n_perline )
103 {
104 if( padding )
105 fprintf( fp, pad_fmt, ' ' );
106 fprintf( fp, "\n" );
107 entries_done = 0;
108 }
109 }
110 }
111 else
112 {
113 for(i=0; i<nval; i++ )
114 {
115 fprintf( fp, fmt, data[i] );
116 entries_done++;
117 if( entries_done == n_perline )
118 {
119 if( padding )
120 fprintf( fp, pad_fmt, ' ' );
121 fprintf( fp, "\n" );
122 entries_done = 0;
123 }
124 }
125 }
126 if( entries_done )
127 {
128 sprintf( pad_fmt, "%%%dc\n", (80 - entries_done*width) );
129 fprintf( fp, pad_fmt, ' ' );
130 }
131 }
132
133 void generate_HB( FILE *fp, dim_t *col_ptr, dim_t *row_ind, double *val )
134 {
135 char buffer[81];
136
137 int val_lines, ind_lines, ptr_lines;
138 int val_perline, ind_perline, ptr_perline;
139 int val_width, ind_width, ptr_width;
140 char ptr_pfmt[7], ind_pfmt[7], val_pfmt[11];
141 char ptr_fmt[10], ind_fmt[10], val_fmt[10];
142
143 /* line 1 */
144 sprintf( buffer, "%-72s%-8s", "Matrix Title", "Key" );
145 buffer[80] = '\0';
146 fprintf( fp, "%s\n", buffer );
147
148 /* line 2 */
149 ptr_width = calc_digits( nz+1 );
150 fmt_str( N+1, 1, &ptr_width, &ptr_lines, &ptr_perline, ptr_pfmt, ptr_fmt );
151 ind_width = calc_digits( N );
152 fmt_str( nz, 1, &ind_width, &ind_lines, &ind_perline, ind_pfmt, ind_fmt );
153 val_width = 13;
154 fmt_str( nz, 0, &val_width, &val_lines, &val_perline, val_pfmt, val_fmt );
155 sprintf( buffer, "%14d%14d%14d%14d%14d%10c", (ptr_lines+ind_lines+val_lines), ptr_lines, ind_lines, val_lines, 0, ' ' );
156 buffer[80] = '\0';
157 fprintf( fp, "%s\n", buffer );
158
159 /* line 3 */
160 sprintf( buffer, "%c%c%c%11c%14d%14d%14d%14d%10c", 'R', 'U', 'A', ' ', M, N, nz, 0, ' ' );
161 buffer[80] = '\0';
162 fprintf( fp, "%s\n", buffer );
163
164 /* line 4 */
165 sprintf( buffer, "%16s%16s%20s%28c", ptr_pfmt, ind_pfmt, val_pfmt, ' ');
166 buffer[80]='\0';
167 fprintf( fp, "%s\n", buffer );
168
169 /* line 5 */
170 /* NOT PRESENT */
171
172 /* write the actual data */
173 print_data( fp, ptr_perline, ptr_width, (N+1), ptr_fmt, col_ptr, 1, 1 );
174 print_data( fp, ind_perline, ind_width, nz, ind_fmt, row_ind, 1, 1 );
175 print_data( fp, val_perline, val_width, nz, val_fmt, val, 0, 0 );
176 }
177
178 void Paso_SystemMatrix_saveHB( Paso_SystemMatrix *A_p, char *filename_p )
179 {
180 #ifdef PASO_MPI
181 Paso_setError(TYPE_ERROR,"Paso_SystemMatrix_saveHB: MPI is currently supported.\n");
182 return;
183 #endif
184 if (A_p->val == NULL) {
185 Paso_setError(TYPE_ERROR,"Paso_SystemMatrix_saveHB: unsupported format detected.\n");
186 return;
187 }
188 int i, curr_col,j ;
189 int iPtr, iCol, ir, ic;
190 index_t index_offset=(A_p->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
191 int nz = A_p->myLen;
192 /* open the file */
193 FILE *fileHandle_p = fopen( filename_p, "w" );
194 if( fileHandle_p == NULL )
195 {
196 Paso_setError(IO_ERROR,"File could not be opened for writing.");
197 return;
198 }
199
200
201 /* check the internal storage type
202 * currently ONLY support for CSC
203 */
204 if ( A_p->type & MATRIX_FORMAT_CSC) {
205 if( A_p->row_block_size == 1 && A_p->col_block_size == 1 )
206 {
207 M = A_p->numRows;
208 N = A_p->numCols;
209 generate_HB( fileHandle_p, A_p->pattern->ptr, A_p->pattern->index, A_p->val );
210 }
211 else
212 {
213
214 // fprintf( fileHandle_p, "NEED unrolling!\n" );
215 M = A_p->numRows*A_p->row_block_size;
216 N = A_p->numCols*A_p->col_block_size;
217 nz = A_p->len;
218
219 dim_t *row_ind = MEMALLOC( nz, dim_t );
220 dim_t *col_ind = MEMALLOC( nz, dim_t );
221
222 i = 0;
223 for( iCol=0; iCol<A_p->pattern->myNumOutput; iCol++ )
224 for( ic=0; ic<A_p->col_block_size; ic++)
225 for( iPtr=A_p->pattern->ptr[iCol]-index_offset; iPtr<A_p->pattern->ptr[iCol+1]-index_offset; iPtr++)
226 for( ir=0; ir<A_p->row_block_size; ir++ )
227 {
228 row_ind[i] = (A_p->pattern->index[iPtr]-index_offset)*A_p->row_block_size+ir+1;
229 col_ind[i] = iCol*A_p->col_block_size+ic+1;
230 i++;
231 }
232
233 /* get the col_ptr */
234 dim_t *col_ptr = MEMALLOC( (N+1), dim_t );
235
236 curr_col = 0;
237 for(j=0; (j<nz && curr_col<N); curr_col++ )
238 {
239 while( col_ind[j] != curr_col )
240 j++;
241 col_ptr[curr_col] = j;
242 }
243 col_ptr[N] = nz;
244
245 /* generate the HB file */
246 generate_HB( fileHandle_p, col_ptr, row_ind, A_p->val );
247
248 /* free the allocated memory */
249 MEMFREE( col_ptr );
250 MEMFREE( col_ind );
251 MEMFREE( row_ind );
252 }
253 } else {
254 Paso_setError(TYPE_ERROR,"Paso_SystemMatrix_saveHB: only CSC is currently supported.\n");
255 }
256
257 /* close the file */
258 fclose( fileHandle_p );
259
260 return;
261 }
262
263 /*
264 * $Log$
265 * Revision 1.2 2005/09/15 03:44:39 jgs
266 * Merge of development branch dev-02 back to main trunk on 2005-09-15
267 *
268 * Revision 1.1.2.1 2005/09/05 06:29:48 gross
269 * These files have been extracted from finley to define a stand alone libray for iterative
270 * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
271 * has not been tested yet.
272 *
273 *
274 */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26