1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "MatLoad test for loading matrices that are created by DMCreateMatrix() and\n\ 3c4762a1bSJed Brown stored in binary via MatView_MPI_DA.MatView_MPI_DA stores the matrix\n\ 4c4762a1bSJed Brown in natural ordering. Hence MatLoad() has to read the matrix first in\n\ 5c4762a1bSJed Brown natural ordering and then permute it back to the application ordering.This\n\ 6c4762a1bSJed Brown example is used for testing the subroutine MatLoad_MPI_DA\n\n"; 7c4762a1bSJed Brown 8c4762a1bSJed Brown #include <petscdm.h> 9c4762a1bSJed Brown #include <petscdmda.h> 10c4762a1bSJed Brown 11c4762a1bSJed Brown int main(int argc,char **argv) 12c4762a1bSJed Brown { 13c4762a1bSJed Brown PetscInt X = 10,Y = 8,Z=8; 14c4762a1bSJed Brown DM da; 15c4762a1bSJed Brown PetscViewer viewer; 16c4762a1bSJed Brown Mat A; 17c4762a1bSJed Brown 18*9566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 19*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp.dat",FILE_MODE_WRITE,&viewer)); 20c4762a1bSJed Brown 21c4762a1bSJed Brown /* Read options */ 22*9566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL,NULL,"-X",&X,NULL)); 23*9566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL,NULL,"-Y",&Y,NULL)); 24*9566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL,NULL,"-Z",&Z,NULL)); 25c4762a1bSJed Brown 26c4762a1bSJed Brown /* Create distributed array and get vectors */ 27*9566063dSJacob Faibussowitsch PetscCall(DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,X,Y,Z,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,1,1,NULL,NULL,NULL,&da)); 28*9566063dSJacob Faibussowitsch PetscCall(DMSetFromOptions(da)); 29*9566063dSJacob Faibussowitsch PetscCall(DMSetUp(da)); 30*9566063dSJacob Faibussowitsch PetscCall(DMSetMatType(da,MATMPIAIJ)); 31*9566063dSJacob Faibussowitsch PetscCall(DMCreateMatrix(da,&A)); 32*9566063dSJacob Faibussowitsch PetscCall(MatShift(A,X)); 33*9566063dSJacob Faibussowitsch PetscCall(MatView(A,viewer)); 34*9566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 35*9566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 36c4762a1bSJed Brown 37*9566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp.dat",FILE_MODE_READ,&viewer)); 38*9566063dSJacob Faibussowitsch PetscCall(DMCreateMatrix(da,&A)); 39*9566063dSJacob Faibussowitsch PetscCall(MatLoad(A,viewer)); 40c4762a1bSJed Brown 41c4762a1bSJed Brown /* Free memory */ 42*9566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 43*9566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 44*9566063dSJacob Faibussowitsch PetscCall(DMDestroy(&da)); 45*9566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 46b122ec5aSJacob Faibussowitsch return 0; 47c4762a1bSJed Brown } 48