Actual source code: ex50.c

  1: static char help[] = "Test DMStagVecSplitToDMDA()\n\n";

  3: #include <petscdm.h>
  4: #include <petscdmstag.h>

  6: int main(int argc, char **argv)
  7: {
  8:   DM                    dm;
  9:   Vec                   x;
 10:   PetscBool             coords;
 11:   PetscInt              dim, dof[4];
 12:   PetscInt              n_loc[4];
 13:   DMStagStencilLocation loc[4][3];

 15:   PetscFunctionBeginUser;
 16:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 17:   dim = 2;
 18:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL));
 19:   coords = PETSC_TRUE;
 20:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-coords", &coords, NULL));

 22:   // Create DMStag and setup set of locations to test
 23:   switch (dim) {
 24:   case 1:
 25:     PetscCall(DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 3, 1, 2, DMSTAG_STENCIL_BOX, 1, NULL, &dm));
 26:     n_loc[0]  = 1;
 27:     loc[0][0] = DMSTAG_LEFT;
 28:     n_loc[1]  = 1;
 29:     loc[1][0] = DMSTAG_ELEMENT;
 30:     n_loc[2]  = 0;
 31:     n_loc[3]  = 0;
 32:     break;
 33:   case 2:
 34:     PetscCall(DMStagCreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 2, PETSC_DECIDE, PETSC_DECIDE, 2, 1, 3, DMSTAG_STENCIL_BOX, 1, NULL, NULL, &dm));
 35:     n_loc[0]  = 1;
 36:     loc[0][0] = DMSTAG_DOWN_LEFT;
 37:     n_loc[1]  = 2;
 38:     loc[1][0] = DMSTAG_LEFT;
 39:     loc[1][1] = DMSTAG_DOWN;
 40:     n_loc[2]  = 1;
 41:     loc[2][0] = DMSTAG_ELEMENT;
 42:     n_loc[3]  = 0;
 43:     break;
 44:   case 3:
 45:     PetscCall(DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 2, 2, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 2, 3, 1, 2, DMSTAG_STENCIL_BOX, 1, NULL, NULL, NULL, &dm));
 46:     n_loc[0]  = 1;
 47:     loc[0][0] = DMSTAG_BACK_DOWN_LEFT;
 48:     n_loc[1]  = 3;
 49:     loc[1][0] = DMSTAG_DOWN_LEFT;
 50:     loc[1][1] = DMSTAG_BACK_LEFT;
 51:     loc[1][2] = DMSTAG_BACK_DOWN;
 52:     n_loc[2]  = 3;
 53:     loc[2][0] = DMSTAG_LEFT;
 54:     loc[2][1] = DMSTAG_DOWN;
 55:     loc[2][2] = DMSTAG_BACK;
 56:     n_loc[3]  = 1;
 57:     loc[3][0] = DMSTAG_ELEMENT;
 58:     break;
 59:   default:
 60:     SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "No support for dimension %" PetscInt_FMT, dim);
 61:   }

 63:   PetscCall(DMSetFromOptions(dm));
 64:   PetscCall(DMSetUp(dm));
 65:   if (coords) PetscCall(DMStagSetUniformCoordinatesProduct(dm, -1.0, 1.0, -2.0, 2.0, -3.0, 3.0));

 67:   PetscCall(DMCreateGlobalVector(dm, &x));
 68:   PetscCall(VecSet(x, 1.2345));

 70:   PetscCall(DMStagGetDOF(dm, &dof[0], &dof[1], &dof[2], &dof[3]));
 71:   for (PetscInt stratum = 0; stratum < dim + 1; ++stratum) {
 72:     for (PetscInt i_loc = 0; i_loc < n_loc[stratum]; ++i_loc) {
 73:       // Extract 3 components, padding or truncating
 74:       {
 75:         DM  da;
 76:         Vec x_da;

 78:         PetscCall(DMStagVecSplitToDMDA(dm, x, loc[stratum][i_loc], -3, &da, &x_da));
 79:         PetscCall(DMView(da, PETSC_VIEWER_STDOUT_WORLD));
 80:         PetscCall(VecView(x_da, PETSC_VIEWER_STDOUT_WORLD));
 81:         PetscCall(DMDestroy(&da));
 82:         PetscCall(VecDestroy(&x_da));
 83:       }

 85:       // Extract individual components
 86:       for (PetscInt c = 0; c < dof[stratum]; ++c) {
 87:         DM  da;
 88:         Vec x_da;

 90:         PetscCall(DMStagVecSplitToDMDA(dm, x, loc[stratum][i_loc], c, &da, &x_da));
 91:         PetscCall(DMView(da, PETSC_VIEWER_STDOUT_WORLD));
 92:         PetscCall(VecView(x_da, PETSC_VIEWER_STDOUT_WORLD));
 93:         PetscCall(DMDestroy(&da));
 94:         PetscCall(VecDestroy(&x_da));
 95:       }
 96:     }
 97:   }

 99:   PetscCall(VecDestroy(&x));
100:   PetscCall(DMDestroy(&dm));
101:   PetscCall(PetscFinalize());
102:   return 0;
103: }

105: /*TEST

107:    test:
108:       requires: !complex
109:       args: -dim {{1 2 3}separate output} -coords {{true false}separate output}

111: TEST*/