1 2 static char help[] = "Tests various DMComposite routines.\n\n"; 3 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 #include <petscdmcomposite.h> 7 8 int main(int argc,char **argv) 9 { 10 PetscMPIInt rank; 11 PetscErrorCode ierr; 12 DM da1,da2,packer; 13 Vec local,global,globals[2],buffer; 14 PetscScalar value; 15 PetscViewer viewer; 16 17 ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 18 19 CHKERRQ(DMCompositeCreate(PETSC_COMM_WORLD,&packer)); 20 CHKERRQ(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,8,1,1,NULL,&da1)); 21 CHKERRQ(DMSetFromOptions(da1)); 22 CHKERRQ(DMSetUp(da1)); 23 CHKERRQ(DMCompositeAddDM(packer,da1)); 24 CHKERRQ(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,6,1,1,NULL,&da2)); 25 CHKERRQ(DMSetFromOptions(da2)); 26 CHKERRQ(DMSetUp(da2)); 27 CHKERRQ(DMCompositeAddDM(packer,da2)); 28 29 CHKERRQ(DMCreateGlobalVector(packer,&global)); 30 CHKERRQ(DMCreateLocalVector(packer,&local)); 31 CHKERRQ(DMCreateLocalVector(packer,&buffer)); 32 33 CHKERRQ(DMCompositeGetAccessArray(packer,global,2,NULL,globals)); 34 value = 1; 35 CHKERRQ(VecSet(globals[0], value)); 36 value = -1; 37 CHKERRQ(VecSet(globals[1], value)); 38 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 39 value = rank + 1; 40 CHKERRQ(VecScale(globals[0], value)); 41 CHKERRQ(VecScale(globals[1], value)); 42 CHKERRQ(DMCompositeRestoreAccessArray(packer,global,2,NULL,globals)); 43 44 /* Test GlobalToLocal in insert mode */ 45 CHKERRQ(DMGlobalToLocalBegin(packer,global,INSERT_VALUES,local)); 46 CHKERRQ(DMGlobalToLocalEnd(packer,global,INSERT_VALUES,local)); 47 48 CHKERRQ(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 49 CHKERRQ(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nLocal Vector: processor %d\n",rank)); 50 CHKERRQ(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&viewer)); 51 CHKERRQ(VecView(local,viewer)); 52 CHKERRQ(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&viewer)); 53 CHKERRQ(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 54 CHKERRQ(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 55 56 /* Test LocalToGlobal in insert mode */ 57 CHKERRQ(DMLocalToGlobalBegin(packer,local,INSERT_VALUES,global)); 58 CHKERRQ(DMLocalToGlobalEnd(packer,local,INSERT_VALUES,global)); 59 60 CHKERRQ(VecView(global,PETSC_VIEWER_STDOUT_WORLD)); 61 62 /* Test LocalToLocal in insert mode */ 63 CHKERRQ(DMLocalToLocalBegin(packer,local,INSERT_VALUES,buffer)); 64 CHKERRQ(DMLocalToLocalEnd(packer,local,INSERT_VALUES,buffer)); 65 66 CHKERRQ(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 67 CHKERRQ(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nLocal Vector: processor %d\n",rank)); 68 CHKERRQ(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&viewer)); 69 CHKERRQ(VecView(buffer,viewer)); 70 CHKERRQ(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&viewer)); 71 CHKERRQ(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 72 CHKERRQ(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 73 74 CHKERRQ(VecDestroy(&buffer)); 75 CHKERRQ(VecDestroy(&local)); 76 CHKERRQ(VecDestroy(&global)); 77 CHKERRQ(DMDestroy(&packer)); 78 CHKERRQ(DMDestroy(&da2)); 79 CHKERRQ(DMDestroy(&da1)); 80 81 ierr = PetscFinalize(); 82 83 return ierr; 84 } 85 86 /*TEST 87 88 test: 89 nsize: 3 90 91 TEST*/ 92