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 PetscFunctionBeginUser; 24 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 25 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 26 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 27 28 PetscCall(PetscOptionsGetBool(NULL, NULL, "-gather_add", &gather_add, NULL)); 29 30 PetscCall(DMCompositeCreate(PETSC_COMM_WORLD, &packer)); 31 32 PetscCall(DMRedundantCreate(PETSC_COMM_WORLD, 0, nredundant1, &dmred1)); 33 PetscCall(DMCreateLocalVector(dmred1, &redundant1)); 34 PetscCall(DMCompositeAddDM(packer, dmred1)); 35 36 PetscCall(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 8, 1, 1, NULL, &da1)); 37 PetscCall(DMSetFromOptions(da1)); 38 PetscCall(DMSetUp(da1)); 39 PetscCall(DMCreateLocalVector(da1, &local1)); 40 PetscCall(DMCompositeAddDM(packer, da1)); 41 42 PetscCall(DMRedundantCreate(PETSC_COMM_WORLD, 1 % size, nredundant2, &dmred2)); 43 PetscCall(DMCreateLocalVector(dmred2, &redundant2)); 44 PetscCall(DMCompositeAddDM(packer, dmred2)); 45 46 PetscCall(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 6, 1, 1, NULL, &da2)); 47 PetscCall(DMSetFromOptions(da2)); 48 PetscCall(DMSetUp(da2)); 49 PetscCall(DMCreateLocalVector(da2, &local2)); 50 PetscCall(DMCompositeAddDM(packer, da2)); 51 52 PetscCall(DMCreateGlobalVector(packer, &global)); 53 PetscCall(PFCreate(PETSC_COMM_WORLD, 1, 1, &pf)); 54 PetscCall(PFSetType(pf, PFIDENTITY, NULL)); 55 PetscCall(PFApplyVec(pf, NULL, global)); 56 PetscCall(PFDestroy(&pf)); 57 PetscCall(VecView(global, PETSC_VIEWER_STDOUT_WORLD)); 58 59 PetscCall(DMCompositeScatter(packer, global, redundant1, local1, redundant2, local2)); 60 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 61 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of redundant1 vector\n", rank)); 62 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 63 PetscCall(VecView(redundant1, sviewer)); 64 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 65 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 66 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of da1 vector\n", rank)); 67 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 68 PetscCall(VecView(local1, sviewer)); 69 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 70 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 71 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of redundant2 vector\n", rank)); 72 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 73 PetscCall(VecView(redundant2, sviewer)); 74 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 75 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 76 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of da2 vector\n", rank)); 77 PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 78 PetscCall(VecView(local2, sviewer)); 79 PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer)); 80 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 81 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 82 83 PetscCall(VecGetArray(redundant1, &redundant1a)); 84 PetscCall(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 PetscCall(VecRestoreArray(redundant1, &redundant1a)); 88 PetscCall(VecRestoreArray(redundant2, &redundant2a)); 89 90 PetscCall(DMCompositeGather(packer, gather_add ? ADD_VALUES : INSERT_VALUES, global, redundant1, local1, redundant2, local2)); 91 PetscCall(VecView(global, PETSC_VIEWER_STDOUT_WORLD)); 92 93 /* get the global numbering for each subvector element */ 94 PetscCall(DMCompositeGetISLocalToGlobalMappings(packer, <og)); 95 96 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of redundant1 vector\n")); 97 PetscCall(ISLocalToGlobalMappingView(ltog[0], PETSC_VIEWER_STDOUT_WORLD)); 98 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of local1 vector\n")); 99 PetscCall(ISLocalToGlobalMappingView(ltog[1], PETSC_VIEWER_STDOUT_WORLD)); 100 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of redundant2 vector\n")); 101 PetscCall(ISLocalToGlobalMappingView(ltog[2], PETSC_VIEWER_STDOUT_WORLD)); 102 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of local2 vector\n")); 103 PetscCall(ISLocalToGlobalMappingView(ltog[3], PETSC_VIEWER_STDOUT_WORLD)); 104 105 for (i = 0; i < 4; i++) PetscCall(ISLocalToGlobalMappingDestroy(<og[i])); 106 PetscCall(PetscFree(ltog)); 107 108 PetscCall(DMDestroy(&da1)); 109 PetscCall(DMDestroy(&dmred1)); 110 PetscCall(DMDestroy(&dmred2)); 111 PetscCall(DMDestroy(&da2)); 112 PetscCall(VecDestroy(&redundant1)); 113 PetscCall(VecDestroy(&redundant2)); 114 PetscCall(VecDestroy(&local1)); 115 PetscCall(VecDestroy(&local2)); 116 PetscCall(VecDestroy(&global)); 117 PetscCall(DMDestroy(&packer)); 118 PetscCall(PetscFinalize()); 119 return 0; 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