1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Reads a PETSc matrix and computes the 2 norm of the columns\n\n"; 3c4762a1bSJed Brown 4c4762a1bSJed Brown /* 5c4762a1bSJed Brown Include "petscmat.h" so that we can use matrices. 6c4762a1bSJed Brown automatically includes: 7c4762a1bSJed Brown petscsys.h - base PETSc routines petscvec.h - vectors 8c4762a1bSJed Brown petscmat.h - matrices 9c4762a1bSJed Brown petscis.h - index sets petscviewer.h - viewers 10c4762a1bSJed Brown */ 11c4762a1bSJed Brown #include <petscmat.h> 12c4762a1bSJed Brown 13c4762a1bSJed Brown int main(int argc,char **args) 14c4762a1bSJed Brown { 15c4762a1bSJed Brown Mat A; /* matrix */ 16c4762a1bSJed Brown PetscViewer fd; /* viewer */ 17c4762a1bSJed Brown char file[PETSC_MAX_PATH_LEN]; /* input file name */ 18c4762a1bSJed Brown PetscReal *norms; 19c4762a1bSJed Brown PetscInt n,cstart,cend; 20c4762a1bSJed Brown PetscBool flg; 21c4762a1bSJed Brown PetscViewerFormat format; 22c4762a1bSJed Brown 23*327415f7SBarry Smith PetscFunctionBeginUser; 249566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 25c4762a1bSJed Brown /* 26c4762a1bSJed Brown Determine files from which we read the matrix 27c4762a1bSJed Brown */ 289566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg)); 2928b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option"); 30c4762a1bSJed Brown 31c4762a1bSJed Brown /* 32c4762a1bSJed Brown Open binary file. Note that we use FILE_MODE_READ to indicate 33c4762a1bSJed Brown reading from this file. 34c4762a1bSJed Brown */ 359566063dSJacob Faibussowitsch PetscCall(PetscViewerCreate(PETSC_COMM_WORLD,&fd)); 369566063dSJacob Faibussowitsch PetscCall(PetscViewerSetType(fd,PETSCVIEWERBINARY)); 379566063dSJacob Faibussowitsch PetscCall(PetscViewerSetFromOptions(fd)); 389566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetEnum(NULL,NULL,"-viewer_format",PetscViewerFormats,(PetscEnum*)&format,&flg)); 399566063dSJacob Faibussowitsch if (flg) PetscCall(PetscViewerPushFormat(fd,format)); 409566063dSJacob Faibussowitsch PetscCall(PetscViewerFileSetMode(fd,FILE_MODE_READ)); 419566063dSJacob Faibussowitsch PetscCall(PetscViewerFileSetName(fd,file)); 42c4762a1bSJed Brown 43c4762a1bSJed Brown /* 44c4762a1bSJed Brown Load the matrix; then destroy the viewer. 45c4762a1bSJed Brown Matrix type is set automatically but you can override it by MatSetType() prior to MatLoad(). 46c4762a1bSJed Brown Do that only if you really insist on the given type. 47c4762a1bSJed Brown */ 489566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 499566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(A,"a_")); 509566063dSJacob Faibussowitsch PetscCall(PetscObjectSetName((PetscObject) A,"A")); 519566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(A)); 529566063dSJacob Faibussowitsch PetscCall(MatLoad(A,fd)); 539566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&fd)); 54c4762a1bSJed Brown 559566063dSJacob Faibussowitsch PetscCall(MatGetSize(A,NULL,&n)); 569566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRangeColumn(A,&cstart,&cend)); 579566063dSJacob Faibussowitsch PetscCall(PetscMalloc1(n,&norms)); 589566063dSJacob Faibussowitsch PetscCall(MatGetColumnNorms(A,NORM_2,norms)); 599566063dSJacob Faibussowitsch PetscCall(PetscRealView(cend-cstart,norms+cstart,PETSC_VIEWER_STDOUT_WORLD)); 609566063dSJacob Faibussowitsch PetscCall(PetscFree(norms)); 61c4762a1bSJed Brown 629566063dSJacob Faibussowitsch PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)A,PETSC_VIEWER_STDOUT_WORLD)); 639566063dSJacob Faibussowitsch PetscCall(MatViewFromOptions(A,NULL,"-mat_view")); 64c4762a1bSJed Brown 659566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 669566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 67b122ec5aSJacob Faibussowitsch return 0; 68c4762a1bSJed Brown } 69c4762a1bSJed Brown 70c4762a1bSJed Brown /*TEST 71c4762a1bSJed Brown 72c4762a1bSJed Brown test: 73c4762a1bSJed Brown suffix: mpiaij 74c4762a1bSJed Brown nsize: 2 75dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 76c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type mpiaij 77c4762a1bSJed Brown args: -a_matload_symmetric 78c4762a1bSJed Brown 79c4762a1bSJed Brown test: 80c4762a1bSJed Brown suffix: mpiaij_hdf5 81c4762a1bSJed Brown nsize: 2 82dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 83c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat 84c4762a1bSJed Brown args: -a_matload_symmetric 85c4762a1bSJed Brown 86c4762a1bSJed Brown test: 87c4762a1bSJed Brown suffix: mpiaij_rect_hdf5 88c4762a1bSJed Brown nsize: 2 89dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 90c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat 91c4762a1bSJed Brown 92c4762a1bSJed Brown test: 93c4762a1bSJed Brown # test for more processes than rows 94c4762a1bSJed Brown suffix: mpiaij_hdf5_tiny 95c4762a1bSJed Brown nsize: 8 96dfd57a17SPierre Jolivet requires: double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 97c4762a1bSJed Brown args: -f ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system_with_x0.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat 98c4762a1bSJed Brown args: -a_matload_symmetric 99c4762a1bSJed Brown 100c4762a1bSJed Brown test: 101c4762a1bSJed Brown # test for more processes than rows, complex 102c4762a1bSJed Brown TODO: not yet implemented for MATLAB complex format 103c4762a1bSJed Brown suffix: mpiaij_hdf5_tiny_complex 104c4762a1bSJed Brown nsize: 8 105dfd57a17SPierre Jolivet requires: double complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 106c4762a1bSJed Brown args: -f ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system_with_x0_complex.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat 107c4762a1bSJed Brown args: -a_matload_symmetric 108c4762a1bSJed Brown 109c4762a1bSJed Brown test: 110c4762a1bSJed Brown TODO: mpibaij not supported yet 111c4762a1bSJed Brown suffix: mpibaij_hdf5 112c4762a1bSJed Brown nsize: 2 113dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 114c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type mpibaij -a_mat_block_size 2 -viewer_type hdf5 -viewer_format hdf5_mat 115c4762a1bSJed Brown args: -a_matload_symmetric 116c4762a1bSJed Brown 117c4762a1bSJed Brown test: 118c4762a1bSJed Brown suffix: mpidense 119c4762a1bSJed Brown nsize: 2 120dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 121c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type mpidense 122c4762a1bSJed Brown args: -a_matload_symmetric 123c4762a1bSJed Brown 124c4762a1bSJed Brown test: 125c4762a1bSJed Brown suffix: seqaij 126dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 127c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type seqaij 128c4762a1bSJed Brown args: -a_matload_symmetric 129c4762a1bSJed Brown 130c4762a1bSJed Brown test: 131c4762a1bSJed Brown suffix: seqaij_hdf5 132dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 133c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type seqaij -viewer_type hdf5 -viewer_format hdf5_mat 134c4762a1bSJed Brown args: -a_matload_symmetric 135c4762a1bSJed Brown 136c4762a1bSJed Brown test: 137c4762a1bSJed Brown suffix: seqaij_rect_hdf5 138dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 139c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect.mat -a_mat_type seqaij -viewer_type hdf5 -viewer_format hdf5_mat 140c4762a1bSJed Brown 141c4762a1bSJed Brown test: 142c4762a1bSJed Brown suffix: seqdense 143dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 144c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type seqdense 145c4762a1bSJed Brown args: -a_matload_symmetric 146c4762a1bSJed Brown 147c4762a1bSJed Brown test: 148c4762a1bSJed Brown suffix: seqdense_hdf5 149dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 150c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_dense.mat -a_mat_type seqdense -viewer_type hdf5 -viewer_format hdf5_mat 151c4762a1bSJed Brown args: -a_matload_symmetric 152c4762a1bSJed Brown 153c4762a1bSJed Brown test: 154c4762a1bSJed Brown suffix: seqdense_rect_hdf5 155dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 156c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect_dense.mat -a_mat_type seqdense -viewer_type hdf5 -viewer_format hdf5_mat 157c4762a1bSJed Brown 158c4762a1bSJed Brown test: 159c4762a1bSJed Brown suffix: mpidense_hdf5 160c4762a1bSJed Brown nsize: 2 161dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 162c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_dense.mat -a_mat_type mpidense -viewer_type hdf5 -viewer_format hdf5_mat 163c4762a1bSJed Brown args: -a_matload_symmetric 164c4762a1bSJed Brown 165c4762a1bSJed Brown test: 166c4762a1bSJed Brown suffix: mpidense_rect_hdf5 167c4762a1bSJed Brown nsize: 2 168dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB) 169c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect_dense.mat -a_mat_type mpidense -viewer_type hdf5 -viewer_format hdf5_mat 170c4762a1bSJed Brown TEST*/ 171