iohb.c

http://gams.cam.nist.gov/acmd/Staff/RPozo/sparselib++.html

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


    FILE *out_file;

    int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
    char *ptrfmt, *indfmt, *valfmt;
    char *Type;
    char *Titlefill, *Keyfill;
    int filllen, i;
    int entry, finfo;

    out_file = fopen( filename, "w");

    ptrcrd = (N+1)/8;
    if ( (N+1)%8 != 0) ptrcrd++;

    indcrd = nz/8;
    if ( nz%8 != 0) indcrd++;

    valcrd = nz/4;
    if ( nz%4 != 0) valcrd++;

    rhscrd = nrhs*M/4; 
    if ( nrhs*M%4 != 0) rhscrd++;

    totcrd = 4+ptrcrd+indcrd+valcrd+rhscrd;

    ptrfmt = "(8I10)          ";
    indfmt = ptrfmt;
    valfmt = "(4E20.16)       ";

    Type = "RUA";

/*  Print header information:  */

    fprintf(out_file,"%-72s%-8s\n%14d%14d%14d%14d%14d\n",Title, Key, totcrd,
            ptrcrd, indcrd, valcrd, rhscrd);
    fprintf(out_file,"%3s%11s%14d%14d%14d\n",Type,"          ", M, N, nz);
    fprintf(out_file,"%16s%16s%16s\n", ptrfmt, indfmt, valfmt);

/*  Print column pointers:   */
    for (i=0;i<N+1;i++)
    {
       entry = colptr[i]+1;
       fprintf(out_file,"%10d",entry);
       if ( (i+1)%8 == 0 ) fprintf(out_file,"\n");
    }

   if ( (N+1) % 8 != 0 ) fprintf(out_file,"\n");

/*  Print row indices:       */
    for (i=0;i<nz;i++)
    {
       entry = rowind[i]+1;
       fprintf(out_file,"%10d",entry);
       if ( (i+1)%8 == 0 ) fprintf(out_file,"\n");
    }

   if ( nz % 8 != 0 ) fprintf(out_file,"\n");

/*  Print values:            */
    for (i=0;i<nz;i++)
    {
       fprintf(out_file,"% 20.12E",val[i]);
       if ( (i+1)%4 == 0 ) fprintf(out_file,"\n");
    }

    if ( nz % 4 != 0 ) fprintf(out_file,"\n");

/*  Print right hand sides:  */
    if ( nrhs > 0 ) {
       for (i=0;i<nrhs*M;i++)
       {
          fprintf(out_file,"% 20.12E",rhs[i]);
          if ( (i+1)%4 == 0 ) fprintf(out_file,"\n");
       }
    }

    finfo = fclose(out_file);
    if (finfo != 0) printf("Error closing file in writeHB_mat_float().\n");

    return;
}

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

              /********************************************/
              /*  Utility functions supporting readHB:    */
              /********************************************/

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

/**************************************************/
/*  Check to ensure header information is present */
/**************************************************/
void nullchk(line_ptr)
char *line_ptr;
{
    if (line_ptr == NULL)
    {
          printf("Cannot complete reading file information.\n ");
          exit(1);
    }
}

/*************************************************/
/*  Parse an *integer* format field to determine */
/*  width and number of elements per line.       */
/*************************************************/
void ParseIfmt(fmt, perline, width)
char *fmt; int *perline; int *width;
{
    char *PPL, *WID;
    upcase(fmt);
    PPL = substr_before(fmt,'I');
    PPL = substr(PPL,1,strlen(PPL)-1);
    *perline = atoi(PPL);
    WID = substr_after(fmt,'I');
    WID = substr(WID,0,strlen(WID)-1);
    *width = atoi(WID);
}

/*************************************************/
/*  Parse an *real* format field to determine    */
/*  width and number of elements per line.       */
/*  Also sets flag indicating 'E' or 'D' format. */
/*************************************************/
void ParseRfmt(fmt, perline, width, flag)
char *fmt; int *perline; int *width; int *flag;
{
    int foundE, foundD, foundF, foundP;
    char *PPL, *WID;

    if (fmt == NULL )
    {
      *perline = 0; *width = 0; *flag = '\0'; return;
    }

    upcase(fmt);
    foundP = my_index(fmt,'P');
    foundE = my_index(fmt,'E');
    foundD = my_index(fmt,'D');
    foundF = my_index(fmt,'F');      /* Fixed format */
    if (foundP != -1 )            /* Remove any scaling factor, which */
    {                             /* affects output only, not input */
      fmt = fmt + foundP + 1;
    }
    if (foundE != -1 )
    {
      *flag = 'E';
      PPL = substr_before(fmt,'E');
      PPL = substr(PPL,1,strlen(PPL)-1);
      *perline = atoi(PPL);
      WID = substr_after(fmt,'E');
      WID = substr_through(WID,'.');
      WID = substr(WID,0,strlen(WID)-1);
      *width = atoi(WID);
    }
    else if (foundD != -1)
    {
      *flag = 'D';
      PPL = substr_before(fmt,'D');
      PPL = substr(PPL,1,strlen(PPL)-1);
      *perline = atoi(PPL);
      WID = substr_after(fmt,'D');
      WID = substr_through(WID,'.');
      WID = substr(WID,0,strlen(WID)-1);
      *width = atoi(WID);
    }
    else if (foundF != -1)
    {
      *flag = 'F';
      PPL = substr_before(fmt,'F');
      PPL = substr(PPL,1,strlen(PPL)-1);
      *perline = atoi(PPL);
      WID = substr_after(fmt,'F');
      WID = substr_through(WID,'.');
      WID = substr(WID,0,strlen(WID)-1);
      *width = atoi(WID);
    }
    else
    {
      printf("Real format in H/B file not supported.\n");
      exit(1);
    }

}


/*****************************************************/
/*  Converts a line with real data in 'D' format     */
/*  to 'E' format by simply changing the character   */
/*  'D' to 'E' by adding 1 to it.                    */
/*****************************************************/
void convertDtoE(line)
char *line;
{
    int len, i;
    len = strlen(line);
    for (i=0;i<len;i++)
       if ( line[i] == 'D' || line[i] == 'd' ) line[i] = line[i]+1;
    return;
}


/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
/*       small string processing library for readhb()            */
/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/

/*****************************************************************/
/*  This file contains some minimal string functions to support  */
/*  support parsing the data formats in a Harwell-Boeing file.   */
/*****************************************************************/

char* substr(S, pos, len)
char *S;
int pos, len;
{
    int i;
    char *SubS;
    SubS = malloc(len+1);
    if ( (pos+len) > strlen(S))
    {
       printf("Error: Substring (%s, %d, %d) will read beyond string boundary.\n",
               S, pos, len);
       exit(1);
    }
    for (i=0;i<len;i++)
       SubS[i] = S[pos+i];
    SubS[len] = '\0';
    return SubS;
}

char* substr_after(S, M)
char *S; char M;
{
/*  Return the substring of S from the char after the first   */
/*  occurrence of M to the end of the string                  */
    int i, pos, SubSlen;
    char * SubS;
    pos = 0;
    for (i=0;i< strlen(S);i++)
    {
       if ( S[i] == M )
       {
          pos = i+1;
          break;
       }
    }
    if (pos != 0)
    {
       SubSlen = strlen(S) - pos;
       SubS = malloc(SubSlen);
       for (i=0;i< SubSlen;i++)
          SubS[i] = S[pos+i];
       SubS[SubSlen] = '\0';
       return SubS;
     }
     else
     {
       printf("Character %s not found in input string.\n", M);
       exit(1);
     }

       return NULL;       
 }



char* substr_before(S, M)
char *S; char M;
{
/*  Return the substring of S from the first char of S to the   */
/*  char before the first occurrence of M                       */
    int i, pos, SubSlen;
    char* SubS;
    pos = 0;
    for (i=0;i< strlen(S);i++)
    {
       if ( S[i] == M )
       {
          pos = i-1;
          break;
       }
    }
    if (pos != 0)
    {
       SubSlen = pos+1;
       SubS = malloc(SubSlen);
       for (i=0;i< SubSlen;i++)
          SubS[i] = S[i];
       SubS[SubSlen] = '\0';
       return SubS;
     }
     else
     {
       printf("Character %s not found in input string.\n", M);
       exit(1);
     }
     return NULL;   
 }


char* substr_through(S, M)
char *S; char M;
{
/*  Similer to substr_before, but include M         */
    int i, pos, SubSlen;
    char *SubS;
    pos = 0;
    for (i=0;i< strlen(S);i++)
    {
       if ( S[i] == M )
       {
          pos = i;
          break;
       }
    }
    if (pos != 0)
    {
       SubSlen = pos+1;
       SubS = malloc(SubSlen);
       for (i=0;i< SubSlen;i++)
          SubS[i] = S[i];
       SubS[SubSlen] = '\0';
       return SubS;
     }
     else
     {
       printf("Character %s not found in input string.\n", M);
       exit(1);
     }
       return NULL;
 }


void upcase(S)
char *S;
{
/*  Convert S to uppercase     */
    int i;
    for (i=0;i< strlen(S);i++)
       if ( S[i] >= 'a' && S[i] <= 'z' ) S[i] = S[i] - 32;
}


int my_index(S, M)
char *S; char M;
{
/*  Return the index of the first occurrence of M in S  */
/*  Return -1 if M does not occur in S                  */
    int i, pos;
    pos = -1;
    for (i=0;i< strlen(S);i++)
    {
       if ( S[i] == M )
       {
          pos = i;
          break;
       }
    }
    return pos;
}


void readHB_newmat_double(filename, M, N, nonzeros, colptr, rowind, val)
char *filename; int *M; int *N; int *nonzeros; int **colptr; 
int **rowind; double **val;
{
    int nrhs;

    readHB_info(filename, M, N, nonzeros, &nrhs);
   *colptr = (int *)malloc((*N+1)*sizeof(int));
   *rowind = (int *)malloc(*nonzeros*sizeof(int));
   *val = (double *)malloc(*nonzeros*sizeof(double));

    readHB_mat_double(filename, *colptr, *rowind, *val);

}

void readHB_newmat_float(filename, M, N, nonzeros, colptr, rowind, val)
char *filename; int *M; int *N; int *nonzeros; int **colptr; 
int **rowind; float **val;
{
    int nrhs;

    readHB_info(filename, M, N, nonzeros, &nrhs);
    *colptr = (int *)malloc((*N+1)*sizeof(int));
    *rowind = (int *)malloc(*nonzeros*sizeof(int));
    *val = (float *)malloc(*nonzeros*sizeof(float));

    readHB_mat_float(filename, *colptr, *rowind, *val);

}