[cig-commits] r19926 - in short/3D/PyLith/trunk/libsrc/pylith: meshio problems topology
knepley at geodynamics.org
knepley at geodynamics.org
Thu Apr 5 19:56:24 PDT 2012
Author: knepley
Date: 2012-04-05 19:56:24 -0700 (Thu, 05 Apr 2012)
New Revision: 19926
Modified:
short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc
short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshIO.cc
short/3D/PyLith/trunk/libsrc/pylith/problems/SolverNonlinear.cc
short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.cc
short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.hh
short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.icc
Log:
Added code to create a DMComplex
Modified: short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc 2012-04-05 23:55:18 UTC (rev 19925)
+++ short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshBuilder.cc 2012-04-06 02:56:24 UTC (rev 19926)
@@ -22,6 +22,7 @@
#include "pylith/topology/Mesh.hh" // USES Mesh
#include "pylith/utils/array.hh" // USES scalar_array, int_array
+#include "pylith/utils/petscerror.h" // USES CHECK_PETSC_ERROR
#include "spatialdata/units/Nondimensional.hh" // USES Nondimensional
@@ -54,6 +55,7 @@
MPI_Comm comm = mesh->comm();
int dim = meshDim;
int rank = 0;
+ PetscErrorCode err;
{ // Check to make sure every vertex is in at least one cell.
// This is required by Sieve
@@ -77,6 +79,10 @@
mesh->createSieveMesh(dim);
const ALE::Obj<SieveMesh>& sieveMesh = mesh->sieveMesh();
assert(!sieveMesh.isNull());
+ /* DMComplex */
+ mesh->createDMMesh(dim);
+ DM complexMesh = mesh->dmMesh();
+ assert(complexMesh);
ALE::Obj<SieveMesh::sieve_type> sieve =
new SieveMesh::sieve_type(mesh->comm());
@@ -117,6 +123,22 @@
delete[] coneO; coneO = 0;
// Symmetrize to fill up supports
sieve->symmetrize();
+ /* DMComplex */
+ err = DMComplexSetChart(complexMesh, 0, numCells+numVertices);CHECK_PETSC_ERROR(err);
+ for(PetscInt c = 0; c < numCells; ++c) {
+ err = DMComplexSetConeSize(complexMesh, c, numCorners);CHECK_PETSC_ERROR(err);
+ }
+ err = DMSetUp(complexMesh);CHECK_PETSC_ERROR(err);
+ PetscInt *cone2 = new PetscInt[numCorners];
+ for(PetscInt c = 0; c < numCells; ++c) {
+ for(PetscInt v = 0; v < numCorners; ++v) {
+ cone2[v] = cells[c*numCorners+v]+numCells;
+ }
+ err = DMComplexSetCone(complexMesh, c, cone2);CHECK_PETSC_ERROR(err);
+ } // for
+ delete [] cone2; cone2 = 0;
+ err = DMComplexSymmetrize(complexMesh);CHECK_PETSC_ERROR(err);
+ err = DMComplexStratify(complexMesh);CHECK_PETSC_ERROR(err);
} else {
// Same old thing
ALE::Obj<SieveFlexMesh::sieve_type> s =
@@ -144,20 +166,6 @@
const ALE::Obj<SieveMesh::label_type>& depth =
sieveMesh->createLabel("depth");
-#ifdef IMESH_NEW_LABELS
- height->setChart(sieveMesh->getSieve()->getChart());
- depth->setChart(sieveMesh->getSieve()->getChart());
- for(int c = 0; c < numCells+numVertices; ++c) {
- height->setConeSize(c, 1);
- depth->setConeSize(c, 1);
- } // for
- if (numCells+numVertices) {
- height->setSupportSize(0, numCells+numVertices);
- depth->setSupportSize(0, numCells+numVertices);
- } // if
- height->allocate();
- depth->allocate();
-#endif
for(int c = 0; c < numCells; ++c) {
height->setCone(0, c);
depth->setCone(1, c);
@@ -166,10 +174,6 @@
height->setCone(1, v);
depth->setCone(0, v);
} // for
-#ifdef IMESH_NEW_LABELS
- height->recalculateLabel();
- depth->recalculateLabel();
-#endif
sieveMesh->setHeight(1);
sieveMesh->setDepth(1);
} else {
@@ -188,6 +192,32 @@
logger.stagePush("MeshCoordinates");
ALE::SieveBuilder<SieveMesh>::buildCoordinates(sieveMesh, spaceDim,
&(*coordinates)[0]);
+ /* DMComplex */
+ PetscSection coordSection;
+ Vec coordVec;
+ PetscScalar *coords;
+ PetscInt coordSize;
+
+ err = DMComplexGetCoordinateSection(complexMesh, &coordSection);CHECK_PETSC_ERROR(err);
+ err = DMComplexGetCoordinateVec(complexMesh, &coordVec);CHECK_PETSC_ERROR(err);
+ err = PetscSectionSetChart(coordSection, numCells, numCells+numVertices);CHECK_PETSC_ERROR(err);
+ for(PetscInt v = numCells; v < numCells+numVertices; ++v) {
+ err = PetscSectionSetDof(coordSection, v, spaceDim);CHECK_PETSC_ERROR(err);
+ }
+ err = PetscSectionSetUp(coordSection);CHECK_PETSC_ERROR(err);
+ err = PetscSectionGetStorageSize(coordSection, &coordSize);CHECK_PETSC_ERROR(err);
+ err = VecSetSizes(coordVec, coordSize, PETSC_DETERMINE);CHECK_PETSC_ERROR(err);
+ err = VecSetFromOptions(coordVec);CHECK_PETSC_ERROR(err);
+ err = VecGetArray(coordVec, &coords);CHECK_PETSC_ERROR(err);
+ for(PetscInt v = 0; v < numVertices; ++v) {
+ PetscInt off;
+
+ err = PetscSectionGetOffset(coordSection, v+numCells, &off);CHECK_PETSC_ERROR(err);
+ for(PetscInt d = 0; d < spaceDim; ++d) {
+ coords[off+d] = (*coordinates)[v*spaceDim+d];
+ }
+ }
+ err = VecRestoreArray(coordVec, &coords);CHECK_PETSC_ERROR(err);
logger.stagePop(); // Coordinates
sieveMesh->getFactory()->clear();
Modified: short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshIO.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshIO.cc 2012-04-05 23:55:18 UTC (rev 19925)
+++ short/3D/PyLith/trunk/libsrc/pylith/meshio/MeshIO.cc 2012-04-06 02:56:24 UTC (rev 19926)
@@ -233,24 +233,11 @@
} // if
int i = 0;
-#ifdef IMESH_NEW_LABELS
- labelMaterials->setChart(sieveMesh->getSieve()->getChart());
for(SieveMesh::label_sequence::iterator e_iter = cellsBegin;
e_iter != cellsEnd;
++e_iter) {
- labelMaterials->setConeSize(*e_iter, 1);
- }
- if (cells->size()) {labelMaterials->setSupportSize(0, cells->size());}
- labelMaterials->allocate();
-#endif
- for(SieveMesh::label_sequence::iterator e_iter = cellsBegin;
- e_iter != cellsEnd;
- ++e_iter) {
sieveMesh->setValue(labelMaterials, *e_iter, materialIds[i++]);
}
-#ifdef IMESH_NEW_LABELS
- labelMaterials->recalculateLabel();
-#endif
} // if
logger.stagePop();
///logger.setDebug(0);
Modified: short/3D/PyLith/trunk/libsrc/pylith/problems/SolverNonlinear.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/problems/SolverNonlinear.cc 2012-04-05 23:55:18 UTC (rev 19925)
+++ short/3D/PyLith/trunk/libsrc/pylith/problems/SolverNonlinear.cc 2012-04-06 02:56:24 UTC (rev 19926)
@@ -114,7 +114,7 @@
PetscSNESLineSearch ls;
err = SNESGetSNESLineSearch(_snes, &ls); CHECK_PETSC_ERROR(err);
- err = SNESLineSearchSetType(ls, SNESSHELL); CHECK_PETSC_ERROR(err);
+ err = SNESLineSearchSetType(ls, SNESLINESEARCHSHELL); CHECK_PETSC_ERROR(err);
err = SNESLineSearchShellSetUserFunc(ls, lineSearch, (void*) formulation); CHECK_PETSC_ERROR(err);
if (formulation->splitFields()) {
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.cc 2012-04-05 23:55:18 UTC (rev 19925)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.cc 2012-04-06 02:56:24 UTC (rev 19926)
@@ -24,6 +24,7 @@
#include "spatialdata/geocoords/CoordSys.hh" // USES CoordSys
#include "spatialdata/units/Nondimensional.hh" // USES Nondimensional
#include "pylith/utils/array.hh" // USES scalar_array
+#include "pylith/utils/petscerror.h" // USES CHECK_PETSC_ERROR
#include <stdexcept> // USES std::runtime_error
#include <sstream> // USES std::ostringstream
@@ -32,6 +33,7 @@
// ----------------------------------------------------------------------
// Default constructor
pylith::topology::Mesh::Mesh(void) :
+ _newMesh(PETSC_NULL),
_coordsys(0),
_comm(PETSC_COMM_WORLD),
_debug(false)
@@ -43,10 +45,12 @@
pylith::topology::Mesh::Mesh(const int dim,
const MPI_Comm& comm) :
_mesh(new SieveMesh(comm, dim)),
+ _newMesh(PETSC_NULL),
_coordsys(0),
_comm(comm),
_debug(false)
{ // constructor
+ createDMMesh(dim);
_mesh->setName("domain");
assert(!_mesh->getFactory().isNull());
_mesh->getFactory()->clear();
@@ -66,6 +70,7 @@
{ // deallocate
delete _coordsys; _coordsys = 0;
_mesh.destroy();
+ PetscErrorCode err = DMDestroy(&_newMesh);CHECK_PETSC_ERROR(err);
} // deallocate
// ----------------------------------------------------------------------
@@ -80,6 +85,19 @@
assert(!_mesh->getFactory().isNull());
_mesh->getFactory()->clear();
} // createSieveMesh
+
+// ----------------------------------------------------------------------
+// Create DMComplex mesh.
+void
+pylith::topology::Mesh::createDMMesh(const int dim)
+{ // createDMMesh
+ PetscErrorCode err;
+ err = DMDestroy(&_newMesh);CHECK_PETSC_ERROR(err);
+ err = DMCreate(_comm, &_newMesh);CHECK_PETSC_ERROR(err);
+ err = DMSetType(_newMesh, DMCOMPLEX);CHECK_PETSC_ERROR(err);
+ err = DMComplexSetDimension(_newMesh, dim);CHECK_PETSC_ERROR(err);
+ err = PetscObjectSetName((PetscObject) _newMesh, "domain");CHECK_PETSC_ERROR(err);
+} // createDMMesh
// ----------------------------------------------------------------------
// Set coordinate system.
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.hh 2012-04-05 23:55:18 UTC (rev 19925)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.hh 2012-04-06 02:56:24 UTC (rev 19926)
@@ -106,7 +106,19 @@
* @returns Sieve mesh.
*/
ALE::Obj<SieveMesh>& sieveMesh(void);
+
+ /** Create DMComplex mesh.
+ *
+ * @param dim Dimension associated with mesh cells.
+ */
+ void createDMMesh(const int dim=3);
+ /** Get DMComplex mesh.
+ *
+ * @returns DMComplex mesh.
+ */
+ DM dmMesh(void);
+
/** Set coordinate system.
*
* @param cs Coordinate system.
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.icc 2012-04-05 23:55:18 UTC (rev 19925)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/Mesh.icc 2012-04-06 02:56:24 UTC (rev 19926)
@@ -34,6 +34,13 @@
return _mesh;
}
+// Get DMComplex mesh.
+inline
+DM
+pylith::topology::Mesh::dmMesh(void) {
+ return _newMesh;
+}
+
// Get coordinate system.
inline
const spatialdata::geocoords::CoordSys*
More information about the CIG-COMMITS
mailing list