1 static char help[] = "Tests MatLoad() for dense matrix with uneven dimensions set in program\n\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc,char **args) 6 { 7 Mat A; 8 PetscViewer fd; 9 PetscErrorCode ierr; 10 PetscMPIInt rank; 11 PetscScalar *Av; 12 PetscInt i; 13 14 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 15 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr); 16 17 ierr = MatCreateDense(PETSC_COMM_WORLD,6,6,12,12,NULL,&A);CHKERRQ(ierr); 18 ierr = MatDenseGetArray(A,&Av);CHKERRQ(ierr); 19 for (i=0; i<6*12; i++) Av[i] = (PetscScalar) i; 20 ierr = MatDenseRestoreArray(A,&Av);CHKERRQ(ierr); 21 22 /* Load matrices */ 23 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"ex191matrix",FILE_MODE_WRITE,&fd);CHKERRQ(ierr); 24 ierr = PetscViewerPushFormat(fd,PETSC_VIEWER_NATIVE);CHKERRQ(ierr); 25 ierr = MatView(A,fd);CHKERRQ(ierr); 26 ierr = MatDestroy(&A);CHKERRQ(ierr); 27 ierr = PetscViewerPopFormat(fd);CHKERRQ(ierr); 28 ierr = PetscViewerDestroy(&fd);CHKERRQ(ierr); 29 30 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 31 ierr = MatSetType(A,MATDENSE);CHKERRQ(ierr); 32 if (!rank) { 33 ierr = MatSetSizes(A, 4, PETSC_DETERMINE, PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr); 34 } else { 35 ierr = MatSetSizes(A, 8, PETSC_DETERMINE, PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr); 36 } 37 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"ex191matrix",FILE_MODE_READ,&fd);CHKERRQ(ierr); 38 ierr = MatLoad(A,fd);CHKERRQ(ierr); 39 ierr = PetscViewerDestroy(&fd);CHKERRQ(ierr); 40 ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 41 ierr = PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL);CHKERRQ(ierr); 42 ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 43 ierr = PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 44 ierr = MatDestroy(&A);CHKERRQ(ierr); 45 ierr = PetscFinalize(); 46 return ierr; 47 } 48 49 /*TEST 50 51 test: 52 nsize: 2 53 filter: grep -v alloced 54 55 TEST*/ 56