Actual source code: ex16.c


  2: static char help[] = "Tests DMComposite routines.\n\n";

  4: #include <petscdmredundant.h>
  5: #include <petscdm.h>
  6: #include <petscdmda.h>
  7: #include <petscdmcomposite.h>
  8: #include <petscpf.h>

 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;

 24:   PetscInitialize(&argc, &argv, (char *)0, help);
 25:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 26:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 28:   PetscOptionsGetBool(NULL, NULL, "-gather_add", &gather_add, NULL);

 30:   DMCompositeCreate(PETSC_COMM_WORLD, &packer);

 32:   DMRedundantCreate(PETSC_COMM_WORLD, 0, nredundant1, &dmred1);
 33:   DMCreateLocalVector(dmred1, &redundant1);
 34:   DMCompositeAddDM(packer, dmred1);

 36:   DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 8, 1, 1, NULL, &da1);
 37:   DMSetFromOptions(da1);
 38:   DMSetUp(da1);
 39:   DMCreateLocalVector(da1, &local1);
 40:   DMCompositeAddDM(packer, da1);

 42:   DMRedundantCreate(PETSC_COMM_WORLD, 1 % size, nredundant2, &dmred2);
 43:   DMCreateLocalVector(dmred2, &redundant2);
 44:   DMCompositeAddDM(packer, dmred2);

 46:   DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 6, 1, 1, NULL, &da2);
 47:   DMSetFromOptions(da2);
 48:   DMSetUp(da2);
 49:   DMCreateLocalVector(da2, &local2);
 50:   DMCompositeAddDM(packer, da2);

 52:   DMCreateGlobalVector(packer, &global);
 53:   PFCreate(PETSC_COMM_WORLD, 1, 1, &pf);
 54:   PFSetType(pf, PFIDENTITY, NULL);
 55:   PFApplyVec(pf, NULL, global);
 56:   PFDestroy(&pf);
 57:   VecView(global, PETSC_VIEWER_STDOUT_WORLD);

 59:   DMCompositeScatter(packer, global, redundant1, local1, redundant2, local2);
 60:   PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD);
 61:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of redundant1 vector\n", rank);
 62:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 63:   VecView(redundant1, sviewer);
 64:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 65:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
 66:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of da1 vector\n", rank);
 67:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 68:   VecView(local1, sviewer);
 69:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 70:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
 71:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of redundant2 vector\n", rank);
 72:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 73:   VecView(redundant2, sviewer);
 74:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 75:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
 76:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of da2 vector\n", rank);
 77:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 78:   VecView(local2, sviewer);
 79:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
 80:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
 81:   PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD);

 83:   VecGetArray(redundant1, &redundant1a);
 84:   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:   VecRestoreArray(redundant1, &redundant1a);
 88:   VecRestoreArray(redundant2, &redundant2a);

 90:   DMCompositeGather(packer, gather_add ? ADD_VALUES : INSERT_VALUES, global, redundant1, local1, redundant2, local2);
 91:   VecView(global, PETSC_VIEWER_STDOUT_WORLD);

 93:   /* get the global numbering for each subvector element */
 94:   DMCompositeGetISLocalToGlobalMappings(packer, &ltog);

 96:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of redundant1 vector\n");
 97:   ISLocalToGlobalMappingView(ltog[0], PETSC_VIEWER_STDOUT_WORLD);
 98:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of local1 vector\n");
 99:   ISLocalToGlobalMappingView(ltog[1], PETSC_VIEWER_STDOUT_WORLD);
100:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of redundant2 vector\n");
101:   ISLocalToGlobalMappingView(ltog[2], PETSC_VIEWER_STDOUT_WORLD);
102:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of local2 vector\n");
103:   ISLocalToGlobalMappingView(ltog[3], PETSC_VIEWER_STDOUT_WORLD);

105:   for (i = 0; i < 4; i++) ISLocalToGlobalMappingDestroy(&ltog[i]);
106:   PetscFree(ltog);

108:   DMDestroy(&da1);
109:   DMDestroy(&dmred1);
110:   DMDestroy(&dmred2);
111:   DMDestroy(&da2);
112:   VecDestroy(&redundant1);
113:   VecDestroy(&redundant2);
114:   VecDestroy(&local1);
115:   VecDestroy(&local2);
116:   VecDestroy(&global);
117:   DMDestroy(&packer);
118:   PetscFinalize();
119:   return 0;
120: }

122: /*TEST

124:    build:
125:       requires: !complex

127:    test:
128:       nsize: 3

130:    test:
131:       suffix: 2
132:       nsize: 3
133:       args: -gather_add

135: TEST*/