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

Contents of /trunk/paso/src/SystemMatrix_saveHB.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 150 - (show annotations)
Thu Sep 15 03:44:45 2005 UTC (14 years, 1 month ago) by jgs
Original Path: trunk/esys2/paso/src/SystemMatrix_saveHB.c
File MIME type: text/plain
File size: 6597 byte(s)
Merge of development branch dev-02 back to main trunk on 2005-09-15

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26