1: #include <petscsnes.h>
  2: #include <petscdmda.h>

  4: int main(int argc, char *argv[]) {
  5:   DM              da, daX, daY;
  6:   DMDALocalInfo     info;
  7:   MPI_Comm        commX, commY;
  8:   Vec             basisX, basisY;
  9:   PetscScalar   **arrayX, **arrayY;
 10:   const PetscInt *lx, *ly;
 11:   PetscInt        M = 3, N = 3;
 12:   PetscInt        p = 1;
 13:   PetscInt        numGP = 3;
 14:   PetscInt        dof = 2*(p+1)*numGP;
 15:   PetscMPIInt     rank, subsize, subrank;
 16:   PetscErrorCode  ierr;

 18:   PetscInitialize(&argc,&argv,0,0);
 19:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 20:   /* Create 2D DMDA */
 21:   DMDACreate2d(PETSC_COMM_WORLD, DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_NONE, DMDA_STENCIL_STAR, M, N, PETSC_DECIDE, PETSC_DECIDE, 1, 1, PETSC_NULL, PETSC_NULL, &da);
 22:   /* Create 1D DMDAs along two directions */
 23:   DMDAGetOwnershipRanges(da, &lx, &ly, PETSC_NULL);
 24:   DMDAGetLocalInfo(da, &info);
 25:   DMDAGetProcessorSubsets(da, DMDA_X, &commX);
 26:   DMDAGetProcessorSubsets(da, DMDA_Y, &commY);
 27:   MPI_Comm_size(commX, &subsize);
 28:   MPI_Comm_rank(commX, &subrank);
 29:   PetscPrintf(PETSC_COMM_SELF, "[%d]X subrank: %d subsize: %d\n", rank, subrank, subsize);
 30:   MPI_Comm_size(commY, &subsize);
 31:   MPI_Comm_rank(commY, &subrank);
 32:   PetscPrintf(PETSC_COMM_SELF, "[%d]Y subrank: %d subsize: %d\n", rank, subrank, subsize);
 33:   DMDACreate1d(commX, DMDA_BOUNDARY_NONE, M, dof, 1, lx, &daX);
 34:   DMDACreate1d(commY, DMDA_BOUNDARY_NONE, N, dof, 1, ly, &daY);
 35:   /* Create 1D vectors for basis functions */
 36:   DMGetGlobalVector(daX, &basisX);
 37:   DMGetGlobalVector(daY, &basisY);
 38:   /* Extract basis functions */
 39:   DMDAVecGetArrayDOF(daX, basisX, &arrayX);
 40:   DMDAVecGetArrayDOF(daY, basisY, &arrayY);
 41:   //arrayX[i][ndof];
 42:   //arrayY[j][ndof];
 43:   DMDAVecRestoreArrayDOF(daX, basisX, &arrayX);
 44:   DMDAVecRestoreArrayDOF(daY, basisY, &arrayY);
 45:   /* Return basis vectors */
 46:   DMRestoreGlobalVector(daX, &basisX);
 47:   DMRestoreGlobalVector(daY, &basisY);
 48:   /* Cleanup */
 49:   DMDestroy(&daX);
 50:   DMDestroy(&daY);
 51:   DMDestroy(&da);
 52:   PetscFinalize();
 53:   return 0;
 54: }