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