[cig-commits] r7328 - in short/3D/PyLith/trunk: libsrc/faults
libsrc/materials libsrc/meshio modulesrc/topology
pylith/meshio pylith/topology
knepley at geodynamics.org
knepley at geodynamics.org
Wed Jun 20 15:49:33 PDT 2007
Author: knepley
Date: 2007-06-20 15:49:32 -0700 (Wed, 20 Jun 2007)
New Revision: 7328
Modified:
short/3D/PyLith/trunk/libsrc/faults/CohesiveTopology.cc
short/3D/PyLith/trunk/libsrc/faults/FaultCohesive.cc
short/3D/PyLith/trunk/libsrc/materials/Material.cc
short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc
short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh
short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc
short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src
short/3D/PyLith/trunk/pylith/meshio/MeshIO.py
short/3D/PyLith/trunk/pylith/topology/Mesh.py
Log:
The first parallel test now runs (correctness not verified)
Modified: short/3D/PyLith/trunk/libsrc/faults/CohesiveTopology.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/CohesiveTopology.cc 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/libsrc/faults/CohesiveTopology.cc 2007-06-20 22:49:32 UTC (rev 7328)
@@ -172,11 +172,17 @@
PointArray origVertices;
PointArray newVertices;
int oppositeVertex;
- const int numCorners = sieve->nCone(*mesh->heightStratum(0)->begin(), mesh->depth())->size();
- const int faceSize = _numFaceVertices(*mesh->heightStratum(0)->begin(), mesh);
- int *indices = new int[faceSize];
- const int firstCohesiveCell = newPoint;
-
+ int numCorners = 0;
+ int faceSize = 0;
+ int *indices = NULL;
+ int firstCohesiveCell;
+
+ if (!(*fault)->commRank()) {
+ numCorners = sieve->nCone(*mesh->heightStratum(0)->begin(), mesh->depth())->size();
+ faceSize = _numFaceVertices(*mesh->heightStratum(0)->begin(), mesh);
+ indices = new int[faceSize];
+ firstCohesiveCell = newPoint;
+ }
for(Mesh::label_sequence::iterator f_iter = faces->begin();
f_iter != faces->end();
++f_iter, ++newPoint) {
@@ -255,7 +261,7 @@
}
mesh->setValue(material, newPoint, materialId);
} // for
- delete [] indices;
+ if (!(*fault)->commRank()) delete [] indices;
mesh->stratify();
const ALE::Obj<Mesh::label_type>& label = mesh->createLabel(std::string("censored depth"));
const ALE::Obj<std::set<Mesh::point_type> > modifiedPoints = new std::set<Mesh::point_type>();
Modified: short/3D/PyLith/trunk/libsrc/faults/FaultCohesive.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/FaultCohesive.cc 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/libsrc/faults/FaultCohesive.cc 2007-06-20 22:49:32 UTC (rev 7328)
@@ -48,7 +48,7 @@
assert(!groupField.isNull());
CohesiveTopology::create(_faultMesh, mesh, groupField, id(),
- _useLagrangeConstraints());
+ _useLagrangeConstraints());
} // adjustTopology
// ----------------------------------------------------------------------
Modified: short/3D/PyLith/trunk/libsrc/materials/Material.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.cc 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.cc 2007-06-20 22:49:32 UTC (rev 7328)
@@ -70,11 +70,13 @@
const ALE::Mesh::label_sequence::iterator cellsEnd = cells->end();
// Check to make sure we have cells
+#if 0
if (0 == cells->size()) {
std::ostringstream msg;
msg << "Could not find any cells for material '" << _label << "'.";
throw std::runtime_error(msg.str());
} // if
+#endif
// Create sections to hold parameters for physical properties
if (0 != _parameters) // Can't delete NULL pointer that holds reference
Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc 2007-06-20 22:49:32 UTC (rev 7328)
@@ -84,25 +84,34 @@
const int meshDim)
{ // _buildMesh
assert(0 != _mesh);
- assert(coordinates.size() == numVertices*spaceDim);
- assert(cells.size() == numCells*numCorners);
+ MPI_Comm comm = PETSC_COMM_WORLD;
+ int dim = meshDim;
+ int rank;
+ MPI_Bcast(&dim, 1, MPI_INT, 0, comm);
// :BUG: This causes a memory leak.
- *_mesh = new Mesh(PETSC_COMM_WORLD, meshDim);
+ *_mesh = new Mesh(PETSC_COMM_WORLD, dim);
_mesh->addRef();
assert(!_mesh->isNull());
(*_mesh)->setDebug(_debug);
ALE::Obj<sieve_type> sieve = new sieve_type((*_mesh)->comm());
+ (*_mesh)->setSieve(sieve);
- ALE::SieveBuilder<Mesh>::buildTopology(sieve, meshDim,
- numCells,
- const_cast<int*>(&cells[0]),
- numVertices,
- _interpolate, numCorners, -1,
- (*_mesh)->getArrowSection("orientation"));
- (*_mesh)->setSieve(sieve);
+ MPI_Comm_rank(comm, &rank);
+ if (!rank) {
+ assert(coordinates.size() == numVertices*spaceDim);
+ assert(cells.size() == numCells*numCorners);
+ ALE::SieveBuilder<Mesh>::buildTopology(sieve, meshDim,
+ numCells,
+ const_cast<int*>(&cells[0]),
+ numVertices,
+ _interpolate, numCorners, -1,
+ (*_mesh)->getArrowSection("orientation"));
+ } else {
+ (*_mesh)->getArrowSection("orientation");
+ }
(*_mesh)->stratify();
ALE::SieveBuilder<Mesh>::buildCoordinates(*_mesh, spaceDim, &coordinates[0]);
} // _buildMesh
@@ -193,27 +202,29 @@
{ // _setMaterials
assert(0 != _mesh);
assert(!_mesh->isNull());
-
- const ALE::Obj<Mesh::label_sequence>& cells = (*_mesh)->heightStratum(0);
- assert(!cells.isNull());
- const int numCells = materialIds.size();
- if (cells->size() != numCells) {
- std::ostringstream msg;
- msg << "Mismatch in size of materials identifier array ("
- << numCells << ") and number of cells in mesh ("
- << cells->size() << ").";
- throw std::runtime_error(msg.str());
- } // if
-
const ALE::Obj<Mesh::label_type>& labelMaterials =
(*_mesh)->createLabel("material-id");
+
+ if (!(*_mesh)->commRank()) {
+ const ALE::Obj<Mesh::label_sequence>& cells = (*_mesh)->heightStratum(0);
+ assert(!cells.isNull());
+
+ const int numCells = materialIds.size();
+ if (cells->size() != numCells) {
+ std::ostringstream msg;
+ msg << "Mismatch in size of materials identifier array ("
+ << numCells << ") and number of cells in mesh ("
+ << cells->size() << ").";
+ throw std::runtime_error(msg.str());
+ } // if
- int i = 0;
- for(Mesh::label_sequence::iterator e_iter = cells->begin();
- e_iter != cells->end();
- ++e_iter)
- (*_mesh)->setValue(labelMaterials, *e_iter, materialIds[i++]);
+ int i = 0;
+ for(Mesh::label_sequence::iterator e_iter = cells->begin();
+ e_iter != cells->end();
+ ++e_iter)
+ (*_mesh)->setValue(labelMaterials, *e_iter, materialIds[i++]);
+ }
} // _setMaterials
// ----------------------------------------------------------------------
@@ -270,6 +281,47 @@
} // _setGroup
// ----------------------------------------------------------------------
+// Create empty groups on other processes
+void
+pylith::meshio::MeshIO::_distributeGroups()
+{ // _distributeGroups
+ assert(0 != _mesh);
+ assert(!_mesh->isNull());
+
+ if (!(*_mesh)->commRank()) {
+ const ALE::Obj<std::set<std::string> >& sectionNames =
+ (*_mesh)->getIntSections();
+ int numGroups = sectionNames->size();
+
+ MPI_Bcast(&numGroups, 1, MPI_INT, 0, (*_mesh)->comm());
+ for (std::set<std::string>::const_iterator name=sectionNames->begin();
+ name != sectionNames->end(); ++name) {
+ int len = name->size();
+
+ MPI_Bcast(&len, 1, MPI_INT, 0, (*_mesh)->comm());
+ MPI_Bcast((void *) name->c_str(), len, MPI_CHAR, 0, (*_mesh)->comm());
+ }
+ } else {
+ int numGroups;
+
+ MPI_Bcast(&numGroups, 1, MPI_INT, 0, (*_mesh)->comm());
+ for(int g = 0; g < numGroups; ++g) {
+ char *name;
+ int len;
+
+ MPI_Bcast(&len, 1, MPI_INT, 0, (*_mesh)->comm());
+ name = new char[len+1];
+ MPI_Bcast(name, len, MPI_CHAR, 0, (*_mesh)->comm());
+ name[len] = 0;
+ const ALE::Obj<int_section_type>& groupField = (*_mesh)->getIntSection(name);
+ assert(!groupField.isNull());
+ (*_mesh)->allocate(groupField);
+ delete [] name;
+ }
+ }
+} // _distributeGroups
+
+// ----------------------------------------------------------------------
// Get names of all groups in mesh.
void
pylith::meshio::MeshIO::_getGroupNames(string_vector* names) const
Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh 2007-06-20 22:49:32 UTC (rev 7328)
@@ -190,6 +190,10 @@
GroupPtType* type,
const char *name) const;
+ /** Create empty groups on other processes
+ */
+ void _distributeGroups();
+
// PRIVATE MEMBERS //////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc 2007-06-20 22:49:32 UTC (rev 7328)
@@ -18,6 +18,8 @@
#include "journal/info.h" // USES journal::info_t
+#include <petsc.h> // USES MPI
+
#include <iomanip> // USES setw(), setiosflags(), resetiosflags()
#include <assert.h> // USES assert()
#include <fstream> // USES std::ifstream, std::ofstream
@@ -46,6 +48,8 @@
void
pylith::meshio::MeshIOAscii::_read(void)
{ // _read
+ MPI_Comm comm = PETSC_COMM_WORLD;
+ int rank;
int meshDim = 0;
int spaceDim = 0;
int numVertices = 0;
@@ -54,81 +58,89 @@
double_array coordinates;
int_array cells;
int_array materialIds;
-
- std::ifstream filein(_filename.c_str());
- if (!filein.is_open() || !filein.good()) {
- std::ostringstream msg;
- msg << "Could not open mesh file '" << _filename
- << "' for reading.\n";
- throw std::runtime_error(msg.str());
- } // if
- std::string token;
- const int maxIgnore = 1024;
- filein >> token;
- if (0 != strcasecmp(token.c_str(), "mesh")) {
- std::ostringstream msg;
- msg << "Expected 'mesh' token but encountered '" << token << "'\n";
- throw std::runtime_error(msg.str());
- }
+ MPI_Comm_rank(comm, &rank);
+ if (!rank) {
+ std::ifstream filein(_filename.c_str());
+ if (!filein.is_open() || !filein.good()) {
+ std::ostringstream msg;
+ msg << "Could not open mesh file '" << _filename
+ << "' for reading.\n";
+ throw std::runtime_error(msg.str());
+ } // if
- bool readDim = false;
- bool readCells = false;
- bool readVertices = false;
- bool builtMesh = false;
+ std::string token;
+ const int maxIgnore = 1024;
+ filein >> token;
+ if (0 != strcasecmp(token.c_str(), "mesh")) {
+ std::ostringstream msg;
+ msg << "Expected 'mesh' token but encountered '" << token << "'\n";
+ throw std::runtime_error(msg.str());
+ }
- filein.ignore(maxIgnore, '{');
- filein >> token;
- while (filein.good() && token != "}") {
- if (0 == strcasecmp(token.c_str(), "dimension")) {
- filein.ignore(maxIgnore, '=');
- filein >> meshDim;
- readDim = true;
- } else if (0 == strcasecmp(token.c_str(), "use-index-zero")) {
- filein.ignore(maxIgnore, '=');
- std::string flag = "";
- filein >> flag;
- if (0 == strcasecmp(flag.c_str(), "true"))
- useIndexZero(true);
- else
- useIndexZero(false);
- } else if (0 == strcasecmp(token.c_str(), "vertices")) {
- filein.ignore(maxIgnore, '{');
- _readVertices(filein, &coordinates, &numVertices, &spaceDim);
- readVertices = true;
- } else if (0 == strcasecmp(token.c_str(), "cells")) {
- filein.ignore(maxIgnore, '{');
- _readCells(filein, &cells, &materialIds, &numCells, &numCorners);
- readCells = true;
- } else if (0 == strcasecmp(token.c_str(), "group")) {
- std::string name;
- GroupPtType type;
- int numPoints = 0;
- int_array points;
+ bool readDim = false;
+ bool readCells = false;
+ bool readVertices = false;
+ bool builtMesh = false;
- if (!builtMesh)
- throw std::runtime_error("Both 'vertices' and 'cells' must "
- "precede any groups in mesh file.");
- filein.ignore(maxIgnore, '{');
- _readGroup(filein, &points, &type, &name);
- _setGroup(name, type, points);
- } else {
- std::ostringstream msg;
- msg << "Could not parse '" << token << "' into a mesh setting.";
- throw std::runtime_error(msg.str());
- } // else
+ filein.ignore(maxIgnore, '{');
+ filein >> token;
+ while (filein.good() && token != "}") {
+ if (0 == strcasecmp(token.c_str(), "dimension")) {
+ filein.ignore(maxIgnore, '=');
+ filein >> meshDim;
+ readDim = true;
+ } else if (0 == strcasecmp(token.c_str(), "use-index-zero")) {
+ filein.ignore(maxIgnore, '=');
+ std::string flag = "";
+ filein >> flag;
+ if (0 == strcasecmp(flag.c_str(), "true"))
+ useIndexZero(true);
+ else
+ useIndexZero(false);
+ } else if (0 == strcasecmp(token.c_str(), "vertices")) {
+ filein.ignore(maxIgnore, '{');
+ _readVertices(filein, &coordinates, &numVertices, &spaceDim);
+ readVertices = true;
+ } else if (0 == strcasecmp(token.c_str(), "cells")) {
+ filein.ignore(maxIgnore, '{');
+ _readCells(filein, &cells, &materialIds, &numCells, &numCorners);
+ readCells = true;
+ } else if (0 == strcasecmp(token.c_str(), "group")) {
+ std::string name;
+ GroupPtType type;
+ int numPoints = 0;
+ int_array points;
- if (readDim && readCells && readVertices && !builtMesh) {
- // Can now build mesh
- _buildMesh(coordinates, numVertices, spaceDim,
- cells, numCells, numCorners, meshDim);
- _setMaterials(materialIds);
- builtMesh = true;
- } // if
+ if (!builtMesh)
+ throw std::runtime_error("Both 'vertices' and 'cells' must "
+ "precede any groups in mesh file.");
+ filein.ignore(maxIgnore, '{');
+ _readGroup(filein, &points, &type, &name);
+ _setGroup(name, type, points);
+ } else {
+ std::ostringstream msg;
+ msg << "Could not parse '" << token << "' into a mesh setting.";
+ throw std::runtime_error(msg.str());
+ } // else
- filein >> token;
- } // while
- filein.close();
+ if (readDim && readCells && readVertices && !builtMesh) {
+ // Can now build mesh
+ _buildMesh(coordinates, numVertices, spaceDim,
+ cells, numCells, numCorners, meshDim);
+ _setMaterials(materialIds);
+ builtMesh = true;
+ } // if
+
+ filein >> token;
+ } // while
+ filein.close();
+ } else {
+ _buildMesh(coordinates, numVertices, spaceDim,
+ cells, numCells, numCorners, meshDim);
+ _setMaterials(materialIds);
+ }
+ _distributeGroups();
} // read
// ----------------------------------------------------------------------
Modified: short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src 2007-06-20 22:49:32 UTC (rev 7328)
@@ -259,8 +259,34 @@
fieldVptr = PyCObject_AsVoidPtr(field)
ptr = Mesh_createMatrix(self.thisptr, fieldVptr)
return PyCObject_FromVoidPtr(ptr, PetscMat_destructor)
+
+ def view(self):
+ """
+ View the mesh.
+ """
+ # create shim for view
+ #embed{ void* Mesh_view(void* objVptr)
+ try {
+ ALE::Obj<ALE::Mesh>* mesh = (ALE::Obj<ALE::Mesh>*) objVptr;
+ assert(0 != mesh);
+ assert(!mesh->isNull());
+ (*mesh)->view("");
+ } catch (const std::exception& err) {
+ PyErr_SetString(PyExc_RuntimeError,
+ const_cast<char*>(err.what()));
+ } catch (const ALE::Exception& err) {
+ PyErr_SetString(PyExc_RuntimeError,
+ const_cast<char*>(err.msg().c_str()));
+ } catch (...) {
+ PyErr_SetString(PyExc_RuntimeError,
+ "Caught unknown C++ exception.");
+ } // try/catch
+ #}embed
+ Mesh_view(self.thisptr)
+ return
+
def _createHandle(self):
"""
Wrap pointer to C++ object in PyCObject.
Modified: short/3D/PyLith/trunk/pylith/meshio/MeshIO.py
===================================================================
--- short/3D/PyLith/trunk/pylith/meshio/MeshIO.py 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/pylith/meshio/MeshIO.py 2007-06-20 22:49:32 UTC (rev 7328)
@@ -59,9 +59,10 @@
mesh.initialize(self.coordsys)
# Read mesh
- import mpi
- if 0 == mpi.MPI_Comm_rank(mesh.comm()):
- self.cppHandle.read(mesh.cppHandle)
+ #import mpi
+ #if 0 == mpi.MPI_Comm_rank(mesh.comm()):
+ self.cppHandle.read(mesh.cppHandle)
+ mesh.view()
return mesh
Modified: short/3D/PyLith/trunk/pylith/topology/Mesh.py
===================================================================
--- short/3D/PyLith/trunk/pylith/topology/Mesh.py 2007-06-20 16:34:18 UTC (rev 7327)
+++ short/3D/PyLith/trunk/pylith/topology/Mesh.py 2007-06-20 22:49:32 UTC (rev 7328)
@@ -63,8 +63,17 @@
if not self.cppHandle is None:
dim = self.cppHandle.dimension
return dim
+
+ def view(self):
+ """
+ View the mesh.
+ """
+ if not self.cppHandle is None:
+ self.cppHandle.view
+ return
+
def comm(self):
"""
Get MPI communicator associated with mesh.
More information about the cig-commits
mailing list