1c4762a1bSJed Brown /* 2c4762a1bSJed Brown Test DMCreateMatrix() for structure_only 3c4762a1bSJed Brown */ 4c4762a1bSJed Brown 5c4762a1bSJed Brown #include <petscdmda.h> 6c4762a1bSJed Brown 7c4762a1bSJed Brown int main(int argc, char *argv[]) 8c4762a1bSJed Brown { 9c4762a1bSJed Brown PetscErrorCode ierr; 10c4762a1bSJed Brown PetscInt nx=6,ny=6,nz=6,dim=1,dof=2; 11c4762a1bSJed Brown DM da; 12c4762a1bSJed Brown Mat A; 13c4762a1bSJed Brown PetscBool struct_only=PETSC_TRUE; 14c4762a1bSJed Brown 15c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,NULL,NULL);if (ierr) return ierr; 16*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-dim",&dim,NULL)); 17*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-dof",&dof,NULL)); 18c4762a1bSJed Brown switch (dim) { 19c4762a1bSJed Brown case 1: 20*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,nx,dof,1,NULL,&da)); 21c4762a1bSJed Brown break; 22c4762a1bSJed Brown case 2: 23*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,nx,ny,PETSC_DECIDE,PETSC_DECIDE,dof,1,NULL,NULL,&da)); 24c4762a1bSJed Brown break; 25c4762a1bSJed Brown default: 26c4762a1bSJed Brown ierr = DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,nx,ny,nz, 27c4762a1bSJed Brown PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,2,NULL,NULL,NULL,&da);CHKERRQ(ierr); 28c4762a1bSJed Brown } 29c4762a1bSJed Brown 30*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetFromOptions(da)); 31*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetUp(da)); 32*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMView(da,PETSC_VIEWER_STDOUT_WORLD)); 33c4762a1bSJed Brown 34*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetBool(NULL,NULL,"-struct_only",&struct_only,NULL)); 35*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetMatrixStructureOnly(da,struct_only)); 36*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMCreateMatrix(da,&A)); 37c4762a1bSJed Brown /* Set da->structure_only to default PETSC_FALSE in case da is being used to create new matrices */ 38*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetMatrixStructureOnly(da,PETSC_FALSE)); 39c4762a1bSJed Brown 40*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatView(A,PETSC_VIEWER_STDOUT_WORLD)); 41*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&A)); 42*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDestroy(&da)); 43c4762a1bSJed Brown ierr = PetscFinalize(); 44c4762a1bSJed Brown return ierr; 45c4762a1bSJed Brown } 46c4762a1bSJed Brown 47c4762a1bSJed Brown /*TEST 48c4762a1bSJed Brown 49c4762a1bSJed Brown test: 50c4762a1bSJed Brown 51c4762a1bSJed Brown test: 52c4762a1bSJed Brown suffix: 2 53c4762a1bSJed Brown args: -dm_mat_type baij -dim 2 54c4762a1bSJed Brown 55c4762a1bSJed Brown TEST*/ 56