1. petsc
  2. PETSc
  3. petsc

Commits

Matt Knepley  committed 22d884e

DMPlex: Now geometry tests can take an input mesh
- Use Exodus format
- Input all mesh data in options
- Made builder runs for 3 hexes

  • Participants
  • Parent commits 7565392
  • Branches master

Comments (0)

Files changed (5)

File config/builder.py

View file
                         'src/dm/impls/plex/examples/tests/ex8': [{'numProcs': 1, 'args': '-dm_view ::ascii_info_detail'},
                                                                  {'numProcs': 1, 'args': '-interpolate -dm_view ::ascii_info_detail'},
                                                                  {'numProcs': 1, 'args': '-transform'},
-                                                                 {'numProcs': 1, 'args': '-interpolate -transform'}],
+                                                                 {'numProcs': 1, 'args': '-interpolate -transform'},
+                                                                 {'numProcs': 1, 'args': '-run_type file -filename %(meshes)s/simpleblock-100.exo -dm_view ::ascii_info_detail -v0 -1.5,-0.5,0.5,-0.5,-0.5,0.5,0.5,-0.5,0.5 -J 0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0 -invJ 0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0 -detJ 0.125,0.125,0.125', 'requires': ['exodusii']},
+                                                                 {'numProcs': 1, 'args': '-interpolate -run_type file -filename %(meshes)s/simpleblock-100.exo -dm_view ::ascii_info_detail -v0 -1.5,-0.5,0.5,-0.5,-0.5,0.5,0.5,-0.5,0.5 -J 0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0 -invJ 0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0 -detJ 0.125,0.125,0.125 -centroid -1.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0 -normal 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 -vol 1.0,1.0,1.0', 'requires': ['exodusii']}],
                         'src/dm/impls/plex/examples/tests/ex1f90': [{'numProcs': 1, 'args': ''}],
                         'src/dm/impls/plex/examples/tests/ex2f90': [{'numProcs': 1, 'args': ''}],
                         'src/dm/impls/plex/examples/tutorials/ex1': [{'numProcs': 1, 'args': ''},

File share/petsc/datafiles/meshes/simpleblock-100.exo

Binary file added.

File src/dm/impls/plex/examples/tests/ex8.c

View file
 static char help[] = "Tests for cell geometry\n\n";
 
-/* TODO
-*/
-
 #include <petscdmplex.h>
+#if defined(PETSC_HAVE_EXODUSII)
+#include <exodusII.h>
+#endif
 
 typedef enum {RUN_REFERENCE, RUN_FILE} RunType;
 
   char      filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */
   PetscBool interpolate;                  /* Interpolate the mesh */
   PetscBool transform;                    /* Use random coordinate transformations */
+  /* Data for input meshes */
+  PetscReal *v0, *J, *invJ, *detJ;        /* FEM data */
+  PetscReal *centroid, *normal, *vol;     /* FVM data */
 } AppCtx;
 
 #undef __FUNCT__
+#define __FUNCT__ "ReadMesh"
+PetscErrorCode ReadMesh(MPI_Comm comm, const char *filename, AppCtx *user, DM *dm)
+{
+  int            CPU_word_size = 0, IO_word_size = 0, exoid = -1;
+  float          version;
+  PetscMPIInt    rank;
+  PetscErrorCode ierr;
+
+  PetscFunctionBegin;
+  ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
+#if defined(PETSC_HAVE_EXODUSII)
+  if (!rank) {
+    exoid = ex_open(filename, EX_READ, &CPU_word_size, &IO_word_size, &version);
+    if (exoid <= 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_LIB, "ex_open(\"%s\",...) did not return a valid file ID", filename);
+  }
+  ierr = DMPlexCreateExodus(comm, exoid, PETSC_FALSE, dm);CHKERRQ(ierr);
+  if (!rank) {ierr = ex_close(exoid);CHKERRQ(ierr);}
+#else
+  SETERRQ(comm, PETSC_ERR_SUP, "Loading meshes requires ExodusII support. Reconfigure using --download-exodusii");
+#endif
+  if (user->interpolate) {
+    DM interpolatedMesh = NULL;
+
+    ierr = DMPlexInterpolate(*dm, &interpolatedMesh);CHKERRQ(ierr);
+    ierr = DMPlexCopyCoordinates(*dm, interpolatedMesh);CHKERRQ(ierr);
+    ierr = DMDestroy(dm);CHKERRQ(ierr);
+    *dm  = interpolatedMesh;
+  }
+  ierr = PetscObjectSetName((PetscObject) *dm, "Input Mesh");CHKERRQ(ierr);
+  ierr = DMSetFromOptions(*dm);CHKERRQ(ierr);
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
 #define __FUNCT__ "ProcessOptions"
 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 {
-  const char    *runTypes[2] = {"full", "test"};
+  const char    *runTypes[2] = {"reference", "file"};
   PetscInt       run;
   PetscErrorCode ierr;
 
   ierr = PetscOptionsString("-filename", "The mesh file", "ex8.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr);
   ierr = PetscOptionsBool("-interpolate", "Interpolate the mesh", "ex8.c", options->interpolate, &options->interpolate, NULL);CHKERRQ(ierr);
   ierr = PetscOptionsBool("-transform", "Use random transforms", "ex8.c", options->transform, &options->transform, NULL);CHKERRQ(ierr);
+
+  if (options->runType == RUN_FILE) {
+    PetscInt  dim, cStart, cEnd, numCells, n;
+    PetscBool flag;
+
+    ierr = ReadMesh(PETSC_COMM_WORLD, options->filename, options, &options->dm);CHKERRQ(ierr);
+    ierr = DMPlexGetDimension(options->dm, &dim);CHKERRQ(ierr);
+    ierr = DMPlexGetHeightStratum(options->dm, 0, &cStart, &cEnd);CHKERRQ(ierr);
+    numCells = cEnd-cStart;
+    ierr = PetscMalloc7(numCells*dim,PetscReal,&options->v0,numCells*dim*dim,PetscReal,&options->J,numCells*dim*dim,PetscReal,&options->invJ,numCells,PetscReal,&options->detJ,
+                        numCells*dim,PetscReal,&options->centroid,numCells*dim,PetscReal,&options->normal,numCells,PetscReal,&options->vol);CHKERRQ(ierr);
+    n    = numCells*dim;
+    ierr = PetscOptionsRealArray("-v0", "Input v0 for each cell", "ex8.c", options->v0, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of v0 %d should be %d", n, numCells*dim);
+    n    = numCells*dim*dim;
+    ierr = PetscOptionsRealArray("-J", "Input Jacobian for each cell", "ex8.c", options->J, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells*dim*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of J %d should be %d", n, numCells*dim*dim);
+    n    = numCells*dim*dim;
+    ierr = PetscOptionsRealArray("-invJ", "Input inverse Jacobian for each cell", "ex8.c", options->invJ, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells*dim*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of invJ %d should be %d", n, numCells*dim*dim);
+    n    = numCells;
+    ierr = PetscOptionsRealArray("-detJ", "Input Jacobian determinant for each cell", "ex8.c", options->detJ, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of detJ %d should be %d", n, numCells);
+    n    = numCells*dim;
+    ierr = PetscOptionsRealArray("-centroid", "Input centroid for each cell", "ex8.c", options->centroid, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of centroid %d should be %d", n, numCells*dim);
+    n    = numCells*dim;
+    ierr = PetscOptionsRealArray("-normal", "Input normal for each cell", "ex8.c", options->normal, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of normal %d should be %d", n, numCells*dim);
+    n    = numCells;
+    ierr = PetscOptionsRealArray("-vol", "Input volume for each cell", "ex8.c", options->vol, &n, &flag);CHKERRQ(ierr);
+    if (flag && n != numCells) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of vol %d should be %d", n, numCells);
+  }
   ierr = PetscOptionsEnd();
 
   if (options->transform) {ierr = PetscPrintf(comm, "Using random transforms");CHKERRQ(ierr);}
     ierr = TestQuadrilateral(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr);
     ierr = TestTetrahedron(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr);
     ierr = TestHexahedron(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr);
+  } else if (user.runType == RUN_FILE) {
+    PetscInt dim, cStart, cEnd, c;
+
+    ierr = DMPlexGetDimension(user.dm, &dim);CHKERRQ(ierr);
+    ierr = DMPlexGetHeightStratum(user.dm, 0, &cStart, &cEnd);CHKERRQ(ierr);
+    for (c = 0; c < cEnd-cStart; ++c) {
+      ierr = CheckFEMGeometry(user.dm, c+cStart, dim, &user.v0[c*dim], &user.J[c*dim*dim], &user.invJ[c*dim*dim], user.detJ[c]);CHKERRQ(ierr);
+      if (user.interpolate) {ierr = CheckFVMGeometry(user.dm, c+cStart, dim, &user.centroid[c*dim], &user.normal[c*dim], user.vol[c]);CHKERRQ(ierr);}
+    }
+    ierr = PetscFree7(user.v0,user.J,user.invJ,user.detJ,user.centroid,user.normal,user.vol);CHKERRQ(ierr);
   }
   ierr = PetscFinalize();
   return 0;

File src/dm/impls/plex/examples/tests/output/ex8_4.out

View file
+Mesh 'Input Mesh':
+Max sizes cone: 8 support: 2
+orientation is missing
+cap --> base:
+[0]: 3 ----> 0
+[0]: 4 ----> 0
+[0]: 5 ----> 0
+[0]: 6 ----> 0
+[0]: 7 ----> 0
+[0]: 7 ----> 1
+[0]: 8 ----> 0
+[0]: 8 ----> 1
+[0]: 9 ----> 0
+[0]: 9 ----> 1
+[0]: 10 ----> 0
+[0]: 10 ----> 1
+[0]: 11 ----> 1
+[0]: 11 ----> 2
+[0]: 12 ----> 1
+[0]: 12 ----> 2
+[0]: 13 ----> 1
+[0]: 13 ----> 2
+[0]: 14 ----> 1
+[0]: 14 ----> 2
+[0]: 15 ----> 2
+[0]: 16 ----> 2
+[0]: 17 ----> 2
+[0]: 18 ----> 2
+base <-- cap:
+[0]: 0 <---- 3 (0)
+[0]: 0 <---- 6 (0)
+[0]: 0 <---- 5 (0)
+[0]: 0 <---- 4 (0)
+[0]: 0 <---- 7 (0)
+[0]: 0 <---- 8 (0)
+[0]: 0 <---- 9 (0)
+[0]: 0 <---- 10 (0)
+[0]: 1 <---- 7 (0)
+[0]: 1 <---- 10 (0)
+[0]: 1 <---- 9 (0)
+[0]: 1 <---- 8 (0)
+[0]: 1 <---- 11 (0)
+[0]: 1 <---- 12 (0)
+[0]: 1 <---- 13 (0)
+[0]: 1 <---- 14 (0)
+[0]: 2 <---- 11 (0)
+[0]: 2 <---- 14 (0)
+[0]: 2 <---- 13 (0)
+[0]: 2 <---- 12 (0)
+[0]: 2 <---- 15 (0)
+[0]: 2 <---- 16 (0)
+[0]: 2 <---- 17 (0)
+[0]: 2 <---- 18 (0)
+coordinates with 1 fields
+  field 0 with 3 components
+Process 0:
+  (   3) dim  3 offset   0 -1.5 -0.5 0.5
+  (   4) dim  3 offset   3 -1.5 -0.5 -0.5
+  (   5) dim  3 offset   6 -1.5 0.5 -0.5
+  (   6) dim  3 offset   9 -1.5 0.5 0.5
+  (   7) dim  3 offset  12 -0.5 -0.5 0.5
+  (   8) dim  3 offset  15 -0.5 -0.5 -0.5
+  (   9) dim  3 offset  18 -0.5 0.5 -0.5
+  (  10) dim  3 offset  21 -0.5 0.5 0.5
+  (  11) dim  3 offset  24 0.5 -0.5 0.5
+  (  12) dim  3 offset  27 0.5 -0.5 -0.5
+  (  13) dim  3 offset  30 0.5 0.5 -0.5
+  (  14) dim  3 offset  33 0.5 0.5 0.5
+  (  15) dim  3 offset  36 1.5 -0.5 0.5
+  (  16) dim  3 offset  39 1.5 -0.5 -0.5
+  (  17) dim  3 offset  42 1.5 0.5 -0.5
+  (  18) dim  3 offset  45 1.5 0.5 0.5

File src/dm/impls/plex/examples/tests/output/ex8_5.out

View file
+Mesh 'Input Mesh':
+Max sizes cone: 6 support: 4
+orientation is missing
+cap --> base:
+[0]: 3 ----> 35
+[0]: 3 ----> 38
+[0]: 3 ----> 44
+[0]: 4 ----> 37
+[0]: 4 ----> 38
+[0]: 4 ----> 43
+[0]: 5 ----> 36
+[0]: 5 ----> 37
+[0]: 5 ----> 46
+[0]: 6 ----> 35
+[0]: 6 ----> 36
+[0]: 6 ----> 45
+[0]: 7 ----> 39
+[0]: 7 ----> 42
+[0]: 7 ----> 44
+[0]: 7 ----> 52
+[0]: 8 ----> 39
+[0]: 8 ----> 40
+[0]: 8 ----> 43
+[0]: 8 ----> 51
+[0]: 9 ----> 40
+[0]: 9 ----> 41
+[0]: 9 ----> 46
+[0]: 9 ----> 54
+[0]: 10 ----> 41
+[0]: 10 ----> 42
+[0]: 10 ----> 45
+[0]: 10 ----> 53
+[0]: 11 ----> 47
+[0]: 11 ----> 50
+[0]: 11 ----> 52
+[0]: 11 ----> 60
+[0]: 12 ----> 47
+[0]: 12 ----> 48
+[0]: 12 ----> 51
+[0]: 12 ----> 59
+[0]: 13 ----> 48
+[0]: 13 ----> 49
+[0]: 13 ----> 54
+[0]: 13 ----> 62
+[0]: 14 ----> 49
+[0]: 14 ----> 50
+[0]: 14 ----> 53
+[0]: 14 ----> 61
+[0]: 15 ----> 55
+[0]: 15 ----> 58
+[0]: 15 ----> 60
+[0]: 16 ----> 55
+[0]: 16 ----> 56
+[0]: 16 ----> 59
+[0]: 17 ----> 56
+[0]: 17 ----> 57
+[0]: 17 ----> 62
+[0]: 18 ----> 57
+[0]: 18 ----> 58
+[0]: 18 ----> 61
+[0]: 19 ----> 0
+[0]: 20 ----> 0
+[0]: 20 ----> 1
+[0]: 21 ----> 0
+[0]: 22 ----> 0
+[0]: 23 ----> 0
+[0]: 24 ----> 0
+[0]: 25 ----> 1
+[0]: 25 ----> 2
+[0]: 26 ----> 1
+[0]: 27 ----> 1
+[0]: 28 ----> 1
+[0]: 29 ----> 1
+[0]: 30 ----> 2
+[0]: 31 ----> 2
+[0]: 32 ----> 2
+[0]: 33 ----> 2
+[0]: 34 ----> 2
+[0]: 35 ----> 19
+[0]: 35 ----> 24
+[0]: 36 ----> 19
+[0]: 36 ----> 22
+[0]: 37 ----> 19
+[0]: 37 ----> 23
+[0]: 38 ----> 19
+[0]: 38 ----> 21
+[0]: 39 ----> 20
+[0]: 39 ----> 21
+[0]: 39 ----> 26
+[0]: 40 ----> 20
+[0]: 40 ----> 23
+[0]: 40 ----> 28
+[0]: 41 ----> 20
+[0]: 41 ----> 22
+[0]: 41 ----> 27
+[0]: 42 ----> 20
+[0]: 42 ----> 24
+[0]: 42 ----> 29
+[0]: 43 ----> 21
+[0]: 43 ----> 23
+[0]: 44 ----> 21
+[0]: 44 ----> 24
+[0]: 45 ----> 22
+[0]: 45 ----> 24
+[0]: 46 ----> 22
+[0]: 46 ----> 23
+[0]: 47 ----> 25
+[0]: 47 ----> 26
+[0]: 47 ----> 31
+[0]: 48 ----> 25
+[0]: 48 ----> 28
+[0]: 48 ----> 33
+[0]: 49 ----> 25
+[0]: 49 ----> 27
+[0]: 49 ----> 32
+[0]: 50 ----> 25
+[0]: 50 ----> 29
+[0]: 50 ----> 34
+[0]: 51 ----> 26
+[0]: 51 ----> 28
+[0]: 52 ----> 26
+[0]: 52 ----> 29
+[0]: 53 ----> 27
+[0]: 53 ----> 29
+[0]: 54 ----> 27
+[0]: 54 ----> 28
+[0]: 55 ----> 30
+[0]: 55 ----> 31
+[0]: 56 ----> 30
+[0]: 56 ----> 33
+[0]: 57 ----> 30
+[0]: 57 ----> 32
+[0]: 58 ----> 30
+[0]: 58 ----> 34
+[0]: 59 ----> 31
+[0]: 59 ----> 33
+[0]: 60 ----> 31
+[0]: 60 ----> 34
+[0]: 61 ----> 32
+[0]: 61 ----> 34
+[0]: 62 ----> 32
+[0]: 62 ----> 33
+base <-- cap:
+[0]: 0 <---- 19 (0)
+[0]: 0 <---- 20 (0)
+[0]: 0 <---- 21 (0)
+[0]: 0 <---- 22 (0)
+[0]: 0 <---- 23 (0)
+[0]: 0 <---- 24 (0)
+[0]: 1 <---- 20 (-4)
+[0]: 1 <---- 25 (0)
+[0]: 1 <---- 26 (0)
+[0]: 1 <---- 27 (0)
+[0]: 1 <---- 28 (0)
+[0]: 1 <---- 29 (0)
+[0]: 2 <---- 25 (-4)
+[0]: 2 <---- 30 (0)
+[0]: 2 <---- 31 (0)
+[0]: 2 <---- 32 (0)
+[0]: 2 <---- 33 (0)
+[0]: 2 <---- 34 (0)
+[0]: 19 <---- 35 (0)
+[0]: 19 <---- 36 (0)
+[0]: 19 <---- 37 (0)
+[0]: 19 <---- 38 (0)
+[0]: 20 <---- 39 (0)
+[0]: 20 <---- 40 (0)
+[0]: 20 <---- 41 (0)
+[0]: 20 <---- 42 (0)
+[0]: 21 <---- 38 (-2)
+[0]: 21 <---- 43 (0)
+[0]: 21 <---- 39 (-2)
+[0]: 21 <---- 44 (0)
+[0]: 22 <---- 36 (-2)
+[0]: 22 <---- 45 (0)
+[0]: 22 <---- 41 (-2)
+[0]: 22 <---- 46 (0)
+[0]: 23 <---- 37 (-2)
+[0]: 23 <---- 46 (-2)
+[0]: 23 <---- 40 (-2)
+[0]: 23 <---- 43 (-2)
+[0]: 24 <---- 44 (-2)
+[0]: 24 <---- 42 (-2)
+[0]: 24 <---- 45 (-2)
+[0]: 24 <---- 35 (-2)
+[0]: 25 <---- 47 (0)
+[0]: 25 <---- 48 (0)
+[0]: 25 <---- 49 (0)
+[0]: 25 <---- 50 (0)
+[0]: 26 <---- 39 (0)
+[0]: 26 <---- 51 (0)
+[0]: 26 <---- 47 (-2)
+[0]: 26 <---- 52 (0)
+[0]: 27 <---- 41 (0)
+[0]: 27 <---- 53 (0)
+[0]: 27 <---- 49 (-2)
+[0]: 27 <---- 54 (0)
+[0]: 28 <---- 40 (0)
+[0]: 28 <---- 54 (-2)
+[0]: 28 <---- 48 (-2)
+[0]: 28 <---- 51 (-2)
+[0]: 29 <---- 52 (-2)
+[0]: 29 <---- 50 (-2)
+[0]: 29 <---- 53 (-2)
+[0]: 29 <---- 42 (0)
+[0]: 30 <---- 55 (0)
+[0]: 30 <---- 56 (0)
+[0]: 30 <---- 57 (0)
+[0]: 30 <---- 58 (0)
+[0]: 31 <---- 47 (0)
+[0]: 31 <---- 59 (0)
+[0]: 31 <---- 55 (-2)
+[0]: 31 <---- 60 (0)
+[0]: 32 <---- 49 (0)
+[0]: 32 <---- 61 (0)
+[0]: 32 <---- 57 (-2)
+[0]: 32 <---- 62 (0)
+[0]: 33 <---- 48 (0)
+[0]: 33 <---- 62 (-2)
+[0]: 33 <---- 56 (-2)
+[0]: 33 <---- 59 (-2)
+[0]: 34 <---- 60 (-2)
+[0]: 34 <---- 58 (-2)
+[0]: 34 <---- 61 (-2)
+[0]: 34 <---- 50 (0)
+[0]: 35 <---- 3 (0)
+[0]: 35 <---- 6 (0)
+[0]: 36 <---- 6 (0)
+[0]: 36 <---- 5 (0)
+[0]: 37 <---- 5 (0)
+[0]: 37 <---- 4 (0)
+[0]: 38 <---- 4 (0)
+[0]: 38 <---- 3 (0)
+[0]: 39 <---- 7 (0)
+[0]: 39 <---- 8 (0)
+[0]: 40 <---- 8 (0)
+[0]: 40 <---- 9 (0)
+[0]: 41 <---- 9 (0)
+[0]: 41 <---- 10 (0)
+[0]: 42 <---- 10 (0)
+[0]: 42 <---- 7 (0)
+[0]: 43 <---- 4 (0)
+[0]: 43 <---- 8 (0)
+[0]: 44 <---- 7 (0)
+[0]: 44 <---- 3 (0)
+[0]: 45 <---- 6 (0)
+[0]: 45 <---- 10 (0)
+[0]: 46 <---- 9 (0)
+[0]: 46 <---- 5 (0)
+[0]: 47 <---- 11 (0)
+[0]: 47 <---- 12 (0)
+[0]: 48 <---- 12 (0)
+[0]: 48 <---- 13 (0)
+[0]: 49 <---- 13 (0)
+[0]: 49 <---- 14 (0)
+[0]: 50 <---- 14 (0)
+[0]: 50 <---- 11 (0)
+[0]: 51 <---- 8 (0)
+[0]: 51 <---- 12 (0)
+[0]: 52 <---- 11 (0)
+[0]: 52 <---- 7 (0)
+[0]: 53 <---- 10 (0)
+[0]: 53 <---- 14 (0)
+[0]: 54 <---- 13 (0)
+[0]: 54 <---- 9 (0)
+[0]: 55 <---- 15 (0)
+[0]: 55 <---- 16 (0)
+[0]: 56 <---- 16 (0)
+[0]: 56 <---- 17 (0)
+[0]: 57 <---- 17 (0)
+[0]: 57 <---- 18 (0)
+[0]: 58 <---- 18 (0)
+[0]: 58 <---- 15 (0)
+[0]: 59 <---- 12 (0)
+[0]: 59 <---- 16 (0)
+[0]: 60 <---- 15 (0)
+[0]: 60 <---- 11 (0)
+[0]: 61 <---- 14 (0)
+[0]: 61 <---- 18 (0)
+[0]: 62 <---- 17 (0)
+[0]: 62 <---- 13 (0)
+coordinates with 1 fields
+  field 0 with 3 components
+Process 0:
+  (   3) dim  3 offset   0 -1.5 -0.5 0.5
+  (   4) dim  3 offset   3 -1.5 -0.5 -0.5
+  (   5) dim  3 offset   6 -1.5 0.5 -0.5
+  (   6) dim  3 offset   9 -1.5 0.5 0.5
+  (   7) dim  3 offset  12 -0.5 -0.5 0.5
+  (   8) dim  3 offset  15 -0.5 -0.5 -0.5
+  (   9) dim  3 offset  18 -0.5 0.5 -0.5
+  (  10) dim  3 offset  21 -0.5 0.5 0.5
+  (  11) dim  3 offset  24 0.5 -0.5 0.5
+  (  12) dim  3 offset  27 0.5 -0.5 -0.5
+  (  13) dim  3 offset  30 0.5 0.5 -0.5
+  (  14) dim  3 offset  33 0.5 0.5 0.5
+  (  15) dim  3 offset  36 1.5 -0.5 0.5
+  (  16) dim  3 offset  39 1.5 -0.5 -0.5
+  (  17) dim  3 offset  42 1.5 0.5 -0.5
+  (  18) dim  3 offset  45 1.5 0.5 0.5