1*c4762a1bSJed Brown 2*c4762a1bSJed Brown static char help[] = "Reads a PETSc matrix and computes the 2 norm of the columns\n\n"; 3*c4762a1bSJed Brown 4*c4762a1bSJed Brown /*T 5*c4762a1bSJed Brown Concepts: Mat^loading a binary matrix; 6*c4762a1bSJed Brown Processors: n 7*c4762a1bSJed Brown T*/ 8*c4762a1bSJed Brown 9*c4762a1bSJed Brown /* 10*c4762a1bSJed Brown Include "petscmat.h" so that we can use matrices. 11*c4762a1bSJed Brown automatically includes: 12*c4762a1bSJed Brown petscsys.h - base PETSc routines petscvec.h - vectors 13*c4762a1bSJed Brown petscmat.h - matrices 14*c4762a1bSJed Brown petscis.h - index sets petscviewer.h - viewers 15*c4762a1bSJed Brown */ 16*c4762a1bSJed Brown #include <petscmat.h> 17*c4762a1bSJed Brown 18*c4762a1bSJed Brown 19*c4762a1bSJed Brown int main(int argc,char **args) 20*c4762a1bSJed Brown { 21*c4762a1bSJed Brown Mat A; /* matrix */ 22*c4762a1bSJed Brown PetscViewer fd; /* viewer */ 23*c4762a1bSJed Brown char file[PETSC_MAX_PATH_LEN]; /* input file name */ 24*c4762a1bSJed Brown PetscErrorCode ierr; 25*c4762a1bSJed Brown PetscReal *norms; 26*c4762a1bSJed Brown PetscInt n,cstart,cend; 27*c4762a1bSJed Brown PetscBool flg; 28*c4762a1bSJed Brown PetscViewerFormat format; 29*c4762a1bSJed Brown 30*c4762a1bSJed Brown ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 31*c4762a1bSJed Brown /* 32*c4762a1bSJed Brown Determine files from which we read the matrix 33*c4762a1bSJed Brown */ 34*c4762a1bSJed Brown ierr = PetscOptionsGetString(NULL,NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 35*c4762a1bSJed Brown if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option"); 36*c4762a1bSJed Brown 37*c4762a1bSJed Brown /* 38*c4762a1bSJed Brown Open binary file. Note that we use FILE_MODE_READ to indicate 39*c4762a1bSJed Brown reading from this file. 40*c4762a1bSJed Brown */ 41*c4762a1bSJed Brown ierr = PetscViewerCreate(PETSC_COMM_WORLD,&fd);CHKERRQ(ierr); 42*c4762a1bSJed Brown ierr = PetscViewerSetType(fd,PETSCVIEWERBINARY);CHKERRQ(ierr); 43*c4762a1bSJed Brown ierr = PetscViewerSetFromOptions(fd);CHKERRQ(ierr); 44*c4762a1bSJed Brown ierr = PetscOptionsGetEnum(NULL,NULL,"-viewer_format",PetscViewerFormats,(PetscEnum*)&format,&flg);CHKERRQ(ierr); 45*c4762a1bSJed Brown if (flg) {ierr = PetscViewerPushFormat(fd,format);CHKERRQ(ierr);} 46*c4762a1bSJed Brown ierr = PetscViewerFileSetMode(fd,FILE_MODE_READ);CHKERRQ(ierr); 47*c4762a1bSJed Brown ierr = PetscViewerFileSetName(fd,file);CHKERRQ(ierr); 48*c4762a1bSJed Brown 49*c4762a1bSJed Brown /* 50*c4762a1bSJed Brown Load the matrix; then destroy the viewer. 51*c4762a1bSJed Brown Matrix type is set automatically but you can override it by MatSetType() prior to MatLoad(). 52*c4762a1bSJed Brown Do that only if you really insist on the given type. 53*c4762a1bSJed Brown */ 54*c4762a1bSJed Brown ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 55*c4762a1bSJed Brown ierr = MatSetOptionsPrefix(A,"a_");CHKERRQ(ierr); 56*c4762a1bSJed Brown ierr = PetscObjectSetName((PetscObject) A,"A");CHKERRQ(ierr); 57*c4762a1bSJed Brown ierr = MatSetFromOptions(A);CHKERRQ(ierr); 58*c4762a1bSJed Brown ierr = MatLoad(A,fd);CHKERRQ(ierr); 59*c4762a1bSJed Brown ierr = PetscViewerDestroy(&fd);CHKERRQ(ierr); 60*c4762a1bSJed Brown 61*c4762a1bSJed Brown ierr = MatGetSize(A,NULL,&n);CHKERRQ(ierr); 62*c4762a1bSJed Brown ierr = MatGetOwnershipRangeColumn(A,&cstart,&cend);CHKERRQ(ierr); 63*c4762a1bSJed Brown ierr = PetscMalloc1(n,&norms);CHKERRQ(ierr); 64*c4762a1bSJed Brown ierr = MatGetColumnNorms(A,NORM_2,norms);CHKERRQ(ierr); 65*c4762a1bSJed Brown ierr = PetscRealView(cend-cstart,norms+cstart,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 66*c4762a1bSJed Brown ierr = PetscFree(norms);CHKERRQ(ierr); 67*c4762a1bSJed Brown 68*c4762a1bSJed Brown ierr = PetscObjectPrintClassNamePrefixType((PetscObject)A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 69*c4762a1bSJed Brown ierr = MatGetOption(A,MAT_SYMMETRIC,&flg);CHKERRQ(ierr); 70*c4762a1bSJed Brown ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"MAT_SYMMETRIC: %D\n",flg);CHKERRQ(ierr); 71*c4762a1bSJed Brown ierr = MatViewFromOptions(A,NULL,"-mat_view");CHKERRQ(ierr); 72*c4762a1bSJed Brown 73*c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 74*c4762a1bSJed Brown ierr = PetscFinalize(); 75*c4762a1bSJed Brown return ierr; 76*c4762a1bSJed Brown } 77*c4762a1bSJed Brown 78*c4762a1bSJed Brown 79*c4762a1bSJed Brown 80*c4762a1bSJed Brown /*TEST 81*c4762a1bSJed Brown 82*c4762a1bSJed Brown test: 83*c4762a1bSJed Brown suffix: mpiaij 84*c4762a1bSJed Brown nsize: 2 85*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 86*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type mpiaij 87*c4762a1bSJed Brown args: -a_matload_symmetric 88*c4762a1bSJed Brown 89*c4762a1bSJed Brown test: 90*c4762a1bSJed Brown suffix: mpiaij_hdf5 91*c4762a1bSJed Brown nsize: 2 92*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 93*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat 94*c4762a1bSJed Brown args: -a_matload_symmetric 95*c4762a1bSJed Brown 96*c4762a1bSJed Brown test: 97*c4762a1bSJed Brown suffix: mpiaij_rect_hdf5 98*c4762a1bSJed Brown nsize: 2 99*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 100*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat 101*c4762a1bSJed Brown 102*c4762a1bSJed Brown test: 103*c4762a1bSJed Brown # test for more processes than rows 104*c4762a1bSJed Brown suffix: mpiaij_hdf5_tiny 105*c4762a1bSJed Brown nsize: 8 106*c4762a1bSJed Brown requires: double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 107*c4762a1bSJed 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 108*c4762a1bSJed Brown args: -a_matload_symmetric 109*c4762a1bSJed Brown 110*c4762a1bSJed Brown test: 111*c4762a1bSJed Brown # test for more processes than rows, complex 112*c4762a1bSJed Brown TODO: not yet implemented for MATLAB complex format 113*c4762a1bSJed Brown suffix: mpiaij_hdf5_tiny_complex 114*c4762a1bSJed Brown nsize: 8 115*c4762a1bSJed Brown requires: double complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 116*c4762a1bSJed 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 117*c4762a1bSJed Brown args: -a_matload_symmetric 118*c4762a1bSJed Brown 119*c4762a1bSJed Brown test: 120*c4762a1bSJed Brown TODO: mpibaij not supported yet 121*c4762a1bSJed Brown suffix: mpibaij_hdf5 122*c4762a1bSJed Brown nsize: 2 123*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 124*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type mpibaij -a_mat_block_size 2 -viewer_type hdf5 -viewer_format hdf5_mat 125*c4762a1bSJed Brown args: -a_matload_symmetric 126*c4762a1bSJed Brown 127*c4762a1bSJed Brown test: 128*c4762a1bSJed Brown suffix: mpidense 129*c4762a1bSJed Brown nsize: 2 130*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 131*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type mpidense 132*c4762a1bSJed Brown args: -a_matload_symmetric 133*c4762a1bSJed Brown 134*c4762a1bSJed Brown test: 135*c4762a1bSJed Brown suffix: seqaij 136*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 137*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type seqaij 138*c4762a1bSJed Brown args: -a_matload_symmetric 139*c4762a1bSJed Brown 140*c4762a1bSJed Brown test: 141*c4762a1bSJed Brown suffix: seqaij_hdf5 142*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 143*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type seqaij -viewer_type hdf5 -viewer_format hdf5_mat 144*c4762a1bSJed Brown args: -a_matload_symmetric 145*c4762a1bSJed Brown 146*c4762a1bSJed Brown test: 147*c4762a1bSJed Brown suffix: seqaij_rect_hdf5 148*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 149*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect.mat -a_mat_type seqaij -viewer_type hdf5 -viewer_format hdf5_mat 150*c4762a1bSJed Brown 151*c4762a1bSJed Brown test: 152*c4762a1bSJed Brown suffix: seqdense 153*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 154*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/small -a_mat_type seqdense 155*c4762a1bSJed Brown args: -a_matload_symmetric 156*c4762a1bSJed Brown 157*c4762a1bSJed Brown test: 158*c4762a1bSJed Brown suffix: seqdense_hdf5 159*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 160*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_dense.mat -a_mat_type seqdense -viewer_type hdf5 -viewer_format hdf5_mat 161*c4762a1bSJed Brown args: -a_matload_symmetric 162*c4762a1bSJed Brown 163*c4762a1bSJed Brown test: 164*c4762a1bSJed Brown suffix: seqdense_rect_hdf5 165*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 166*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect_dense.mat -a_mat_type seqdense -viewer_type hdf5 -viewer_format hdf5_mat 167*c4762a1bSJed Brown 168*c4762a1bSJed Brown test: 169*c4762a1bSJed Brown suffix: mpidense_hdf5 170*c4762a1bSJed Brown nsize: 2 171*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 172*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_dense.mat -a_mat_type mpidense -viewer_type hdf5 -viewer_format hdf5_mat 173*c4762a1bSJed Brown args: -a_matload_symmetric 174*c4762a1bSJed Brown 175*c4762a1bSJed Brown test: 176*c4762a1bSJed Brown suffix: mpidense_rect_hdf5 177*c4762a1bSJed Brown nsize: 2 178*c4762a1bSJed Brown requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) hdf5 define(PETSC_HDF5_HAVE_ZLIB) 179*c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/matlab/small_rect_dense.mat -a_mat_type mpidense -viewer_type hdf5 -viewer_format hdf5_mat 180*c4762a1bSJed Brown TEST*/ 181