1: static char help[] = "Create a Plex sphere from quads and create a P1 section\n\n";

  3:  #include <petscdmplex.h>

  5: typedef struct {
  6:   PetscInt  dim;     /* Topological problem dimension */
  7:   PetscBool simplex; /* Mesh with simplices */
  8: } AppCtx;

 10: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 11: {

 15:   options->dim     = 2;
 16:   options->simplex = PETSC_FALSE;

 18:   PetscOptionsBegin(comm, "", "Sphere Mesh Options", "DMPLEX");
 19:   PetscOptionsRangeInt("-dim", "Problem dimension", "ex7.c", options->dim, &options->dim, NULL,1,3);
 20:   PetscOptionsBool("-simplex", "Use simplices, or tensor product cells", "ex7.c", options->simplex, &options->simplex, NULL);
 21:   PetscOptionsEnd();
 22:   return(0);
 23: }

 25: static PetscErrorCode ProjectToUnitSphere(DM dm)
 26: {
 27:   Vec            coordinates;
 28:   PetscScalar   *coords;
 29:   PetscInt       Nv, v, bs, dim, d;

 33:   DMGetCoordinatesLocal(dm, &coordinates);
 34:   VecGetLocalSize(coordinates, &Nv);
 35:   VecGetBlockSize(coordinates, &bs);
 36:   DMGetCoordinateDim(dm, &dim);
 37:   if (dim != bs) SETERRQ2(PetscObjectComm((PetscObject)dm),PETSC_ERR_PLIB,"Coordinate bs %D does not match dim %D",bs,dim);
 38:   Nv  /= dim;
 39:   VecGetArray(coordinates, &coords);
 40:   for (v = 0; v < Nv; ++v) {
 41:     PetscReal r = 0.0;

 43:     for (d = 0; d < dim; ++d) r += PetscSqr(PetscRealPart(coords[v*dim+d]));
 44:     r = PetscSqrtReal(r);
 45:     for (d = 0; d < dim; ++d) coords[v*dim+d] /= r;
 46:   }
 47:   VecRestoreArray(coordinates, &coords);
 48:   return(0);
 49: }

 51: static PetscErrorCode SetupSection(DM dm)
 52: {
 53:   PetscSection   s;
 54:   PetscInt       vStart, vEnd, v;

 58:   DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);
 59:   PetscSectionCreate(PetscObjectComm((PetscObject) dm), &s);
 60:   PetscSectionSetNumFields(s, 1);
 61:   PetscSectionSetFieldComponents(s, 0, 1);
 62:   PetscSectionSetChart(s, vStart, vEnd);
 63:   for (v = vStart; v < vEnd; ++v) {
 64:     PetscSectionSetDof(s, v, 1);
 65:     PetscSectionSetFieldDof(s, v, 0, 1);
 66:   }
 67:   PetscSectionSetUp(s);
 68:   DMSetLocalSection(dm, s);
 69:   PetscSectionDestroy(&s);
 70:   return(0);
 71: }

 73: int main(int argc, char **argv)
 74: {
 75:   DM             dm;
 76:   Vec            u;
 77:   AppCtx         ctx;

 80:   PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
 81:   ProcessOptions(PETSC_COMM_WORLD, &ctx);
 82:   DMPlexCreateSphereMesh(PETSC_COMM_WORLD, ctx.dim, ctx.simplex, &dm);
 83:   PetscObjectSetName((PetscObject) dm, "Sphere");
 84:   /* Distribute mesh over processes */
 85:   {
 86:      DM dmDist = NULL;
 87:      PetscPartitioner part;
 88:      DMPlexGetPartitioner(dm, &part);
 89:      PetscPartitionerSetFromOptions(part);
 90:      DMPlexDistribute(dm, 0, NULL, &dmDist);
 91:      if (dmDist) {
 92:        DMDestroy(&dm);
 93:        dm  = dmDist;
 94:      }
 95:   }
 96:   DMSetFromOptions(dm);
 97:   ProjectToUnitSphere(dm);
 98:   DMViewFromOptions(dm, NULL, "-dm_view");
 99:   SetupSection(dm);
100:   DMGetGlobalVector(dm, &u);
101:   VecSet(u, 2);
102:   VecViewFromOptions(u, NULL, "-vec_view");
103:   DMRestoreGlobalVector(dm, &u);
104:   DMDestroy(&dm);
105:   PetscFinalize();
106:   return ierr;
107: }

109: /*TEST

111:   test:
112:     suffix: 2d_quad
113:     requires: !__float128
114:     args: -dm_view

116:   test:
117:     suffix: 2d_quad_parallel
118:     requires: !__float128
119:     args: -dm_view -petscpartitioner_type simple
120:     nsize: 2

122:   test:
123:     suffix: 2d_tri
124:     requires: !__float128
125:     args: -simplex -dm_view

127:   test:
128:     suffix: 2d_tri_parallel
129:     requires: !__float128
130:     args: -simplex -dm_view -petscpartitioner_type simple
131:     nsize: 2

133:   test:
134:     suffix: 3d_tri
135:     requires: !__float128
136:     args: -dim 3 -simplex -dm_view

138:   test:
139:     suffix: 3d_tri_parallel
140:     requires: !__float128
141:     args: -dim 3 -simplex -dm_view -petscpartitioner_type simple
142:     nsize: 2

144: TEST*/