1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Reads a PETSc matrix and vector from a file and reorders it.\n\ 3c4762a1bSJed Brown -f0 <input_file> : first file to load (small system)\n\ 4c4762a1bSJed Brown -f1 <input_file> : second file to load (larger system)\n\n"; 5c4762a1bSJed Brown 6c4762a1bSJed Brown /* 7c4762a1bSJed Brown Include "petscmat.h" so that we can use matrices. 8c4762a1bSJed Brown automatically includes: 9c4762a1bSJed Brown petscsys.h - base PETSc routines petscvec.h - vectors 10c4762a1bSJed Brown petscmat.h - matrices 11c4762a1bSJed Brown petscis.h - index sets petscviewer.h - viewers 12c4762a1bSJed Brown */ 13c4762a1bSJed Brown #include <petscmat.h> 14c4762a1bSJed Brown 15c4762a1bSJed Brown int main(int argc,char **args) 16c4762a1bSJed Brown { 17c4762a1bSJed Brown Mat A; /* matrix */ 18c4762a1bSJed Brown PetscViewer fd; /* viewer */ 19c4762a1bSJed Brown char file[PETSC_MAX_PATH_LEN]; /* input file name */ 20c4762a1bSJed Brown IS isrow,iscol; /* row and column permutations */ 21c4762a1bSJed Brown MatOrderingType rtype = MATORDERINGRCM; 22c4762a1bSJed Brown PetscBool flg; 23c4762a1bSJed Brown 24*327415f7SBarry Smith PetscFunctionBeginUser; 259566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 26c4762a1bSJed Brown /* 27c4762a1bSJed Brown Determine files from which we read the two linear systems 28c4762a1bSJed Brown (matrix and right-hand-side vector). 29c4762a1bSJed Brown */ 309566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg)); 3128b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option"); 32c4762a1bSJed Brown 33c4762a1bSJed Brown /* ----------------------------------------------------------- 34c4762a1bSJed Brown Beginning of loop 35c4762a1bSJed Brown ----------------------------------------------------------- */ 36c4762a1bSJed Brown /* 37c4762a1bSJed Brown Loop through the reordering 2 times. 38c4762a1bSJed Brown - The intention here is to preload and solve a small system; 39c4762a1bSJed Brown then load another (larger) system and solve it as well. 40c4762a1bSJed Brown This process preloads the instructions with the smaller 41c4762a1bSJed Brown system so that more accurate performance monitoring (via 42c4762a1bSJed Brown -log_view) can be done with the larger one (that actually 43c4762a1bSJed Brown is the system of interest). 44c4762a1bSJed Brown */ 45c4762a1bSJed Brown 46c4762a1bSJed Brown /* - - - - - - - - - - - New Stage - - - - - - - - - - - - - 47c4762a1bSJed Brown Load system i 48c4762a1bSJed Brown - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 49c4762a1bSJed Brown 50c4762a1bSJed Brown /* 51c4762a1bSJed Brown Open binary file. Note that we use FILE_MODE_READ to indicate 52c4762a1bSJed Brown reading from this file. 53c4762a1bSJed Brown */ 549566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd)); 55c4762a1bSJed Brown 56c4762a1bSJed Brown /* 57c4762a1bSJed Brown Load the matrix; then destroy the viewer. 58c4762a1bSJed Brown */ 599566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 609566063dSJacob Faibussowitsch PetscCall(MatSetType(A,MATSEQAIJ)); 619566063dSJacob Faibussowitsch PetscCall(MatLoad(A,fd)); 629566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&fd)); 63c4762a1bSJed Brown 649566063dSJacob Faibussowitsch PetscCall(MatGetOrdering(A,rtype,&isrow,&iscol)); 659566063dSJacob Faibussowitsch PetscCall(ISView(isrow,PETSC_VIEWER_STDOUT_WORLD)); 66c4762a1bSJed Brown 67c4762a1bSJed Brown /* 68c4762a1bSJed Brown Free work space. All PETSc objects should be destroyed when they 69c4762a1bSJed Brown are no longer needed. 70c4762a1bSJed Brown */ 719566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 729566063dSJacob Faibussowitsch PetscCall(ISDestroy(&isrow)); 739566063dSJacob Faibussowitsch PetscCall(ISDestroy(&iscol)); 74c4762a1bSJed Brown 759566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 76b122ec5aSJacob Faibussowitsch return 0; 77c4762a1bSJed Brown } 78c4762a1bSJed Brown 79c4762a1bSJed Brown /*TEST 80c4762a1bSJed Brown 81c4762a1bSJed Brown test: 82dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 83c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/medium -viewer_binary_skip_info 84c4762a1bSJed Brown 85c4762a1bSJed Brown TEST*/ 86