[cig-commits] r21312 - in short/3D/PyLith/trunk/libsrc/pylith: faults meshio topology
knepley at geodynamics.org
knepley at geodynamics.org
Thu Jan 31 20:55:00 PST 2013
Author: knepley
Date: 2013-01-31 20:55:00 -0800 (Thu, 31 Jan 2013)
New Revision: 21312
Modified:
short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.cc
short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.hh
short/3D/PyLith/trunk/libsrc/pylith/meshio/DataWriterHDF5.cc
short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc
short/3D/PyLith/trunk/libsrc/pylith/topology/SubMesh.cc
Log:
Just checks and more cohesivee work
Modified: short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.cc 2013-01-31 00:58:01 UTC (rev 21311)
+++ short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.cc 2013-02-01 04:55:00 UTC (rev 21312)
@@ -56,9 +56,16 @@
// Memory logging
ALE::MemoryLogger& logger = ALE::MemoryLogger::singleton();
logger.stagePush("SerialFaultCreation");
+ PetscErrorCode err;
faultMesh->coordsys(mesh);
+ DM dmMesh = mesh.dmMesh();
+ PetscInt dim, depth;
+ assert(dmMesh);
+ err = DMPlexGetDimension(dmMesh, &dim);CHECK_PETSC_ERROR(err);
+ err = DMPlexGetDepth(dmMesh, &depth);CHECK_PETSC_ERROR(err);
+
const ALE::Obj<SieveMesh>& sieveMesh = mesh.sieveMesh();
assert(!sieveMesh.isNull());
const ALE::Obj<SieveMesh::sieve_type>& sieve = sieveMesh->getSieve();
@@ -119,7 +126,26 @@
renumbering.clear();
// Convert fault to a DM
- {
+ if (depth == dim) {
+ DM subdm;
+ DMLabel label;
+ const char *labelName = groupField->getName().c_str();
+
+ // Put fault vertices in a label
+ err = DMPlexCreateLabel(dmMesh, labelName);CHECK_PETSC_ERROR(err);
+ err = DMPlexGetLabel(dmMesh, labelName, &label);CHECK_PETSC_ERROR(err);
+ const IntSection::chart_type& chart = groupField->getChart();
+ const IntSection::chart_type::const_iterator chartEnd = chart.end();
+ for(IntSection::chart_type::const_iterator c_iter = chart.begin(); c_iter != chartEnd; ++c_iter) {
+ if (groupField->getFiberDimension(*c_iter)) {
+ err = DMLabelSetValue(label, *c_iter, 0);CHECK_PETSC_ERROR(err);
+ }
+ }
+
+ err = DMPlexCreateSubmesh(dmMesh, labelName, "faultSurface", &subdm);CHECK_PETSC_ERROR(err);
+ faultMesh->setDMMesh(subdm);
+ } else {
+ // TODO: This leg will be unnecessary
DM dm;
IS subpointMap;
PetscInt *renum;
@@ -888,7 +914,7 @@
pylith::faults::CohesiveTopology::createInterpolated(topology::Mesh* mesh,
const topology::SubMesh& faultMesh,
const ALE::Obj<SieveFlexMesh>& faultBoundary,
- const char groupLabel[],
+ const ALE::Obj<topology::Mesh::IntSection>& groupField,
const int materialId,
int& firstFaultVertex,
int& firstLagrangeVertex,
@@ -897,26 +923,17 @@
{ // createInterpolated
assert(0 != mesh);
assert(!faultBoundary.isNull());
- assert(groupLabel);
+ assert(!groupField.isNull());
+ DM dm = mesh->dmMesh(), sdm;
+ DMLabel label;
+ const char *labelName = "faultSurface";
PetscErrorCode err;
- /* Recreate old mesh in new mesh
- Should be able to use the insertion routines from DMPlex here
- */
- DM dm = mesh->dmMesh();
- DM sdm;
-
- err = DMPlexConstructCohesiveCells(dm, groupLabel, &sdm);CHECK_PETSC_ERROR(err);
- // Renumber labels
- // Add fault vertices to groups and construct renumberings
- // Split the mesh along the fault sieve and create cohesive elements
+ err = DMPlexGetLabel(dm, labelName, &label);CHECK_PETSC_ERROR(err);
// Completes the set of cells scheduled to be replaced
// Have to do internal fault vertices before fault boundary vertices, and this is the only thing I use faultBoundary for
- err = DMPlexSymmetrize(sdm);CHECK_PETSC_ERROR(err);
- err = DMPlexStratify(sdm);CHECK_PETSC_ERROR(err);
- // Various checks for replaced cells
- // Fix coordinates
-
+ err = DMLabelCohesiveComplete(dm, label);CHECK_PETSC_ERROR(err);
+ err = DMPlexConstructCohesiveCells(dm, labelName, &sdm);CHECK_PETSC_ERROR(err);
mesh->setDMMesh(sdm);
} // createInterpolated
Modified: short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.hh 2013-01-31 00:58:01 UTC (rev 21311)
+++ short/3D/PyLith/trunk/libsrc/pylith/faults/CohesiveTopology.hh 2013-02-01 04:55:00 UTC (rev 21312)
@@ -101,7 +101,7 @@
void createInterpolated(topology::Mesh* mesh,
const topology::SubMesh& faultMesh,
const ALE::Obj<SieveFlexMesh>& faultBoundary,
- const char groupLabel[],
+ const ALE::Obj<topology::Mesh::IntSection>& groupField,
const int materialId,
int& firstFaultVertex,
int& firstLagrangeVertex,
Modified: short/3D/PyLith/trunk/libsrc/pylith/meshio/DataWriterHDF5.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/meshio/DataWriterHDF5.cc 2013-01-31 00:58:01 UTC (rev 21311)
+++ short/3D/PyLith/trunk/libsrc/pylith/meshio/DataWriterHDF5.cc 2013-02-01 04:55:00 UTC (rev 21312)
@@ -398,7 +398,7 @@
assert(_viewer);
try {
- const char* context = DataWriter<mesh_type, field_type>::_context.c_str();
+ const char* context = DataWriter<mesh_type, field_type>::_context.c_str();
field.createScatterWithBC(mesh, "", 0, context);
field.scatterSectionToVector(context);
@@ -416,6 +416,18 @@
if (_tstampIndex == istep)
_writeTimeStamp(t, commRank);
+ const int spaceDim = mesh.coordsys()->spaceDim();
+ PetscInt bs;
+ err = VecGetBlockSize(vector, &bs);CHECK_PETSC_ERROR(err);
+ switch (field.vectorFieldType()) {
+ case pylith::topology::FieldBase::VECTOR:
+ if (bs%spaceDim) CHECK_PETSC_ERROR(PETSC_ERR_ARG_WRONG); break;
+ case pylith::topology::FieldBase::TENSOR:
+ if (bs%spaceDim) CHECK_PETSC_ERROR(PETSC_ERR_ARG_WRONG); break;
+ //default:
+ //if (bs > 1) CHECK_PETSC_ERROR(PETSC_ERR_ARG_WRONG); break;
+ }
+
PetscErrorCode err = 0;
err = PetscViewerHDF5PushGroup(_viewer, "/vertex_fields");CHECK_PETSC_ERROR(err);
err = PetscViewerHDF5SetTimestep(_viewer, istep);CHECK_PETSC_ERROR(err);
Modified: short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc 2013-01-31 00:58:01 UTC (rev 21311)
+++ short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc 2013-02-01 04:55:00 UTC (rev 21312)
@@ -138,6 +138,7 @@
delete [] cone2; cone2 = 0;
err = DMPlexSymmetrize(complexMesh);CHECK_PETSC_ERROR(err);
err = DMPlexStratify(complexMesh);CHECK_PETSC_ERROR(err);
+ /* TODO Interpolate mesh if necessary */
} else {
// Same old thing
ALE::Obj<SieveFlexMesh::sieve_type> s =
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/SubMesh.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/SubMesh.cc 2013-01-31 00:58:01 UTC (rev 21311)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/SubMesh.cc 2013-02-01 04:55:00 UTC (rev 21312)
@@ -112,7 +112,7 @@
meshSieveMesh->getRealSection("coordinates_dimensioned"));
/* TODO: Add creation of pointSF for submesh */
- err = DMPlexCreateSubmesh(dmMesh, label, &_newMesh);CHECK_PETSC_ERROR(err);
+ err = DMPlexCreateSubmesh(dmMesh, label, PETSC_NULL, &_newMesh);CHECK_PETSC_ERROR(err);
// Create the parallel overlap
const ALE::Obj<SieveMesh::sieve_type>& sieve = _mesh->getSieve();
More information about the CIG-COMMITS
mailing list