1 2 static char help[] = "Tests DMComposite routines.\n\n"; 3 4 #include <petscdmredundant.h> 5 #include <petscdm.h> 6 #include <petscdmda.h> 7 #include <petscdmcomposite.h> 8 #include <petscpf.h> 9 10 int main(int argc,char **argv) 11 { 12 PetscInt nredundant1 = 5,nredundant2 = 2,i; 13 ISLocalToGlobalMapping *ltog; 14 PetscMPIInt rank,size; 15 DM packer; 16 Vec global,local1,local2,redundant1,redundant2; 17 PF pf; 18 DM da1,da2,dmred1,dmred2; 19 PetscScalar *redundant1a,*redundant2a; 20 PetscViewer sviewer; 21 PetscBool gather_add = PETSC_FALSE; 22 23 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 24 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 25 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 26 27 PetscCall(PetscOptionsGetBool(NULL,NULL,"-gather_add",&gather_add,NULL)); 28 29 PetscCall(DMCompositeCreate(PETSC_COMM_WORLD,&packer)); 30 31 PetscCall(DMRedundantCreate(PETSC_COMM_WORLD,0,nredundant1,&dmred1)); 32 PetscCall(DMCreateLocalVector(dmred1,&redundant1)); 33 PetscCall(DMCompositeAddDM(packer,dmred1)); 34 35 PetscCall(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,8,1,1,NULL,&da1)); 36 PetscCall(DMSetFromOptions(da1)); 37 PetscCall(DMSetUp(da1)); 38 PetscCall(DMCreateLocalVector(da1,&local1)); 39 PetscCall(DMCompositeAddDM(packer,da1)); 40 41 PetscCall(DMRedundantCreate(PETSC_COMM_WORLD,1%size,nredundant2,&dmred2)); 42 PetscCall(DMCreateLocalVector(dmred2,&redundant2)); 43 PetscCall(DMCompositeAddDM(packer,dmred2)); 44 45 PetscCall(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,6,1,1,NULL,&da2)); 46 PetscCall(DMSetFromOptions(da2)); 47 PetscCall(DMSetUp(da2)); 48 PetscCall(DMCreateLocalVector(da2,&local2)); 49 PetscCall(DMCompositeAddDM(packer,da2)); 50 51 PetscCall(DMCreateGlobalVector(packer,&global)); 52 PetscCall(PFCreate(PETSC_COMM_WORLD,1,1,&pf)); 53 PetscCall(PFSetType(pf,PFIDENTITY,NULL)); 54 PetscCall(PFApplyVec(pf,NULL,global)); 55 PetscCall(PFDestroy(&pf)); 56 PetscCall(VecView(global,PETSC_VIEWER_STDOUT_WORLD)); 57 58 PetscCall(DMCompositeScatter(packer,global,redundant1,local1,redundant2,local2)); 59 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 60 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant1 vector\n",rank)); 61 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 62 PetscCall(VecView(redundant1,sviewer)); 63 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 64 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 65 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da1 vector\n",rank)); 66 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 67 PetscCall(VecView(local1,sviewer)); 68 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 69 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 70 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of redundant2 vector\n",rank)); 71 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 72 PetscCall(VecView(redundant2,sviewer)); 73 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 74 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 75 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"[%d] My part of da2 vector\n",rank)); 76 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 77 PetscCall(VecView(local2,sviewer)); 78 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 79 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 80 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 81 82 PetscCall(VecGetArray(redundant1,&redundant1a)); 83 PetscCall(VecGetArray(redundant2,&redundant2a)); 84 for (i=0; i<nredundant1; i++) redundant1a[i] = (rank+2)*i; 85 for (i=0; i<nredundant2; i++) redundant2a[i] = (rank+10)*i; 86 PetscCall(VecRestoreArray(redundant1,&redundant1a)); 87 PetscCall(VecRestoreArray(redundant2,&redundant2a)); 88 89 PetscCall(DMCompositeGather(packer,gather_add ? ADD_VALUES : INSERT_VALUES,global,redundant1,local1,redundant2,local2)); 90 PetscCall(VecView(global,PETSC_VIEWER_STDOUT_WORLD)); 91 92 /* get the global numbering for each subvector element */ 93 PetscCall(DMCompositeGetISLocalToGlobalMappings(packer,<og)); 94 95 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of redundant1 vector\n")); 96 PetscCall(ISLocalToGlobalMappingView(ltog[0],PETSC_VIEWER_STDOUT_WORLD)); 97 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of local1 vector\n")); 98 PetscCall(ISLocalToGlobalMappingView(ltog[1],PETSC_VIEWER_STDOUT_WORLD)); 99 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of redundant2 vector\n")); 100 PetscCall(ISLocalToGlobalMappingView(ltog[2],PETSC_VIEWER_STDOUT_WORLD)); 101 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Local to global mapping of local2 vector\n")); 102 PetscCall(ISLocalToGlobalMappingView(ltog[3],PETSC_VIEWER_STDOUT_WORLD)); 103 104 for (i=0; i<4; i++) PetscCall(ISLocalToGlobalMappingDestroy(<og[i])); 105 PetscCall(PetscFree(ltog)); 106 107 PetscCall(DMDestroy(&da1)); 108 PetscCall(DMDestroy(&dmred1)); 109 PetscCall(DMDestroy(&dmred2)); 110 PetscCall(DMDestroy(&da2)); 111 PetscCall(VecDestroy(&redundant1)); 112 PetscCall(VecDestroy(&redundant2)); 113 PetscCall(VecDestroy(&local1)); 114 PetscCall(VecDestroy(&local2)); 115 PetscCall(VecDestroy(&global)); 116 PetscCall(DMDestroy(&packer)); 117 PetscCall(PetscFinalize()); 118 return 0; 119 } 120 121 /*TEST 122 123 build: 124 requires: !complex 125 126 test: 127 nsize: 3 128 129 test: 130 suffix: 2 131 nsize: 3 132 args: -gather_add 133 134 TEST*/ 135