[cig-commits] r14829 - in short/3D/PyLith/branches/pylith-swig: . libsrc libsrc/faults libsrc/feassemble libsrc/materials libsrc/meshio libsrc/topology modulesrc/topology pylith/topology unittests/libtests/faults unittests/libtests/meshio
brad at geodynamics.org
brad at geodynamics.org
Thu Apr 30 17:34:25 PDT 2009
Author: brad
Date: 2009-04-30 17:34:23 -0700 (Thu, 30 Apr 2009)
New Revision: 14829
Added:
short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Distributor.i
Modified:
short/3D/PyLith/branches/pylith-swig/TODO
short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am
short/3D/PyLith/branches/pylith-swig/libsrc/faults/CohesiveTopology.cc
short/3D/PyLith/branches/pylith-swig/libsrc/faults/FaultCohesiveKin.cc
short/3D/PyLith/branches/pylith-swig/libsrc/feassemble/ElasticityImplicit.cc
short/3D/PyLith/branches/pylith-swig/libsrc/materials/GenMaxwellIsotropic3D.cc
short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileAscii.cc
short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileBinary.cc
short/3D/PyLith/branches/pylith-swig/libsrc/meshio/MeshIOCubit.cc
short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.cc
short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.hh
short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Makefile.am
short/3D/PyLith/branches/pylith-swig/modulesrc/topology/topology.i
short/3D/PyLith/branches/pylith-swig/pylith/topology/Distributor.py
short/3D/PyLith/branches/pylith-swig/unittests/libtests/faults/TestFaultCohesiveKin.cc
short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/Makefile.am
short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.cc
short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.hh
Log:
Updated DataWriterVTK tests for fault mesh. Added missing interface file for Distributor.i.
Modified: short/3D/PyLith/branches/pylith-swig/TODO
===================================================================
--- short/3D/PyLith/branches/pylith-swig/TODO 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/TODO 2009-05-01 00:34:23 UTC (rev 14829)
@@ -6,11 +6,6 @@
0. SWIG conversion [Brad]
- (1) Faults
-
- (2) Output
- TestDataWriterVTKFaultMesh
-
(3) Full-scale tests (few 1-D, 2-D, and 3-D)
Analytical solutions, time stepping, multiple faults
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am 2009-05-01 00:34:23 UTC (rev 14829)
@@ -99,11 +99,11 @@
topology/MeshOps.cc \
topology/SubMesh.cc \
topology/SolutionFields.cc \
+ topology/Distributor.cc \
utils/EventLogger.cc \
utils/TestArray.cc
# materials/GenMaxwellIsotropic3D.cc \
-# topology/Distributor.cc \
# topology/MeshRefiner.cc \
# topology/RefineUniform.cc
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/faults/CohesiveTopology.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/faults/CohesiveTopology.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/faults/CohesiveTopology.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -200,7 +200,7 @@
++name) {
sieveMesh->reallocate(sieveMesh->getIntSection(*name));
} // for
-#if 0
+#if 0 // TEST OF OPTIMIZATION?? [MATT: WHY IS THIS COMMENTED OUT?]
for(SieveSubMesh::label_sequence::iterator v_iter = fVertices->begin();
v_iter != fVerticesEnd;
++v_iter, ++newPoint) {
@@ -372,7 +372,7 @@
} // if
// TODO: Need to reform the material label when sieve is reallocated
sieveMesh->setValue(material, newPoint, materialId);
-#if 0
+#if 0 // TEST OF OPTIMIZATION?? [MATT: WHY IS THIS COMMENTED OUT?]
// OPTIMIZATION
sieveMesh->setHeight(newPoint, 0);
sieveMesh->setDepth(newPoint, 1);
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/faults/FaultCohesiveKin.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/faults/FaultCohesiveKin.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/faults/FaultCohesiveKin.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -475,7 +475,7 @@
stiffnessCell.size(),
&stiffnessCell[0]);
-#if 0
+#if 0 // DEBUGGING
// Check that fault cells match cohesive cells
ALE::ISieveVisitor::PointRetriever<sieve_type> cV(std::max(1, mesh->getSieve()->getMaxConeSize()));
ALE::ISieveVisitor::PointRetriever<sieve_type> cV2(std::max(1, _faultMesh->getSieve()->getMaxConeSize()));
@@ -1130,7 +1130,7 @@
const SieveSubMesh::renumbering_type::const_iterator renumberingEnd =
renumbering.end();
-#if 0 // MOVE TO SEPARATE CHECK METHOD
+#if 0 // DEBUGGING, MOVE TO SEPARATE CHECK METHOD
// Check fault mesh and volume mesh coordinates
const ALE::Obj<RealSection>& coordinates = mesh->getRealSection("coordinates");
const ALE::Obj<RealSection>& fCoordinates = _faultMesh->getRealSection("coordinates");
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/feassemble/ElasticityImplicit.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/feassemble/ElasticityImplicit.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/feassemble/ElasticityImplicit.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -257,7 +257,7 @@
CALL_MEMBER_FN(*this, elasticityResidualFn)(stressCell);
_logger->eventEnd(computeEvent);
-#if 0
+#if 0 // DEBUGGING
std::cout << "Updating residual for cell " << *c_iter << std::endl;
for(int i = 0; i < _quadrature->spaceDim() * _quadrature->numBasis(); ++i) {
std::cout << " v["<<i<<"]: " << _cellVector[i] << std::endl;
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/materials/GenMaxwellIsotropic3D.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/materials/GenMaxwellIsotropic3D.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/materials/GenMaxwellIsotropic3D.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -659,7 +659,7 @@
elasticConsts[19] = 0; // C2313
elasticConsts[20] = elasticConsts[15]; // C1313
-#if 0
+#if 0 // DEBUGGING
std::cout << "_calcElasticConstsViscoelastic" << std::endl;
std::cout << elasticConsts[0] << " " << elasticConsts[1] << " " << elasticConsts[2] << std::endl;
std::cout << elasticConsts[6] << " " << elasticConsts[7] << std::endl;
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileAscii.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileAscii.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileAscii.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -113,7 +113,7 @@
assert(cells.size() == numCells * numCorners);
assert(materialIds.size() == numCells);
-#if 0
+#if 0 // NOT YET IMPLEMENTED
_writeHeader();
_writeVertices(coordinates);
_writeCells(cells);
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileBinary.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileBinary.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/GMVFileBinary.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -118,7 +118,7 @@
assert(cells.size() == numCells * numCorners);
assert(materialIds.size() == numCells);
-#if 0
+#if 0 // NOT YET IMPLEMENTED
_writeHeader();
_writeVertices(coordinates);
_writeCells(cells);
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/MeshIOCubit.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/MeshIOCubit.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/MeshIOCubit.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -347,7 +347,7 @@
assert(cells->size() == numCells*numCorners);
if (2 == meshDim && 4 == numCorners) // QUAD
-#if 0
+#if 0 // OLD
// 0 1 2 3 -> 0 1 3 2
for (int iCell=0; iCell < numCells; ++iCell) {
const int i2 = iCell*numCorners+2;
@@ -360,7 +360,7 @@
; // do nothing
#endif
else if (3 == meshDim && 8 == numCorners) // HEX
-#if 0
+#if 0 // OLD
// 0 1 2 3 4 5 6 7 -> 0 1 3 2 4 5 7 6
for (int iCell=0; iCell < numCells; ++iCell) {
const int i2 = iCell*numCorners+2;
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -14,8 +14,8 @@
#include "Distributor.hh" // implementation of class methods
-#include "pylith/utils/sievetypes.hh" // USES PETSc Mesh
-#include "pylith/utils/vectorfields.hh" // USES SCALAR_FIELD
+#include "pylith/topology/Mesh.hh" // USES Mesh
+#include "pylith/topology/Field.hh" // USES Field<Mesh>
#include "pylith/meshio/DataWriter.hh" // USES DataWriter
#include <cstring> // USES strlen()
@@ -25,6 +25,11 @@
#include <cassert> // USES assert()
// ----------------------------------------------------------------------
+typedef pylith::topology::Mesh::SieveMesh SieveMesh;
+typedef pylith::topology::Mesh::RealSection RealSection;
+typedef pylith::topology::Mesh::IntSection IntSection;
+
+// ----------------------------------------------------------------------
// Constructor
pylith::topology::Distributor::Distributor(void)
{ // constructor
@@ -39,99 +44,199 @@
// ----------------------------------------------------------------------
// Distribute mesh among processors.
void
-pylith::topology::Distributor::distribute(ALE::Obj<Mesh>* const newMesh,
- const ALE::Obj<Mesh>& origMesh,
+pylith::topology::Distributor::distribute(topology::Mesh* const newMesh,
+ const topology::Mesh& origMesh,
const char* partitioner)
{ // distribute
- if (0 == strcasecmp(partitioner, "")) {
- distribute_private<ALE::DistributionNew<Mesh> >(newMesh, origMesh);
+ assert(0 != newMesh);
+
+ newMesh->coordsys(origMesh.coordsys());
+
+ if (0 == strcasecmp(partitioner, ""))
+ _distribute<ALE::DistributionNew<SieveMesh> >(newMesh, origMesh);
#if defined(PETSC_HAVE_CHACO)
- } else if (0 == strcasecmp(partitioner, "chaco")) {
- distribute_private<ALE::DistributionNew<Mesh, ALE::Partitioner<ALE::Chaco::Partitioner<> > > >(newMesh, origMesh);
+ else if (0 == strcasecmp(partitioner, "chaco"))
+ _distribute<ALE::DistributionNew<SieveMesh, ALE::Partitioner<ALE::Chaco::Partitioner<> > > >(newMesh, origMesh);
#endif
#if defined(PETSC_HAVE_PARMETIS)
- } else if (0 == strcasecmp(partitioner, "parmetis")) {
- distribute_private<ALE::DistributionNew<Mesh, ALE::Partitioner<ALE::ParMetis::Partitioner<> > > >(newMesh, origMesh);
+ else if (0 == strcasecmp(partitioner, "parmetis"))
+ _distribute<ALE::DistributionNew<SieveMesh, ALE::Partitioner<ALE::ParMetis::Partitioner<> > > >(newMesh, origMesh);
#endif
- } else {
- std::cout << "ERROR: Using default partitioner instead of unknown partitioner " << partitioner << std::endl;
- distribute_private<ALE::DistributionNew<Mesh> >(newMesh, origMesh);
- }
-}
+ else {
+ std::cerr << "ERROR: Using default partitioner instead of unknown "
+ "partitioner '" << partitioner << "'." << std::endl;
+ _distribute<ALE::DistributionNew<SieveMesh> >(newMesh, origMesh);
+ } // else
+} // distribute
+// ----------------------------------------------------------------------
+// Write partitioning info for distributed mesh.
+void
+pylith::topology::Distributor::write(meshio::DataWriter<topology::Mesh, topology::Field<topology::Mesh> >* const writer,
+ const topology::Mesh& mesh)
+{ // write
+
+ // Setup and allocate field
+ const int fiberDim = 1;
+ topology::Field<topology::Mesh> partition(mesh);
+ partition.scale(1.0);
+ partition.label("partition");
+ partition.vectorFieldType(topology::FieldBase::SCALAR);
+ partition.newSection(topology::FieldBase::CELLS_FIELD, fiberDim);
+ partition.allocate();
+ const ALE::Obj<RealSection>& partitionSection = partition.section();
+ assert(!partitionSection.isNull());
+
+ const ALE::Obj<SieveMesh> sieveMesh = mesh.sieveMesh();
+ assert(!sieveMesh.isNull());
+ double rankReal = double(sieveMesh->commRank());
+ const ALE::Obj<SieveMesh::label_sequence>& cells =
+ sieveMesh->heightStratum(0);
+ assert(!cells.isNull());
+ const SieveMesh::label_sequence::iterator cellsBegin = cells->begin();
+ const SieveMesh::label_sequence::iterator cellsEnd = cells->end();
+ for (SieveMesh::label_sequence::iterator c_iter=cellsBegin;
+ c_iter != cellsEnd;
+ ++c_iter) {
+ partitionSection->updatePoint(*c_iter, &rankReal);
+ } // for
+
+ //partition->view("PARTITION");
+ const double t = 0.0;
+ const int numTimeSteps = 0;
+ writer->open(mesh, numTimeSteps);
+ writer->openTimeStep(t, mesh);
+ writer->writeCellField(t, partition);
+ writer->closeTimeStep();
+ writer->close();
+} // write
+
+// ----------------------------------------------------------------------
template<typename DistributionType>
void
-pylith::topology::Distributor::distribute_private(ALE::Obj<Mesh>* const newMesh,
- const ALE::Obj<Mesh>& origMesh)
+pylith::topology::Distributor::_distribute(topology::Mesh* const newMesh,
+ const topology::Mesh& origMesh)
{ // distribute
- typedef typename Mesh::point_type point_type;
+ typedef typename SieveMesh::point_type point_type;
typedef typename DistributionType::partitioner_type partitioner_type;
typedef typename DistributionType::partition_type partition_type;
- const Obj<Mesh::sieve_type> newSieve = new Mesh::sieve_type(origMesh->comm(), origMesh->debug());
- const Obj<Mesh::send_overlap_type> sendMeshOverlap = new Mesh::send_overlap_type(origMesh->comm(), origMesh->debug());
- const Obj<Mesh::recv_overlap_type> recvMeshOverlap = new Mesh::recv_overlap_type(origMesh->comm(), origMesh->debug());
+ ALE::Obj<SieveMesh>& newSieveMesh = newMesh->sieveMesh();
+ assert(!newSieveMesh.isNull());
+ const ALE::Obj<SieveMesh>& origSieveMesh = origMesh.sieveMesh();
+ assert(!origSieveMesh.isNull());
- *newMesh = new Mesh(origMesh->comm(), origMesh->getDimension(), origMesh->debug());
- (*newMesh)->setSieve(newSieve);
- // IMESH_TODO
- // This might be unnecessary, since the overlap for submeshes is just the restriction of the overlaps
- // std::map<point_type,point_type> renumbering;
- Mesh::renumbering_type& renumbering = (*newMesh)->getRenumbering();
+ const ALE::Obj<SieveMesh::sieve_type> newSieve =
+ new SieveMesh::sieve_type(origSieveMesh->comm(), origSieveMesh->debug());
+ assert(!newSieve.isNull());
+ const ALE::Obj<SieveMesh::send_overlap_type> sendMeshOverlap =
+ new SieveMesh::send_overlap_type(origSieveMesh->comm(),
+ origSieveMesh->debug());
+ assert(!sendMeshOverlap.isNull());
+ const ALE::Obj<SieveMesh::recv_overlap_type> recvMeshOverlap =
+ new SieveMesh::recv_overlap_type(origSieveMesh->comm(),
+ origSieveMesh->debug());
+ assert(!recvMeshOverlap.isNull());
+
+ newSieveMesh = new SieveMesh(origSieveMesh->comm(),
+ origSieveMesh->getDimension(),
+ origSieveMesh->debug());
+ assert(!newSieveMesh.isNull());
+ newSieveMesh->setSieve(newSieve);
+ // IMESH_TODO This might be unnecessary, since the overlap for
+ // submeshes is just the restriction of the overlaps
+ // std::map<point_type,point_type> renumbering;
+ SieveMesh::renumbering_type& renumbering = newSieveMesh->getRenumbering();
// Distribute the mesh
- Obj<partition_type> partition = DistributionType::distributeMeshV(origMesh, (*newMesh), renumbering, sendMeshOverlap, recvMeshOverlap);
- if (origMesh->debug()) {
- std::cout << "["<<origMesh->commRank()<<"]: Mesh Renumbering:" << std::endl;
- for(Mesh::renumbering_type::const_iterator r_iter = renumbering.begin(); r_iter != renumbering.end(); ++r_iter) {
- std::cout << "["<<origMesh->commRank()<<"]: global point " << r_iter->first << " --> " << " local point " << r_iter->second << std::endl;
- }
- }
+ ALE::Obj<partition_type> partition =
+ DistributionType::distributeMeshV(origSieveMesh, newSieveMesh,
+ renumbering,
+ sendMeshOverlap, recvMeshOverlap);
+ if (origSieveMesh->debug()) {
+ std::cout << "["<<origSieveMesh->commRank()<<"]: Mesh Renumbering:"
+ << std::endl;
+ for (SieveMesh::renumbering_type::const_iterator r_iter = renumbering.begin();
+ r_iter != renumbering.end();
+ ++r_iter) {
+ std::cout << "["<<origSieveMesh->commRank()<<"]: global point "
+ << r_iter->first << " --> " << " local point "
+ << r_iter->second << std::endl;
+ } // for
+ } // if
// Check overlap
int localSendOverlapSize = 0, sendOverlapSize;
int localRecvOverlapSize = 0, recvOverlapSize;
- for(int p = 0; p < sendMeshOverlap->commSize(); ++p) {
+ const int commSize = sendMeshOverlap->commSize();
+ for (int p = 0; p < commSize; ++p) {
localSendOverlapSize += sendMeshOverlap->cone(p)->size();
localRecvOverlapSize += recvMeshOverlap->support(p)->size();
- }
- MPI_Allreduce(&localSendOverlapSize, &sendOverlapSize, 1, MPI_INT, MPI_SUM, sendMeshOverlap->comm());
- MPI_Allreduce(&localRecvOverlapSize, &recvOverlapSize, 1, MPI_INT, MPI_SUM, recvMeshOverlap->comm());
- if(sendOverlapSize != recvOverlapSize) {
- std::cout <<"["<<sendMeshOverlap->commRank()<<"]: Size mismatch " << sendOverlapSize << " != " << recvOverlapSize << std::endl;
+ } // for
+ MPI_Allreduce(&localSendOverlapSize, &sendOverlapSize, 1, MPI_INT, MPI_SUM,
+ sendMeshOverlap->comm());
+ MPI_Allreduce(&localRecvOverlapSize, &recvOverlapSize, 1, MPI_INT, MPI_SUM,
+ recvMeshOverlap->comm());
+ if (sendOverlapSize != recvOverlapSize) {
+ std::cout <<"["<<sendMeshOverlap->commRank()<<"]: Size mismatch " <<
+ sendOverlapSize << " != " << recvOverlapSize << std::endl;
sendMeshOverlap->view("Send Overlap");
recvMeshOverlap->view("Recv Overlap");
throw ALE::Exception("Invalid Overlap");
- }
+ } // if
// Distribute the coordinates
- const Obj<real_section_type>& coordinates = origMesh->getRealSection("coordinates");
- const Obj<real_section_type>& parallelCoordinates = (*newMesh)->getRealSection("coordinates");
+ const ALE::Obj<RealSection>& coordinates =
+ origSieveMesh->getRealSection("coordinates");
+ assert(!coordinates.isNull());
+ const ALE::Obj<RealSection>& parallelCoordinates =
+ newSieveMesh->getRealSection("coordinates");
+ assert(!parallelCoordinates.isNull());
- (*newMesh)->setupCoordinates(parallelCoordinates);
- DistributionType::distributeSection(coordinates, partition, renumbering, sendMeshOverlap, recvMeshOverlap, parallelCoordinates);
+ newSieveMesh->setupCoordinates(parallelCoordinates);
+ DistributionType::distributeSection(coordinates, partition, renumbering,
+ sendMeshOverlap, recvMeshOverlap,
+ parallelCoordinates);
// Distribute other sections
- if (origMesh->getRealSections()->size() > 1) {
- Obj<std::set<std::string> > names = origMesh->getRealSections();
- int n = 0;
+ if (origSieveMesh->getRealSections()->size() > 1) {
+ ALE::Obj<std::set<std::string> > names = origSieveMesh->getRealSections();
+ assert(!names.isNull());
+ int n = 0;
- for(std::set<std::string>::const_iterator n_iter = names->begin(); n_iter != names->end(); ++n_iter) {
- if (*n_iter == "coordinates") continue;
- if (*n_iter == "replaced_cells") continue;
- std::cout << "ERROR: Did not distribute real section " << *n_iter << std::endl;
+ const std::set<std::string>::const_iterator namesBegin = names->begin();
+ const std::set<std::string>::const_iterator namesEnd = names->end();
+ for (std::set<std::string>::const_iterator n_iter = namesBegin;
+ n_iter != namesEnd;
+ ++n_iter) {
+ if (*n_iter == "coordinates") continue; // already copied
+ if (*n_iter == "replaced_cells") continue; // ignore
+ std::cerr << "ERROR: Did not distribute real section '" << *n_iter
+ << "'." << std::endl;
++n;
- }
- if (n) {throw ALE::Exception("Need to distribute more real sections");}
+ } // if
+ if (n)
+ throw std::logic_error("Need to distribute more real sections");
}
- if (origMesh->getIntSections()->size() > 0) {
- Obj<std::set<std::string> > names = origMesh->getIntSections();
+ if (origSieveMesh->getIntSections()->size() > 0) {
+ ALE::Obj<std::set<std::string> > names = origSieveMesh->getIntSections();
+ assert(!names.isNull());
- for(std::set<std::string>::const_iterator n_iter = names->begin(); n_iter != names->end(); ++n_iter) {
- const Obj<Mesh::int_section_type>& origSection = origMesh->getIntSection(*n_iter);
- const Obj<Mesh::int_section_type>& newSection = (*newMesh)->getIntSection(*n_iter);
+ std::set<std::string>::const_iterator namesBegin = names->begin();
+ std::set<std::string>::const_iterator namesEnd = names->end();
+ for (std::set<std::string>::const_iterator n_iter = namesBegin;
+ n_iter != namesEnd;
+ ++n_iter) {
+ const ALE::Obj<IntSection>& origSection =
+ origSieveMesh->getIntSection(*n_iter);
+ assert(!origSection.isNull());
+ const ALE::Obj<IntSection>& newSection =
+ newSieveMesh->getIntSection(*n_iter);
+ assert(!newSection.isNull());
// We assume all integer sections are complete sections
- newSection->setChart((*newMesh)->getSieve()->getChart());
- DistributionType::distributeSection(origSection, partition, renumbering, sendMeshOverlap, recvMeshOverlap, newSection);
-#if 0
+ newSection->setChart(newSieveMesh->getSieve()->getChart());
+ DistributionType::distributeSection(origSection, partition, renumbering,
+ sendMeshOverlap, recvMeshOverlap,
+ newSection);
+#if 0 // DEBUGGING
std::string serialName("Serial ");
std::string parallelName("Parallel ");
serialName += *n_iter;
@@ -139,36 +244,44 @@
origSection->view(serialName.c_str());
newSection->view(parallelName.c_str());
#endif
- }
- }
- if (origMesh->getArrowSections()->size() > 1) {
- throw ALE::Exception("Need to distribute more arrow sections");
- }
+ } // for
+ } // if
+ if (origSieveMesh->getArrowSections()->size() > 1)
+ throw std::logic_error("Need to distribute more arrow sections");
+
// Distribute labels
- const Mesh::labels_type& labels = origMesh->getLabels();
+ const SieveMesh::labels_type& labels = origSieveMesh->getLabels();
+ const SieveMesh::labels_type::const_iterator labelsBegin = labels.begin();
+ const SieveMesh::labels_type::const_iterator labelsEnd = labels.end();
- for(Mesh::labels_type::const_iterator l_iter = labels.begin(); l_iter != labels.end(); ++l_iter) {
- if ((*newMesh)->hasLabel(l_iter->first)) continue;
+ for (SieveMesh::labels_type::const_iterator l_iter = labelsBegin;
+ l_iter != labelsEnd;
+ ++l_iter) {
+ if (newSieveMesh->hasLabel(l_iter->first)) continue;
+ const ALE::Obj<SieveMesh::label_type>& origLabel = l_iter->second;
+ assert(!origLabel.isNull());
+ const ALE::Obj<SieveMesh::label_type>& newLabel =
+ newSieveMesh->createLabel(l_iter->first);
+ assert(!newLabel.isNull());
+
#ifdef IMESH_NEW_LABELS
- const Obj<Mesh::label_type>& origLabel = l_iter->second;
- const Obj<Mesh::label_type>& newLabel = (*newMesh)->createLabel(l_iter->first);
-
newLabel->setChart(newSieve->getChart());
// Size the local mesh
- partitioner_type::sizeLocalSieveV(origLabel, partition, renumbering, newLabel);
+ partitioner_type::sizeLocalSieveV(origLabel, partition, renumbering,
+ newLabel);
// Create the remote meshes
- DistributionType::completeConesV(origLabel, newLabel, renumbering, sendMeshOverlap, recvMeshOverlap);
+ DistributionType::completeConesV(origLabel, newLabel, renumbering,
+ sendMeshOverlap, recvMeshOverlap);
// Create the local mesh
- partitioner_type::createLocalSieveV(origLabel, partition, renumbering, newLabel);
+ partitioner_type::createLocalSieveV(origLabel, partition, renumbering,
+ newLabel);
newLabel->symmetrize();
#else
- const Obj<Mesh::label_type>& origLabel = l_iter->second;
- const Obj<Mesh::label_type>& newLabel = (*newMesh)->createLabel(l_iter->first);
// Get remote labels
- ALE::New::Completion<Mesh,Mesh::point_type>::scatterCones(origLabel, newLabel, sendMeshOverlap, recvMeshOverlap, renumbering);
+ ALE::New::Completion<SieveMesh,SieveMesh::point_type>::scatterCones(origLabel, newLabel, sendMeshOverlap, recvMeshOverlap, renumbering);
// Create local label
- newLabel->add(origLabel, (*newMesh)->getSieve(), renumbering);
-#if 0
+ newLabel->add(origLabel, newSieveMesh->getSieve(), renumbering);
+#if 0 // DEBUGGING
std::string serialName("Serial ");
std::string parallelName("Parallel ");
serialName += l_iter->first;
@@ -177,54 +290,23 @@
newLabel->view(parallelName.c_str());
#endif
#endif
- }
+ } // for
+
// Create the parallel overlap
- Obj<Mesh::send_overlap_type> sendParallelMeshOverlap = (*newMesh)->getSendOverlap();
- Obj<Mesh::recv_overlap_type> recvParallelMeshOverlap = (*newMesh)->getRecvOverlap();
+ ALE::Obj<SieveMesh::send_overlap_type> sendParallelMeshOverlap =
+ newSieveMesh->getSendOverlap();
+ assert(!sendParallelMeshOverlap.isNull());
+ ALE::Obj<SieveMesh::recv_overlap_type> recvParallelMeshOverlap =
+ newSieveMesh->getRecvOverlap();
+ assert(!recvParallelMeshOverlap.isNull());
+
// Can I figure this out in a nicer way?
ALE::SetFromMap<std::map<point_type,point_type> > globalPoints(renumbering);
-
- ALE::OverlapBuilder<>::constructOverlap(globalPoints, renumbering, sendParallelMeshOverlap, recvParallelMeshOverlap);
- (*newMesh)->setCalculatedOverlap(true);
+ ALE::OverlapBuilder<>::constructOverlap(globalPoints, renumbering,
+ sendParallelMeshOverlap,
+ recvParallelMeshOverlap);
+ newSieveMesh->setCalculatedOverlap(true);
} // distribute
-// ----------------------------------------------------------------------
-// Write partitioning info for distributed mesh.
-void
-pylith::topology::Distributor::write(meshio::DataWriter* const writer,
- const ALE::Obj<Mesh>& mesh,
- const spatialdata::geocoords::CoordSys* cs)
-{ // write
-
- // Setup and allocate field
- const int fiberDim = 1;
- ALE::Obj<real_section_type> partition =
- new real_section_type(mesh->comm(), mesh->debug());
- const ALE::Obj<Mesh::label_sequence>& cells = mesh->heightStratum(0);
- assert(!cells.isNull());
- partition->setChart(real_section_type::chart_type(*std::min_element(cells->begin(), cells->end()),
- *std::max_element(cells->begin(), cells->end())+1));
- partition->setFiberDimension(cells, fiberDim);
- mesh->allocate(partition);
- const int rank = mesh->commRank();
- double rankReal = double(rank);
- const Mesh::label_sequence::iterator cellsEnd = cells->end();
- for (Mesh::label_sequence::iterator c_iter=cells->begin();
- c_iter != cellsEnd;
- ++c_iter) {
- partition->updatePoint(*c_iter, &rankReal);
- } // for
-
- //partition->view("PARTITION");
- const double t = 0.0;
- const int numTimeSteps = 0;
- writer->open(mesh, cs, numTimeSteps);
- writer->openTimeStep(t, mesh, cs);
- writer->writeCellField(t, "partition", partition, SCALAR_FIELD, mesh);
- writer->closeTimeStep();
- writer->close();
-} // write
-
-
// End of file
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.hh 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/topology/Distributor.hh 2009-05-01 00:34:23 UTC (rev 14829)
@@ -19,25 +19,11 @@
#if !defined(pylith_topology_distributor_hh)
#define pylith_topology_distributor_hh
-#include "pylith/utils/sievetypes.hh" // USES PETSc Mesh
+// Include directives ---------------------------------------------------
+#include "topologyfwd.hh" // forward declarations
-namespace pylith {
- namespace topology {
- class Distributor;
- class TestDistributor;
- } // topology
+#include "pylith/meshio/meshiofwd.hh" // USES DataWriter<Mesh>
- namespace meshio {
- class DataWriter;
- } // meshio
-} // pylith
-
-namespace spatialdata {
- namespace geocoords {
- class CoordSys;
- } // geocoords
-} // spatialdata
-
class pylith::topology::Distributor
{ // Distributor
friend class TestDistributor; // unit testing
@@ -54,12 +40,12 @@
/** Distribute mesh among processors.
*
* @param newMesh Distributed mesh (result).
- * @param mesh Mesh to distribute.
+ * @param origMesh Mesh to distribute.
* @param partitioner Name of partitioner to use in distributing mesh.
*/
static
- void distribute(ALE::Obj<Mesh>* const newMesh,
- const ALE::Obj<Mesh>& mesh,
+ void distribute(topology::Mesh* const newMesh,
+ const topology::Mesh& origMesh,
const char* partitioner);
/** Write partitioning info for distributed mesh.
@@ -69,25 +55,29 @@
* @param cs Coordinate system for mesh.
*/
static
- void write(meshio::DataWriter* const writer,
- const ALE::Obj<Mesh>& mesh,
- const spatialdata::geocoords::CoordSys* cs);
+ void write(meshio::DataWriter<topology::Mesh, topology::Field<topology::Mesh> >* const writer,
+ const topology::Mesh& mesh);
-// NOT IMPLEMENTED //////////////////////////////////////////////////////
+// PRIVATE //////////////////////////////////////////////////////////////
private :
- /// Not implemented
- Distributor(const Distributor&);
-
- /// Not implemented
- const Distributor& operator=(const Distributor&);
-
+ /** Distribute mesh among processors.
+ *
+ * @param newMesh Distributed mesh (result).
+ * @param origMesh Mesh to distribute.
+ */
template<typename DistributionType>
static
void
- distribute_private(ALE::Obj<Mesh>* const newMesh,
- const ALE::Obj<Mesh>& origMesh);
+ _distribute(topology::Mesh* const newMesh,
+ const topology::Mesh& origMesh);
+// NOT IMPLEMENTED //////////////////////////////////////////////////////
+private :
+
+ Distributor(const Distributor&); ///< Not implemented
+ const Distributor& operator=(const Distributor&); ///< Not implemented
+
}; // Distributor
#endif // pylith_topology_distributor_hh
Added: short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Distributor.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Distributor.i (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Distributor.i 2009-05-01 00:34:23 UTC (rev 14829)
@@ -0,0 +1,61 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+/**
+ * @file modulesrc/topology/Distributor.i
+ *
+ * @brief Python interface to C++ Distributor object.
+ */
+
+namespace pylith {
+ namespace topology {
+
+ class Distributor
+ { // Distributor
+
+ // PUBLIC MEMBERS /////////////////////////////////////////////////
+ public :
+
+ /// Constructor
+ Distributor(void);
+
+ /// Destructor
+ ~Distributor(void);
+
+ /** Distribute mesh among processors.
+ *
+ * @param newMesh Distributed mesh (result).
+ * @param origMesh Mesh to distribute.
+ * @param partitioner Name of partitioner to use in distributing mesh.
+ */
+ static
+ void distribute(pylith::topology::Mesh* const newMesh,
+ const pylith::topology::Mesh& origMesh,
+ const char* partitioner);
+
+ /** Write partitioning info for distributed mesh.
+ *
+ * @param writer Data writer for partition information.
+ * @param mesh Distributed mesh.
+ * @param cs Coordinate system for mesh.
+ */
+ static
+ void write(pylith::meshio::DataWriter<pylith::topology::Mesh, pylith::topology::Field<pylith::topology::Mesh> >* const writer,
+ const pylith::topology::Mesh& mesh);
+
+ }; // Distributor
+
+ } // topology
+} // pylith
+
+
+// End of file
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Makefile.am 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Makefile.am 2009-05-01 00:34:23 UTC (rev 14829)
@@ -26,7 +26,8 @@
Field.i \
Fields.i \
SolutionFields.i \
- Jacobian.i
+ Jacobian.i \
+ Distributor.i
swig_generated = \
topology_wrap.cxx \
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/topology/topology.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/topology/topology.i 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/topology/topology.i 2009-05-01 00:34:23 UTC (rev 14829)
@@ -23,6 +23,7 @@
#include "pylith/topology/Fields.hh"
#include "pylith/topology/SolutionFields.hh"
#include "pylith/topology/Jacobian.hh"
+#include "pylith/topology/Distributor.hh"
%}
%include "exception.i"
@@ -56,6 +57,7 @@
%include "Fields.i"
%include "SolutionFields.i"
%include "Jacobian.i"
+%include "Distributor.i"
// Template instatiation
%template(MeshField) pylith::topology::Field<pylith::topology::Mesh>;
Modified: short/3D/PyLith/branches/pylith-swig/pylith/topology/Distributor.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/topology/Distributor.py 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/pylith/topology/Distributor.py 2009-05-01 00:34:23 UTC (rev 14829)
@@ -16,61 +16,55 @@
##
## Factory: mesh_distributor.
-from pyre.components.Component import Component
+from pylith.utils.PetscComponent import PetscComponent
+from topology import Distributor as ModuleDistributor
# Distributor class
-class Distributor(Component):
+class Distributor(PetscComponent, ModuleDistributor):
"""
Python manager for distributing mesh among processors.
+ Inventory
+
+ \b Properties
+ @li \b partitioner Name of mesh partitioner {"parmetis", "chaco"}.
+ @li \b debug Write partition information to file.
+
+ \b Facilities
+ @li \b writer Data writer for for partition information.
+
Factory: mesh_distributor
"""
# INVENTORY //////////////////////////////////////////////////////////
- class Inventory(Component.Inventory):
- """
- Python object for managing Distributor facilities and properties.
- """
+ import pyre.inventory
+
+ partitioner = pyre.inventory.str("partitioner", default="chaco",
+ validator=pyre.inventory.choice(["chaco",
+ "parmetis"]))
+ partitioner.meta['tip'] = "Name of mesh partitioner."
+
+ debug = pyre.inventory.bool("debug", default=False)
+ debug.meta['tip'] = "Write partition information to file."
+
+ from pylith.meshio.DataWriterVTKMesh import DataWriterVTKMesh
+ dataWriter = pyre.inventory.facility("data_writer", factory=DataWriterVTKMesh,
+ family="output_data_writer")
+ dataWriter.meta['tip'] = "Data writer for partition information."
- ## @class Inventory
- ## Python object for managing Distributor facilities and properties.
- ##
- ## \b Properties
- ## @li \b partitioner Name of mesh partitioner {"parmetis", "chaco"}.
- ## @li \b debug Write partition information to file.
- ##
- ## \b Facilities
- ## @li \b writer Data writer for for partition information.
-
- import pyre.inventory
-
- partitioner = pyre.inventory.str("partitioner", default="chaco",
- validator=pyre.inventory.choice(["chaco",
- "parmetis"]))
- partitioner.meta['tip'] = "Name of mesh partitioner."
-
- debug = pyre.inventory.bool("debug", default=False)
- debug.meta['tip'] = "Write partition information to file."
-
- from pylith.meshio.DataWriterVTK import DataWriterVTK
- dataWriter = pyre.inventory.facility("data_writer", factory=DataWriterVTK,
- family="output_data_writer")
- dataWriter.meta['tip'] = "Data writer for partition information."
-
# PUBLIC METHODS /////////////////////////////////////////////////////
def __init__(self, name="partitioner"):
"""
Constructor.
"""
- Component.__init__(self, name, facility="partitioner")
- self.cppHandle = None
- self.debug = False
+ PetscComponent.__init__(self, name, facility="partitioner")
+ ModuleDistributor.__init__(self)
return
- def distribute(self, mesh):
+ def distribute(self, mesh, normalizer):
"""
Distribute a Mesh
"""
@@ -78,18 +72,13 @@
logEvent = "%sdistribute" % self._loggingPrefix
self._logger.eventBegin(logEvent)
- self._createCppHandle()
-
- from Mesh import Mesh
+ from pylith.topology.Mesh import Mesh
newMesh = Mesh()
- newMesh.cppHandle = self.cppHandle.distribute(mesh.cppHandle,
- self.partitioner)
- newMesh.coordsys = mesh.coordsys
+ ModuleDistributor.distribute(newMesh, mesh, self.partitioner)
if self.debug:
- self.dataWriter.initialize()
- self.cppHandle.write(self.dataWriter.cppHandle,
- newMesh.cppHandle, newMesh.coordsys.cppHandle)
+ self.dataWriter.initialize(normalizer)
+ ModuleDistributor.write(self.dataWriter, newMesh)
self._logger.eventEnd(logEvent)
return newMesh
@@ -101,23 +90,13 @@
"""
Set members based using inventory.
"""
- Component._configure(self)
+ PetscComponent._configure(self)
self.partitioner = self.inventory.partitioner
self.debug = self.inventory.debug
self.dataWriter = self.inventory.dataWriter
return
- def _createCppHandle(self):
- """
- Create handle to C++ object.
- """
- if None == self.cppHandle:
- import pylith.topology.topology as bindings
- self.cppHandle = bindings.Distributor()
- return
-
-
def _setupLogging(self):
"""
Setup event logging.
Modified: short/3D/PyLith/branches/pylith-swig/unittests/libtests/faults/TestFaultCohesiveKin.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/libtests/faults/TestFaultCohesiveKin.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/faults/TestFaultCohesiveKin.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -95,7 +95,7 @@
const char* names[] = {"one", "two"};
sources[0] = &eqsrcA;
sources[1] = &eqsrcB;
- fault.eqsrcs(names, sources, nsrcs);
+ fault.eqsrcs(names, nsrcs, sources, nsrcs);
CPPUNIT_ASSERT(&eqsrcA == fault._eqSrcs["one"]);
CPPUNIT_ASSERT(&eqsrcB == fault._eqSrcs["two"]);
delete[] sources; sources = 0;
@@ -878,7 +878,7 @@
fault->label(_data->label);
fault->quadrature(_quadrature);
- fault->eqsrcs(const_cast<const char**>(names), sources, nsrcs);
+ fault->eqsrcs(const_cast<const char**>(names), nsrcs, sources, nsrcs);
fault->adjustTopology(mesh, _flipFault);
const double upDir[] = { 0.0, 0.0, 1.0 };
Modified: short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/Makefile.am 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/Makefile.am 2009-05-01 00:34:23 UTC (rev 14829)
@@ -46,16 +46,14 @@
TestDataWriterVTKSubMeshHex8.cc \
TestOutputManager.cc \
TestOutputSolnSubset.cc \
+ TestDataWriterVTKFaultMesh.cc \
+ TestDataWriterVTKFaultMeshTri3.cc \
+ TestDataWriterVTKFaultMeshQuad4.cc \
+ TestDataWriterVTKFaultMeshTet4.cc \
+ TestDataWriterVTKFaultMeshHex8.cc \
test_meshio.cc
-# TestDataWriterVTKFaultMesh.cc \
-# TestDataWriterVTKFaultMeshTri3.cc \
-# TestDataWriterVTKFaultMeshQuad4.cc \
-# TestDataWriterVTKFaultMeshTet4.cc \
-# TestDataWriterVTKFaultMeshHex8.cc
-
-
noinst_HEADERS = \
TestCellFilterAvg.hh \
TestDataWriterVTK.hh \
Modified: short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.cc 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.cc 2009-05-01 00:34:23 UTC (rev 14829)
@@ -16,34 +16,308 @@
#include "data/DataWriterVTKData.hh" // USES DataWriterVTKData
+#include "pylith/topology/Mesh.hh" // USES Mesh
+#include "pylith/topology/Field.hh" // USES Field
+#include "pylith/topology/Fields.hh" // USES Fields
#include "pylith/meshio/MeshIOAscii.hh" // USES MeshIOAscii
+#include "pylith/meshio/DataWriterVTK.hh" // USES DataWriterVTK
#include "pylith/faults/FaultCohesiveKin.hh" // USES FaultCohesiveKin
#include "pylith/faults/CohesiveTopology.hh" // USES CohesiveTopology
#include <map> // USES std::map
// ----------------------------------------------------------------------
+CPPUNIT_TEST_SUITE_REGISTRATION( pylith::meshio::TestDataWriterVTKFaultMesh );
+
+// ----------------------------------------------------------------------
+typedef pylith::topology::Field<pylith::topology::SubMesh> MeshField;
+
+// ----------------------------------------------------------------------
+// Setup testing data.
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::setUp(void)
+{ // setUp
+ TestDataWriterVTK::setUp();
+ _mesh = new topology::Mesh();
+ _faultMesh = new topology::SubMesh();
+ _flipFault = false;
+} // setUp
+
+// ----------------------------------------------------------------------
+// Tear down testing data.
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::tearDown(void)
+{ // tearDown
+ TestDataWriterVTK::tearDown();
+ delete _mesh; _mesh = 0;
+ delete _faultMesh; _faultMesh = 0;
+} // tearDown
+
+// ----------------------------------------------------------------------
+// Test constructor
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::testConstructor(void)
+{ // testConstructor
+ DataWriterVTK<topology::SubMesh, MeshField> writer;
+
+ CPPUNIT_ASSERT(0 == writer._viewer);
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+} // testConstructor
+
+// ----------------------------------------------------------------------
+// Test openTimeStep() and closeTimeStep()
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::testTimeStep(void)
+{ // testTimeStep
+ CPPUNIT_ASSERT(0 != _mesh);
+ CPPUNIT_ASSERT(0 != _data);
+
+ DataWriterVTK<topology::SubMesh, MeshField> writer;
+
+ writer.filename(_data->timestepFilename);
+ writer.timeFormat(_data->timeFormat);
+
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+
+ const double t = _data->time;
+ const int numTimeSteps = 1;
+ if (0 == _data->cellsLabel) {
+ writer.open(*_faultMesh, numTimeSteps);
+ writer.openTimeStep(t, *_faultMesh);
+ } else {
+ const char* label = _data->cellsLabel;
+ const int id = _data->labelId;
+ writer.open(*_faultMesh, numTimeSteps, label, id);
+ writer.openTimeStep(t, *_faultMesh, label, id);
+ } // else
+
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+
+ writer.closeTimeStep();
+ writer.close();
+
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+
+ checkFile(_data->timestepFilename, t, _data->timeFormat);
+} // testTimeStep
+
+// ----------------------------------------------------------------------
+// Test writeVertexField.
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::testWriteVertexField(void)
+{ // testWriteVertexField
+ CPPUNIT_ASSERT(0 != _mesh);
+ CPPUNIT_ASSERT(0 != _data);
+
+ DataWriterVTK<topology::SubMesh, MeshField> writer;
+
+ topology::Fields<MeshField> vertexFields(*_faultMesh);
+ _createVertexFields(&vertexFields);
+
+ writer.filename(_data->vertexFilename);
+ writer.timeFormat(_data->timeFormat);
+
+ const int nfields = _data->numVertexFields;
+
+ const double t = _data->time;
+ const int numTimeSteps = 1;
+ if (0 == _data->cellsLabel) {
+ writer.open(*_faultMesh, numTimeSteps);
+ writer.openTimeStep(t, *_faultMesh);
+ } else {
+ const char* label = _data->cellsLabel;
+ const int id = _data->labelId;
+ writer.open(*_faultMesh, numTimeSteps, label, id);
+ writer.openTimeStep(t, *_faultMesh, label, id);
+ } // else
+ for (int i=0; i < nfields; ++i) {
+ const MeshField& field = vertexFields.get(_data->vertexFieldsInfo[i].name);
+ writer.writeVertexField(t, field);
+ CPPUNIT_ASSERT(writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+ } // for
+ writer.closeTimeStep();
+ writer.close();
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+
+ checkFile(_data->vertexFilename, t, _data->timeFormat);
+} // testWriteVertexField
+
+// ----------------------------------------------------------------------
+// Test writeCellField.
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::testWriteCellField(void)
+{ // testWriteCellField
+ CPPUNIT_ASSERT(0 != _mesh);
+ CPPUNIT_ASSERT(0 != _data);
+
+ DataWriterVTK<topology::SubMesh, MeshField> writer;
+
+ topology::Fields<MeshField> cellFields(*_faultMesh);
+ _createCellFields(&cellFields);
+
+ writer.filename(_data->cellFilename);
+ writer.timeFormat(_data->timeFormat);
+
+ const int nfields = _data->numCellFields;
+
+ const double t = _data->time;
+ const int numTimeSteps = 1;
+ if (0 == _data->cellsLabel) {
+ writer.open(*_faultMesh, numTimeSteps);
+ writer.openTimeStep(t, *_faultMesh);
+ for (int i=0; i < nfields; ++i) {
+ const MeshField& field = cellFields.get(_data->cellFieldsInfo[i].name);
+ writer.writeCellField(t, field);
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(writer._wroteCellHeader);
+ } // for
+ } else {
+ const char* label = _data->cellsLabel;
+ const int id = _data->labelId;
+ writer.open(*_faultMesh, numTimeSteps, label, id);
+ writer.openTimeStep(t, *_faultMesh, label, id);
+ for (int i=0; i < nfields; ++i) {
+ const MeshField& field = cellFields.get(_data->cellFieldsInfo[i].name);
+ writer.writeCellField(t, field, label, id);
+ CPPUNIT_ASSERT(false == writer._wroteVertexHeader);
+ CPPUNIT_ASSERT(writer._wroteCellHeader);
+ } // for
+ } // else
+ writer.closeTimeStep();
+ writer.close();
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+ CPPUNIT_ASSERT(false == writer._wroteCellHeader);
+
+ checkFile(_data->cellFilename, t, _data->timeFormat);
+} // testWriteCellField
+
+// ----------------------------------------------------------------------
// Initialize mesh.
void
pylith::meshio::TestDataWriterVTKFaultMesh::_initialize(void)
{ // _initialize
CPPUNIT_ASSERT(0 != _data);
+ CPPUNIT_ASSERT(0 != _mesh);
+ CPPUNIT_ASSERT(0 != _faultMesh);
MeshIOAscii iohandler;
iohandler.filename(_data->meshFilename);
- iohandler.read(&_meshDomain);
- CPPUNIT_ASSERT(!_meshDomain.isNull());
+ iohandler.read(_mesh);
faults::FaultCohesiveKin fault;
fault.label(_data->faultLabel);
fault.id(_data->faultId);
- fault.adjustTopology(_meshDomain, _flipFault);
+ fault.adjustTopology(_mesh, _flipFault);
const bool constraintCell = true;
std::map<Mesh::point_type, Mesh::point_type> cohesiveToFault;
- faults::CohesiveTopology::createParallel(&_mesh, &cohesiveToFault,
- _meshDomain, _data->faultId,
- constraintCell);
+ faults::CohesiveTopology::createFaultParallel(_faultMesh, &cohesiveToFault,
+ *_mesh, _data->faultId,
+ constraintCell);
} // _initialize
+// ----------------------------------------------------------------------
+// Create vertex fields.
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::_createVertexFields(
+ topology::Fields<MeshField>* fields) const
+{ // _createVertexFields
+ CPPUNIT_ASSERT(0 != fields);
+ CPPUNIT_ASSERT(0 != _faultMesh);
+ CPPUNIT_ASSERT(0 != _data);
+ try {
+ const int nfields = _data->numVertexFields;
+
+ const ALE::Obj<topology::SubMesh::SieveMesh>& sieveFaultMesh =
+ _faultMesh->sieveMesh();
+ CPPUNIT_ASSERT(!sieveFaultMesh.isNull());
+ const ALE::Obj<topology::SubMesh::SieveMesh::label_sequence>& vertices =
+ sieveFaultMesh->depthStratum(0);
+ CPPUNIT_ASSERT(!vertices.isNull());
+ const topology::SubMesh::SieveMesh::label_sequence::iterator verticesEnd =
+ vertices->end();
+
+ // Set vertex fields
+ for (int i=0; i < nfields; ++i) {
+ const char* name = _data->vertexFieldsInfo[i].name;
+ const int fiberDim = _data->vertexFieldsInfo[i].fiber_dim;
+ fields->add(name, name);
+ MeshField& field = fields->get(name);
+ field.newSection(topology::FieldBase::VERTICES_FIELD, fiberDim);
+ field.allocate();
+ field.vectorFieldType(_data->vertexFieldsInfo[i].field_type);
+
+ const ALE::Obj<topology::SubMesh::RealSection>& section = field.section();
+ CPPUNIT_ASSERT(!section.isNull());
+ int ipt = 0;
+ for (topology::SubMesh::SieveMesh::label_sequence::iterator v_iter=vertices->begin();
+ v_iter != verticesEnd;
+ ++v_iter, ++ipt) {
+ const double* values = &_data->vertexFields[i][ipt*fiberDim];
+ section->updatePoint(*v_iter, values);
+ } // for
+ CPPUNIT_ASSERT_EQUAL(_data->numVertices, ipt);
+ } // for
+ } catch (const ALE::Exception& err) {
+ throw std::runtime_error(err.msg());
+ } // catch
+} // _createVertexFields
+
+// ----------------------------------------------------------------------
+// Create cell fields.
+void
+pylith::meshio::TestDataWriterVTKFaultMesh::_createCellFields(
+ topology::Fields<MeshField>* fields) const
+{ // _createCellFields
+ CPPUNIT_ASSERT(0 != fields);
+ CPPUNIT_ASSERT(0 != _mesh);
+ CPPUNIT_ASSERT(0 != _data);
+
+ try {
+ const int nfields = _data->numCellFields;
+
+ const ALE::Obj<topology::SubMesh::SieveMesh>& sieveFaultMesh =
+ _faultMesh->sieveMesh();
+ CPPUNIT_ASSERT(!sieveFaultMesh.isNull());
+ const ALE::Obj<topology::SubMesh::SieveMesh::label_sequence>& cells =
+ sieveFaultMesh->heightStratum(0);
+ assert(!cells.isNull());
+ const topology::SubMesh::SieveMesh::label_sequence::iterator cellsBegin =
+ cells->begin();
+ const topology::SubMesh::SieveMesh::label_sequence::iterator cellsEnd =
+ cells->end();
+
+ // Set cell fields
+ for (int i=0; i < nfields; ++i) {
+ const char* name = _data->cellFieldsInfo[i].name;
+ const int fiberDim = _data->cellFieldsInfo[i].fiber_dim;
+ fields->add(name, name);
+ MeshField& field = fields->get(name);
+ field.newSection(topology::FieldBase::CELLS_FIELD, fiberDim);
+ field.allocate();
+ field.vectorFieldType(_data->cellFieldsInfo[i].field_type);
+
+ const ALE::Obj<topology::SubMesh::RealSection>& section = field.section();
+ CPPUNIT_ASSERT(!section.isNull());
+ int icell = 0;
+ for (topology::SubMesh::SieveMesh::label_sequence::iterator c_iter=cellsBegin;
+ c_iter != cellsEnd;
+ ++c_iter, ++icell) {
+ const double* values = &_data->cellFields[i][icell*fiberDim];
+ section->updatePoint(*c_iter, values);
+ } // for
+ CPPUNIT_ASSERT_EQUAL(_data->numCells, icell);
+ } // for
+ } catch (const ALE::Exception& err) {
+ throw std::runtime_error(err.msg());
+ } // catch
+} // _createCellFields
+
+
// End of file
Modified: short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.hh 2009-04-30 23:46:54 UTC (rev 14828)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestDataWriterVTKFaultMesh.hh 2009-05-01 00:34:23 UTC (rev 14829)
@@ -23,6 +23,8 @@
#include "TestDataWriterVTK.hh"
+#include "pylith/topology/topologyfwd.hh" // USES Mesh, SubMesh, Field
+
/// Namespace for pylith package
namespace pylith {
namespace meshio {
@@ -34,6 +36,34 @@
class pylith::meshio::TestDataWriterVTKFaultMesh : public TestDataWriterVTK
{ // class TestDataWriterVTKFaultMesh
+ // CPPUNIT TEST SUITE /////////////////////////////////////////////////
+ CPPUNIT_TEST_SUITE( TestDataWriterVTKFaultMesh );
+
+ CPPUNIT_TEST( testConstructor );
+
+ CPPUNIT_TEST_SUITE_END();
+
+ // PUBLIC METHODS /////////////////////////////////////////////////////
+public :
+
+ /// Setup testing data.
+ void setUp(void);
+
+ /// Tear down testing data.
+ void tearDown(void);
+
+ /// Test constructor
+ void testConstructor(void);
+
+ /// Test openTimeStep() and closeTimeStep()
+ void testTimeStep(void);
+
+ /// Test writeVertexField.
+ void testWriteVertexField(void);
+
+ /// Test writeCellField.
+ void testWriteCellField(void);
+
// PROTECTED MEMBERS //////////////////////////////////////////////////
protected :
@@ -43,8 +73,27 @@
// PROTECTED MEMBERS //////////////////////////////////////////////////
protected :
- ALE::Obj<Mesh> _meshDomain; ///< Mesh for domain.
+ topology::Mesh* _mesh; ///< Mesh for domain
+ topology::SubMesh* _faultMesh; ///< Fault mesh.
+ bool _flipFault; ///< If true, flip fault orientation.
+ // PRIVATE MEMBERS ////////////////////////////////////////////////////
+private :
+
+ /** Create vertex fields.
+ *
+ * @param fields Vertex fields.
+ */
+ void
+ _createVertexFields(topology::Fields<topology::Field<topology::SubMesh> >* fields) const;
+
+ /** Create cell fields.
+ *
+ * @param fields Cell fields.
+ */
+ void
+ _createCellFields(topology::Fields<topology::Field<topology::SubMesh> >* fields) const;
+
}; // class TestDataWriterVTKFaultMesh
#endif // pylith_meshio_testdatawritervtkfaultmesh_hh
More information about the CIG-COMMITS
mailing list