Setting Dirichlet BCs from field: Input array needs to be sorted

[baagaard-usgs]

I am reading Cubit and Gmsh files with boundaries specified by faces and when I try to set boundary values for Dirichlet BCs I get the error below:

[0]PETSC ERROR: Invalid argument
[0]PETSC ERROR: Input array needs to be sorted
-- snip –
[0]PETSC ERROR: #1 PetscFindInt() at /Users/baagaard/software/unix/petsc-dev/src/sys/utils/sorti.c:513
[0]PETSC ERROR: #2 DMGetEnclosurePoint() at /Users/baagaard/software/unix/petsc-dev/src/dm/impls/plex/plexsubmesh.c:4476
[0]PETSC ERROR: #3 DMProjectPoint_BdField_Private() at /Users/baagaard/software/unix/petsc-dev/src/dm/impls/plex/plexproject.c:375
[0]PETSC ERROR: #4 DMProjectPoint_Private() at /Users/baagaard/software/unix/petsc-dev/src/dm/impls/plex/plexproject.c:488
[0]PETSC ERROR: #5 DMProjectLocal_Generic_Plex() at /Users/baagaard/software/unix/petsc-dev/src/dm/impls/plex/plexproject.c:970
[0]PETSC ERROR: #6 DMProjectBdFieldLabelLocal_Plex() at /Users/baagaard/software/unix/petsc-dev/src/dm/impls/plex/plexproject.c:1060
[0]PETSC ERROR: #7 DMProjectBdFieldLabelLocal() at /Users/baagaard/software/unix/petsc-dev/src/dm/interface/dm.c:8508
[0]PETSC ERROR: #8 DMPlexInsertBoundaryValuesEssentialBdField() at /Users/baagaard/software/unix/petsc-dev/src/dm/impls/plex/plexfem.c:994

The MMS tests are working, so there is something different between using a function for setting values as opposed to a field (projection).

This works:

    err = DMPlexLabelAddCells(dmField, dmLabel);PYLITH_CHECK_ERROR(err);
    err = DMPlexInsertBoundaryValuesEssential(dmField, t, fieldIndex, numConstrained, &constraint._constrainedDOF[0], dmLabel, 1,
                                              &constraint._labelValue, fn, context, field->getLocalVector());PYLITH_CHECK_ERROR(err);
    err = DMPlexLabelClearCells(dmField, dmLabel);PYLITH_CHECK_ERROR(err);

This generates the error:

    err = DMPlexLabelAddFaceCells(dmSoln, dmLabel);PYLITH_CHECK_ERROR(err);
    err = DMPlexInsertBoundaryValuesEssentialBdField(dmSoln, t, solution->getLocalVector(), fieldIndex,
                                                     numConstrained, &_constrainedDOF[0], dmLabel, 1, &_labelValue,
                                                     _kernelConstraint, context, solution->getLocalVector());PYLITH_CHECK_ERROR(err);
    err = DMPlexLabelClearCells(dmSoln, dmLabel);PYLITH_CHECK_ERROR(err);

Reproduce

• PETSc branch: knepley/pylith
• PyLith fork: https://github.com/baagaard-usgs/pylith.git
• PyLith branch: feature-meshio-hdf5

# BUILD and then run the full-scale test
cd tests/fullscale/linearelasticity/nofaults-2d
make export-data
pylith axialdisp.cfg axialdisp_quad.cfg

[baagaard-usgs]

Matt says DMPlexFilter() does not create a submesh with the domain cells hanging off it. As a result, when we later try to project from the submesh to the domain mesh to set boundary conditions, the projection fails.

Matt says we should go back to using DMPlexCreateSubmesh() like we have been doing when we use vertices to specify the BCs.

[baagaard-usgs]

Resolved by using DMPlexCreateSubmesh().