Actual source code: ex6.c

  1: static char help[] = "Spot test DMStag->DMDA routines in 3d\n\n";
  2: #include <petscdm.h>
  3: #include <petscdmstag.h>

  5: int main(int argc, char **argv)
  6: {
  7:   DM  dm;
  8:   Vec vec;

 11:   PetscInitialize(&argc, &argv, (char *)0, help);
 12:   DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 4, 4, 4, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 3, 3, 3, 3, DMSTAG_STENCIL_STAR, 1, NULL, NULL, NULL, &dm);
 13:   DMSetFromOptions(dm);
 14:   DMSetUp(dm);
 15:   DMStagSetUniformCoordinatesProduct(dm, 0.0, 10.0, 0.0, 10.0, 0.0, 10.0);

 17:   DMCreateGlobalVector(dm, &vec);
 18:   VecSet(vec, 1.234);

 20:   /* All element values */
 21:   {
 22:     DM  da;
 23:     Vec vecda;
 24:     DMStagVecSplitToDMDA(dm, vec, DMSTAG_ELEMENT, -3, &da, &vecda);
 25:     DMDestroy(&da);
 26:     VecDestroy(&vecda);
 27:   }

 29:   /* Pad element values */
 30:   {
 31:     DM  da;
 32:     Vec vecda;
 33:     DMStagVecSplitToDMDA(dm, vec, DMSTAG_ELEMENT, -5, &da, &vecda);
 34:     DMDestroy(&da);
 35:     VecDestroy(&vecda);
 36:   }

 38:   /* 2 element values */
 39:   {
 40:     DM  da;
 41:     Vec vecda;
 42:     DMStagVecSplitToDMDA(dm, vec, DMSTAG_ELEMENT, -2, &da, &vecda);
 43:     DMDestroy(&da);
 44:     VecDestroy(&vecda);
 45:   }

 47:   /* One corner value */
 48:   {
 49:     DM  da;
 50:     Vec vecda;
 51:     DMStagVecSplitToDMDA(dm, vec, DMSTAG_FRONT_DOWN_LEFT, 2, &da, &vecda);
 52:     DMDestroy(&da);
 53:     VecDestroy(&vecda);
 54:   }

 56:   /* One edge value */
 57:   {
 58:     DM  da;
 59:     Vec vecda;
 60:     DMStagVecSplitToDMDA(dm, vec, DMSTAG_BACK_RIGHT, 1, &da, &vecda);
 61:     DMDestroy(&da);
 62:     VecDestroy(&vecda);
 63:   }

 65:   /* One face value */
 66:   {
 67:     DM  da;
 68:     Vec vecda;
 69:     DMStagVecSplitToDMDA(dm, vec, DMSTAG_DOWN, 0, &da, &vecda);
 70:     DMDestroy(&da);
 71:     VecDestroy(&vecda);
 72:   }

 74:   VecDestroy(&vec);
 75:   DMDestroy(&dm);
 76:   PetscFinalize();
 77:   return 0;
 78: }

 80: /*TEST

 82:    test:
 83:       suffix: 1
 84:       nsize: 12
 85:       args: -stag_ranks_x 2 -stag_ranks_y 3 -stag_ranks_z 2

 87: TEST*/