iohb.c

Go to the documentation of this file.
00001 /*
00002 Fri Aug 15 16:29:47 EDT 1997
00003 
00004                        Harwell-Boeing File I/O in C
00005                                 V. 1.0
00006 
00007            National Institute of Standards and Technology, MD.
00008                              K.A. Remington
00009 
00010 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00011                                 NOTICE
00012 
00013  Permission to use, copy, modify, and distribute this software and
00014  its documentation for any purpose and without fee is hereby granted
00015  provided that the above copyright notice appear in all copies and
00016  that both the copyright notice and this permission notice appear in
00017  supporting documentation.
00018 
00019  Neither the Author nor the Institution (National Institute of Standards
00020  and Technology) make any representations about the suitability of this 
00021  software for any purpose. This software is provided "as is" without 
00022  expressed or implied warranty.
00023 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00024 
00025                          ---------------------
00026                          INTERFACE DESCRIPTION
00027                          ---------------------
00028   ---------------
00029   QUERY FUNCTIONS
00030   ---------------
00031 
00032   FUNCTION:
00033 
00034   int readHB_info(const char *filename, int *M, int *N, int *nz,
00035   char **Type, int *Nrhs)
00036 
00037   DESCRIPTION:
00038 
00039   The readHB_info function opens and reads the header information from
00040   the specified Harwell-Boeing file, and reports back the number of rows
00041   and columns in the stored matrix (M and N), the number of nonzeros in
00042   the matrix (nz), the 3-character matrix type(Type), and the number of
00043   right-hand-sides stored along with the matrix (Nrhs).  This function
00044   is designed to retrieve basic size information which can be used to 
00045   allocate arrays.
00046 
00047   FUNCTION:
00048 
00049   int  readHB_header(FILE* in_file, char* Title, char* Key, char* Type, 
00050                     int* Nrow, int* Ncol, int* Nnzero, int* Nrhs,
00051                     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt, 
00052                     int* Ptrcrd, int* Indcrd, int* Valcrd, int* Rhscrd, 
00053                     char *Rhstype)
00054 
00055   DESCRIPTION:
00056 
00057   More detailed than the readHB_info function, readHB_header() reads from 
00058   the specified Harwell-Boeing file all of the header information.  
00059 
00060 
00061   ------------------------------
00062   DOUBLE PRECISION I/O FUNCTIONS
00063   ------------------------------
00064   FUNCTION:
00065 
00066   int readHB_newmat_double(const char *filename, int *M, int *N, *int nz,
00067   int **colptr, int **rowind,  double**val)
00068 
00069   int readHB_mat_double(const char *filename, int *colptr, int *rowind,
00070   double*val)
00071 
00072 
00073   DESCRIPTION:
00074 
00075   This function opens and reads the specified file, interpreting its
00076   contents as a sparse matrix stored in the Harwell/Boeing standard
00077   format.  (See readHB_aux_double to read auxillary vectors.)
00078         -- Values are interpreted as double precision numbers. --
00079 
00080   The "mat" function uses _pre-allocated_ vectors to hold the index and 
00081   nonzero value information.
00082 
00083   The "newmat" function allocates vectors to hold the index and nonzero
00084   value information, and returns pointers to these vectors along with
00085   matrix dimension and number of nonzeros.
00086 
00087   FUNCTION:
00088 
00089   int readHB_aux_double(const char* filename, const char AuxType, double b[])
00090 
00091   int readHB_newaux_double(const char* filename, const char AuxType, double** b)
00092 
00093   DESCRIPTION:
00094 
00095   This function opens and reads from the specified file auxillary vector(s).
00096   The char argument Auxtype determines which type of auxillary vector(s)
00097   will be read (if present in the file).
00098 
00099                   AuxType = 'F'   right-hand-side 
00100                   AuxType = 'G'   initial estimate (Guess)
00101                   AuxType = 'X'   eXact solution
00102 
00103   If Nrhs > 1, all of the Nrhs vectors of the given type are read and 
00104   stored in column-major order in the vector b.
00105 
00106   The "newaux" function allocates a vector to hold the values retrieved.
00107   The "mat" function uses a _pre-allocated_ vector to hold the values.
00108 
00109   FUNCTION:
00110 
00111   int writeHB_mat_double(const char* filename, int M, int N, 
00112                         int nz, const int colptr[], const int rowind[], 
00113                         const double val[], int Nrhs, const double rhs[], 
00114                         const double guess[], const double exact[],
00115                         const char* Title, const char* Key, const char* Type, 
00116                         char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
00117                         const char* Rhstype)
00118 
00119   DESCRIPTION:
00120 
00121   The writeHB_mat_double function opens the named file and writes the specified
00122   matrix and optional auxillary vector(s) to that file in Harwell-Boeing
00123   format.  The format arguments (Ptrfmt,Indfmt,Valfmt, and Rhsfmt) are
00124   character strings specifying "Fortran-style" output formats -- as they
00125   would appear in a Harwell-Boeing file.  They are used to produce output
00126   which is as close as possible to what would be produced by Fortran code,
00127   but note that "D" and "P" edit descriptors are not supported.
00128   If NULL, the following defaults will be used:
00129                     Ptrfmt = Indfmt = "(8I10)"
00130                     Valfmt = Rhsfmt = "(4E20.13)"
00131 
00132   -----------------------
00133   CHARACTER I/O FUNCTIONS 
00134   -----------------------
00135   FUNCTION: 
00136 
00137   int readHB_mat_char(const char* filename, int colptr[], int rowind[], 
00138                                            char val[], char* Valfmt)
00139   int readHB_newmat_char(const char* filename, int* M, int* N, int* nonzeros, 
00140                           int** colptr, int** rowind, char** val, char** Valfmt)
00141 
00142   DESCRIPTION: 
00143 
00144   This function opens and reads the specified file, interpreting its
00145   contents as a sparse matrix stored in the Harwell/Boeing standard
00146   format.  (See readHB_aux_char to read auxillary vectors.)
00147               -- Values are interpreted as char strings.     --
00148   (Used to translate exact values from the file into a new storage format.)
00149 
00150   The "mat" function uses _pre-allocated_ arrays to hold the index and 
00151   nonzero value information.
00152 
00153   The "newmat" function allocates char arrays to hold the index 
00154   and nonzero value information, and returns pointers to these arrays 
00155   along with matrix dimension and number of nonzeros.
00156 
00157   FUNCTION:
00158 
00159   int readHB_aux_char(const char* filename, const char AuxType, char b[])
00160   int readHB_newaux_char(const char* filename, const char AuxType, char** b, 
00161                          char** Rhsfmt)
00162 
00163   DESCRIPTION:
00164 
00165   This function opens and reads from the specified file auxillary vector(s).
00166   The char argument Auxtype determines which type of auxillary vector(s)
00167   will be read (if present in the file).
00168 
00169                   AuxType = 'F'   right-hand-side 
00170                   AuxType = 'G'   initial estimate (Guess)
00171                   AuxType = 'X'   eXact solution
00172 
00173   If Nrhs > 1, all of the Nrhs vectors of the given type are read and 
00174   stored in column-major order in the vector b.
00175 
00176   The "newaux" function allocates a character array to hold the values 
00177                 retrieved.
00178   The "mat" function uses a _pre-allocated_ array to hold the values.
00179 
00180   FUNCTION:
00181 
00182   int writeHB_mat_char(const char* filename, int M, int N, 
00183                         int nz, const int colptr[], const int rowind[], 
00184                         const char val[], int Nrhs, const char rhs[], 
00185                         const char guess[], const char exact[], 
00186                         const char* Title, const char* Key, const char* Type, 
00187                         char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
00188                         const char* Rhstype)
00189 
00190   DESCRIPTION:
00191 
00192   The writeHB_mat_char function opens the named file and writes the specified
00193   matrix and optional auxillary vector(s) to that file in Harwell-Boeing
00194   format.  The format arguments (Ptrfmt,Indfmt,Valfmt, and Rhsfmt) are
00195   character strings specifying "Fortran-style" output formats -- as they
00196   would appear in a Harwell-Boeing file.  Valfmt and Rhsfmt must accurately 
00197   represent the character representation of the values stored in val[] 
00198   and rhs[]. 
00199 
00200   If NULL, the following defaults will be used for the integer vectors:
00201                     Ptrfmt = Indfmt = "(8I10)"
00202                     Valfmt = Rhsfmt = "(4E20.13)"
00203 
00204 
00205 */
00206 
00207 /*---------------------------------------------------------------------*/
00208 /* If zero-based indexing is desired, _SP_base should be set to 0      */
00209 /* This will cause indices read from H-B files to be decremented by 1  */
00210 /*             and indices written to H-B files to be incremented by 1 */
00211 /*            <<<  Standard usage is _SP_base = 1  >>>                 */
00212 #ifndef _SP_base
00213 #define _SP_base 1
00214 #endif
00215 /*---------------------------------------------------------------------*/
00216 
00217 #include "iohb.h"
00218 #include<stdio.h>
00219 #include<stdlib.h>
00220 #include<string.h>
00221 #include<math.h>
00222 #include<malloc.h>
00223 
00224 char* substr(const char* S, const int pos, const int len);
00225 void upcase(char* S);
00226 void IOHBTerminate(char* message);
00227 
00228 int readHB_info(const char* filename, int* M, int* N, int* nz, char** Type, 
00229                                                       int* Nrhs)
00230 {
00231 /****************************************************************************/
00232 /*  The readHB_info function opens and reads the header information from    */
00233 /*  the specified Harwell-Boeing file, and reports back the number of rows  */
00234 /*  and columns in the stored matrix (M and N), the number of nonzeros in   */
00235 /*  the matrix (nz), and the number of right-hand-sides stored along with   */
00236 /*  the matrix (Nrhs).                                                      */
00237 /*                                                                          */
00238 /*  For a description of the Harwell Boeing standard, see:                  */
00239 /*            Duff, et al.,  ACM TOMS Vol.15, No.1, March 1989              */
00240 /*                                                                          */
00241 /*    ----------                                                            */
00242 /*    **CAVEAT**                                                            */
00243 /*    ----------                                                            */
00244 /*  **  If the input file does not adhere to the H/B format, the  **        */
00245 /*  **             results will be unpredictable.                 **        */
00246 /*                                                                          */
00247 /****************************************************************************/
00248     FILE *in_file;
00249     int Ptrcrd, Indcrd, Valcrd, Rhscrd; 
00250     int Nrow, Ncol, Nnzero;
00251     char *mat_type;
00252     char Title[73], Key[9], Rhstype[4];
00253     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
00254 
00255     mat_type = (char *) malloc(4);
00256     if ( mat_type == NULL ) IOHBTerminate("Insufficient memory for mat_typen");
00257     
00258     if ( (in_file = fopen( filename, "r")) == NULL ) {
00259        fprintf(stderr,"Error: Cannot open file: %s\n",filename);
00260        return 0;
00261     }
00262 
00263     readHB_header(in_file, Title, Key, mat_type, &Nrow, &Ncol, &Nnzero, Nrhs,
00264                   Ptrfmt, Indfmt, Valfmt, Rhsfmt, 
00265                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
00266     fclose(in_file);
00267     *Type = mat_type;
00268     *(*Type+3) = (char) NULL;
00269     *M    = Nrow;
00270     *N    = Ncol;
00271     *nz   = Nnzero;
00272     if (Rhscrd == 0) {*Nrhs = 0;}
00273 
00274 /*  In verbose mode, print some of the header information:   */
00275 /*
00276     if (verbose == 1)
00277     {
00278         printf("Reading from Harwell-Boeing file %s (verbose on)...\n",filename);
00279         printf("  Title: %s\n",Title);
00280         printf("  Key:   %s\n",Key);
00281         printf("  The stored matrix is %i by %i with %i nonzeros.\n", 
00282                 *M, *N, *nz );
00283         printf("  %i right-hand--side(s) stored.\n",*Nrhs);
00284     }
00285 */
00286  
00287     return 1;
00288 
00289 }
00290 
00291 
00292 
00293 int readHB_header(FILE* in_file, char* Title, char* Key, char* Type, 
00294                     int* Nrow, int* Ncol, int* Nnzero, int* Nrhs,
00295                     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt, 
00296                     int* Ptrcrd, int* Indcrd, int* Valcrd, int* Rhscrd, 
00297                     char *Rhstype)
00298 {
00299 /*************************************************************************/
00300 /*  Read header information from the named H/B file...                   */
00301 /*************************************************************************/
00302     int Totcrd,Neltvl,Nrhsix;
00303     char line[BUFSIZ];
00304 
00305 /*  First line:   */
00306     fgets(line, BUFSIZ, in_file);
00307     if ( sscanf(line,"%*s") < 0 ) 
00308         IOHBTerminate("iohb.c: Null (or blank) first line of HB file.\n");
00309     (void) sscanf(line, "%72c%8[^\n]", Title, Key);
00310     *(Key+8) = (char) NULL;
00311     *(Title+72) = (char) NULL;
00312 
00313 /*  Second line:  */
00314     fgets(line, BUFSIZ, in_file);
00315     if ( sscanf(line,"%*s") < 0 ) 
00316         IOHBTerminate("iohb.c: Null (or blank) second line of HB file.\n");
00317     if ( sscanf(line,"%i",&Totcrd) != 1) Totcrd = 0;
00318     if ( sscanf(line,"%*i%i",Ptrcrd) != 1) *Ptrcrd = 0;
00319     if ( sscanf(line,"%*i%*i%i",Indcrd) != 1) *Indcrd = 0;
00320     if ( sscanf(line,"%*i%*i%*i%i",Valcrd) != 1) *Valcrd = 0;
00321     if ( sscanf(line,"%*i%*i%*i%*i%i",Rhscrd) != 1) *Rhscrd = 0;
00322 
00323 /*  Third line:   */
00324     fgets(line, BUFSIZ, in_file);
00325     if ( sscanf(line,"%*s") < 0 ) 
00326         IOHBTerminate("iohb.c: Null (or blank) third line of HB file.\n");
00327     if ( sscanf(line, "%3c", Type) != 1) 
00328         IOHBTerminate("iohb.c: Invalid Type info, line 3 of Harwell-Boeing file.\n");
00329     upcase(Type);
00330     if ( sscanf(line,"%*3c%i",Nrow) != 1) *Nrow = 0 ;
00331     if ( sscanf(line,"%*3c%*i%i",Ncol) != 1) *Ncol = 0 ;
00332     if ( sscanf(line,"%*3c%*i%*i%i",Nnzero) != 1) *Nnzero = 0 ;
00333     if ( sscanf(line,"%*3c%*i%*i%*i%i",&Neltvl) != 1) Neltvl = 0 ;
00334 
00335 /*  Fourth line:  */
00336     fgets(line, BUFSIZ, in_file);
00337     if ( sscanf(line,"%*s") < 0 ) 
00338         IOHBTerminate("iohb.c: Null (or blank) fourth line of HB file.\n");
00339     if ( sscanf(line, "%16c",Ptrfmt) != 1)
00340         IOHBTerminate("iohb.c: Invalid format info, line 4 of Harwell-Boeing file.\n"); 
00341     if ( sscanf(line, "%*16c%16c",Indfmt) != 1)
00342         IOHBTerminate("iohb.c: Invalid format info, line 4 of Harwell-Boeing file.\n"); 
00343     if ( sscanf(line, "%*16c%*16c%20c",Valfmt) != 1) 
00344         IOHBTerminate("iohb.c: Invalid format info, line 4 of Harwell-Boeing file.\n"); 
00345     sscanf(line, "%*16c%*16c%*20c%20c",Rhsfmt);
00346     *(Ptrfmt+16) = (char) NULL;
00347     *(Indfmt+16) = (char) NULL;
00348     *(Valfmt+20) = (char) NULL;
00349     *(Rhsfmt+20) = (char) NULL;
00350    
00351 /*  (Optional) Fifth line: */
00352     if (*Rhscrd != 0 )
00353     { 
00354        fgets(line, BUFSIZ, in_file);
00355        if ( sscanf(line,"%*s") < 0 ) 
00356            IOHBTerminate("iohb.c: Null (or blank) fifth line of HB file.\n");
00357        if ( sscanf(line, "%3c", Rhstype) != 1) 
00358          IOHBTerminate("iohb.c: Invalid RHS type information, line 5 of Harwell-Boeing file.\n");
00359        if ( sscanf(line, "%*3c%i", Nrhs) != 1) *Nrhs = 0;
00360        if ( sscanf(line, "%*3c%*i%i", &Nrhsix) != 1) Nrhsix = 0;
00361     }
00362     return 1;
00363 }
00364 
00365 
00366 int readHB_mat_double(const char* filename, int colptr[], int rowind[], 
00367                                                                  double val[])
00368 {
00369 /****************************************************************************/
00370 /*  This function opens and reads the specified file, interpreting its      */
00371 /*  contents as a sparse matrix stored in the Harwell/Boeing standard       */
00372 /*  format and creating compressed column storage scheme vectors to hold    */
00373 /*  the index and nonzero value information.                                */
00374 /*                                                                          */
00375 /*    ----------                                                            */
00376 /*    **CAVEAT**                                                            */
00377 /*    ----------                                                            */
00378 /*  Parsing real formats from Fortran is tricky, and this file reader       */
00379 /*  does not claim to be foolproof.   It has been tested for cases when     */
00380 /*  the real values are printed consistently and evenly spaced on each      */
00381 /*  line, with Fixed (F), and Exponential (E or D) formats.                 */
00382 /*                                                                          */
00383 /*  **  If the input file does not adhere to the H/B format, the  **        */
00384 /*  **             results will be unpredictable.                 **        */
00385 /*                                                                          */
00386 /****************************************************************************/
00387     FILE *in_file;
00388     int i,j,ind,col,offset,count,last,Nrhs;
00389     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
00390     int Nrow, Ncol, Nnzero, Nentries;
00391     int Ptrperline, Ptrwidth, Indperline, Indwidth;
00392     int Valperline, Valwidth, Valprec;
00393     int Valflag;           /* Indicates 'E','D', or 'F' float format */
00394     char* ThisElement;
00395     char Title[73], Key[8], Type[4] = "XXX\0", Rhstype[4];
00396     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
00397     char line[BUFSIZ];
00398 
00399     if ( (in_file = fopen( filename, "r")) == NULL ) {
00400        fprintf(stderr,"Error: Cannot open file: %s\n",filename);
00401        return 0;
00402     }
00403 
00404     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
00405                   Ptrfmt, Indfmt, Valfmt, Rhsfmt,
00406                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
00407 
00408 /*  Parse the array input formats from Line 3 of HB file  */
00409     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
00410     ParseIfmt(Indfmt,&Indperline,&Indwidth);
00411     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
00412     ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
00413     }
00414 
00415 /*  Read column pointer array:   */
00416 
00417     offset = 1-_SP_base;  /* if base 0 storage is declared (via macro definition), */
00418                           /* then storage entries are offset by 1                  */
00419 
00420     ThisElement = (char *) malloc(Ptrwidth+1);
00421     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
00422     *(ThisElement+Ptrwidth) = (char) NULL;
00423     count=0;
00424     for (i=0;i<Ptrcrd;i++)
00425     {
00426        fgets(line, BUFSIZ, in_file);
00427        if ( sscanf(line,"%*s") < 0 ) 
00428          IOHBTerminate("iohb.c: Null (or blank) line in pointer data region of HB file.\n");
00429        col =  0;
00430        for (ind = 0;ind<Ptrperline;ind++)
00431        {
00432           if (count > Ncol) break;
00433           strncpy(ThisElement,line+col,Ptrwidth);
00434   /* ThisElement = substr(line,col,Ptrwidth); */
00435           colptr[count] = atoi(ThisElement)-offset;
00436           count++; col += Ptrwidth;
00437        }
00438     }
00439     free(ThisElement);
00440 
00441 /*  Read row index array:  */
00442 
00443     ThisElement = (char *) malloc(Indwidth+1);
00444     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
00445     *(ThisElement+Indwidth) = (char) NULL;
00446     count = 0;
00447     for (i=0;i<Indcrd;i++)
00448     {
00449        fgets(line, BUFSIZ, in_file);
00450        if ( sscanf(line,"%*s") < 0 ) 
00451          IOHBTerminate("iohb.c: Null (or blank) line in index data region of HB file.\n");
00452        col =  0;
00453        for (ind = 0;ind<Indperline;ind++)
00454        {
00455           if (count == Nnzero) break;
00456           strncpy(ThisElement,line+col,Indwidth);
00457 /*        ThisElement = substr(line,col,Indwidth); */
00458           rowind[count] = atoi(ThisElement)-offset;
00459           count++; col += Indwidth;
00460        }
00461     }
00462     free(ThisElement);
00463 
00464 /*  Read array of values:  */
00465 
00466     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
00467 
00468        if ( Type[0] == 'C' ) Nentries = 2*Nnzero;
00469            else Nentries = Nnzero;
00470 
00471     ThisElement = (char *) malloc(Valwidth+1);
00472     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
00473     *(ThisElement+Valwidth) = (char) NULL;
00474     count = 0;
00475     for (i=0;i<Valcrd;i++)
00476     {
00477        fgets(line, BUFSIZ, in_file);
00478        if ( sscanf(line,"%*s") < 0 ) 
00479          IOHBTerminate("iohb.c: Null (or blank) line in value data region of HB file.\n");
00480        if (Valflag == 'D')  {
00481           while( strchr(line,'D') ) *strchr(line,'D') = 'E';
00482 /*           *strchr(Valfmt,'D') = 'E'; */
00483        }
00484        col =  0;
00485        for (ind = 0;ind<Valperline;ind++)
00486        {
00487           if (count == Nentries) break;
00488           strncpy(ThisElement,line+col,Valwidth);
00489           /*ThisElement = substr(line,col,Valwidth);*/
00490           if ( Valflag != 'F' && strchr(ThisElement,'E') == NULL ) { 
00491              /* insert a char prefix for exp */
00492              last = strlen(ThisElement);
00493              for (j=last+1;j>=0;j--) {
00494                 ThisElement[j] = ThisElement[j-1];
00495                 if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
00496                    ThisElement[j-1] = Valflag;                    
00497                    break;
00498                 }
00499              }
00500           }
00501           val[count] = atof(ThisElement);
00502           count++; col += Valwidth;
00503        }
00504     }
00505     free(ThisElement);
00506     }
00507 
00508     fclose(in_file);
00509     return 1;
00510 }
00511 
00512 int readHB_newmat_double(const char* filename, int* M, int* N, int* nonzeros, 
00513                          int** colptr, int** rowind, double** val)
00514 {
00515   int Nrhs;
00516         char *Type;
00517 
00518   readHB_info(filename, M, N, nonzeros, &Type, &Nrhs);
00519 
00520         *colptr = (int *)malloc((*N+1)*sizeof(int));
00521         if ( *colptr == NULL ) IOHBTerminate("Insufficient memory for colptr.\n");
00522         *rowind = (int *)malloc(*nonzeros*sizeof(int));
00523         if ( *rowind == NULL ) IOHBTerminate("Insufficient memory for rowind.\n");
00524         if ( Type[0] == 'C' ) {
00525 /*
00526    fprintf(stderr, "Warning: Reading complex data from HB file %s.\n",filename);
00527    fprintf(stderr, "         Real and imaginary parts will be interlaced in val[].\n");
00528 */
00529            /* Malloc enough space for real AND imaginary parts of val[] */
00530            *val = (double *)malloc(*nonzeros*sizeof(double)*2);
00531            if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
00532         } else {
00533            if ( Type[0] != 'P' ) {   
00534              /* Malloc enough space for real array val[] */
00535              *val = (double *)malloc(*nonzeros*sizeof(double));
00536              if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
00537            }
00538         }  /* No val[] space needed if pattern only */
00539   return readHB_mat_double(filename, *colptr, *rowind, *val);
00540 
00541 }
00542 
00543 int readHB_aux_double(const char* filename, const char AuxType, double b[])
00544 {
00545 /****************************************************************************/
00546 /*  This function opens and reads the specified file, placing auxillary     */
00547 /*  vector(s) of the given type (if available) in b.                        */
00548 /*  Return value is the number of vectors successfully read.                */
00549 /*                                                                          */
00550 /*                AuxType = 'F'   full right-hand-side vector(s)            */
00551 /*                AuxType = 'G'   initial Guess vector(s)                   */
00552 /*                AuxType = 'X'   eXact solution vector(s)                  */
00553 /*                                                                          */
00554 /*    ----------                                                            */
00555 /*    **CAVEAT**                                                            */
00556 /*    ----------                                                            */
00557 /*  Parsing real formats from Fortran is tricky, and this file reader       */
00558 /*  does not claim to be foolproof.   It has been tested for cases when     */
00559 /*  the real values are printed consistently and evenly spaced on each      */
00560 /*  line, with Fixed (F), and Exponential (E or D) formats.                 */
00561 /*                                                                          */
00562 /*  **  If the input file does not adhere to the H/B format, the  **        */
00563 /*  **             results will be unpredictable.                 **        */
00564 /*                                                                          */
00565 /****************************************************************************/
00566     FILE *in_file;
00567     int i,j,n,maxcol,start,stride,col,last,linel;
00568     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
00569     int Nrow, Ncol, Nnzero, Nentries;
00570     int Nrhs, nvecs, rhsi;
00571     int Rhsperline, Rhswidth, Rhsprec;
00572     int Rhsflag;
00573     char *ThisElement;
00574     char Title[73], Key[9], Type[4] = "XXX\0", Rhstype[4];
00575     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
00576     char line[BUFSIZ];
00577 
00578     if ((in_file = fopen( filename, "r")) == NULL) {
00579       fprintf(stderr,"Error: Cannot open file: %s\n",filename);
00580       return 0;
00581      }
00582 
00583     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
00584                   Ptrfmt, Indfmt, Valfmt, Rhsfmt,
00585                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
00586 
00587     if (Nrhs <= 0)
00588     {
00589       fprintf(stderr, "Warn: Attempt to read auxillary vector(s) when none are present.\n");
00590       return 0;
00591     }
00592     if (Rhstype[0] != 'F' )
00593     {
00594       fprintf(stderr,"Warn: Attempt to read auxillary vector(s) which are not stored in Full form.\n");
00595       fprintf(stderr,"       Rhs must be specified as full. \n");
00596       return 0;
00597     }
00598 
00599 /* If reading complex data, allow for interleaved real and imaginary values. */ 
00600     if ( Type[0] == 'C' ) {
00601        Nentries = 2*Nrow;
00602      } else {
00603        Nentries = Nrow;
00604     }
00605 
00606     nvecs = 1;
00607     
00608     if ( Rhstype[1] == 'G' ) nvecs++;
00609     if ( Rhstype[2] == 'X' ) nvecs++;
00610 
00611     if ( AuxType == 'G' && Rhstype[1] != 'G' ) {
00612       fprintf(stderr, "Warn: Attempt to read auxillary Guess vector(s) when none are present.\n");
00613       return 0;
00614     }
00615     if ( AuxType == 'X' && Rhstype[2] != 'X' ) {
00616       fprintf(stderr, "Warn: Attempt to read auxillary eXact solution vector(s) when none are present.\n");
00617       return 0;
00618     }
00619 
00620     ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec,&Rhsflag);
00621     maxcol = Rhsperline*Rhswidth;
00622 
00623 /*  Lines to skip before starting to read RHS values... */
00624     n = Ptrcrd + Indcrd + Valcrd;
00625 
00626     for (i = 0; i < n; i++)
00627       fgets(line, BUFSIZ, in_file);
00628 
00629 /*  start  - number of initial aux vector entries to skip   */
00630 /*           to reach first  vector requested               */
00631 /*  stride - number of aux vector entries to skip between   */
00632 /*           requested vectors                              */
00633     if ( AuxType == 'F' ) start = 0;
00634     else if ( AuxType == 'G' ) start = Nentries;
00635     else start = (nvecs-1)*Nentries;
00636     stride = (nvecs-1)*Nentries;
00637 
00638     fgets(line, BUFSIZ, in_file);
00639     linel= strchr(line,'\n')-line;
00640     col = 0;
00641 /*  Skip to initial offset */
00642 
00643     for (i=0;i<start;i++) {
00644        if ( col >=  ( maxcol<linel?maxcol:linel ) ) {
00645            fgets(line, BUFSIZ, in_file);
00646            linel= strchr(line,'\n')-line;
00647            col = 0;
00648        }
00649        col += Rhswidth;
00650     }
00651     if (Rhsflag == 'D')  {
00652        while( strchr(line,'D') ) *strchr(line,'D') = 'E';
00653     }
00654 
00655 /*  Read a vector of desired type, then skip to next */
00656 /*  repeating to fill Nrhs vectors                   */
00657 
00658   ThisElement = (char *) malloc(Rhswidth+1);
00659   if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
00660   *(ThisElement+Rhswidth) = (char) NULL;
00661   for (rhsi=0;rhsi<Nrhs;rhsi++) {
00662 
00663     for (i=0;i<Nentries;i++) {
00664        if ( col >= ( maxcol<linel?maxcol:linel ) ) {
00665            fgets(line, BUFSIZ, in_file);
00666            linel= strchr(line,'\n')-line;
00667            if (Rhsflag == 'D')  {
00668               while( strchr(line,'D') ) *strchr(line,'D') = 'E';
00669            }
00670            col = 0;
00671        }
00672        strncpy(ThisElement,line+col,Rhswidth);
00673        /*ThisElement = substr(line, col, Rhswidth);*/
00674           if ( Rhsflag != 'F' && strchr(ThisElement,'E') == NULL ) { 
00675              /* insert a char prefix for exp */
00676              last = strlen(ThisElement);
00677              for (j=last+1;j>=0;j--) {
00678                 ThisElement[j] = ThisElement[j-1];
00679                 if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
00680                    ThisElement[j-1] = Rhsflag;                    
00681                    break;
00682                 }
00683              }
00684           }
00685        b[i] = atof(ThisElement);
00686        col += Rhswidth;
00687     }
00688  
00689 /*  Skip any interleaved Guess/eXact vectors */
00690 
00691     for (i=0;i<stride;i++) {
00692        if ( col >= ( maxcol<linel?maxcol:linel ) ) {
00693            fgets(line, BUFSIZ, in_file);
00694            linel= strchr(line,'\n')-line;
00695            col = 0;
00696        }
00697        col += Rhswidth;
00698     }
00699 
00700   }
00701   free(ThisElement);
00702     
00703 
00704     fclose(in_file);
00705     return Nrhs;
00706 }
00707 
00708 int readHB_newaux_double(const char* filename, const char AuxType, double** b)
00709 {
00710         int Nrhs,M,N,nonzeros;
00711         char *Type;
00712 
00713   readHB_info(filename, &M, &N, &nonzeros, &Type, &Nrhs);
00714         if ( Nrhs <= 0 ) {
00715           fprintf(stderr,"Warn: Requested read of aux vector(s) when none are present.\n");
00716           return 0;
00717         } else { 
00718           if ( Type[0] == 'C' ) {
00719             fprintf(stderr, "Warning: Reading complex aux vector(s) from HB file %s.",filename);
00720             fprintf(stderr, "         Real and imaginary parts will be interlaced in b[].");
00721             *b = (double *)malloc(M*Nrhs*sizeof(double)*2);
00722             if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
00723             return readHB_aux_double(filename, AuxType, *b);
00724           } else {
00725             *b = (double *)malloc(M*Nrhs*sizeof(double));
00726             if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
00727       return readHB_aux_double(filename, AuxType, *b);
00728           }
00729         }
00730 }
00731 
00732 int writeHB_mat_double(const char* filename, int M, int N, 
00733                         int nz, const int colptr[], const int rowind[], 
00734                         const double val[], int Nrhs, const double rhs[], 
00735                         const double guess[], const double exact[],
00736                         const char* Title, const char* Key, const char* Type, 
00737                         char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
00738                         const char* Rhstype)
00739 {
00740 /****************************************************************************/
00741 /*  The writeHB function opens the named file and writes the specified      */
00742 /*  matrix and optional right-hand-side(s) to that file in Harwell-Boeing   */
00743 /*  format.                                                                 */
00744 /*                                                                          */
00745 /*  For a description of the Harwell Boeing standard, see:                  */
00746 /*            Duff, et al.,  ACM TOMS Vol.15, No.1, March 1989              */
00747 /*                                                                          */
00748 /****************************************************************************/
00749     FILE *out_file;
00750     int i,j,entry,offset,acount,linemod;
00751     int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
00752     int nvalentries, nrhsentries;
00753     int Ptrperline, Ptrwidth, Indperline, Indwidth;
00754     int Rhsperline, Rhswidth, Rhsprec;
00755     int Rhsflag;
00756     int Valperline, Valwidth, Valprec;
00757     int Valflag;           /* Indicates 'E','D', or 'F' float format */
00758     char pformat[16],iformat[16],vformat[19],rformat[19];
00759 
00760     if ( Type[0] == 'C' ) {
00761          nvalentries = 2*nz;
00762          nrhsentries = 2*M;
00763     } else {
00764          nvalentries = nz;
00765          nrhsentries = M;
00766     }
00767 
00768     if ( filename != NULL ) {
00769        if ( (out_file = fopen( filename, "w")) == NULL ) {
00770          fprintf(stderr,"Error: Cannot open file: %s\n",filename);
00771          return 0;
00772        }
00773     } else out_file = stdout;
00774 
00775     if ( Ptrfmt == NULL ) Ptrfmt = "(8I10)";
00776     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
00777     sprintf(pformat,"%%%dd",Ptrwidth);
00778     ptrcrd = (N+1)/Ptrperline;
00779     if ( (N+1)%Ptrperline != 0) ptrcrd++;
00780    
00781     if ( Indfmt == NULL ) Indfmt =  Ptrfmt;
00782     ParseIfmt(Indfmt,&Indperline,&Indwidth);
00783     sprintf(iformat,"%%%dd",Indwidth);
00784     indcrd = nz/Indperline;
00785     if ( nz%Indperline != 0) indcrd++;
00786 
00787     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
00788       if ( Valfmt == NULL ) Valfmt = "(4E20.13)";
00789       ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
00790       if (Valflag == 'D') *strchr(Valfmt,'D') = 'E';
00791       if (Valflag == 'F')
00792          sprintf(vformat,"%% %d.%df",Valwidth,Valprec);
00793       else
00794          sprintf(vformat,"%% %d.%dE",Valwidth,Valprec);
00795       valcrd = nvalentries/Valperline;
00796       if ( nvalentries%Valperline != 0) valcrd++;
00797     } else valcrd = 0;
00798 
00799     if ( Nrhs > 0 ) {
00800        if ( Rhsfmt == NULL ) Rhsfmt = Valfmt;
00801        ParseRfmt(Rhsfmt,&Rhsperline,&Rhswidth,&Rhsprec, &Rhsflag);
00802        if (Rhsflag == 'F')
00803           sprintf(rformat,"%% %d.%df",Rhswidth,Rhsprec);
00804        else
00805           sprintf(rformat,"%% %d.%dE",Rhswidth,Rhsprec);
00806        if (Rhsflag == 'D') *strchr(Rhsfmt,'D') = 'E';
00807        rhscrd = nrhsentries/Rhsperline; 
00808        if ( nrhsentries%Rhsperline != 0) rhscrd++;
00809        if ( Rhstype[1] == 'G' ) rhscrd+=rhscrd;
00810        if ( Rhstype[2] == 'X' ) rhscrd+=rhscrd;
00811        rhscrd*=Nrhs;
00812     } else rhscrd = 0;
00813 
00814     totcrd = 4+ptrcrd+indcrd+valcrd+rhscrd;
00815 
00816 
00817 /*  Print header information:  */
00818 
00819     fprintf(out_file,"%-72s%-8s\n%14d%14d%14d%14d%14d\n",Title, Key, totcrd,
00820             ptrcrd, indcrd, valcrd, rhscrd);
00821     fprintf(out_file,"%3s%11s%14d%14d%14d\n",Type,"          ", M, N, nz);
00822     fprintf(out_file,"%-16s%-16s%-20s", Ptrfmt, Indfmt, Valfmt);
00823     if ( Nrhs != 0 ) {
00824 /*     Print Rhsfmt on fourth line and                                    */
00825 /*           optional fifth header line for auxillary vector information: */
00826        fprintf(out_file,"%-20s\n%-14s%d\n",Rhsfmt,Rhstype,Nrhs);
00827     } else fprintf(out_file,"\n");
00828 
00829     offset = 1-_SP_base;  /* if base 0 storage is declared (via macro definition), */
00830                           /* then storage entries are offset by 1                  */
00831 
00832 /*  Print column pointers:   */
00833     for (i=0;i<N+1;i++)
00834     {
00835        entry = colptr[i]+offset;
00836        fprintf(out_file,pformat,entry);
00837        if ( (i+1)%Ptrperline == 0 ) fprintf(out_file,"\n");
00838     }
00839 
00840    if ( (N+1) % Ptrperline != 0 ) fprintf(out_file,"\n");
00841 
00842 /*  Print row indices:       */
00843     for (i=0;i<nz;i++)
00844     {
00845        entry = rowind[i]+offset;
00846        fprintf(out_file,iformat,entry);
00847        if ( (i+1)%Indperline == 0 ) fprintf(out_file,"\n");
00848     }
00849 
00850    if ( nz % Indperline != 0 ) fprintf(out_file,"\n");
00851 
00852 /*  Print values:            */
00853 
00854     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
00855 
00856     for (i=0;i<nvalentries;i++)
00857     {
00858        fprintf(out_file,vformat,val[i]);
00859        if ( (i+1)%Valperline == 0 ) fprintf(out_file,"\n");
00860     }
00861 
00862     if ( nvalentries % Valperline != 0 ) fprintf(out_file,"\n");
00863 
00864 /*  If available,  print right hand sides, 
00865            guess vectors and exact solution vectors:  */
00866     acount = 1;
00867     linemod = 0;
00868     if ( Nrhs > 0 ) {
00869        for (i=0;i<Nrhs;i++)
00870        {
00871           for ( j=0;j<nrhsentries;j++ ) {
00872             fprintf(out_file,rformat,rhs[j]);
00873             if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
00874           }
00875           if ( (acount-1)%Rhsperline != linemod ) {
00876             fprintf(out_file,"\n");
00877             linemod = (acount-1)%Rhsperline;
00878           }
00879           rhs += nrhsentries;
00880           if ( Rhstype[1] == 'G' ) {
00881             for ( j=0;j<nrhsentries;j++ ) {
00882               fprintf(out_file,rformat,guess[j]);
00883               if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
00884             }
00885             if ( (acount-1)%Rhsperline != linemod ) {
00886               fprintf(out_file,"\n");
00887               linemod = (acount-1)%Rhsperline;
00888             }
00889             guess += nrhsentries;
00890           }
00891           if ( Rhstype[2] == 'X' ) {
00892             for ( j=0;j<nrhsentries;j++ ) {
00893               fprintf(out_file,rformat,exact[j]);
00894               if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
00895             }
00896             if ( (acount-1)%Rhsperline != linemod ) {
00897               fprintf(out_file,"\n");
00898               linemod = (acount-1)%Rhsperline;
00899             }
00900             exact += nrhsentries;
00901           }
00902        }
00903     }
00904 
00905     }
00906 
00907     if ( fclose(out_file) != 0){
00908       fprintf(stderr,"Error closing file in writeHB_mat_double().\n");
00909       return 0;
00910     } else return 1;
00911     
00912 }
00913 
00914 int readHB_mat_char(const char* filename, int colptr[], int rowind[], 
00915                                            char val[], char* Valfmt)
00916 {
00917 /****************************************************************************/
00918 /*  This function opens and reads the specified file, interpreting its      */
00919 /*  contents as a sparse matrix stored in the Harwell/Boeing standard       */
00920 /*  format and creating compressed column storage scheme vectors to hold    */
00921 /*  the index and nonzero value information.                                */
00922 /*                                                                          */
00923 /*    ----------                                                            */
00924 /*    **CAVEAT**                                                            */
00925 /*    ----------                                                            */
00926 /*  Parsing real formats from Fortran is tricky, and this file reader       */
00927 /*  does not claim to be foolproof.   It has been tested for cases when     */
00928 /*  the real values are printed consistently and evenly spaced on each      */
00929 /*  line, with Fixed (F), and Exponential (E or D) formats.                 */
00930 /*                                                                          */
00931 /*  **  If the input file does not adhere to the H/B format, the  **        */
00932 /*  **             results will be unpredictable.                 **        */
00933 /*                                                                          */
00934 /****************************************************************************/
00935     FILE *in_file;
00936     int i,j,ind,col,offset,count,last;
00937     int Nrow,Ncol,Nnzero,Nentries,Nrhs;
00938     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
00939     int Ptrperline, Ptrwidth, Indperline, Indwidth;
00940     int Valperline, Valwidth, Valprec;
00941     int Valflag;           /* Indicates 'E','D', or 'F' float format */
00942     char* ThisElement;
00943     char line[BUFSIZ];
00944     char Title[73], Key[8], Type[4] = "XXX\0", Rhstype[4];
00945     char Ptrfmt[17], Indfmt[17], Rhsfmt[21];
00946 
00947     if ( (in_file = fopen( filename, "r")) == NULL ) {
00948        fprintf(stderr,"Error: Cannot open file: %s\n",filename);
00949        return 0;
00950     }
00951 
00952     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
00953                   Ptrfmt, Indfmt, Valfmt, Rhsfmt,
00954                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
00955 
00956 /*  Parse the array input formats from Line 3 of HB file  */
00957     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
00958     ParseIfmt(Indfmt,&Indperline,&Indwidth);
00959     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
00960        ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
00961        if (Valflag == 'D') {
00962           *strchr(Valfmt,'D') = 'E';
00963        }
00964     }
00965 
00966 /*  Read column pointer array:   */
00967 
00968     offset = 1-_SP_base;  /* if base 0 storage is declared (via macro definition), */
00969                           /* then storage entries are offset by 1                  */
00970 
00971     ThisElement = (char *) malloc(Ptrwidth+1);
00972     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
00973     *(ThisElement+Ptrwidth) = (char) NULL;
00974     count=0; 
00975     for (i=0;i<Ptrcrd;i++)
00976     {
00977        fgets(line, BUFSIZ, in_file);
00978        if ( sscanf(line,"%*s") < 0 ) 
00979          IOHBTerminate("iohb.c: Null (or blank) line in pointer data region of HB file.\n");
00980        col =  0;
00981        for (ind = 0;ind<Ptrperline;ind++)
00982        {
00983           if (count > Ncol) break;
00984           strncpy(ThisElement,line+col,Ptrwidth);
00985           /*ThisElement = substr(line,col,Ptrwidth);*/
00986           colptr[count] = atoi(ThisElement)-offset;
00987           count++; col += Ptrwidth;
00988        }
00989     }
00990     free(ThisElement);
00991 
00992 /*  Read row index array:  */
00993 
00994     ThisElement = (char *) malloc(Indwidth+1);
00995     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
00996     *(ThisElement+Indwidth) = (char) NULL;
00997     count = 0;
00998     for (i=0;i<Indcrd;i++)
00999     {
01000        fgets(line, BUFSIZ, in_file);
01001        if ( sscanf(line,"%*s") < 0 ) 
01002          IOHBTerminate("iohb.c: Null (or blank) line in index data region of HB file.\n");
01003        col =  0;
01004        for (ind = 0;ind<Indperline;ind++)
01005        {
01006           if (count == Nnzero) break;
01007           strncpy(ThisElement,line+col,Indwidth);
01008           /*ThisElement = substr(line,col,Indwidth);*/
01009           rowind[count] = atoi(ThisElement)-offset;
01010           count++; col += Indwidth;
01011        }
01012     }
01013     free(ThisElement);
01014 
01015 /*  Read array of values:  AS CHARACTERS*/
01016 
01017     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
01018 
01019        if ( Type[0] == 'C' ) Nentries = 2*Nnzero;
01020            else Nentries = Nnzero;
01021 
01022     ThisElement = (char *) malloc(Valwidth+1);
01023     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
01024     *(ThisElement+Valwidth) = (char) NULL;
01025     count = 0;
01026     for (i=0;i<Valcrd;i++)
01027     {
01028        fgets(line, BUFSIZ, in_file);
01029        if ( sscanf(line,"%*s") < 0 ) 
01030          IOHBTerminate("iohb.c: Null (or blank) line in value data region of HB file.\n");
01031        if (Valflag == 'D') {
01032           while( strchr(line,'D') ) *strchr(line,'D') = 'E';
01033        }
01034        col =  0;
01035        for (ind = 0;ind<Valperline;ind++)
01036        {
01037           if (count == Nentries) break;
01038           ThisElement = &val[count*Valwidth];
01039           strncpy(ThisElement,line+col,Valwidth);
01040           /*strncpy(ThisElement,substr(line,col,Valwidth),Valwidth);*/
01041           if ( Valflag != 'F' && strchr(ThisElement,'E') == NULL ) { 
01042              /* insert a char prefix for exp */
01043              last = strlen(ThisElement);
01044              for (j=last+1;j>=0;j--) {
01045                 ThisElement[j] = ThisElement[j-1];
01046                 if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
01047                    ThisElement[j-1] = Valflag;                    
01048                    break;
01049                 }
01050              }
01051           }
01052           count++; col += Valwidth;
01053        }
01054     }
01055     }
01056 
01057     return 1;
01058 }
01059 
01060 int readHB_newmat_char(const char* filename, int* M, int* N, int* nonzeros, int** colptr, 
01061                           int** rowind, char** val, char** Valfmt)
01062 {
01063     FILE *in_file;
01064     int Nrhs;
01065     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
01066     int Valperline, Valwidth, Valprec;
01067     int Valflag;           /* Indicates 'E','D', or 'F' float format */
01068     char Title[73], Key[9], Type[4] = "XXX\0", Rhstype[4];
01069     char Ptrfmt[17], Indfmt[17], Rhsfmt[21];
01070 
01071     if ((in_file = fopen( filename, "r")) == NULL) {
01072       fprintf(stderr,"Error: Cannot open file: %s\n",filename);
01073       return 0;
01074      }
01075     
01076     *Valfmt = (char *)malloc(21*sizeof(char));
01077     if ( *Valfmt == NULL ) IOHBTerminate("Insufficient memory for Valfmt.");
01078     readHB_header(in_file, Title, Key, Type, M, N, nonzeros, &Nrhs,
01079                   Ptrfmt, Indfmt, (*Valfmt), Rhsfmt,
01080                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
01081     fclose(in_file);
01082     ParseRfmt(*Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
01083 
01084         *colptr = (int *)malloc((*N+1)*sizeof(int));
01085         if ( *colptr == NULL ) IOHBTerminate("Insufficient memory for colptr.\n");
01086         *rowind = (int *)malloc(*nonzeros*sizeof(int));
01087         if ( *rowind == NULL ) IOHBTerminate("Insufficient memory for rowind.\n");
01088         if ( Type[0] == 'C' ) {
01089 /*
01090    fprintf(stderr, "Warning: Reading complex data from HB file %s.\n",filename);
01091    fprintf(stderr, "         Real and imaginary parts will be interlaced in val[].\n");
01092 */
01093            /* Malloc enough space for real AND imaginary parts of val[] */
01094            *val = (char *)malloc(*nonzeros*Valwidth*sizeof(char)*2);
01095            if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
01096         } else {
01097            if ( Type[0] != 'P' ) {   
01098              /* Malloc enough space for real array val[] */
01099              *val = (char *)malloc(*nonzeros*Valwidth*sizeof(char));
01100              if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
01101            }
01102         }  /* No val[] space needed if pattern only */
01103   return readHB_mat_char(filename, *colptr, *rowind, *val, *Valfmt);
01104 
01105 }
01106 
01107 int readHB_aux_char(const char* filename, const char AuxType, char b[])
01108 {
01109 /****************************************************************************/
01110 /*  This function opens and reads the specified file, placing auxilary      */
01111 /*  vector(s) of the given type (if available) in b :                       */
01112 /*  Return value is the number of vectors successfully read.                */
01113 /*                                                                          */
01114 /*                AuxType = 'F'   full right-hand-side vector(s)            */
01115 /*                AuxType = 'G'   initial Guess vector(s)                   */
01116 /*                AuxType = 'X'   eXact solution vector(s)                  */
01117 /*                                                                          */
01118 /*    ----------                                                            */
01119 /*    **CAVEAT**                                                            */
01120 /*    ----------                                                            */
01121 /*  Parsing real formats from Fortran is tricky, and this file reader       */
01122 /*  does not claim to be foolproof.   It has been tested for cases when     */
01123 /*  the real values are printed consistently and evenly spaced on each      */
01124 /*  line, with Fixed (F), and Exponential (E or D) formats.                 */
01125 /*                                                                          */
01126 /*  **  If the input file does not adhere to the H/B format, the  **        */
01127 /*  **             results will be unpredictable.                 **        */
01128 /*                                                                          */
01129 /****************************************************************************/
01130     FILE *in_file;
01131     int i,j,n,maxcol,start,stride,col,last,linel,nvecs,rhsi;
01132     int Nrow, Ncol, Nnzero, Nentries,Nrhs;
01133     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
01134     int Rhsperline, Rhswidth, Rhsprec;
01135     int Rhsflag;
01136     char Title[73], Key[9], Type[4] = "XXX\0", Rhstype[4];
01137     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
01138     char line[BUFSIZ];
01139     char *ThisElement;
01140 
01141     if ((in_file = fopen( filename, "r")) == NULL) {
01142       fprintf(stderr,"Error: Cannot open file: %s\n",filename);
01143       return 0;
01144      }
01145 
01146     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
01147                   Ptrfmt, Indfmt, Valfmt, Rhsfmt,
01148                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
01149 
01150     if (Nrhs <= 0)
01151     {
01152       fprintf(stderr, "Warn: Attempt to read auxillary vector(s) when none are present.\n");
01153       return 0;
01154     }
01155     if (Rhstype[0] != 'F' )
01156     {
01157       fprintf(stderr,"Warn: Attempt to read auxillary vector(s) which are not stored in Full form.\n");
01158       fprintf(stderr,"       Rhs must be specified as full. \n");
01159       return 0;
01160     }
01161 
01162 /* If reading complex data, allow for interleaved real and imaginary values. */ 
01163     if ( Type[0] == 'C' ) {
01164        Nentries = 2*Nrow;
01165      } else {
01166        Nentries = Nrow;
01167     }
01168 
01169     nvecs = 1;
01170     
01171     if ( Rhstype[1] == 'G' ) nvecs++;
01172     if ( Rhstype[2] == 'X' ) nvecs++;
01173 
01174     if ( AuxType == 'G' && Rhstype[1] != 'G' ) {
01175       fprintf(stderr, "Warn: Attempt to read auxillary Guess vector(s) when none are present.\n");
01176       return 0;
01177     }
01178     if ( AuxType == 'X' && Rhstype[2] != 'X' ) {
01179       fprintf(stderr, "Warn: Attempt to read auxillary eXact solution vector(s) when none are present.\n");
01180       return 0;
01181     }
01182 
01183     ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec,&Rhsflag);
01184     maxcol = Rhsperline*Rhswidth;
01185 
01186 /*  Lines to skip before starting to read RHS values... */
01187     n = Ptrcrd + Indcrd + Valcrd;
01188 
01189     for (i = 0; i < n; i++)
01190       fgets(line, BUFSIZ, in_file);
01191 
01192 /*  start  - number of initial aux vector entries to skip   */
01193 /*           to reach first  vector requested               */
01194 /*  stride - number of aux vector entries to skip between   */
01195 /*           requested vectors                              */
01196     if ( AuxType == 'F' ) start = 0;
01197     else if ( AuxType == 'G' ) start = Nentries;
01198     else start = (nvecs-1)*Nentries;
01199     stride = (nvecs-1)*Nentries;
01200 
01201     fgets(line, BUFSIZ, in_file);
01202     linel= strchr(line,'\n')-line;
01203     if ( sscanf(line,"%*s") < 0 ) 
01204        IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
01205     col = 0;
01206 /*  Skip to initial offset */
01207 
01208     for (i=0;i<start;i++) {
01209        col += Rhswidth;
01210        if ( col >= ( maxcol<linel?maxcol:linel ) ) {
01211            fgets(line, BUFSIZ, in_file);
01212            linel= strchr(line,'\n')-line;
01213        if ( sscanf(line,"%*s") < 0 ) 
01214        IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
01215            col = 0;
01216        }
01217     }
01218 
01219     if (Rhsflag == 'D')  {
01220       while( strchr(line,'D') ) *strchr(line,'D') = 'E';
01221     }
01222 /*  Read a vector of desired type, then skip to next */
01223 /*  repeating to fill Nrhs vectors                   */
01224 
01225   for (rhsi=0;rhsi<Nrhs;rhsi++) {
01226 
01227     for (i=0;i<Nentries;i++) {
01228        if ( col >= ( maxcol<linel?maxcol:linel ) ) {
01229            fgets(line, BUFSIZ, in_file);
01230            linel= strchr(line,'\n')-line;
01231        if ( sscanf(line,"%*s") < 0 ) 
01232        IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
01233            if (Rhsflag == 'D')  {
01234               while( strchr(line,'D') ) *strchr(line,'D') = 'E';
01235            }
01236            col = 0;
01237        }
01238        ThisElement = &b[i*Rhswidth]; 
01239        strncpy(ThisElement,line+col,Rhswidth);
01240           if ( Rhsflag != 'F' && strchr(ThisElement,'E') == NULL ) { 
01241              /* insert a char prefix for exp */
01242              last = strlen(ThisElement);
01243              for (j=last+1;j>=0;j--) {
01244                 ThisElement[j] = ThisElement[j-1];
01245                 if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
01246                    ThisElement[j-1] = Rhsflag;                    
01247                    break;
01248                 }
01249              }
01250           }
01251        col += Rhswidth;
01252     }
01253     b+=Nentries*Rhswidth;
01254  
01255 /*  Skip any interleaved Guess/eXact vectors */
01256 
01257     for (i=0;i<stride;i++) {
01258        col += Rhswidth;
01259        if ( col >= ( maxcol<linel?maxcol:linel ) ) {
01260            fgets(line, BUFSIZ, in_file);
01261            linel= strchr(line,'\n')-line;
01262        if ( sscanf(line,"%*s") < 0 ) 
01263        IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
01264            col = 0;
01265        }
01266     }
01267 
01268   }
01269     
01270 
01271     fclose(in_file);
01272     return Nrhs;
01273 }
01274 
01275 int readHB_newaux_char(const char* filename, const char AuxType, char** b, char** Rhsfmt)
01276 {
01277     FILE *in_file;
01278     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
01279     int Nrow,Ncol,Nnzero,Nrhs;
01280     int Rhsperline, Rhswidth, Rhsprec;
01281     int Rhsflag;
01282     char Title[73], Key[9], Type[4] = "XXX\0", Rhstype[4];
01283     char Ptrfmt[17], Indfmt[17], Valfmt[21];
01284 
01285     if ((in_file = fopen( filename, "r")) == NULL) {
01286       fprintf(stderr,"Error: Cannot open file: %s\n",filename);
01287       return 0;
01288      }
01289 
01290     *Rhsfmt = (char *)malloc(21*sizeof(char));
01291     if ( *Rhsfmt == NULL ) IOHBTerminate("Insufficient memory for Rhsfmt.");
01292     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
01293                   Ptrfmt, Indfmt, Valfmt, (*Rhsfmt),
01294                   &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
01295      fclose(in_file);
01296         if ( Nrhs == 0 ) {
01297           fprintf(stderr,"Warn: Requested read of aux vector(s) when none are present.\n");
01298           return 0;
01299         } else {
01300           ParseRfmt(*Rhsfmt,&Rhsperline,&Rhswidth,&Rhsprec,&Rhsflag);
01301           if ( Type[0] == 'C' ) {
01302             fprintf(stderr, "Warning: Reading complex aux vector(s) from HB file %s.",filename);
01303             fprintf(stderr, "         Real and imaginary parts will be interlaced in b[].");
01304             *b = (char *)malloc(Nrow*Nrhs*Rhswidth*sizeof(char)*2);
01305             if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
01306       return readHB_aux_char(filename, AuxType, *b);
01307           } else {
01308             *b = (char *)malloc(Nrow*Nrhs*Rhswidth*sizeof(char));
01309             if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
01310       return readHB_aux_char(filename, AuxType, *b);
01311           }
01312         } 
01313 }
01314 
01315 int writeHB_mat_char(const char* filename, int M, int N, 
01316                         int nz, const int colptr[], const int rowind[], 
01317                         const char val[], int Nrhs, const char rhs[], 
01318                         const char guess[], const char exact[], 
01319                         const char* Title, const char* Key, const char* Type, 
01320                         char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
01321                         const char* Rhstype)
01322 {
01323 /****************************************************************************/
01324 /*  The writeHB function opens the named file and writes the specified      */
01325 /*  matrix and optional right-hand-side(s) to that file in Harwell-Boeing   */
01326 /*  format.                                                                 */
01327 /*                                                                          */
01328 /*  For a description of the Harwell Boeing standard, see:                  */
01329 /*            Duff, et al.,  ACM TOMS Vol.15, No.1, March 1989              */
01330 /*                                                                          */
01331 /****************************************************************************/
01332     FILE *out_file;
01333     int i,j,acount,linemod,entry,offset;
01334     int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
01335     int nvalentries, nrhsentries;
01336     int Ptrperline, Ptrwidth, Indperline, Indwidth;
01337     int Rhsperline, Rhswidth, Rhsprec;
01338     int Rhsflag;
01339     int Valperline, Valwidth, Valprec;
01340     int Valflag;           /* Indicates 'E','D', or 'F' float format */
01341     char pformat[16],iformat[16],vformat[19],rformat[19];
01342 
01343     if ( Type[0] == 'C' ) {
01344          nvalentries = 2*nz;
01345          nrhsentries = 2*M;
01346     } else {
01347          nvalentries = nz;
01348          nrhsentries = M;
01349     }
01350 
01351     if ( filename != NULL ) {
01352        if ( (out_file = fopen( filename, "w")) == NULL ) {
01353          fprintf(stderr,"Error: Cannot open file: %s\n",filename);
01354          return 0;
01355        }
01356     } else out_file = stdout;
01357 
01358     if ( Ptrfmt == NULL ) Ptrfmt = "(8I10)";
01359     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
01360     sprintf(pformat,"%%%dd",Ptrwidth);
01361    
01362     if ( Indfmt == NULL ) Indfmt =  Ptrfmt;
01363     ParseIfmt(Indfmt,&Indperline,&Indwidth);
01364     sprintf(iformat,"%%%dd",Indwidth);
01365 
01366     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
01367       if ( Valfmt == NULL ) Valfmt = "(4E20.13)";
01368       ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
01369       sprintf(vformat,"%%%ds",Valwidth);
01370     }
01371 
01372     ptrcrd = (N+1)/Ptrperline;
01373     if ( (N+1)%Ptrperline != 0) ptrcrd++;
01374 
01375     indcrd = nz/Indperline;
01376     if ( nz%Indperline != 0) indcrd++;
01377 
01378     valcrd = nvalentries/Valperline;
01379     if ( nvalentries%Valperline != 0) valcrd++;
01380 
01381     if ( Nrhs > 0 ) {
01382        if ( Rhsfmt == NULL ) Rhsfmt = Valfmt;
01383        ParseRfmt(Rhsfmt,&Rhsperline,&Rhswidth,&Rhsprec, &Rhsflag);
01384        sprintf(rformat,"%%%ds",Rhswidth);
01385        rhscrd = nrhsentries/Rhsperline; 
01386        if ( nrhsentries%Rhsperline != 0) rhscrd++;
01387        if ( Rhstype[1] == 'G' ) rhscrd+=rhscrd;
01388        if ( Rhstype[2] == 'X' ) rhscrd+=rhscrd;
01389        rhscrd*=Nrhs;
01390     } else rhscrd = 0;
01391 
01392     totcrd = 4+ptrcrd+indcrd+valcrd+rhscrd;
01393 
01394 
01395 /*  Print header information:  */
01396 
01397     fprintf(out_file,"%-72s%-8s\n%14d%14d%14d%14d%14d\n",Title, Key, totcrd,
01398             ptrcrd, indcrd, valcrd, rhscrd);
01399     fprintf(out_file,"%3s%11s%14d%14d%14d\n",Type,"          ", M, N, nz);
01400     fprintf(out_file,"%-16s%-16s%-20s", Ptrfmt, Indfmt, Valfmt);
01401     if ( Nrhs != 0 ) {
01402 /*     Print Rhsfmt on fourth line and                                    */
01403 /*           optional fifth header line for auxillary vector information: */
01404        fprintf(out_file,"%-20s\n%-14s%d\n",Rhsfmt,Rhstype,Nrhs);
01405     } else fprintf(out_file,"\n");
01406 
01407     offset = 1-_SP_base;  /* if base 0 storage is declared (via macro definition), */
01408                           /* then storage entries are offset by 1                  */
01409 
01410 /*  Print column pointers:   */
01411     for (i=0;i<N+1;i++)
01412     {
01413        entry = colptr[i]+offset;
01414        fprintf(out_file,pformat,entry);
01415        if ( (i+1)%Ptrperline == 0 ) fprintf(out_file,"\n");
01416     }
01417 
01418    if ( (N+1) % Ptrperline != 0 ) fprintf(out_file,"\n");
01419 
01420 /*  Print row indices:       */
01421     for (i=0;i<nz;i++)
01422     {
01423        entry = rowind[i]+offset;
01424        fprintf(out_file,iformat,entry);
01425        if ( (i+1)%Indperline == 0 ) fprintf(out_file,"\n");
01426     }
01427 
01428    if ( nz % Indperline != 0 ) fprintf(out_file,"\n");
01429 
01430 /*  Print values:            */
01431 
01432     if ( Type[0] != 'P' ) {          /* Skip if pattern only  */
01433     for (i=0;i<nvalentries;i++)
01434     {
01435        fprintf(out_file,vformat,val+i*Valwidth);
01436        if ( (i+1)%Valperline == 0 ) fprintf(out_file,"\n");
01437     }
01438 
01439     if ( nvalentries % Valperline != 0 ) fprintf(out_file,"\n");
01440 
01441 /*  Print right hand sides:  */
01442     acount = 1;
01443     linemod=0;
01444     if ( Nrhs > 0 ) {
01445       for (j=0;j<Nrhs;j++) {
01446        for (i=0;i<nrhsentries;i++)
01447        {
01448           fprintf(out_file,rformat,rhs+i*Rhswidth);
01449           if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
01450        }
01451        if ( acount%Rhsperline != linemod ) {
01452           fprintf(out_file,"\n");
01453           linemod = (acount-1)%Rhsperline;
01454        }
01455        if ( Rhstype[1] == 'G' ) {
01456          for (i=0;i<nrhsentries;i++)
01457          {
01458            fprintf(out_file,rformat,guess+i*Rhswidth);
01459            if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
01460          }
01461          if ( acount%Rhsperline != linemod ) {
01462             fprintf(out_file,"\n");
01463             linemod = (acount-1)%Rhsperline;
01464          }
01465        }
01466        if ( Rhstype[2] == 'X' ) {
01467          for (i=0;i<nrhsentries;i++)
01468          {
01469            fprintf(out_file,rformat,exact+i*Rhswidth);
01470            if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
01471          }
01472          if ( acount%Rhsperline != linemod ) {
01473             fprintf(out_file,"\n");
01474             linemod = (acount-1)%Rhsperline;
01475          }
01476        }
01477       }
01478     }
01479 
01480     }
01481 
01482     if ( fclose(out_file) != 0){
01483       fprintf(stderr,"Error closing file in writeHB_mat_char().\n");
01484       return 0;
01485     } else return 1;
01486     
01487 }
01488 
01489 int ParseIfmt(char* fmt, int* perline, int* width)
01490 {
01491 /*************************************************/
01492 /*  Parse an *integer* format field to determine */
01493 /*  width and number of elements per line.       */
01494 /*************************************************/
01495     char *tmp;
01496     if (fmt == NULL ) {
01497       *perline = 0; *width = 0; return 0;
01498     }
01499     upcase(fmt);
01500     tmp = strchr(fmt,'(');
01501     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,'I') - tmp - 1);
01502     *perline = atoi(tmp);
01503     tmp = strchr(fmt,'I');
01504     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,')') - tmp - 1);
01505     return *width = atoi(tmp);
01506 }
01507 
01508 int ParseRfmt(char* fmt, int* perline, int* width, int* prec, int* flag)
01509 {
01510 /*************************************************/
01511 /*  Parse a *real* format field to determine     */
01512 /*  width and number of elements per line.       */
01513 /*  Also sets flag indicating 'E' 'F' 'P' or 'D' */
01514 /*  format.                                      */
01515 /*************************************************/
01516     char* tmp;
01517     char* tmp2;
01518     char* tmp3;
01519     int len;
01520 
01521     if (fmt == NULL ) {
01522       *perline = 0; 
01523       *width = 0; 
01524       flag = NULL;  
01525       return 0;
01526     }
01527 
01528     upcase(fmt);
01529     if (strchr(fmt,'(') != NULL)  fmt = strchr(fmt,'(');
01530     if (strchr(fmt,')') != NULL)  {
01531        tmp2 = strchr(fmt,')');
01532        while ( strchr(tmp2+1,')') != NULL ) {
01533           tmp2 = strchr(tmp2+1,')');
01534        }
01535        *(tmp2+1) = (int) NULL;
01536     }
01537     if (strchr(fmt,'P') != NULL)  /* Remove any scaling factor, which */
01538     {                             /* affects output only, not input */
01539       if (strchr(fmt,'(') != NULL) {
01540         tmp = strchr(fmt,'P');
01541         if ( *(++tmp) == ',' ) tmp++;
01542         tmp3 = strchr(fmt,'(')+1;
01543         len = tmp-tmp3;
01544         tmp2 = tmp3;
01545         while ( *(tmp2+len) != (int) NULL ) {
01546            *tmp2=*(tmp2+len);
01547            tmp2++; 
01548         }
01549         *(strchr(fmt,')')+1) = (int) NULL;
01550       }
01551     }
01552     if (strchr(fmt,'E') != NULL) { 
01553        *flag = 'E';
01554     } else if (strchr(fmt,'D') != NULL) { 
01555        *flag = 'D';
01556     } else if (strchr(fmt,'F') != NULL) { 
01557        *flag = 'F';
01558     } else {
01559       fprintf(stderr,"Real format %s in H/B file not supported.\n",fmt);
01560       return 0;
01561     }
01562     tmp = strchr(fmt,'(');
01563     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,*flag) - tmp - 1);
01564     *perline = atoi(tmp);
01565     tmp = strchr(fmt,*flag);
01566     if ( strchr(fmt,'.') ) {
01567       *prec = atoi( substr( fmt, strchr(fmt,'.') - fmt + 1, strchr(fmt,')') - strchr(fmt,'.')-1) );
01568       tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,'.') - tmp - 1);
01569     } else {
01570       tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,')') - tmp - 1);
01571     }
01572     return *width = atoi(tmp);
01573 }
01574 
01575 char* substr(const char* S, const int pos, const int len)
01576 {
01577     int i;
01578     char *SubS;
01579     if ( pos+len <= strlen(S)) {
01580     SubS = (char *)malloc(len+1);
01581     if ( SubS == NULL ) IOHBTerminate("Insufficient memory for SubS.");
01582     for (i=0;i<len;i++) SubS[i] = S[pos+i];
01583     SubS[len] = (char) NULL;
01584     } else {
01585       SubS = NULL;
01586     }
01587     return SubS;
01588 }
01589 
01590 #include<ctype.h>
01591 void upcase(char* S)
01592 {
01593 /*  Convert S to uppercase     */
01594     int i,len;
01595     len = strlen(S);
01596     for (i=0;i< len;i++)
01597        S[i] = toupper(S[i]);
01598 }
01599 
01600 void IOHBTerminate(char* message) 
01601 {
01602    fprintf(stderr,message);
01603    exit(1);
01604 }
01605 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Wed Apr 13 10:05:35 2011 for Ifpack Package Browser (Single Doxygen Collection) by  doxygen 1.6.3