1c4762a1bSJed Brown static char help[] = "Tests MatView()/MatLoad() with binary viewers for BAIJ matrices.\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown #include <petscviewer.h> 5c4762a1bSJed Brown 6c4762a1bSJed Brown #include <petsc/private/hashtable.h> 7c4762a1bSJed Brown static PetscReal MakeValue(PetscInt i,PetscInt j,PetscInt M) 8c4762a1bSJed Brown { 9c4762a1bSJed Brown PetscHash_t h = PetscHashCombine(PetscHashInt(i),PetscHashInt(j)); 10c4762a1bSJed Brown return (PetscReal) ((h % 5 == 0) ? (1 + i + j*M) : 0); 11c4762a1bSJed Brown } 12c4762a1bSJed Brown 13c4762a1bSJed Brown static PetscErrorCode CheckValuesAIJ(Mat A) 14c4762a1bSJed Brown { 15c4762a1bSJed Brown PetscInt M,N,rstart,rend,i,j; 16c4762a1bSJed Brown PetscReal v,w; 17c4762a1bSJed Brown PetscScalar val; 18c4762a1bSJed Brown PetscErrorCode ierr; 19c4762a1bSJed Brown 20c4762a1bSJed Brown PetscFunctionBegin; 21c4762a1bSJed Brown ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 22c4762a1bSJed Brown ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 23c4762a1bSJed Brown for (i=rstart; i<rend; i++) { 24c4762a1bSJed Brown for (j=0; j<N; j++) { 25c4762a1bSJed Brown ierr = MatGetValue(A,i,j,&val);CHKERRQ(ierr); 26c4762a1bSJed Brown v = MakeValue(i,j,M); w = PetscRealPart(val); 27*2c71b3e2SJacob Faibussowitsch PetscCheckFalse(PetscAbsReal(v-w) > 0,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Matrix entry (%" PetscInt_FMT ",%" PetscInt_FMT ") should be %g, got %g",i,j,(double)v,(double)w); 28c4762a1bSJed Brown } 29c4762a1bSJed Brown } 30c4762a1bSJed Brown PetscFunctionReturn(0); 31c4762a1bSJed Brown } 32c4762a1bSJed Brown 33c4762a1bSJed Brown int main(int argc,char **args) 34c4762a1bSJed Brown { 35c4762a1bSJed Brown Mat A; 36c4762a1bSJed Brown PetscInt M = 24,N = 48,bs = 2; 37c4762a1bSJed Brown PetscInt rstart,rend,i,j; 38c4762a1bSJed Brown PetscErrorCode ierr; 39c4762a1bSJed Brown PetscViewer view; 40c4762a1bSJed Brown 41c4762a1bSJed Brown ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr; 42c4762a1bSJed Brown /* 43c4762a1bSJed Brown Create a parallel BAIJ matrix shared by all processors 44c4762a1bSJed Brown */ 45c4762a1bSJed Brown ierr = MatCreateBAIJ(PETSC_COMM_WORLD, 46c4762a1bSJed Brown bs, 47c4762a1bSJed Brown PETSC_DECIDE,PETSC_DECIDE, 48c4762a1bSJed Brown M,N, 49c4762a1bSJed Brown PETSC_DECIDE,NULL, 50c4762a1bSJed Brown PETSC_DECIDE,NULL, 51c4762a1bSJed Brown &A);CHKERRQ(ierr); 52c4762a1bSJed Brown 53c4762a1bSJed Brown /* 54c4762a1bSJed Brown Set values into the matrix 55c4762a1bSJed Brown */ 56c4762a1bSJed Brown ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 57c4762a1bSJed Brown ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr); 58c4762a1bSJed Brown for (i=rstart; i<rend; i++) { 59c4762a1bSJed Brown for (j=0; j<N; j++) { 60c4762a1bSJed Brown PetscReal v = MakeValue(i,j,M); 61c4762a1bSJed Brown if (PetscAbsReal(v) > 0) { 62c4762a1bSJed Brown ierr = MatSetValue(A,i,j,v,INSERT_VALUES);CHKERRQ(ierr); 63c4762a1bSJed Brown } 64c4762a1bSJed Brown } 65c4762a1bSJed Brown } 66c4762a1bSJed Brown ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 67c4762a1bSJed Brown ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 68c4762a1bSJed Brown ierr = MatViewFromOptions(A,NULL,"-mat_base_view");CHKERRQ(ierr); 69c4762a1bSJed Brown 70c4762a1bSJed Brown /* 71c4762a1bSJed Brown Store the binary matrix to a file 72c4762a1bSJed Brown */ 73c4762a1bSJed Brown ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr); 74c4762a1bSJed Brown for (i=0; i<3; i++) { 75c4762a1bSJed Brown ierr = MatView(A,view);CHKERRQ(ierr); 76c4762a1bSJed Brown } 77c4762a1bSJed Brown ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 78c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 79c4762a1bSJed Brown 80c4762a1bSJed Brown /* 81c4762a1bSJed Brown Now reload the matrix and check its values 82c4762a1bSJed Brown */ 83c4762a1bSJed Brown ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 84c4762a1bSJed Brown ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 85c4762a1bSJed Brown ierr = MatSetType(A,MATBAIJ);CHKERRQ(ierr); 86c4762a1bSJed Brown for (i=0; i<3; i++) { 87c4762a1bSJed Brown if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 88c4762a1bSJed Brown ierr = MatLoad(A,view);CHKERRQ(ierr); 89c4762a1bSJed Brown ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 90c4762a1bSJed Brown } 91c4762a1bSJed Brown ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 92c4762a1bSJed Brown ierr = MatViewFromOptions(A,NULL,"-mat_load_view");CHKERRQ(ierr); 93c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 94c4762a1bSJed Brown 95c4762a1bSJed Brown /* 96c4762a1bSJed Brown Reload in SEQBAIJ matrix and check its values 97c4762a1bSJed Brown */ 98c4762a1bSJed Brown ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 99c4762a1bSJed Brown ierr = MatCreate(PETSC_COMM_SELF,&A);CHKERRQ(ierr); 100c4762a1bSJed Brown ierr = MatSetType(A,MATSEQBAIJ);CHKERRQ(ierr); 101c4762a1bSJed Brown for (i=0; i<3; i++) { 102c4762a1bSJed Brown if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 103c4762a1bSJed Brown ierr = MatLoad(A,view);CHKERRQ(ierr); 104c4762a1bSJed Brown ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 105c4762a1bSJed Brown } 106c4762a1bSJed Brown ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 107c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 108c4762a1bSJed Brown 109c4762a1bSJed Brown /* 110c4762a1bSJed Brown Reload in MPIBAIJ matrix and check its values 111c4762a1bSJed Brown */ 112c4762a1bSJed Brown ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr); 113c4762a1bSJed Brown ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 114c4762a1bSJed Brown ierr = MatSetType(A,MATMPIBAIJ);CHKERRQ(ierr); 115c4762a1bSJed Brown for (i=0; i<3; i++) { 116c4762a1bSJed Brown if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);} 117c4762a1bSJed Brown ierr = MatLoad(A,view);CHKERRQ(ierr); 118c4762a1bSJed Brown ierr = CheckValuesAIJ(A);CHKERRQ(ierr); 119c4762a1bSJed Brown } 120c4762a1bSJed Brown ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 121c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 122c4762a1bSJed Brown 123c4762a1bSJed Brown ierr = PetscFinalize(); 124c4762a1bSJed Brown return ierr; 125c4762a1bSJed Brown } 126c4762a1bSJed Brown 127c4762a1bSJed Brown /*TEST 128c4762a1bSJed Brown 129c4762a1bSJed Brown testset: 130c4762a1bSJed Brown args: -viewer_binary_mpiio 0 131c4762a1bSJed Brown output_file: output/ex45.out 132c4762a1bSJed Brown test: 133c4762a1bSJed Brown suffix: stdio_1 134c4762a1bSJed Brown nsize: 1 135c4762a1bSJed Brown test: 136c4762a1bSJed Brown suffix: stdio_2 137c4762a1bSJed Brown nsize: 2 138c4762a1bSJed Brown test: 139c4762a1bSJed Brown suffix: stdio_3 140c4762a1bSJed Brown nsize: 3 141c4762a1bSJed Brown test: 142c4762a1bSJed Brown suffix: stdio_4 143c4762a1bSJed Brown nsize: 4 144c4762a1bSJed Brown test: 145c4762a1bSJed Brown suffix: stdio_5 146c4762a1bSJed Brown nsize: 4 147c4762a1bSJed Brown 148c4762a1bSJed Brown testset: 149c4762a1bSJed Brown requires: mpiio 150c4762a1bSJed Brown args: -viewer_binary_mpiio 1 151c4762a1bSJed Brown output_file: output/ex45.out 152c4762a1bSJed Brown test: 153c4762a1bSJed Brown suffix: mpiio_1 154c4762a1bSJed Brown nsize: 1 155c4762a1bSJed Brown test: 156c4762a1bSJed Brown suffix: mpiio_2 157c4762a1bSJed Brown nsize: 2 158c4762a1bSJed Brown test: 159c4762a1bSJed Brown suffix: mpiio_3 160c4762a1bSJed Brown nsize: 3 161c4762a1bSJed Brown test: 162c4762a1bSJed Brown suffix: mpiio_4 163c4762a1bSJed Brown nsize: 4 164c4762a1bSJed Brown test: 165c4762a1bSJed Brown suffix: mpiio_5 166c4762a1bSJed Brown nsize: 5 167c4762a1bSJed Brown 168c4762a1bSJed Brown TEST*/ 169