Actual source code: ex50.c

  1: /*$Id: ex50.c,v 1.26 2001/04/10 19:35:44 bsmith Exp $*/

  3: static char help[] = "Reads in a matrix and vector in ASCII format. Writesn
  4: them using the PETSc sparse format. Input parameters are:n
  5:   -fin <filename> : input filen
  6:   -fout <filename> : output filenn";

  8: #include "petscmat.h"

 10: int main(int argc,char **args)
 11: {
 12:   Mat         A;
 13:   Vec         b;
 14:   char        filein[256],finname[256],fileout[256];
 15:   int         n,ierr,col,row;
 16:   int         rowin;
 17:   PetscTruth  flg;
 18:   Scalar      val,*array;
 19:   FILE*       file;
 20:   PetscViewer view;

 22:   PetscInitialize(&argc,&args,(char *)0,help);

 24:   /* Read in matrix and RHS */
 25:   PetscOptionsGetString(PETSC_NULL,"-fin",filein,255,&flg);
 26:   if (!flg) SETERRQ(1,"Must indicate file for reading");
 27:   PetscOptionsGetString(PETSC_NULL,"-fout",fileout,255,&flg);
 28:   if (!flg) SETERRQ(1,"Must indicate file for writing");

 30:   PetscFixFilename(filein,finname);
 31:   if (!(file = fopen(finname,"r"))) {
 32:     SETERRQ(1,"cannot open input filen");
 33:   }
 34:   fscanf(file,"%dn",&n);

 36:   MatCreate(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,n,n,&A);
 37:   MatSetFromOptions(A);
 38:   VecCreate(PETSC_COMM_WORLD,PETSC_DECIDE,n,&b);
 39:   VecSetFromOptions(b);

 41:   for (row=0; row<n; row++) {
 42:     fscanf(file,"row %d:",&rowin);
 43:     if (rowin != row) SETERRQ(1,"Bad file");
 44:     while (fscanf(file," %d %le",&col,&val)) {
 45:       MatSetValues(A,1,&row,1,&col,&val,INSERT_VALUES);
 46:     }
 47:   }
 48:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 49:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 50:   VecGetArray(b,&array);
 51:   for (row=0; row<n; row++) {
 52:     fscanf(file," ii= %d %le",&col,array+row);
 53:   }
 54:   VecRestoreArray(b,&array);

 56:   fclose(file);

 58:   PetscPrintf(PETSC_COMM_SELF,"Reading matrix complete.n");
 59:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,fileout,PETSC_BINARY_CREATE,&view);
 60:   MatView(A,view);
 61:   VecView(b,view);
 62:   PetscViewerDestroy(view);

 64:   VecDestroy(b);
 65:   MatDestroy(A);

 67:   PetscFinalize();
 68:   return 0;
 69: }