[cig-commits] r14584 - in short/3D/PyLith/branches/pylith-swig: . libsrc libsrc/meshio libsrc/topology modulesrc/topology unittests/libtests/meshio

brad at geodynamics.org brad at geodynamics.org
Fri Apr 3 17:02:42 PDT 2009


Author: brad
Date: 2009-04-03 17:02:41 -0700 (Fri, 03 Apr 2009)
New Revision: 14584

Modified:
   short/3D/PyLith/branches/pylith-swig/TODO
   short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/DataWriter.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/meshio/meshiofwd.hh
   short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.cc
   short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.hh
   short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Field.i
   short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/Makefile.am
   short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestCellFilterAvg.cc
   short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestVertexFilterVecNorm.cc
Log:
Finished updating output. Started updating output unit tests.

Modified: short/3D/PyLith/branches/pylith-swig/TODO
===================================================================
--- short/3D/PyLith/branches/pylith-swig/TODO	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/TODO	2009-04-04 00:02:41 UTC (rev 14584)
@@ -7,8 +7,6 @@
 0. SWIG conversion [Brad]
 
   (1) Output
-    OutputSolnSubset
-    VertexFilterNorm
 
   (2) Full-scale tests (few 1-D, 2-D, and 3-D)
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am	2009-04-04 00:02:41 UTC (rev 14584)
@@ -75,6 +75,7 @@
 	meshio/PsetFile.cc \
 	meshio/PsetFileAscii.cc \
 	meshio/PsetFileBinary.cc \
+	meshio/OutputSolnSubset.cc \
 	problems/Formulation.cc \
 	problems/Solver.cc \
 	problems/SolverLinear.cc \
@@ -105,7 +106,6 @@
 #	meshio/DataWriter.cc \
 # 	meshio/DataWriterVTK.cc \
 # 	meshio/OutputManager.cc \
-# 	meshio/OutputSolnSubset.cc \
 # 	meshio/VertexFilter.cc \
 # 	meshio/VertexFilterVecNorm.cc \
 # 	meshio/UCDFaultFile.cc \

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -37,16 +37,17 @@
   _quadrature(0)
 { // copy constructor
   if (0 != f._quadrature)
-    _quadrature = f._quadrature->clone();
+    _quadrature = new feassemble::Quadrature<mesh_type>(*f._quadrature);
 } // copy constructor
 
 // ----------------------------------------------------------------------
 // Set quadrature associated with cells.
 template<typename mesh_type>
 void
-pylith::meshio::CellFilter<mesh_type>::quadrature(const feassemble::Quadrature* q)
+pylith::meshio::CellFilter<mesh_type>::quadrature(const feassemble::Quadrature<mesh_type>* q)
 { // quadrature
-    delete _quadrature; _quadrature = (0 != q) ? q->clone() : 0;
+  delete _quadrature; 
+  _quadrature = (0 != q) ? new feassemble::Quadrature<mesh_type>(*q) : 0;
 } // quadrature
 
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilter.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -23,7 +23,7 @@
 // Include directives ---------------------------------------------------
 #include "meshiofwd.hh" // forward declarations
 
-#include "pylith/feassemble/topologyfwd.hh" // HOLDSA Quadrature<Mesh>
+#include "pylith/topology/topologyfwd.hh" // HOLDSA Quadrature<Mesh>
 #include "pylith/feassemble/feassemblefwd.hh" // HOLDSA Quadrature<Mesh>
 
 // CellFilter -----------------------------------------------------------

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -17,7 +17,8 @@
 // ----------------------------------------------------------------------
 // Constructor
 template<typename mesh_type>
-pylith::meshio::CellFilterAvg<mesh_type>::CellFilterAvg(void)
+pylith::meshio::CellFilterAvg<mesh_type>::CellFilterAvg(void) :
+  _fieldAvg(0)
 { // constructor
 } // constructor
 
@@ -26,13 +27,15 @@
 template<typename mesh_type>
 pylith::meshio::CellFilterAvg<mesh_type>::~CellFilterAvg(void)
 { // destructor
+  delete _fieldAvg; _fieldAvg = 0;
 } // destructor  
 
 // ----------------------------------------------------------------------
 // Copy constructor.
 template<typename mesh_type>
 pylith::meshio::CellFilterAvg<mesh_type>::CellFilterAvg(const CellFilterAvg& f) :
-  CellFilter(f)
+  CellFilter<mesh_type>(f),
+  _fieldAvg(0)
 { // copy constructor
 } // copy constructor
 
@@ -54,19 +57,23 @@
 				  const char* label,
 				  const int labelId)
 { // filter
-  assert(0 != _quadrature);
+  typedef typename mesh_type::RealSection RealSection;
+  typedef typename mesh_type::SieveMesh SieveMesh;
+  typedef typename SieveMesh::label_sequence label_sequence;
 
-  const int numQuadPts = _quadrature->numQuadPts();
-  const double_array& wts = _quadrature->quadWts();
+  assert(0 != CellFilter<mesh_type>::_quadrature);
+
+  const int numQuadPts = CellFilter<mesh_type>::_quadrature->numQuadPts();
+  const double_array& wts = CellFilter<mesh_type>::_quadrature->quadWts();
   
   const ALE::Obj<SieveMesh>& sieveMesh = fieldIn.mesh().sieveMesh();
   assert(!sieveMesh.isNull());
 
-  const ALE::Obj<SieveMesh::label_sequence>& cells = (0 == label) ?
-    mesh->heightStratum(0) :
-    mesh->getLabelStratum(label, labelId);
+  const ALE::Obj<label_sequence>& cells = (0 == label) ?
+    sieveMesh->heightStratum(0) :
+    sieveMesh->getLabelStratum(label, labelId);
   assert(!cells.isNull());
-  const SieveMesh::label_sequence::iterator cellsEnd = cells->end();
+  const typename label_sequence::iterator cellsEnd = cells->end();
 
   // Only processors with cells for output get the correct fiber dimension.
   const ALE::Obj<RealSection>& sectionIn = fieldIn.section();
@@ -79,28 +86,28 @@
   if (0 == _fieldAvg) {
     _fieldAvg = new topology::Field<mesh_type>(fieldIn.mesh());
     assert(0 != _fieldAvg);
-    _fieldAvg->newSection(fieldIn->getChart(), fiberDim);
+    _fieldAvg->newSection(sectionIn->getChart(), fiberDim);
     _fieldAvg->allocate();
 
     _fieldAvg->label(fieldIn.label());
     switch (fieldIn.vectorFieldType())
       { // switch
-      case FieldBase::MULTI_SCALAR:
-	_fieldAvg->vectorFieldType(FieldBase::SCALAR);
+      case topology::FieldBase::MULTI_SCALAR:
+	_fieldAvg->vectorFieldType(topology::FieldBase::SCALAR);
 	break;
-      case FieldBase::MULTI_VECTOR:
-	_fieldAvg->vectorFieldType(FieldBase::VECTOR);
+      case topology::FieldBase::MULTI_VECTOR:
+	_fieldAvg->vectorFieldType(topology::FieldBase::VECTOR);
 	break;
-      case FieldBase::MULTI_TENSOR:
-	_fieldAvg->vectorFieldType(FieldBase::TENSOR);
+      case topology::FieldBase::MULTI_TENSOR:
+	_fieldAvg->vectorFieldType(topology::FieldBase::TENSOR);
 	break;
-      case FieldBase::MULTI_OTHER:
-	_fieldAvg->vectorFieldType(FieldBase::OTHER);
+      case topology::FieldBase::MULTI_OTHER:
+	_fieldAvg->vectorFieldType(topology::FieldBase::OTHER);
 	break;
-      case FieldBase::SCALAR:
-      case FieldBase::VECTOR:
-      case FieldBase::TENSOR:
-      case FieldBase::OTHER:
+      case topology::FieldBase::SCALAR:
+      case topology::FieldBase::VECTOR:
+      case topology::FieldBase::TENSOR:
+      case topology::FieldBase::OTHER:
       default :
 	std::cerr << "Bad vector field type for CellFilterAvg." << std::endl;
 	assert(0);
@@ -115,7 +122,7 @@
     scalar += wts[iQuad];
 
   // Loop over cells
-  for (SieveMesh::label_sequence::iterator c_iter=cells->begin();
+  for (typename label_sequence::iterator c_iter=cells->begin();
        c_iter != cellsEnd;
        ++c_iter) {
     const double* values = sectionIn->restrictPoint(*c_iter);
@@ -125,7 +132,7 @@
       for (int i=0; i < fiberDim; ++i)
 	fieldAvgCell[i] += wts[iQuad] / scalar * values[iQuad*fiberDim+i];
 
-    _sectionAvg->updatePoint(*c_iter, &fieldAvgCell[0]);
+    sectionAvg->updatePoint(*c_iter, &fieldAvgCell[0]);
     PetscLogFlops( numQuadPts*fiberDim*3 );
   } // for
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/CellFilterAvg.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -19,12 +19,13 @@
 
 #if !defined(pylith_meshio_cellfilteravg_hh)
 #define pylith_meshio_cellfilteravg_hh
+
 // Include directives ---------------------------------------------------
 #include "CellFilter.hh" // ISA CellFilter
 
 // CellFilter -----------------------------------------------------------
 template<typename mesh_type>
-class pylith::meshio::CellFilterAvg : public CellFilter
+class pylith::meshio::CellFilterAvg : public CellFilter<mesh_type>
 { // CellFilterAvg
 
 // PUBLIC METHODS ///////////////////////////////////////////////////////
@@ -40,7 +41,7 @@
    *
    * @returns Copy of filter.
    */
-  CellFilter* clone(void) const;
+  CellFilter<mesh_type>* clone(void) const;
 
   /** Filter field over cells.
    *
@@ -69,7 +70,7 @@
 private :
 
   /// Not implemented.
-  const CellFilter& operator=(const CellFilter&);
+  const CellFilterAvg& operator=(const CellFilterAvg&);
 
 // PRIVATE MEMBERS //////////////////////////////////////////////////////
 private :

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/DataWriter.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/DataWriter.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/DataWriter.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -22,6 +22,8 @@
 // Include directives ---------------------------------------------------
 #include "meshiofwd.hh" // forward declarations
 
+#include "pylith/topology/topologyfwd.hh" // USES Field
+
 // DataWriter -----------------------------------------------------------
 template<typename mesh_type>
 class pylith::meshio::DataWriter

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -12,6 +12,7 @@
 
 #include <portinfo>
 
+#include "DataWriter.hh" // USES DataWriter
 #include "VertexFilter.hh" // USES VertexFilter
 #include "CellFilter.hh" // USES CellFilter
 
@@ -43,7 +44,8 @@
 // Set coordinate system in output. The vertex fields in the output
 template<typename mesh_type>
 void
-pylith::meshio::OutputManager<mesh_type>::coordsys(const spatialdata::geocoords::CoordSys* cs)
+pylith::meshio::OutputManager<mesh_type>::coordsys(
+				  const spatialdata::geocoords::CoordSys* cs)
 { // coordsys
   delete _coordsys; _coordsys = (0 != cs) ? cs->clone() : 0;
 } // coordsys
@@ -52,7 +54,8 @@
 // Set writer to write data to file.
 template<typename mesh_type>
 void
-pylith::meshio::OutputManager<mesh_type>::writer(const DataWriter* datawriter)
+pylith::meshio::OutputManager<mesh_type>::writer(
+				       const DataWriter<mesh_type>* datawriter)
 { // writer
   delete _writer; _writer = (0 != datawriter) ? datawriter->clone() : 0;
 } // writer
@@ -61,7 +64,8 @@
 // Set filter for vertex data.
 template<typename mesh_type>
 void
-pylith::meshio::OutputManager<mesh_type>::vertexFilter(const VertexFilter* filter)
+pylith::meshio::OutputManager<mesh_type>::vertexFilter(
+					const VertexFilter<mesh_type>* filter)
 { // vertexFilter
   delete _vertexFilter; _vertexFilter = (0 != filter) ? filter->clone() : 0;
 } // vertexFilter
@@ -70,7 +74,8 @@
 // Set filter for cell data.
 template<typename mesh_type>
 void
-pylith::meshio::OutputManager<mesh_type>::cellFilter(const CellFilter* filter)
+pylith::meshio::OutputManager<mesh_type>::cellFilter(
+				         const CellFilter<mesh_type>* filter)
 { // cellFilter
   delete _cellFilter; _cellFilter = (0 != filter) ? filter->clone() : 0;
 } // cellFilter

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputManager.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -19,17 +19,14 @@
 #if !defined(pylith_meshio_outputmanager_hh)
 #define pylith_meshio_outputmanager_hh
 
-
-
-#include "pylith/utils/sievetypes.hh" // USES PETSc Mesh, real_section_type
-#include "pylith/utils/vectorfields.hh" // USES VectorFieldEnum
-
 // Include directives ---------------------------------------------------
 #include "meshiofwd.hh" // forward declarations
 
-#include "DataWriter.hh" // USES DataWriter in templated methods
+#include "pylith/topology/topologyfwd.hh" // USES Field
+#include "spatialdata/geocoords/geocoordsfwd.hh" // USES CoordSys
 
 // OutputManager --------------------------------------------------------
+template<typename mesh_type>
 class pylith::meshio::OutputManager
 { // OutputManager
   friend class TestOutputManager; // unit testing
@@ -61,13 +58,13 @@
    *
    * @param filter Filter to apply to vertex data before writing.
    */
-  void vertexFilter(const VertexFilter* filter);
+  void vertexFilter(const VertexFilter<mesh_type>* filter);
 
   /** Set filter for cell data.
    *
    * @param filter Filter to apply to cell data before writing.
    */
-  void cellFilter(const CellFilter* filter);
+  void cellFilter(const CellFilter<mesh_type>* filter);
 
   /** Prepare for output.
    *
@@ -77,7 +74,6 @@
    *   (=0 means use all cells in mesh).
    * @param labelId Value of label defining which cells to include.
    */
-  template<typename mesh_type>
   void open(const mesh_type& mesh,
 	    const int numTimeSteps,
 	    const char* label =0,
@@ -94,7 +90,6 @@
    *   (=0 means use all cells in mesh).
    * @param labelId Value of label defining which cells to include.
    */
-  template<typename mesh_type>
   void openTimeStep(const double t,
 		    const mesh_type& mesh,
 		    const char* label =0,
@@ -108,24 +103,18 @@
    * @param t Time associated with field.
    * @param field Vertex field.
    */
-  template<typename mesh_type>
   void appendVertexField(const double t,
 			 const topology::Field<mesh_type>& field);
 
   /** Append finite-element cell field to file.
    *
    * @param t Time associated with field.
-   * @param name Name of field.
    * @param field Cell field.
-   * @param fieldType Type of field.
-   * @param mesh PETSc mesh object.
    * @param label Name of label defining cells to include in output
    *   (=0 means use all cells in mesh).
    * @param labelId Value of label defining which cells to include.
    */
-  template<typename mesh_type>
   void appendCellField(const double t,
-		       const char* name,
 		       const topology::Field<mesh_type>& field,
 		       const char* label =0,
 		       const int labelId =0);
@@ -139,14 +128,16 @@
 // PRIVATE MEMBERS //////////////////////////////////////////////////////
 private :
 
-  spatialdata::geocoords::CoordSys* _coordsys; ///< Coordinate system for output.
+  /// Coordinate system for output.
+  spatialdata::geocoords::CoordSys* _coordsys;
+
   DataWriter<mesh_type>* _writer; ///< Writer for data.
-  VertexFilter* _vertexFilter; ///< Filter applied to vertex data.
-  CellFilter* _cellFilter; ///< Filter applied to cell data.
+  VertexFilter<mesh_type>* _vertexFilter; ///< Filter applied to vertex data.
+  CellFilter<mesh_type>* _cellFilter; ///< Filter applied to cell data.
 
 }; // OutputManager
 
-#include "OutputManager.icc" // template methods
+#include "OutputManager.cc" // template methods
 
 #endif // pylith_meshio_outputmanager_hh
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -14,12 +14,13 @@
 
 #include "OutputSolnSubset.hh" // implementation of class methods
 
-#include <Selection.hh> // USES submesh algorithms
+#include "pylith/topology/Mesh.hh" // USES Mesh
 
 // ----------------------------------------------------------------------
 // Constructor
 pylith::meshio::OutputSolnSubset::OutputSolnSubset(void) :
-  _label("")
+  _label(""),
+  _submesh(0)
 { // constructor
 } // constructor
 
@@ -27,6 +28,7 @@
 // Destructor
 pylith::meshio::OutputSolnSubset::~OutputSolnSubset(void)
 { // destructor
+  delete _submesh; _submesh = 0;
 } // destructor  
 
 // ----------------------------------------------------------------------
@@ -40,11 +42,12 @@
 // ----------------------------------------------------------------------
 // Verify configuration is acceptable.
 void
-pylith::meshio::OutputSolnSubset::verifyConfiguration(const ALE::Obj<Mesh>& mesh) const
+pylith::meshio::OutputSolnSubset::verifyConfiguration(const topology::Mesh& mesh) const
 { // verifyConfiguration
-  assert(!mesh.isNull());
+  const ALE::Obj<topology::Mesh::SieveMesh>& sieveMesh = mesh.sieveMesh();
+  assert(!sieveMesh.isNull());
 
-  if (!mesh->hasIntSection(_label)) {
+  if (!sieveMesh->hasIntSection(_label)) {
     std::ostringstream msg;
     msg << "Mesh missing group of vertices '" << _label
 	<< " for subdomain output.";
@@ -54,36 +57,12 @@
 
 // ----------------------------------------------------------------------
 // Get mesh associated with subdomain.
-const ALE::Obj<pylith::Mesh>&
-pylith::meshio::OutputSolnSubset::subdomainMesh(const ALE::Obj<Mesh>& mesh)
+const pylith::topology::SubMesh&
+pylith::meshio::OutputSolnSubset::subdomainMesh(const topology::Mesh& mesh)
 { // subdomainMesh
-  _mesh =
-    ALE::Selection<Mesh>::submeshV<SubMesh>(mesh, mesh->getIntSection(_label));
-  if (_mesh.isNull()) {
-    std::ostringstream msg;
-    msg << "Could not construct mesh of subdomain " << _label << "'.";
-    throw std::runtime_error(msg.str());
-  } // if
-  _mesh->setRealSection("coordinates", 
-			mesh->getRealSection("coordinates"));
-  // Create the parallel overlap
-  const Obj<Mesh::sieve_type>& sieve                   = _mesh->getSieve();
-  Obj<Mesh::send_overlap_type> sendParallelMeshOverlap = _mesh->getSendOverlap();
-  Obj<Mesh::recv_overlap_type> recvParallelMeshOverlap = _mesh->getRecvOverlap();
-  Mesh::renumbering_type&      renumbering             = _mesh->getRenumbering();
-  Mesh::renumbering_type&      oldRenumbering          = mesh->getRenumbering();
-  for(Mesh::renumbering_type::const_iterator r_iter = oldRenumbering.begin(); r_iter != oldRenumbering.end(); ++r_iter) {
-    if (sieve->getChart().hasPoint(r_iter->second) && (sieve->getConeSize(r_iter->second) || sieve->getSupportSize(r_iter->second))) {
-      renumbering[r_iter->first] = r_iter->second;
-    }
-  }
-  //   Can I figure this out in a nicer way?
-  ALE::SetFromMap<std::map<Mesh::point_type,Mesh::point_type> > globalPoints(renumbering);
-
-  ALE::OverlapBuilder<>::constructOverlap(globalPoints, renumbering, sendParallelMeshOverlap, recvParallelMeshOverlap);
-  _mesh->setCalculatedOverlap(true);
-
-  return _mesh;
+  delete _submesh; _submesh = new topology::SubMesh(mesh, _label.c_str());
+  assert(0 != _submesh);
+  return *_submesh;
 } // subdomainMesh
 
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/OutputSolnSubset.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -20,20 +20,16 @@
 #if !defined(pylith_meshio_outputsolnsubset_hh)
 #define pylith_meshio_outputsolnsubset_hh
 
+// Include directives ---------------------------------------------------
+#include "meshiofwd.hh" // forward declarations
+
+#include "pylith/topology/SubMesh.hh" // ISA OutputManager<SubMesh>
 #include "OutputManager.hh" // ISA OutputManager
 
-#include "pylith/utils/sievetypes.hh" // HASA PETSc Mesh
 #include <string> // HASA std::string
 
-namespace pylith {
-  namespace meshio {
-    class OutputSolnSubset;
-
-    class TestOutputSolnSubset; // unit testing
-  } // meshio
-} // pylith
-
-class pylith::meshio::OutputSolnSubset : public OutputManager
+// OutputSolnSubset -----------------------------------------------------
+class pylith::meshio::OutputSolnSubset : public OutputManager<topology::SubMesh>
 { // OutputSolnSubset
   friend class TestOutputSolnSubset; // unit testing
 
@@ -56,13 +52,13 @@
    *
    * @param mesh PETSc mesh
    */
-  void verifyConfiguration(const ALE::Obj<Mesh>& mesh) const;
+  void verifyConfiguration(const topology::Mesh& mesh) const;
 
   /** Get mesh associated with subdomain.
    *
    * @returns Mesh associated with subdomain.
    */
-  const ALE::Obj<Mesh>& subdomainMesh(const ALE::Obj<Mesh>& mesh);
+  const topology::SubMesh& subdomainMesh(const topology::Mesh& mesh);
   
 
 // NOT IMPLEMENTED //////////////////////////////////////////////////////
@@ -75,7 +71,7 @@
 private :
 
   std::string _label; ///< Label of subdomain.
-  ALE::Obj<SubMesh> _mesh; ///< Mesh of subdomain.
+  topology::SubMesh* _submesh; ///< Mesh of subdomain.
 
 }; // OutputSolnSubset
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -36,7 +36,7 @@
 // ----------------------------------------------------------------------
 // operator=.
 template<typename mesh_type>
-const pylith::meshio::VertexFilter&
+const pylith::meshio::VertexFilter<mesh_type>&
 pylith::meshio::VertexFilter<mesh_type>::operator=(const VertexFilter& f)
 { // operator=
 } // operator=

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilter.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -23,8 +23,7 @@
 // Include directives ---------------------------------------------------
 #include "meshiofwd.hh" // forward declarations
 
-#include "pylith/feassemble/topologyfwd.hh" // HOLDSA Quadrature<Mesh>
-#include "pylith/feassemble/feassemblefwd.hh" // HOLDSA Quadrature<Mesh>
+#include "pylith/topology/topologyfwd.hh" // USES Field
 
 // VertexFilter ---------------------------------------------------------
 template<typename mesh_type>
@@ -65,6 +64,7 @@
    */
   VertexFilter(const VertexFilter& f);
 
+private :
   /** operator=.
   *
   * @param f Filter to copy.

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -12,80 +12,110 @@
 
 #include <portinfo>
 
-#include "VertexFilterVecNorm.hh" // implementation of class methods
+#include "pylith/topology/Field.hh" // USES Field
 
 // ----------------------------------------------------------------------
 // Constructor
-pylith::meshio::VertexFilterVecNorm::VertexFilterVecNorm(void)
+template<typename mesh_type>
+pylith::meshio::VertexFilterVecNorm<mesh_type>::VertexFilterVecNorm(void) :
+  _fieldVecNorm(0)
 { // constructor
 } // constructor
 
 // ----------------------------------------------------------------------
 // Destructor
-pylith::meshio::VertexFilterVecNorm::~VertexFilterVecNorm(void)
+template<typename mesh_type>
+pylith::meshio::VertexFilterVecNorm<mesh_type>::~VertexFilterVecNorm(void)
 { // destructor
+  delete _fieldVecNorm; _fieldVecNorm = 0;
 } // destructor  
 
 // ----------------------------------------------------------------------
 // Copy constructor.
-pylith::meshio::VertexFilterVecNorm::VertexFilterVecNorm(const VertexFilterVecNorm& f) :
-  VertexFilter(f)
+template<typename mesh_type>
+pylith::meshio::VertexFilterVecNorm<mesh_type>::VertexFilterVecNorm(const VertexFilterVecNorm& f) :
+  VertexFilter<mesh_type>(f)
 { // copy constructor
 } // copy constructor
 
 // ----------------------------------------------------------------------
 // Create copy of filter.
-pylith::meshio::VertexFilter*
-pylith::meshio::VertexFilterVecNorm::clone(void) const
+template<typename mesh_type>
+pylith::meshio::VertexFilter<mesh_type>*
+pylith::meshio::VertexFilterVecNorm<mesh_type>::clone(void) const
 { // clone
   return new VertexFilterVecNorm(*this);
 } // clone
 
 // ----------------------------------------------------------------------
 // Filter field.
-const ALE::Obj<pylith::real_section_type>&
-pylith::meshio::VertexFilterVecNorm::filter(
-				  VectorFieldEnum* fieldType,
-				  const ALE::Obj<real_section_type>& fieldIn,
-				  const ALE::Obj<Mesh>& mesh)
+template<typename mesh_type>
+const pylith::topology::Field<mesh_type>&
+pylith::meshio::VertexFilterVecNorm<mesh_type>::filter(
+				   const topology::Field<mesh_type>& fieldIn)
 { // filter
-  assert(0 != fieldType);
-  *fieldType = SCALAR_FIELD;
+  typedef typename mesh_type::RealSection RealSection;
+  typedef typename mesh_type::SieveMesh SieveMesh;
+  typedef typename SieveMesh::label_sequence label_sequence;
 
-  const ALE::Obj<Mesh::label_sequence>& vertices = mesh->depthStratum(0);
+  const ALE::Obj<SieveMesh>& sieveMesh = fieldIn.mesh().sieveMesh();
+  assert(!sieveMesh.isNull());
+
+  const ALE::Obj<label_sequence>& vertices = sieveMesh->depthStratum(0);
   assert(!vertices.isNull());
-  const Mesh::label_sequence::iterator verticesEnd = vertices->end();
+  const typename label_sequence::iterator verticesEnd = vertices->end();
 
-  const int fiberDim = fieldIn->getFiberDimension(*vertices->begin());
+  const ALE::Obj<RealSection>& sectionIn = fieldIn.section();
+  assert(!sectionIn.isNull());
+  const int fiberDimIn = sectionIn->getFiberDimension(*vertices->begin());
+  const int fiberDimNorm = 1;
 
   // Allocation field if necessary
-  if (_fieldVecNorm.isNull() ||
-      1 != _fieldVecNorm->getFiberDimension(*vertices->begin())) {
-    _fieldVecNorm = new real_section_type(mesh->comm(), mesh->debug());
-    _fieldVecNorm->setChart(real_section_type::chart_type(*std::min_element(vertices->begin(), vertices->end()), *std::max_element(vertices->begin(), vertices->end())+1));
-    _fieldVecNorm->setFiberDimension(vertices, 1);
-    mesh->allocate(_fieldVecNorm);
+  if (0 == _fieldVecNorm) {
+    _fieldVecNorm = new topology::Field<mesh_type>(fieldIn.mesh());
+    _fieldVecNorm->newSection(sectionIn->getChart(), fiberDimNorm);
+    _fieldVecNorm->allocate();
+
+    _fieldVecNorm->label(fieldIn.label());    
+    switch (fieldIn.vectorFieldType())
+      { // switch
+      case topology::FieldBase::SCALAR:
+      case topology::FieldBase::VECTOR:
+	_fieldVecNorm->vectorFieldType(topology::FieldBase::SCALAR);
+	break;
+      case topology::FieldBase::MULTI_SCALAR:
+      case topology::FieldBase::MULTI_VECTOR:
+      case topology::FieldBase::MULTI_TENSOR:
+      case topology::FieldBase::MULTI_OTHER:
+      case topology::FieldBase::TENSOR:
+      case topology::FieldBase::OTHER:
+      default :
+	std::cerr << "Bad vector field type for VertexFilterVecNorm." << std::endl;
+	assert(0);
+      } // switch
   } // if
 
+  const ALE::Obj<RealSection>& sectionNorm = 
+    _fieldVecNorm->section();
+  assert(!sectionNorm.isNull());
+
   double norm = 0.0;
-
   // Loop over vertices
-  for (Mesh::label_sequence::iterator v_iter=vertices->begin();
+  for (typename label_sequence::iterator v_iter=vertices->begin();
        v_iter != verticesEnd;
        ++v_iter) {
-    const real_section_type::value_type* values = 
-      fieldIn->restrictPoint(*v_iter);
+    const double* values = sectionIn->restrictPoint(*v_iter);
     
     norm = 0.0;
-    for (int i=0; i < fiberDim; ++i)
+    for (int i=0; i < fiberDimIn; ++i)
       norm += values[i]*values[i];
     norm = sqrt(norm);
 
-    _fieldVecNorm->updatePoint(*v_iter, &norm);
-    PetscLogFlops( 1 + fiberDim*2 );
+    sectionNorm->updatePoint(*v_iter, &norm);
   } // for
+  PetscLogFlops(vertices->size() * (1 + 2*fiberDimIn) );
 
-  return _fieldVecNorm;
+  return *_fieldVecNorm;
 } // filter
 
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/VertexFilterVecNorm.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -20,15 +20,12 @@
 #if !defined(pylith_meshio_cellfiltervecnorm_hh)
 #define pylith_meshio_cellfiltervecnorm_hh
 
+// Include directives ---------------------------------------------------
 #include "VertexFilter.hh" // ISA VertexFilter
 
-namespace pylith {
-  namespace meshio {
-    class VertexFilterVecNorm;
-  } // meshio
-} // pylith
-
-class pylith::meshio::VertexFilterVecNorm : public VertexFilter
+// CellFilter -----------------------------------------------------------
+template<typename mesh_type>
+class pylith::meshio::VertexFilterVecNorm : public VertexFilter<mesh_type>
 { // VertexFilterVecNorm
 
 // PUBLIC METHODS ///////////////////////////////////////////////////////
@@ -44,18 +41,14 @@
    *
    * @returns Copy of filter.
    */
-  VertexFilter* clone(void) const;
+  VertexFilter<mesh_type>* clone(void) const;
 
-  /** Filter field. Field type of filtered field is returned via an argument.
+  /** Filter vertex field.
    *
-   * @param fieldType Field type of filtered field.
    * @param fieldIn Field to filter.
-   * @param mesh PETSc mesh.
    */
-  const ALE::Obj<real_section_type>&
-  filter(VectorFieldEnum* fieldType,
-	 const ALE::Obj<real_section_type>& fieldIn,
-	 const ALE::Obj<Mesh>& mesh);
+  const topology::Field<mesh_type>&
+  filter(const topology::Field<mesh_type>& fieldIn);
 
 // PROTECTED METHODS ////////////////////////////////////////////////////
 protected :
@@ -71,15 +64,17 @@
 private :
 
   /// Not implemented.
-  const VertexFilter& operator=(const VertexFilter&);
+  const VertexFilterVecNorm& operator=(const VertexFilterVecNorm&);
 
 // PRIVATE MEMBERS //////////////////////////////////////////////////////
 private :
 
-  ALE::Obj<real_section_type> _fieldVecNorm; ///< Filtered vertex field
+  topology::Field<mesh_type>* _fieldVecNorm; ///< Filtered vertex field
 
 }; // VertexFilterVecNorm
 
+#include "VertexFilterVecNorm.cc" // template definitions
+
 #endif // pylith_meshio_cellfiltervecnorm_hh
 
 

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/meshio/meshiofwd.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/meshio/meshiofwd.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/meshio/meshiofwd.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -41,12 +41,12 @@
     template<typename mesh_type> class OutputManager;
     template<typename mesh_type> class DataWriter;
     template<typename mesh_type> class DataWriterVTK;
+    template<typename mesh_type> class CellFilter;
+    template<typename mesh_type> class CellFilterAvg;
+    template<typename mesh_type> class VertexFilter;
+    template<typename mesh_type> class VertexFilterVecNorm;
+    class OutputSolnSubset;
 
-    class CellFilter;
-    class CellFilterAvg;
-    class VertexFilter;
-    class VertexFilterVecNorm;
-
     class UCDFaultFile;
 
   } // meshio

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -122,16 +122,17 @@
 template<typename mesh_type>
 void
 pylith::topology::Field<mesh_type>::newSection(const DomainEnum domain,
-				    const int fiberDim)
+					       const int fiberDim,
+					       const int stratum)
 { // newSection
   const ALE::Obj<SieveMesh>& sieveMesh = _mesh.sieveMesh();
   assert(!sieveMesh.isNull());
 
   ALE::Obj<label_sequence> points;
   if (VERTICES_FIELD == domain)
-    points = sieveMesh->depthStratum(0);
+    points = sieveMesh->depthStratum(stratum);
   else if (CELLS_FIELD == domain)
-    points = sieveMesh->heightStratum(1);
+    points = sieveMesh->heightStratum(stratum);
   else {
     std::cerr << "Unknown value for DomainEnum: " << domain << std::endl;
     assert(0);

Modified: short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.hh	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/topology/Field.hh	2009-04-04 00:02:41 UTC (rev 14584)
@@ -147,9 +147,11 @@
    *
    * @param domain Type of points over which to define section.
    * @param dim Fiber dimension for section.
+   * @param stratum Stratum depth (for vertices) and height (for cells).
    */
   void newSection(const DomainEnum domain,
-		  const int fiberDim);
+		  const int fiberDim,
+		  const int stratum =0);
 
   /** Create section given chart. This allows a chart to be reused
    * across multiple fields, reducing memory usage.

Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Field.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Field.i	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/topology/Field.i	2009-04-04 00:02:41 UTC (rev 14584)
@@ -114,9 +114,11 @@
        *
        * @param domain Type of points over which to define section.
        * @param dim Fiber dimension for section.
+       * @param stratum Stratum depth (for vertices) and height (for cells).
        */
       void newSection(const pylith::topology::FieldBase::DomainEnum domain,
-		      const int fiberDim);
+		      const int fiberDim,
+		      const int stratum =0);
 
       /** Create section with same layout (fiber dimension and
        * constraints) as another section. This allows the layout data

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-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/Makefile.am	2009-04-04 00:02:41 UTC (rev 14584)
@@ -24,10 +24,11 @@
 	TestMeshIO.cc \
 	TestMeshIOAscii.cc \
 	TestMeshIOLagrit.cc \
+	TestCellFilterAvg.cc \
+	TestVertexFilterVecNorm.cc \
 	test_meshio.cc
 
 
-#	TestCellFilterAvg.cc \
 # 	TestDataWriterVTK.cc \
 # 	TestDataWriterVTKMesh.cc \
 # 	TestDataWriterVTKMeshLine2.cc \
@@ -52,8 +53,7 @@
 # 	TestDataWriterVTKBCMeshHex8.cc
 
 # 	TestOutputManager.cc \
-# 	TestOutputSolnSubset.cc \
-# 	TestVertexFilterVecNorm.cc
+# 	TestOutputSolnSubset.cc
 
 
 noinst_HEADERS = \

Modified: short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestCellFilterAvg.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestCellFilterAvg.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestCellFilterAvg.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -14,11 +14,11 @@
 
 #include "TestCellFilterAvg.hh" // Implementation of class methods
 
+#include "pylith/topology/Mesh.hh" // USES Mesh
 #include "pylith/meshio/CellFilterAvg.hh"
 
 #include "pylith/meshio/MeshIOAscii.hh" // USES MeshIOAscii
-#include "pylith/topology/Mesh.hh" // USES Mesh
-#include "pylith/feassemble/Quadrature2D.hh" // USES Quadrature
+#include "pylith/feassemble/Quadrature.hh" // USES Quadrature
 
 // ----------------------------------------------------------------------
 CPPUNIT_TEST_SUITE_REGISTRATION( pylith::meshio::TestCellFilterAvg );
@@ -28,7 +28,7 @@
 void
 pylith::meshio::TestCellFilterAvg::testConstructor(void)
 { // testConstructor
-  CellFilterAvg filter;
+  CellFilterAvg<topology::Mesh> filter;
 } // testConstructor
 
 // ----------------------------------------------------------------------
@@ -36,11 +36,15 @@
 void
 pylith::meshio::TestCellFilterAvg::testFilter(void)
 { // testFilter
+  typedef pylith::topology::Mesh::SieveMesh SieveMesh;
+  typedef pylith::topology::Mesh::RealSection RealSection;
+
   const char* filename = "data/quad4.mesh";
   const int fiberDim = 2;
   const int ncells = 2;
-  const char* fieldName = "field data";
-  const VectorFieldEnum fieldType = OTHER_FIELD;
+  const std::string label = "field data";
+  const topology::FieldBase::VectorFieldEnum fieldType = 
+    topology::FieldBase::MULTI_SCALAR;
   const double fieldValues[] = {
     1.1, 1.2,
     2.1, 2.2,
@@ -67,9 +71,8 @@
   const double quadWts[] = { 1.5, 0.5 };
   const double minJacobian = 1.0;
 
-
-
-  const VectorFieldEnum fieldTypeE = OTHER_FIELD;
+  const topology::FieldBase::VectorFieldEnum fieldTypeE = 
+    topology::FieldBase::SCALAR;
   const int fiberDimE = 1;
   const double fieldValuesE[] = {
     (1.5*1.1 + 0.5*1.2)/2.0,
@@ -77,53 +80,58 @@
   };
 
   MeshIOAscii iohandler;
-  topology::Mesh(PETSC_COMM_WORLD, cellDim);
+  topology::Mesh mesh;
   iohandler.filename(filename);
   iohandler.read(&mesh);
-  CPPUNIT_ASSERT(!mesh.isNull());
 
+  // Set cell field
+  topology::Field<topology::Mesh> field(mesh);
+  field.newSection(topology::FieldBase::CELLS_FIELD, fiberDim);
+  field.allocate();
+  field.vectorFieldType(fieldType);
+  field.label(label.c_str());
+
   const ALE::Obj<SieveMesh>& sieveMesh = mesh.sieveMesh();
   CPPUNIT_ASSERT(!sieveMesh.isNull());
-
-  // Set cell field
   const ALE::Obj<SieveMesh::label_sequence>& cells = 
     sieveMesh->heightStratum(0);
+  CPPUNIT_ASSERT(!cells.isNull());
   const SieveMesh::label_sequence::iterator cellsEnd = cells->end();
 
-  ALE::Obj<SieveMesh::real_section_type> field = 
-    new SieveMesh::real_section_type(mesh.comm(), mesh.debug());
-  field->setChart(SieveMesh::real_section_type::chart_type(0, cells->size()));
-  field->setFiberDimension(cells, fiberDim);
-  sieveMesh->allocate(field);
-
+  const ALE::Obj<RealSection>& section = field.section();
+  CPPUNIT_ASSERT(!section.isNull());
   CPPUNIT_ASSERT_EQUAL(ncells, int(cells->size()));
-
   int ipt = 0;
   for (SieveMesh::label_sequence::iterator c_iter=cells->begin();
        c_iter != cellsEnd;
        ++c_iter, ++ipt) {
     const double* values = &fieldValues[ipt*fiberDim];
-    field->updatePoint(*c_iter, values);
+    section->updatePoint(*c_iter, values);
   } // for
 
-  feassemble::Quadrature2D quadrature;
-  quadrature.initialize(basis, basisDerivRef, quadPtsRef, quadWts,
-			cellDim, numBasis, numQuadPts, spaceDim);
+  feassemble::Quadrature<topology::Mesh> quadrature;
+  quadrature.initialize(basis, numQuadPts, numBasis,
+			basisDerivRef, numQuadPts, numBasis, cellDim,
+			quadPtsRef, numQuadPts, cellDim,
+			quadWts, numQuadPts,
+			spaceDim);
 
-  CellFilterAvg filter;
+  CellFilterAvg<topology::Mesh> filter;
   filter.quadrature(&quadrature);
 
-  VectorFieldEnum fieldTypeF = SCALAR_FIELD;
-  const ALE::Obj<real_section_type>& fieldF =
-    filter.filter(&fieldTypeF, field, meshMesh);
+  const topology::Field<topology::Mesh>& fieldF = filter.filter(field);
+  const ALE::Obj<RealSection>& sectionF = fieldF.section();
+  CPPUNIT_ASSERT(!sectionF.isNull());
 
-  CPPUNIT_ASSERT_EQUAL(fieldTypeE, fieldTypeF);
+  CPPUNIT_ASSERT_EQUAL(fieldTypeE, fieldF.vectorFieldType());
+  CPPUNIT_ASSERT_EQUAL(label, std::string(fieldF.label()));
+
   ipt = 0;
   for (SieveMesh::label_sequence::iterator c_iter=cells->begin();
        c_iter != cellsEnd;
        ++c_iter, ++ipt) {
-    CPPUNIT_ASSERT_EQUAL(fiberDimE, fieldF->getFiberDimension(*c_iter));
-    const double* values = fieldF->restrictPoint(*c_iter);
+    CPPUNIT_ASSERT_EQUAL(fiberDimE, sectionF->getFiberDimension(*c_iter));
+    const double* values = sectionF->restrictPoint(*c_iter);
     CPPUNIT_ASSERT(0 != values);
     const double tolerance = 1.0e-06;
     for (int i=0; i < fiberDimE; ++i)

Modified: short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestVertexFilterVecNorm.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestVertexFilterVecNorm.cc	2009-04-03 22:01:53 UTC (rev 14583)
+++ short/3D/PyLith/branches/pylith-swig/unittests/libtests/meshio/TestVertexFilterVecNorm.cc	2009-04-04 00:02:41 UTC (rev 14584)
@@ -14,9 +14,11 @@
 
 #include "TestVertexFilterVecNorm.hh" // Implementation of class methods
 
+#include "pylith/topology/Mesh.hh" // USES Mesh
 #include "pylith/meshio/VertexFilterVecNorm.hh"
 
 #include "pylith/meshio/MeshIOAscii.hh" // USES MeshIOAscii
+#include "pylith/topology/Field.hh" // USES Field
 
 // ----------------------------------------------------------------------
 CPPUNIT_TEST_SUITE_REGISTRATION( pylith::meshio::TestVertexFilterVecNorm );
@@ -26,7 +28,7 @@
 void
 pylith::meshio::TestVertexFilterVecNorm::testConstructor(void)
 { // testConstructor
-  VertexFilterVecNorm filter;
+  VertexFilterVecNorm<topology::Mesh> filter;
 } // testConstructor
 
 // ----------------------------------------------------------------------
@@ -34,18 +36,23 @@
 void
 pylith::meshio::TestVertexFilterVecNorm::testFilter(void)
 { // testFilter
+  typedef pylith::topology::Mesh::SieveMesh SieveMesh;
+  typedef pylith::topology::Mesh::RealSection RealSection;
+
   const char* filename = "data/tri3.mesh";
   const int fiberDim = 2;
   const int nvertices = 4;
-  const char* fieldName = "field data";
-  const VectorFieldEnum fieldType = VECTOR_FIELD;
+  const std::string label = "field data";
+  const topology::FieldBase::VectorFieldEnum fieldType = 
+    topology::FieldBase::VECTOR;
   const double fieldValues[] = {
     1.1, 1.2,
     2.1, 2.2,
     3.1, 3.2,
     4.1, 4.2
   };
-  const VectorFieldEnum fieldTypeE = SCALAR_FIELD;
+  const topology::FieldBase::VectorFieldEnum fieldTypeE = 
+    topology::FieldBase::SCALAR;
   const int fiberDimE = 1;
   const double fieldValuesE[] = {
     sqrt(pow(1.1, 2) + pow(1.2, 2)),
@@ -55,43 +62,49 @@
   };
 
   MeshIOAscii iohandler;
-  ALE::Obj<Mesh> mesh;
+  topology::Mesh mesh;
   iohandler.filename(filename);
   iohandler.read(&mesh);
-  CPPUNIT_ASSERT(!mesh.isNull());
 
   // Set vertex field
-  const ALE::Obj<Mesh::label_sequence>& vertices = mesh->depthStratum(0);
-  const Mesh::label_sequence::iterator verticesEnd = vertices->end();
+  topology::Field<topology::Mesh> field(mesh);
+  field.newSection(topology::FieldBase::VERTICES_FIELD, fiberDim);
+  field.allocate();
+  field.vectorFieldType(fieldType);
+  field.label(label.c_str());
 
-  ALE::Obj<real_section_type> field = 
-    new real_section_type(mesh->comm(), mesh->debug());
-  field->setChart(mesh->getSieve()->getChart());
-  field->setFiberDimension(vertices, fiberDim);
-  mesh->allocate(field);
-
+  const ALE::Obj<SieveMesh>& sieveMesh = mesh.sieveMesh();
+  CPPUNIT_ASSERT(!sieveMesh.isNull());
+  const ALE::Obj<SieveMesh::label_sequence>& vertices = 
+    sieveMesh->depthStratum(0);
+  CPPUNIT_ASSERT(!vertices.isNull());
+  const SieveMesh::label_sequence::iterator verticesEnd = vertices->end();
+  
+  const ALE::Obj<RealSection>& section = field.section();
+  CPPUNIT_ASSERT(!section.isNull());
   CPPUNIT_ASSERT_EQUAL(nvertices, int(vertices->size()));
-
   int ipt = 0;
-  for (Mesh::label_sequence::iterator v_iter=vertices->begin();
+  for (SieveMesh::label_sequence::iterator v_iter=vertices->begin();
        v_iter != verticesEnd;
        ++v_iter, ++ipt) {
     const double* values = &fieldValues[ipt*fiberDim];
-    field->updatePoint(*v_iter, values);
+    section->updatePoint(*v_iter, values);
   } // for
 
-  VertexFilterVecNorm filter;
-  VectorFieldEnum fieldTypeF = OTHER_FIELD;
-  const ALE::Obj<real_section_type>& fieldF =
-    filter.filter(&fieldTypeF, field, mesh);
+  VertexFilterVecNorm<topology::Mesh> filter;
+  const topology::Field<topology::Mesh>& fieldF = filter.filter(field);
+  const ALE::Obj<RealSection>& sectionF = fieldF.section();
+  CPPUNIT_ASSERT(!sectionF.isNull());
 
-  CPPUNIT_ASSERT_EQUAL(fieldTypeE, fieldTypeF);
+  CPPUNIT_ASSERT_EQUAL(fieldTypeE, fieldF.vectorFieldType());
+  CPPUNIT_ASSERT_EQUAL(label, std::string(fieldF.label()));
+
   ipt = 0;
-  for (Mesh::label_sequence::iterator v_iter=vertices->begin();
+  for (SieveMesh::label_sequence::iterator v_iter=vertices->begin();
        v_iter != verticesEnd;
        ++v_iter, ++ipt) {
-    CPPUNIT_ASSERT_EQUAL(fiberDimE, fieldF->getFiberDimension(*v_iter));
-    const double* values = fieldF->restrictPoint(*v_iter);
+    CPPUNIT_ASSERT_EQUAL(fiberDimE, sectionF->getFiberDimension(*v_iter));
+    const double* values = sectionF->restrictPoint(*v_iter);
     CPPUNIT_ASSERT(0 != values);
     const double tolerance = 1.0e-06;
     for (int i=0; i < fiberDimE; ++i)



More information about the CIG-COMMITS mailing list