paso/src/SystemMatrix_saveHB.c

/* $Id$ */

/*
********************************************************************************
*               Copyright   2006 by ACcESS MNRF                                *
*                                                                              * 
*                 http://www.access.edu.au                                     *
*           Primary Business: Queensland, Australia                            *
*     Licensed under the Open Software License version 3.0             *
*        http://www.opensource.org/licenses/osl-3.0.php                        *
********************************************************************************
*/

/**************************************************************/

/* Paso: SystemMatrix is saved to Harwell-Boeing format     */

/**************************************************************/

/* Copyright: ACcESS Australia 2005                           */
/* Author: imran@esscc.uq.edu.au                              */

/**************************************************************/

#include "Paso.h"
#include "SystemMatrix.h"

/* TODO: Refactor the stuff in here into hbio, just like mmio! */

static dim_t M, N, nz;

static int calc_digits( int );
static void fmt_str( int, int, int*, int*, int*, char*, char* );
static void print_data( FILE*, int, int, int, char*, void*, int, int );
static void generate_HB( FILE*, dim_t*, dim_t*, double* );

/* function to get number of digits in an integer */
int calc_digits( int var )
{
    int digits = 1;
    while( (var/=10) )
        digits++;

    return digits;
}

/* function to generate the format string.
 *
 * use maxlen to determine no. of entries per line
 * use nvalues to determine no. of lines
 */
void fmt_str( int nvalues, int integer, int *width, int *nlines, int *nperline, char *pfmt, char *fmt )
{
    int per_line;
    int maxlen = *width;

    if( integer && maxlen < 10 )
        maxlen = 10;
    else
        maxlen = 13;

    per_line = 80 / maxlen;
    *nlines = nvalues / per_line;
    if( nvalues % per_line )
        (*nlines)++;

    *nperline = per_line;
    if( integer )
    {
        sprintf( pfmt, "(%dI%d)", per_line, maxlen );
        sprintf( fmt, "%%%dd", maxlen );
    }
    else
    {
        sprintf( pfmt, "(1P%dE%d.6)", per_line, maxlen );
        sprintf( fmt, "%%%d.6E", maxlen );
    }
    *width = maxlen;
}

/* function to print the actual data in the right format */
void print_data( FILE *fp, int n_perline, int width, int nval, char *fmt, void *ptr, int integer, int adjust )
{
    double *data = ptr;
    int entries_done = 0;
    int padding, i, j;
    char pad_fmt[10];

    padding = 80 - n_perline*width;
    sprintf( pad_fmt, "%%%dc", padding );

    if( adjust != 1 )
        adjust = 0;

    if( integer )
    {
        dim_t *data = ptr;
        for(i=0; i<nval; i++ )
        {
            fprintf( fp, fmt, data[i]+adjust );
            entries_done++;
            if( entries_done == n_perline )
            {
                if( padding )
                    fprintf( fp, pad_fmt, ' ' );
                fprintf( fp, "\n" );
                entries_done = 0;
            }
        }
    }
    else
    {
        for(i=0; i<nval; i++ )
        {
            fprintf( fp, fmt, data[i] );
            entries_done++;
            if( entries_done == n_perline )
            {
                if( padding )
                    fprintf( fp, pad_fmt, ' ' );
                fprintf( fp, "\n" );
                entries_done = 0;
            }
        }
    }
    if( entries_done )
    {
        sprintf( pad_fmt, "%%%dc\n", (80 - entries_done*width) );
        fprintf( fp, pad_fmt, ' ' );
    }
}

void generate_HB( FILE *fp, dim_t *col_ptr, dim_t *row_ind, double *val )
{
    char buffer[81];

    int val_lines, ind_lines, ptr_lines;
    int val_perline, ind_perline, ptr_perline;
    int val_width, ind_width, ptr_width;
    char ptr_pfmt[7], ind_pfmt[7], val_pfmt[11];
    char ptr_fmt[10], ind_fmt[10], val_fmt[10];

    /* line 1 */
    sprintf( buffer, "%-72s%-8s", "Matrix Title", "Key" );
    buffer[80] = '\0';
    fprintf( fp, "%s\n", buffer );

    /* line 2 */
    ptr_width = calc_digits( nz+1 );
    fmt_str( N+1, 1, &ptr_width, &ptr_lines, &ptr_perline, ptr_pfmt, ptr_fmt );
    ind_width = calc_digits( N );
    fmt_str( nz, 1, &ind_width, &ind_lines, &ind_perline, ind_pfmt, ind_fmt );
    val_width = 13;
    fmt_str( nz, 0, &val_width, &val_lines, &val_perline, val_pfmt, val_fmt );
    sprintf( buffer, "%14d%14d%14d%14d%14d%10c", (ptr_lines+ind_lines+val_lines), ptr_lines, ind_lines, val_lines, 0, ' ' );
    buffer[80] = '\0';
    fprintf( fp, "%s\n", buffer );

    /* line 3 */
    sprintf( buffer, "%c%c%c%11c%14d%14d%14d%14d%10c", 'R', 'U', 'A', ' ', M, N, nz, 0, ' ' );
    buffer[80] = '\0';
    fprintf( fp, "%s\n", buffer );

    /* line 4 */
    sprintf( buffer, "%16s%16s%20s%28c", ptr_pfmt, ind_pfmt, val_pfmt, ' ');
    buffer[80]='\0';
    fprintf( fp, "%s\n", buffer );

    /* line 5 */
    /* NOT PRESENT */

    /* write the actual data */
    print_data( fp, ptr_perline, ptr_width, (N+1), ptr_fmt, col_ptr, 1, 1 );
    print_data( fp, ind_perline, ind_width, nz, ind_fmt, row_ind, 1, 1 );
    print_data( fp, val_perline, val_width, nz, val_fmt, val, 0, 0 );
}

void Paso_SystemMatrix_saveHB( Paso_SystemMatrix *A_p, char *filename_p )
{
        int i, curr_col,j ;
    int iPtr, iCol, ir, ic;
        index_t index_offset=(A_p->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
    /* open the file */
    FILE *fileHandle_p = fopen( filename_p, "w" );
    if( fileHandle_p == NULL )
    {
        Paso_setError(IO_ERROR,"File could not be opened for writing.");
        return;
    }

        
    /* check the internal storage type
     * currently ONLY support for CSC
     */
    if ( A_p->type & MATRIX_FORMAT_CSC) {
                if( A_p->row_block_size == 1 && A_p->col_block_size == 1 )
                {
                    M = A_p->num_rows;
                    N = A_p->num_cols;
                    nz = A_p->len;
                    generate_HB( fileHandle_p, A_p->pattern->ptr, A_p->pattern->index, A_p->val );
                }
                else
                {

//                  fprintf( fileHandle_p, "NEED unrolling!\n" );
                    M = A_p->num_rows*A_p->row_block_size;
                    N = A_p->num_cols*A_p->col_block_size;
                    nz = A_p->len;

                    dim_t *row_ind = MEMALLOC( nz, dim_t );
                    dim_t *col_ind = MEMALLOC( nz, dim_t );

                    i = 0;
                    for( iCol=0; iCol<A_p->pattern->n_ptr; iCol++ )
                        for( ic=0; ic<A_p->col_block_size; ic++)
                            for( iPtr=A_p->pattern->ptr[iCol]-index_offset; iPtr<A_p->pattern->ptr[iCol+1]-index_offset; iPtr++)
                                for( ir=0; ir<A_p->row_block_size; ir++ )
                                {
                                    row_ind[i] = (A_p->pattern->index[iPtr]-index_offset)*A_p->row_block_size+ir+1;
                                    col_ind[i] = iCol*A_p->col_block_size+ic+1;
                                    i++;
                                }

                    /* get the col_ptr */
                    dim_t *col_ptr = MEMALLOC( (N+1), dim_t );

                    curr_col = 0;
                    for(j=0; (j<nz && curr_col<N); curr_col++ )
                    {
                        while( col_ind[j] != curr_col )
                            j++;
                        col_ptr[curr_col] = j;
                    }
                    col_ptr[N] = nz;

                    /* generate the HB file */
                    generate_HB( fileHandle_p, col_ptr, row_ind, A_p->val );

                    /* free the allocated memory */
                    MEMFREE( col_ptr );
                    MEMFREE( col_ind );
                    MEMFREE( row_ind );
                }
        } else {
                Paso_setError(TYPE_ERROR,"Paso_SystemMatrix_saveHB: only CSC is currently supported.\n");
    }

    /* close the file */
    fclose( fileHandle_p );

    return;
}

/*
 * $Log$
 * Revision 1.2  2005/09/15 03:44:39  jgs
 * Merge of development branch dev-02 back to main trunk on 2005-09-15
 *
 * Revision 1.1.2.1  2005/09/05 06:29:48  gross
 * These files have been extracted from finley to define a stand alone libray for iterative
 * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
 * has not been tested yet.
 *
 *
 */
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	int i, curr_col,j ;
181	int iPtr, iCol, ir, ic;
182	index_t index_offset=(A_p->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
183	/* open the file */
184	FILE *fileHandle_p = fopen( filename_p, "w" );
185	if( fileHandle_p == NULL )
186	{
187	Paso_setError(IO_ERROR,"File could not be opened for writing.");
188	return;
189	}
190
191
192	/* check the internal storage type
193	* currently ONLY support for CSC
194	*/
195	if ( A_p->type & MATRIX_FORMAT_CSC) {
196	if( A_p->row_block_size == 1 && A_p->col_block_size == 1 )
197	{
198	M = A_p->num_rows;
199	N = A_p->num_cols;
200	nz = A_p->len;
201	generate_HB( fileHandle_p, A_p->pattern->ptr, A_p->pattern->index, A_p->val );
202	}
203	else
204	{
205
206	// fprintf( fileHandle_p, "NEED unrolling!\n" );
207	M = A_p->num_rows*A_p->row_block_size;
208	N = A_p->num_cols*A_p->col_block_size;
209	nz = A_p->len;
210
211	dim_t *row_ind = MEMALLOC( nz, dim_t );
212	dim_t *col_ind = MEMALLOC( nz, dim_t );
213
214	i = 0;
215	for( iCol=0; iCol<A_p->pattern->n_ptr; iCol++ )
216	for( ic=0; ic<A_p->col_block_size; ic++)
217	for( iPtr=A_p->pattern->ptr[iCol]-index_offset; iPtr<A_p->pattern->ptr[iCol+1]-index_offset; iPtr++)
218	for( ir=0; ir<A_p->row_block_size; ir++ )
219	{
220	row_ind[i] = (A_p->pattern->index[iPtr]-index_offset)*A_p->row_block_size+ir+1;
221	col_ind[i] = iCol*A_p->col_block_size+ic+1;
222	i++;
223	}
224
225	/* get the col_ptr */
226	dim_t *col_ptr = MEMALLOC( (N+1), dim_t );
227
228	curr_col = 0;
229	for(j=0; (j<nz && curr_col<N); curr_col++ )
230	{
231	while( col_ind[j] != curr_col )
232	j++;
233	col_ptr[curr_col] = j;
234	}
235	col_ptr[N] = nz;
236
237	/* generate the HB file */
238	generate_HB( fileHandle_p, col_ptr, row_ind, A_p->val );
239
240	/* free the allocated memory */
241	MEMFREE( col_ptr );
242	MEMFREE( col_ind );
243	MEMFREE( row_ind );
244	}
245	} else {
246	Paso_setError(TYPE_ERROR,"Paso_SystemMatrix_saveHB: only CSC is currently supported.\n");
247	}
248
249	/* close the file */
250	fclose( fileHandle_p );
251
252	return;
253	}
254
255	/*
256	* $Log$
257	* Revision 1.2 2005/09/15 03:44:39 jgs
258	* Merge of development branch dev-02 back to main trunk on 2005-09-15
259	*
260	* Revision 1.1.2.1 2005/09/05 06:29:48 gross
261	* These files have been extracted from finley to define a stand alone libray for iterative
262	* linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
263	* has not been tested yet.
264	*
265	*
266	*/
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision