xref: /petsc/src/mat/tests/ex173.c (revision 28b400f66ebc7ae0049166a2294dfcd3df27e64b)
1a5b23f4aSJose E. Roman static char help[] = "Test MatrixMarket outputting.\n\n";
2c4762a1bSJed Brown 
3c4762a1bSJed Brown /*
4c4762a1bSJed Brown   Include "petscmat.h" so that we can use matrices.
5c4762a1bSJed Brown   automatically includes:
6c4762a1bSJed Brown      petscsys.h    - base PETSc routines   petscvec.h    - vectors
7c4762a1bSJed Brown      petscmat.h    - matrices
8c4762a1bSJed Brown      petscis.h     - index sets            petscviewer.h - viewers
9c4762a1bSJed Brown */
10c4762a1bSJed Brown 
11c4762a1bSJed Brown #include <petscmat.h>
12c4762a1bSJed Brown 
13c4762a1bSJed Brown int main(int argc,char **args)
14c4762a1bSJed Brown {
15c4762a1bSJed Brown   Mat            A;
16c4762a1bSJed Brown   PetscViewer    fd;                        /* viewer */
17c4762a1bSJed Brown   char           file[PETSC_MAX_PATH_LEN];  /* input file name */
18c4762a1bSJed Brown   PetscErrorCode ierr;
19c4762a1bSJed Brown   PetscBool      flg;
20c4762a1bSJed Brown 
21c4762a1bSJed Brown   ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
225f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsGetString(NULL,NULL,"-f0",file,sizeof(file),&flg));
23*28b400f6SJacob Faibussowitsch   PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f0 option");
245f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd));
255f80ce2aSJacob Faibussowitsch   CHKERRQ(MatCreate(PETSC_COMM_WORLD,&A));
265f80ce2aSJacob Faibussowitsch   CHKERRQ(MatLoad(A,fd));
275f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscViewerDestroy(&fd));
28c4762a1bSJed Brown 
295f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_MATRIXMARKET));
305f80ce2aSJacob Faibussowitsch   CHKERRQ(MatView(A,PETSC_VIEWER_STDOUT_WORLD));
315f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
325f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&A));
33c4762a1bSJed Brown   ierr = PetscFinalize();
34c4762a1bSJed Brown   return ierr;
35c4762a1bSJed Brown }
36c4762a1bSJed Brown 
37c4762a1bSJed Brown /*TEST
38c4762a1bSJed Brown 
39c4762a1bSJed Brown    test:
40c4762a1bSJed Brown       args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/ns-real-int32-float64
41dfd57a17SPierre Jolivet       requires: !complex double !defined(PETSC_USE_64BIT_INDICES)
42c4762a1bSJed Brown 
43c4762a1bSJed Brown TEST*/
44