[cig-commits] r9327 - in short/3D/PyLith/trunk: examples/twocells/twoquad4 libsrc libsrc/meshio modulesrc/meshio pylith pylith/feassemble pylith/meshio

brad at geodynamics.org brad at geodynamics.org
Wed Feb 13 17:14:20 PST 2008


Author: brad
Date: 2008-02-13 17:14:20 -0800 (Wed, 13 Feb 2008)
New Revision: 9327

Added:
   short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc
   short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.hh
   short/3D/PyLith/trunk/pylith/meshio/CellFilterAvg.py
Modified:
   short/3D/PyLith/trunk/examples/twocells/twoquad4/axialdisp.cfg
   short/3D/PyLith/trunk/libsrc/Makefile.am
   short/3D/PyLith/trunk/libsrc/meshio/CellFilter.cc
   short/3D/PyLith/trunk/libsrc/meshio/CellFilter.hh
   short/3D/PyLith/trunk/libsrc/meshio/Makefile.am
   short/3D/PyLith/trunk/libsrc/meshio/OutputManager.cc
   short/3D/PyLith/trunk/libsrc/meshio/OutputManager.hh
   short/3D/PyLith/trunk/modulesrc/meshio/meshio.pyxe.src
   short/3D/PyLith/trunk/pylith/Makefile.am
   short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py
   short/3D/PyLith/trunk/pylith/meshio/OutputManager.py
   short/3D/PyLith/trunk/pylith/meshio/__init__.py
Log:
Added ability to average cell fields over quadrature points in output.

Modified: short/3D/PyLith/trunk/examples/twocells/twoquad4/axialdisp.cfg
===================================================================
--- short/3D/PyLith/trunk/examples/twocells/twoquad4/axialdisp.cfg	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/examples/twocells/twoquad4/axialdisp.cfg	2008-02-14 01:14:20 UTC (rev 9327)
@@ -110,5 +110,6 @@
 filename = axialdisp.vtk
 
 # Give basename for VTK output of state variables.
-[pylithapp.timedependent.materials.material.output.writer]
-filename = axialdisp-statevars.vtk
+[pylithapp.timedependent.materials.material.output]
+cell_filter = pylith.meshio.CellFilterAvg
+writer.filename = axialdisp-statevars.vtk

Modified: short/3D/PyLith/trunk/libsrc/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/libsrc/Makefile.am	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/libsrc/Makefile.am	2008-02-14 01:14:20 UTC (rev 9327)
@@ -72,6 +72,7 @@
 	materials/ViscoelasticMaxwell.cc \
 	meshio/BinaryIO.cc \
 	meshio/CellFilter.cc \
+	meshio/CellFilterAvg.cc \
 	meshio/DataWriter.cc \
 	meshio/DataWriterVTK.cc \
 	meshio/GMVFile.cc \

Modified: short/3D/PyLith/trunk/libsrc/meshio/CellFilter.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/CellFilter.cc	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/libsrc/meshio/CellFilter.cc	2008-02-14 01:14:20 UTC (rev 9327)
@@ -40,17 +40,6 @@
 } // copy constructor
 
 // ----------------------------------------------------------------------
-// operator=.
-const pylith::meshio::CellFilter&
-pylith::meshio::CellFilter::operator=(const CellFilter& f)
-{ // operator=
-  if (&f != this) {
-    delete _quadrature; 
-    _quadrature = (0 != f._quadrature) ? f._quadrature->clone() : 0;
-  } // if
-} // operator=
-
-// ----------------------------------------------------------------------
 // Set quadrature associated with cells.
 void
 pylith::meshio::CellFilter::quadrature(const feassemble::Quadrature* q)

Modified: short/3D/PyLith/trunk/libsrc/meshio/CellFilter.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/CellFilter.hh	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/libsrc/meshio/CellFilter.hh	2008-02-14 01:14:20 UTC (rev 9327)
@@ -61,11 +61,17 @@
    *
    * @param fieldIn Field to filter.
    * @param mesh PETSc mesh.
+   * @param label Label identifying cells.
+   * @param Value of label of cells to filter.
+   *
+   * @returns Averaged field.
    */
   virtual
   const ALE::Obj<real_section_type>&
   filter(const ALE::Obj<real_section_type>& fieldIn,
-	 const ALE::Obj<ALE::Mesh>& mesh) = 0;
+	 const ALE::Obj<ALE::Mesh>& mesh,
+	 const char* label,
+	 const int labelId) = 0;
 
 // PROTECTED METHODS ////////////////////////////////////////////////////
 protected :
@@ -77,13 +83,12 @@
    */
   CellFilter(const CellFilter& f);
 
-  /** operator=.
-  *
-  * @param f Filter to copy.
-  * @returns Copy of filter.
-  */
-  const CellFilter& operator=(const CellFilter& f);
+// NOT IMPLEMENTED //////////////////////////////////////////////////////
+private :
 
+  /// Not implemented.
+  const CellFilter& operator=(const CellFilter&);
+
 // PROTECTED MEMBERS ////////////////////////////////////////////////////
 protected :
 

Added: short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc	                        (rev 0)
+++ short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc	2008-02-14 01:14:20 UTC (rev 9327)
@@ -0,0 +1,102 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+//                           Brad T. Aagaard
+//                        U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+#include <portinfo>
+
+#include "CellFilterAvg.hh" // implementation of class methods
+
+#include "pylith/feassemble/Quadrature.hh" // USES Quadrature
+
+// ----------------------------------------------------------------------
+// Constructor
+pylith::meshio::CellFilterAvg::CellFilterAvg(void)
+{ // constructor
+} // constructor
+
+// ----------------------------------------------------------------------
+// Destructor
+pylith::meshio::CellFilterAvg::~CellFilterAvg(void)
+{ // destructor
+} // destructor  
+
+// ----------------------------------------------------------------------
+// Copy constructor.
+pylith::meshio::CellFilterAvg::CellFilterAvg(const CellFilterAvg& f) :
+  CellFilter(f)
+{ // copy constructor
+} // copy constructor
+
+// ----------------------------------------------------------------------
+// Create copy of filter.
+pylith::meshio::CellFilter*
+pylith::meshio::CellFilterAvg::clone(void) const
+{ // clone
+  return new CellFilterAvg(*this);
+} // clone
+
+// ----------------------------------------------------------------------
+// Filter field.
+const ALE::Obj<pylith::real_section_type>&
+pylith::meshio::CellFilterAvg::filter(
+				  const ALE::Obj<real_section_type>& fieldIn,
+				  const ALE::Obj<ALE::Mesh>& mesh,
+				  const char* label,
+				  const int labelId)
+{ // filter
+  assert(0 != _quadrature);
+
+  const int numQuadPts = _quadrature->numQuadPts();
+  const double_array& wts = _quadrature->quadWts();
+  
+  const ALE::Obj<Mesh::label_sequence>& cells = (0 == label) ?
+    mesh->heightStratum(0) :
+    mesh->getLabelStratum(label, labelId);
+  assert(!cells.isNull());
+  const Mesh::label_sequence::iterator cellsEnd = cells->end();
+
+  const int totalFiberDim = fieldIn->getFiberDimension(*cells->begin());
+  const int fiberDim = totalFiberDim / numQuadPts;
+  assert(fiberDim * numQuadPts == totalFiberDim);
+
+  // Allocation field if necessary
+  if (_fieldAvg.isNull() ||
+      fiberDim != _fieldAvg->getFiberDimension(*cells->begin())) {
+    _fieldAvg = new real_section_type(mesh->comm(), mesh->debug());
+    _fieldAvg->setFiberDimension(cells, fiberDim);
+    mesh->allocate(_fieldAvg);
+  } // if
+
+  double_array fieldAvgCell(fiberDim);
+  double scalar = 0.0;
+  for (int iQuad=0; iQuad < numQuadPts; ++iQuad)
+    scalar += wts[iQuad];
+
+  // Loop over cells
+  for (Mesh::label_sequence::iterator c_iter=cells->begin();
+       c_iter != cellsEnd;
+       ++c_iter) {
+    const real_section_type::value_type* values = 
+      fieldIn->restrictPoint(*c_iter);
+    
+    fieldAvgCell = 0.0;
+    for (int iQuad=0; iQuad < numQuadPts; ++iQuad)
+      for (int i=0; i < fiberDim; ++i)
+	fieldAvgCell[i] += wts[iQuad] / scalar * values[iQuad*fiberDim+i];
+
+    _fieldAvg->updatePoint(*c_iter, &fieldAvgCell[0]);
+  } // for
+
+  return _fieldAvg;
+} // filter
+
+
+// End of file

Added: short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.hh	                        (rev 0)
+++ short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.hh	2008-02-14 01:14:20 UTC (rev 9327)
@@ -0,0 +1,90 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+//                           Brad T. Aagaard
+//                        U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+/**
+ * @file pylith/meshio/CellFilterAvg.hh
+ *
+ * @brief C++ object for averaging cell fields over quadrature points
+ * when outputing finite-element data.
+ */
+
+#if !defined(pylith_meshio_cellfilteravg_hh)
+#define pylith_meshio_cellfilteravg_hh
+
+#include "CellFilter.hh" // ISA CellFilter
+
+namespace pylith {
+  namespace meshio {
+    class CellFilterAvg;
+  } // meshio
+} // pylith
+
+class pylith::meshio::CellFilterAvg : public CellFilter
+{ // CellFilterAvg
+
+// PUBLIC METHODS ///////////////////////////////////////////////////////
+public :
+
+  /// Constructor
+  CellFilterAvg(void);
+
+  /// Destructor
+  ~CellFilterAvg(void);
+
+  /** Create copy of filter.
+   *
+   * @returns Copy of filter.
+   */
+  CellFilter* clone(void) const;
+
+  /** Filter field.
+   *
+   * @param fieldIn Field to filter.
+   * @param mesh PETSc mesh.
+   * @param label Label identifying cells.
+   * @param Value of label of cells to filter.
+   *
+   * @returns Averaged field.
+   */
+  const ALE::Obj<real_section_type>&
+  filter(const ALE::Obj<real_section_type>& fieldIn,
+	 const ALE::Obj<ALE::Mesh>& mesh,
+	 const char* label,
+	 const int labelId);
+
+// PROTECTED METHODS ////////////////////////////////////////////////////
+protected :
+
+  /** Copy constructor.
+   *
+   * @param f Filter to copy.
+   * @returns Pointer to this.
+   */
+  CellFilterAvg(const CellFilterAvg& f);
+
+// NOT IMPLEMENTED //////////////////////////////////////////////////////
+private :
+
+  /// Not implemented.
+  const CellFilter& operator=(const CellFilter&);
+
+// PRIVATE MEMBERS //////////////////////////////////////////////////////
+private :
+
+  ALE::Obj<real_section_type> _fieldAvg; ///< Averaged cell field
+
+}; // CellFilterAvg
+
+#endif // pylith_meshio_cellfilteravg_hh
+
+
+// End of file 

Modified: short/3D/PyLith/trunk/libsrc/meshio/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/Makefile.am	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/libsrc/meshio/Makefile.am	2008-02-14 01:14:20 UTC (rev 9327)
@@ -15,6 +15,7 @@
 
 subpkginclude_HEADERS = \
 	CellFilter.hh \
+	CellFilterAvg.hh \
 	DataWriter.hh \
 	DataWriterVTK.hh \
 	DataWriterVTK.icc \

Modified: short/3D/PyLith/trunk/libsrc/meshio/OutputManager.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/OutputManager.cc	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/libsrc/meshio/OutputManager.cc	2008-02-14 01:14:20 UTC (rev 9327)
@@ -145,12 +145,14 @@
 				const char* name,
 				const ALE::Obj<real_section_type>& field,
 				const VectorFieldEnum fieldType,
-				const ALE::Obj<ALE::Mesh>& mesh)
+				const ALE::Obj<ALE::Mesh>& mesh,
+				const char* label,
+				const int labelId)
 { // appendCellField
   assert(0 != name);
 
   const ALE::Obj<real_section_type>& fieldFiltered = 
-    (0 == _cellFilter) ? field : _cellFilter->filter(field, mesh);
+    (0 == _cellFilter) ? field : _cellFilter->filter(field, mesh, label, labelId);
 
   _writer->writeCellField(t, name, fieldFiltered, fieldType, mesh);
 } // appendCellField

Modified: short/3D/PyLith/trunk/libsrc/meshio/OutputManager.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/OutputManager.hh	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/libsrc/meshio/OutputManager.hh	2008-02-14 01:14:20 UTC (rev 9327)
@@ -132,12 +132,17 @@
    * @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.
    */
   void appendCellField(const double t,
 		       const char* name,
 		       const ALE::Obj<real_section_type>& field,
 		       const VectorFieldEnum fieldType,
-		       const ALE::Obj<ALE::Mesh>& mesh);
+		       const ALE::Obj<ALE::Mesh>& mesh,
+		       const char* label =0,
+		       const int labelId =0);
 
 // PRIVATE MEMBERS //////////////////////////////////////////////////////
 private :

Modified: short/3D/PyLith/trunk/modulesrc/meshio/meshio.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/meshio/meshio.pyxe.src	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/modulesrc/meshio/meshio.pyxe.src	2008-02-14 01:14:20 UTC (rev 9327)
@@ -12,6 +12,7 @@
 
 #header{
 #include "pylith/meshio/CellFilter.hh"
+#include "pylith/meshio/CellFilterAvg.hh"
 #include "pylith/meshio/DataWriter.hh"
 #include "pylith/meshio/DataWriterVTK.hh"
 #include "pylith/meshio/MeshIO.hh"
@@ -767,12 +768,12 @@
     return
 
 
-  def appendCellField(self, t, name, field, fieldType, mesh):
+  def appendCellField(self, t, name, field, fieldType, mesh, label, labelId):
     """
     Append field over cells to file.
     """
     # create shim for method 'appendCellField'
-    #embed{ void OutputManager_appendCellField(void* objVptr, double t, char* name, void* fieldVptr, int fieldTypeInt, void* meshVptr)
+    #embed{ void OutputManager_appendCellField(void* objVptr, double t, char* name, void* fieldVptr, int fieldTypeInt, void* meshVptr, char* label, int labelId)
     try {
       assert(0 != objVptr);
       assert(0 != fieldVptr);
@@ -786,7 +787,7 @@
       ((pylith::meshio::OutputManager*) objVptr)->appendCellField(t, name,
                                                                   *field,
                                                                   fieldType,
-                                                                  *mesh);
+                                                                  *mesh, label, labelId);
     } catch (const std::exception& err) {
       PyErr_SetString(PyExc_RuntimeError,
                       const_cast<char*>(err.what()));
@@ -803,8 +804,12 @@
             "Argument must be extension module type 'Mesh'."
     cdef void* fieldVptr
     fieldVptr = PyCObject_AsVoidPtr(field)
-    OutputManager_appendCellField(self.thisptr, t, name, fieldVptr, fieldType,
-                                  ptrFromHandle(mesh))
+    if None == label or None == labelId:
+      OutputManager_appendCellField(self.thisptr, t, name, fieldVptr, fieldType,
+                                    ptrFromHandle(mesh), NULL, 0)
+    else:
+      OutputManager_appendCellField(self.thisptr, t, name, fieldVptr, fieldType,
+                                    ptrFromHandle(mesh), label, labelId)
     return
 
 
@@ -1031,4 +1036,98 @@
       DataWriterVTK_timeFormat_set(self.thisptr, format)
 
 
+# ----------------------------------------------------------------------
+cdef void CellFilter_destructor(void* obj):
+  """
+  Destroy CellFilter object.
+  """
+  # create shim for destructor
+  #embed{ void CellFilter_destructor_cpp(void* pObj)
+  pylith::meshio::CellFilter* f = (pylith::meshio::CellFilter*) pObj;
+  delete f;
+  #}embed
+  CellFilter_destructor_cpp(obj)
+  return
+
+cdef class CellFilter:
+
+  cdef void* thisptr # Pointer to C++ object
+  cdef readonly object handle # PyCObject holding pointer to C++ object
+  cdef readonly object name # Identifier for object base type
+
+  def __init__(self):
+    """
+    Constructor.
+    """
+    self.handle = None
+    self.thisptr = NULL
+    self.name = "pylith_meshio_CellFilter"
+    return
+
+
+  property quadrature:
+    def __set__(self, quadrature):
+      """
+      Set quadrature.
+      """
+      # create shim for method 'quadrature'
+      #embed{ void CellFilter_quadrature_set(void* objVptr, void* qVptr)
+      try {
+        assert(0 != objVptr);
+        assert(0 != qVptr);
+        pylith::feassemble::Quadrature* quadrature =
+          (pylith::feassemble::Quadrature*) qVptr;
+        ((pylith::meshio::CellFilter*) objVptr)->quadrature(quadrature);
+      } catch (const std::exception& err) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        const_cast<char*>(err.what()));
+      } catch (const ALE::Exception& err) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        const_cast<char*>(err.msg().c_str()));
+      } catch (...) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        "Caught unknown C++ exception.");
+      } // try/catch
+      #}embed
+      if not quadrature.name == "pylith_feassemble_Quadrature":
+        raise TypeError, \
+              "Argument must be extension module type 'Quadrature'."
+      CellFilter_quadrature_set(self.thisptr, ptrFromHandle(quadrature))
+
+
+  def _createHandle(self):
+    """Wrap pointer to C++ object in PyCObject."""
+    return PyCObject_FromVoidPtr(self.thisptr, CellFilter_destructor)
+
+
+# ----------------------------------------------------------------------
+cdef class CellFilterAvg(CellFilter):
+
+  def __init__(self):
+    """Constructor."""
+    # create shim for constructor
+    #embed{ void* CellFilterAvg_constructor()
+    void* result = 0;
+    try {
+      result = (void*)(new pylith::meshio::CellFilterAvg);
+      assert(0 != result);
+    } catch (const std::exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.what()));
+    } catch (const ALE::Exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.msg().c_str()));
+    } catch (...) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      "Caught unknown C++ exception.");
+    } // try/catch
+    return result;
+    #}embed
+
+    CellFilter.__init__(self)
+    self.thisptr = CellFilterAvg_constructor()
+    self.handle = self._createHandle()
+    return
+
+
 # End of file

Modified: short/3D/PyLith/trunk/pylith/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/pylith/Makefile.am	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/pylith/Makefile.am	2008-02-14 01:14:20 UTC (rev 9327)
@@ -78,6 +78,7 @@
 	materials/MaxwellIsotropic3D.py \
 	meshio/__init__.py \
 	meshio/CellFilter.py \
+	meshio/CellFilterAvg.py \
 	meshio/DataWriter.py \
 	meshio/DataWriterVTK.py \
 	meshio/MeshIO.py \

Modified: short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py	2008-02-14 01:14:20 UTC (rev 9327)
@@ -84,7 +84,7 @@
                    self.material.label)
 
     self.material.initialize(self.mesh, totalTime, numTimeSteps)
-    self.output.initialize(self.quadrature.cppHandle)
+    self.output.initialize(self.quadrature)
     self.output.writeInfo()
     self.output.open(totalTime, numTimeSteps)
 

Added: short/3D/PyLith/trunk/pylith/meshio/CellFilterAvg.py
===================================================================
--- short/3D/PyLith/trunk/pylith/meshio/CellFilterAvg.py	                        (rev 0)
+++ short/3D/PyLith/trunk/pylith/meshio/CellFilterAvg.py	2008-02-14 01:14:20 UTC (rev 9327)
@@ -0,0 +1,62 @@
+#!/usr/bin/env python
+#
+# ----------------------------------------------------------------------
+#
+#                           Brad T. Aagaard
+#                        U.S. Geological Survey
+#
+# <LicenseText>
+#
+# ----------------------------------------------------------------------
+#
+
+## @file pyre/meshio/CellFilterAvg.py
+##
+## @brief Python class for averageing cell fields over each cell's
+## quadrature points when writing finite-element data.
+##
+## Factory: output_cell_filter
+
+from CellFilter import CellFilter
+
+# CellFilterAvg class
+class CellFilterAvg(CellFilter):
+  """
+  Python class for average cell fields over each cell's quadrature
+  points when writing finite-element data.
+
+  Factory: output_cell_filter
+  """
+
+  # PUBLIC METHODS /////////////////////////////////////////////////////
+
+  def __init__(self, name="cellfilteravg"):
+    """
+    Constructor.
+    """
+    CellFilter.__init__(self, name)
+    return
+
+
+  # PRIVATE METHODS ////////////////////////////////////////////////////
+
+  def _createCppHandle(self):
+    """
+    Create handle to C++ object.
+    """
+    if None == self.cppHandle:
+      import meshio as bindings
+      self.cppHandle = bindings.CellFilterAvg()
+    return
+  
+
+# FACTORIES ////////////////////////////////////////////////////////////
+
+def output_cell_filter():
+  """
+  Factory associated with CellFilter.
+  """
+  return CellFilterAvg()
+
+
+# End of file 

Modified: short/3D/PyLith/trunk/pylith/meshio/OutputManager.py
===================================================================
--- short/3D/PyLith/trunk/pylith/meshio/OutputManager.py	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/pylith/meshio/OutputManager.py	2008-02-14 01:14:20 UTC (rev 9327)
@@ -210,7 +210,7 @@
       for name in self.cellInfoFields:
         (field, fieldType) = self.dataProvider.getCellField(name)
         self.cppHandle.appendCellField(t.value, name, field, fieldType, 
-                                       mesh.cppHandle)
+                                       mesh.cppHandle, label, labelId)
 
       self.cppHandle.closeTimeStep()
       self.close()
@@ -241,7 +241,7 @@
       for name in self.cellDataFields:
         (field, fieldType) = self.dataProvider.getCellField(name, fields)
         self.cppHandle.appendCellField(t.value, name, field, fieldType, 
-                                       mesh.cppHandle)
+                                       mesh.cppHandle, label, labelId)
 
       self.cppHandle.closeTimeStep()
 

Modified: short/3D/PyLith/trunk/pylith/meshio/__init__.py
===================================================================
--- short/3D/PyLith/trunk/pylith/meshio/__init__.py	2008-02-14 00:04:44 UTC (rev 9326)
+++ short/3D/PyLith/trunk/pylith/meshio/__init__.py	2008-02-14 01:14:20 UTC (rev 9327)
@@ -15,6 +15,7 @@
 ## @brief Python PyLith meshio module initialization
 
 __all__ = ['CellFilter',
+           'CellFilterAvg',
            'DataWriter',
            'DataWriterVTK',
            'MeshIO',



More information about the cig-commits mailing list