[cig-commits] r5179 - in short/3D/PyLith/trunk: libsrc/meshio pylith/feassemble

brad at geodynamics.org brad at geodynamics.org
Sun Nov 5 19:42:15 PST 2006


Author: brad
Date: 2006-11-05 19:42:14 -0800 (Sun, 05 Nov 2006)
New Revision: 5179

Modified:
   short/3D/PyLith/trunk/libsrc/meshio/Makefile.am
   short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc
   short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh
   short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc
   short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.hh
   short/3D/PyLith/trunk/pylith/feassemble/Field.py
   short/3D/PyLith/trunk/pylith/feassemble/Integrator.py
Log:
Cleaned up MeshIO interface. Moved all Sieve access to MeshIO. Still need to create unit tests.

Modified: short/3D/PyLith/trunk/libsrc/meshio/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/Makefile.am	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/libsrc/meshio/Makefile.am	2006-11-06 03:42:14 UTC (rev 5179)
@@ -16,12 +16,12 @@
 lib_LTLIBRARIES = libpylithmeshio.la
 
 libpylithmeshio_la_SOURCES = \
-	MeshIO.cc \
+	MeshIO.cc 
 	MeshIOAscii.cc
 
 subpkginclude_HEADERS = \
 	MeshIO.hh \
-	MeshIO.icc \
+	MeshIO.icc
 	MeshIOAscii.hh \
 	MeshIOAscii.icc
 

Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIO.cc	2006-11-06 03:42:14 UTC (rev 5179)
@@ -17,7 +17,8 @@
 // ----------------------------------------------------------------------
 // Constructor
 pylith::meshio::MeshIO::MeshIO(void) :
-  _useIndexZero(true)
+  _useIndexZero(true),
+  _mesh(0)
 { // constructor
 } // constructor
 
@@ -27,4 +28,159 @@
 { // destructor
 } // destructor
   
+// ----------------------------------------------------------------------
+// Read mesh from file.
+void 
+pylith::meshio::MeshIO::read(ALE::Obj<Mesh>* mesh)
+{ // read
+  _mesh = mesh;
+  _read();
+  _mesh = 0;
+} // read
+
+// ----------------------------------------------------------------------
+// Write mesh to file.
+void 
+pylith::meshio::MeshIO::write(ALE::Obj<Mesh>* mesh)
+{ // write
+  _mesh = mesh;
+  _write();
+  _mesh = 0;
+} // write
+
+// ----------------------------------------------------------------------
+// Set vertices in mesh.
+void
+pylith::meshio::MeshIO::_buildMesh(const double* coordinates,
+				   const int numVertices,
+				   const int spaceDim,
+				   const int* cells,
+				   const int numCells,
+				   const int numCorners,
+				   const int meshDim)
+{ // _buildMesh
+  assert(0 != _mesh);
+  *_mesh = Mesh(PETSC_COMM_WORLD, meshDim);
+  ALE::Obj<Mesh>& mesh = *_mesh;
+  
+  ALE::Obj<sieve_type> sieve = new sieve_type(mesh->comm());
+  ALE::Obj<topology_type> topology = new topology_type(mesh->comm());
+
+  const bool interpolate = false;
+  ALE::New::SieveBuilder<sieve_type>::buildTopology(sieve, meshDim, 
+						    numCells, 
+						    const_cast<int*>(cells), 
+						    numVertices, 
+						    interpolate, numCorners);
+  sieve->stratify();
+  topology->setPatch(0, sieve);
+  topology->stratify();
+  mesh->setTopology(topology);
+  ALE::New::SieveBuilder<sieve_type>::buildCoordinates(
+					  mesh->getRealSection("coordinates"), 
+					  spaceDim, coordinates);
+} // _buildMesh
+
+// ----------------------------------------------------------------------
+// Get coordinates of vertices in mesh.
+void
+pylith::meshio::MeshIO::_getVertices(double** pCoordinates,
+				     int* pNumVertices,
+				     int* pSpaceDim) const
+{ // _getVertices
+  assert(0 != _mesh);
+  ALE::Obj<Mesh>& mesh = *_mesh;
+
+  const Mesh::real_section_type::patch_type patch = 0;
+  const ALE::Obj<topology_type>& topology = mesh->getTopology();
+
+  const ALE::Obj<Mesh::topology_type::label_sequence>& vertices = 
+    topology->depthStratum(patch, 0);
+  const ALE::Obj<Mesh::real_section_type>& coordsField =
+    mesh->getRealSection("coordinates");
+
+  const int numVertices = vertices->size();
+  const int spaceDim = 
+    coordsField->getFiberDimension(patch, *vertices->begin());
+
+  double* coordinates = 0;
+  const int size = numVertices * spaceDim;
+  if (0 != *pCoordinates && size > 0) {
+    coordinates = new double[size];
+
+    int i = 0;
+    for(Mesh::topology_type::label_sequence::iterator v_iter = 
+	  vertices->begin();
+	v_iter != vertices->end();
+	++v_iter) {
+      const Mesh::real_section_type::value_type *vertexCoords = 
+	coordsField->restrict(patch, *v_iter);
+      for (int iDim=0; iDim < spaceDim; ++iDim)
+	coordinates[i++] = coordinates[iDim];
+    } // for
+  } // if
+
+  if (0 != pCoordinates)
+    *pCoordinates = coordinates;
+  if (0 != pNumVertices)
+    *pNumVertices = numVertices;
+  if (0 != *pSpaceDim)
+    *pSpaceDim = spaceDim;
+} // _getVertices
+
+// ----------------------------------------------------------------------
+// Get cells in mesh.
+void
+pylith::meshio::MeshIO::_getCells(int** pCells,
+				  int* pNumCells,
+				  int* pNumCorners,
+				  int* pMeshDim) const
+{ // _getCells
+  assert(0 != _mesh);
+  ALE::Obj<Mesh>& mesh = *_mesh;
+
+  const topology_type::patch_type patch = 0;
+  const ALE::Obj<topology_type>& topology = mesh->getTopology();
+
+  const ALE::Obj<sieve_type>& sieve = topology->getPatch(patch);
+  const ALE::Obj<Mesh::topology_type::label_sequence>& cells = 
+    topology->heightStratum(patch, 0);
+
+  const int meshDim = mesh->getDimension();
+  const int numCells = cells->size();
+  const int numCorners = sieve->nCone(*cells->begin(), 
+				      topology->depth())->size();
+
+  const ALE::Obj<Mesh::numbering_type>& vNumbering = 
+    mesh->getFactory()->getLocalNumbering(topology, patch, 0);
+
+  int* cellsArray = 0;
+  const int size = numCells * numCorners;
+  if (0 != pCells && size > 0) {
+    cellsArray = new int[size];
+    
+    const int offset = (useIndexZero()) ? 0 : 1;
+    int i = 0;
+    for(Mesh::topology_type::label_sequence::iterator e_iter = cells->begin();
+	e_iter != cells->end();
+	++e_iter) {
+      const ALE::Obj<sieve_type::traits::coneSequence>& cone = 
+	sieve->cone(*e_iter);
+      for(sieve_type::traits::coneSequence::iterator c_iter = cone->begin();
+	  c_iter != cone->end();
+	  ++c_iter)
+	cellsArray[i++] = vNumbering->getIndex(*c_iter) + offset;
+    } // for
+  } // if  
+
+  if (0 != pCells)
+    *pCells = cellsArray;
+  if (0 != pNumCells)
+    *pNumCells = numCells;
+  if (0 != pNumCorners)
+    *pNumCorners = numCorners;
+  if (0 != pMeshDim)
+    *pMeshDim = meshDim;
+} // _getCells
+
 // End of file 

Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIO.hh	2006-11-06 03:42:14 UTC (rev 5179)
@@ -23,6 +23,11 @@
 
 class pylith::meshio::MeshIO
 { // MeshIO
+// PUBLIC TYPEDEFS //////////////////////////////////////////////////////
+public :
+  typedef ALE::Mesh Mesh;
+  typedef ALE::Mesh::sieve_type sieve_type;
+  typedef ALE::Mesh::topology_type topology_type;
   
 // PUBLIC MEMBERS ///////////////////////////////////////////////////////
 public :
@@ -35,19 +40,27 @@
 
   /** Read mesh from file.
    *
-   * @param pMesh Pointer to PETSc mesh object
+   * @param mesh Pointer to PETSc mesh object
    */
-  virtual void read(ALE::Obj<ALE::Mesh>& pMesh, const bool interpolate) = 0;
+  void read(ALE::Obj<Mesh>* mesh);
 
   /** Write mesh to file.
    *
-   * @param mesh PETSc mesh object
+   * @param mesh Pointer to PETSc mesh object
    */
-  virtual void write(const ALE::Obj<ALE::Mesh>& mesh) const = 0;
+  void write(ALE::Obj<Mesh>* mesh);
 
 // PROTECTED MEMBERS ////////////////////////////////////////////////////
 protected :
 
+  /// Write mesh
+  virtual
+  void _write(void) const = 0;
+
+  /// Read mesh
+  virtual
+  void _read(void) = 0;
+
   /** Get flag indicating whether indices start at 0 (True) or 1 (False).
    *
    * @returns True if indices start at 0, false if 1.
@@ -60,11 +73,57 @@
    */
   void useIndexZero(const bool flag);
 
+  /** Build mesh topology and set vertex coordinates.
+   *
+   * @param coordinates Array of coordinates of vertices
+   * @param numVertices Number of vertices
+   * @param spaceDim Dimension of vector space for vertex coordinates
+   * @param cells Array of indices of vertices in cells (first index is 0 or 1)
+   * @param numCells Number of cells
+   * @param numCorners Number of vertices per cell
+   * @param meshDim Dimension of cells in mesh
+   */
+  void _buildMesh(const double* coordinates,
+		  const int numVertices,
+		  const int spaceDim,
+		  const int* cells,
+		  const int numCells,
+		  const int numCorners,
+		  const int meshDim);
+
+  /** Get information about vertices in mesh.
+   *
+   * Method caller is responsible for memory management.
+   *
+   * @param pCoordinates Pointer to array of vertex coordinates
+   * @param pNumVertices Pointer to number of vertices
+   * @param pSpaceDim Poiner to dimension of vector space for coordinates
+   */
+  void _getVertices(double** pCoordinates,
+		    int* pNumVertices,
+		    int* pSpaceDim) const;
+
+  /** Get information about cells in mesh.
+   *
+   * Method caller is responsible for memory management.
+   *
+   * @param pCells Pointer to array of indicates of vertices in each cell
+   * @param pNumCells Pointer to number of cells in mesh
+   * @param pNumCorners Pointer to number of vertices in each cell
+   * @param pMeshDim Pointer to number of dimensions associated with cell
+   */
+  void _getCells(int** pCells,
+		 int* pNumCells,
+		 int* pNumCorners,
+		 int* pMeshDim) const;
+
 // PRIVATE MEMBERS //////////////////////////////////////////////////////
 private :
 
   bool _useIndexZero; ///< Flag indicating if indicates start at 0 (T) or 1 (F)
 
+  ALE::Obj<Mesh>* _mesh; ///< Pointer to PETSc mesh object
+
 }; // MeshIO
 
 #include "MeshIO.icc" // inline methods

Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.cc	2006-11-06 03:42:14 UTC (rev 5179)
@@ -36,11 +36,10 @@
 // ----------------------------------------------------------------------
 // Unpickle mesh
 void
-pylith::meshio::MeshIOAscii::read(ALE::Obj<Mesh>& mesh, 
-				  const bool interpolate)
-{ // read
+pylith::meshio::MeshIOAscii::_read(void)
+{ // _read
   int meshDim = 0;
-  int numDims = 0;
+  int spaceDim = 0;
   int numVertices = 0;
   int numCells = 0;
   int numCorners = 0;
@@ -67,7 +66,7 @@
   bool readDim = false;
   bool readCells = false;
   bool readVertices = false;
-  bool builtTopology = false;
+  bool builtMesh = false;
 
   filein.ignore(maxIgnore, '{');
   filein >> token;
@@ -105,22 +104,11 @@
     } // else
 
     if (readDim && readCells && readVertices && !builtTopology) {
-      // Can now build topology
-      mesh = Mesh(PETSC_COMM_WORLD, meshDim);
-
-      // allow mesh to have different dimension than coordinates
-      ALE::Obj<sieve_type> sieve = new sieve_type(mesh->comm());
-      ALE::Obj<topology_type> topology = new topology_type(mesh->comm());
-
-      ALE::New::SieveBuilder<sieve_type>::buildTopology(sieve, meshDim, numCells, cells, numVertices, interpolate, numCorners);
-      sieve->stratify();
-      topology->setPatch(0, sieve);
-      topology->stratify();
-      mesh->setTopology(topology);
-      ALE::New::SieveBuilder<sieve_type>::buildCoordinates(mesh->getRealSection("coordinates"), 
-						  meshDim, coordinates);
-      delete[] coordinates; coordinates = NULL;
-      delete[] cells; cells = NULL;
+      // Can now build mesh
+      _buildMesh(coordinates, numVertices, spaceDim,
+		 cells, numCells, numCorners, meshDim);
+      delete[] coordinates; coordinates = 0;
+      delete[] cells; cells = 0;
       builtTopology = true;
     } // if
 
@@ -132,7 +120,7 @@
 // ----------------------------------------------------------------------
 // Write mesh to file.
 void
-pylith::meshio::MeshIOAscii::write(const ALE::Obj<Mesh>& mesh) const
+pylith::meshio::MeshIOAscii::_write(void)
 { // write
   std::ofstream fileout(_filename.c_str());
   if (!fileout.is_open() || !fileout.good()) {
@@ -142,15 +130,24 @@
     throw std::runtime_error(msg.str());
   } // if
 
-  const int dimension = mesh->getDimension();
+  int meshDim = 0;
+  int spaceDim = 0;
+  int numVertices = 0;
+  int numCells = 0;
+  int numCorners = 0;
+  double* coordinates = 0;
+  int* cells = 0;
 
+  _getVertices(&coordinates, &numVertices, &spaceDim);
+  _getCells(&cells, &numCells, &numCorners, &meshDim);
+  
   fileout
     << "mesh = {\n"
-    << "  dimension = " << dimension << "\n"
+    << "  dimension = " << meshDim << "\n"
     << "  use-index-zero = " << (useIndexZero() ? "true" : "false") << "\n";
 
-  _writeVertices(fileout, mesh);
-  _writeCells(fileout, mesh);
+  _writeVertices(fileout, coordinates, numVertices, spaceDim);
+  _writeCells(fileout, cells, numCells, numCorners);
 
   // LOOP OVER GROUPS
   // _writeGroup(fileout, mesh, nameIter->c_str());
@@ -217,35 +214,22 @@
 // Write mesh vertices.
 void
 pylith::meshio::MeshIOAscii::_writeVertices(std::ostream& fileout,
-					    const ALE::Obj<Mesh>& mesh) const
+					    const double* coordinates,
+					    const int numVertices,
+					    const int spaceDim) const
 { // _writeVertices
-  const ALE::Obj<Mesh::real_section_type>& 
-    coords_field = mesh->getRealSection("coordinates");
-  const ALE::Obj<Mesh::topology_type>& topology = mesh->getTopology();
-  const Mesh::real_section_type::patch_type patch = 0;
-  const ALE::Obj<Mesh::topology_type::label_sequence>& vertices = 
-    topology->depthStratum(patch, 0);
-  const int embedDim = 
-    coords_field->getFiberDimension(patch, *vertices->begin());
-
   fileout
     << "  vertices = {\n"
-    << "    dimension = " << embedDim << "\n"
-    << "    count = " << vertices->size() << "\n"
+    << "    dimension = " << spaceDim << "\n"
+    << "    count = " << numVertices << "\n"
     << "    coordinates = {\n"
     << std::resetiosflags(std::ios::fixed)
     << std::setiosflags(std::ios::scientific)
     << std::setprecision(6);
-  for(Mesh::topology_type::label_sequence::iterator v_iter = vertices->begin();
-      v_iter != vertices->end();
-      ++v_iter) {
-    const Mesh::real_section_type::value_type *coordinates = 
-      coords_field->restrict(patch, *v_iter);
-
+  for(int iVertex=0, i=0; iVertex < numVertices; ++iVertex) {
     fileout << "      ";
-    for(int d = 0; d < embedDim; ++d) {
-      fileout << std::setw(18) << coordinates[d];
-    }
+    for(int iDim=0; iDim < spaceDim; ++iDim)
+      fileout << std::setw(18) << coordinates[i++];
     fileout << "\n";
   } // for
   fileout
@@ -319,36 +303,21 @@
 // Write mesh cells.
 void
 pylith::meshio::MeshIOAscii::_writeCells(std::ostream& fileout,
-					 const ALE::Obj<Mesh>& mesh) const
+					 const int* cells,
+					 const int numCells,
+					 const int numCorners) const
 { // _writeCells
-  const ALE::Obj<topology_type>& topology = mesh->getTopology();
-  const topology_type::patch_type patch = 0;
-  const ALE::Obj<sieve_type>& sieve = topology->getPatch(patch);
-  const ALE::Obj<Mesh::topology_type::label_sequence>& cells = 
-    topology->heightStratum(patch, 0);
-  const ALE::Obj<Mesh::numbering_type>& vNumbering = 
-    mesh->getFactory()->getLocalNumbering(topology, patch, 0);
-  const int numCorners = sieve->nCone(*cells->begin(), 
-				      topology->depth())->size();
-
   fileout
     << "  cells = {\n"
-    << "    count = " << cells->size() << "\n"
+    << "    count = " << numCells << "\n"
     << "    num-corners = " << numCorners << "\n"
     << "    simplices = {\n";
 
-  const int offset = (useIndexZero()) ? 0 : 1;
-  for(Mesh::topology_type::label_sequence::iterator e_iter = cells->begin();
-      e_iter != cells->end();
-      ++e_iter) {
+  for(int iCell=0, i=0; iCell < numCells; ++iCell) {
     fileout << "      ";
-    const ALE::Obj<sieve_type::traits::coneSequence>& cone = 
-      sieve->cone(*e_iter);
-    for(sieve_type::traits::coneSequence::iterator c_iter = cone->begin();
-	c_iter != cone->end();
-	++c_iter)
+    for (int iCorner=0; iCorner < numCorners; ++iCorner)
       fileout << std::setw(8)
-	      << vNumbering->getIndex(*c_iter) + offset;
+	      << cells[i++];
     fileout << "\n";
   } // for
   fileout
@@ -356,6 +325,7 @@
     << "  }\n";
 } // _writeCells
 
+#if 0
 // ----------------------------------------------------------------------
 // Read mesh group.
 void
@@ -403,7 +373,6 @@
   if (!filein.good())
     throw std::runtime_error("I/O error while parsing group settings.");
 
-#if 0
   assert(!mesh.isNull());
   ALE::Obj<Mesh::field_type> groupField = mesh->getField(name);
   const int meshDim = mesh->getDimension();
@@ -423,7 +392,6 @@
   const double zero = 0;
   for (int i=0; i < count; ++i)
     groupField->update(patch, Mesh::point_type(0, i), &zero);
-#endif
 } // _readGroup
 
 // ----------------------------------------------------------------------
@@ -449,5 +417,6 @@
     << "    }\n"
     << "  }\n";
 } // _writeGroup
+#endif
   
 // End of file 

Modified: short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.hh	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/libsrc/meshio/MeshIOAscii.hh	2006-11-06 03:42:14 UTC (rev 5179)
@@ -15,7 +15,6 @@
 
 #include <iosfwd> // USES std::istream, std::ostream
 #include <string> // HASA std::string
-#include "MeshIO.hh"
 
 namespace pylith {
   namespace meshio {
@@ -25,11 +24,6 @@
 
 class pylith::meshio::MeshIOAscii : public pylith::meshio::MeshIO
 { // MeshIOAscii
-// PUBLIC TYPEDEFS -------------------------------------------------------
-public :
-  typedef ALE::Mesh                Mesh;
-  typedef ALE::Mesh::sieve_type    sieve_type;
-  typedef ALE::Mesh::topology_type topology_type;
 // PUBLIC METHODS -------------------------------------------------------
 public :
 
@@ -51,20 +45,15 @@
    */
   const char* filename(void) const;
 
-  /** Read mesh from file.
-   *
-   * @param pMesh Pointer to PETSc mesh object
-   * @param interpolate Flag indicating whether to build intermediate topology
-   */
-  void read(ALE::Obj<Mesh>& mesh, 
-	    const bool interpolate =false);
+// PROTECTED METHODS ----------------------------------------------------
+protected :
 
-  /** Write mesh to file.
-   *
-   * @param mesh PETSc mesh object
-   */
-  void write(const ALE::Obj<Mesh>& mesh) const;
+  /// Write mesh
+  void _write(void) const;
 
+  /// Read mesh
+  void _read(void);
+
 // PRIVATE METHODS ------------------------------------------------------
 private :
 
@@ -73,60 +62,50 @@
    * @param filein Input stream
    * @param pCoordinates Pointer to array of vertex coordinates
    * @param pNumVertices Pointer to number of vertices
-   * @param pNumDims Pointer to number of dimensions
+   * @param pSpaceDim Pointer to dimension of coordinates vector space
    */
   void _readVertices(std::istream& filein,
 		     double** pCoordinates,
-		     int* pNumVertices, 
-		     int* pNumDims) const;
+		     int* pNumVertices,
+		     int* pSpaceDim) const;
   
   /** Write mesh vertices.
    *
    * @param fileout Output stream
-   * @param mesh PETSc mesh
+   * @param coordinates Array of vertex coordinates
+   * @param numVertices Number of vertices
+   * @param spaceDim Dimension of coordinates vector space
    */
   void _writeVertices(std::ostream& fileout,
-		      const ALE::Obj<Mesh>& mesh) const;
+		      const double* coordinates,
+		      const int numVertices,
+		      const int spaceDim) const;
   
   /** Read mesh cells.
    *
    * @param filein Input stream
-   * @param pCells Pointer to array of cells (vertices in each element)
+   * @param pCells Pointer to array of indices of cell vertices
    * @param pNumCells Pointer to number of cells
-   * @param pNumCorners Pointer to number of corners in each element
+   * @param pNumCorners Pointer to number of corners
    */
   void _readCells(std::istream& filein,
 		  int** pCells,
-		  int* pNumCells, 
+		  int* pNumCells,
 		  int* pNumCorners) const;
   
   /** Write mesh cells.
    *
    * @param fileout Output stream
-   * @param mesh PETSc mesh
+   * @param cells Array of indices of cell vertices
+   * @param numCells Number of cells
+   * @param numCorners Number of corners
    */
   void _writeCells(std::ostream& fileout,
-		   const ALE::Obj<Mesh>& mesh) const;
+		   const int* cells,
+		   const int numCells,
+		   const int numCorners) const;
 
-  /** Read mesh group.
-   *
-   * @param filein Output stream
-   * @param pMesh Pointer to PETSc mesh
-   */
-  void _readGroup(std::istream& filein,
-		  const ALE::Obj<Mesh>& pMesh) const;
-
-  /** Write mesh group.
-   *
-   * @param fileout Output stream
-   * @param mesh PETSc mesh
-   * @param name Name of group
-   */
-  void _writeGroup(std::ostream& fileout,
-		   const ALE::Obj<Mesh>& mesh,
-		   const char* name) const;
-
-// PRIVATE MEMBERS ------------------------------------------------------
+  // PRIVATE MEMBERS ----------------------------------------------------
 private :
 
   std::string _filename; ///< Name of file
@@ -137,7 +116,4 @@
 
 #endif // pylith_meshio_meshioascii_hh
 
-// version
-// $Id$
-
 // End of file 

Modified: short/3D/PyLith/trunk/pylith/feassemble/Field.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/Field.py	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/pylith/feassemble/Field.py	2006-11-06 03:42:14 UTC (rev 5179)
@@ -10,31 +10,15 @@
 # ----------------------------------------------------------------------
 #
 
-## @file pyre/feassemble/Field.py
+## @file pylith/feassemble/Field.py
+
 ## @brief Python PyLith field.
-import FIAT.shapes
 
-elementShapes = {'LINE':        FIAT.shapes.LINE,
-                 'TRIANGLE':    FIAT.shapes.TRIANGLE,
-                 'TETRAHEDRON': FIAT.shapes.TETRAHEDRON}
-
-from pyre.components.Component import Component
-
-def validateFamily(value):
-  try:
-    __import__('FIAT.'+str(value))
-  except ImportError:
-    raise ValueError, 'Invalid element family: '+str(value)
-  return value
-
-def validateShape(value):
-  if not str(value).upper() in elementShapes:
-    raise ValueError, 'Invalid element shape: '+str(value)
-  return value.upper()
-
 # Field class
 class Field(Component):
-  """Python finite-element assembler."""
+  """
+  Python PyLith field.
+  """
 
   # INVENTORY //////////////////////////////////////////////////////////
 
@@ -45,47 +29,18 @@
     ## Python object for managing Field facilities and properties.
     ##
     ## \b Properties
-    ## @li \b family Element family for field
-    ## @li \b family_order Order of element family
-    ## @li \b quadrature_order Order for quadrature
-    ## @li \b shape Element shape
+    ## @li None
     ##
     ## \b Facilities
     ## @li None
 
     import pyre.inventory
 
-    family = pyre.inventory.str("family", default="Lagrange",
-                                validator=validateFamily)
-    family.meta['tip'] = "Element family for field"
-
-    familyOrder = pyre.inventory.int("family_order", default=1,
-                                     validator=pyre.inventory.greaterEqual(0))
-    familyOrder.meta['tip'] = "Order of element family"
-
-    quadratureOrder = pyre.inventory.int("quadrature_order", default=1,
-                                         validator=pyre.inventory.greater(0))
-    quadratureOrder.meta['tip'] = "Order for quadrature."
-
-    shape = pyre.inventory.str("shape", default="tetrahedron",
-                               validator=validateShape)
-    shape.meta['tip'] = "Element shape."
-
   # PUBLIC METHODS /////////////////////////////////////////////////////
 
   def initialize(self):
-    """Setup basis fns and quadrature info."""
-    import FIAT.quadrature
-    from pylith.utils import importing
-
-    self._info.log('Creating the '+self.family+' element of order ' +
-                   str(self.familyOrder))
-    self.element = getattr(importing.importModule('FIAT.'+self.family),
-                           self.family)(self.shape, self.familyOrder)
-    self.quadrature = FIAT.quadrature.make_quadrature_by_degree(self.shape,
-                                                        2*self.element.n - 1)
-    self._info.log('Created quadrature of order '+str(self.quadrature.degree))
-    self._evaluateBasis()
+    """
+    """
     return
 
 
@@ -93,41 +48,17 @@
     """Constructor."""
     Component.__init__(self, name, facility="field")
     self.sieveField = None
-    self.element = None
-    self.quadrature = None
-    self.basis = None
-    self.basisder = None
     return
 
 
   # PRIVATE METHODS /////////////////////////////////////////////////////
 
-  def _evaluateBasis(self):
-    '''Return Numeric arrays with the basis functions and their
-    derivatives evalauted at the quadrature points
-       - FIAT uses a reference element of (-1,-1):(1,-1):(-1,1)'''
-    import Numeric
 
-    points = self.quadrature.get_points()
-    basis = self.element.function_space()
-    dim = FIAT.shapes.dimension(basis.base.shape)
-    self.basis = Numeric.transpose(basis.tabulate(points))
-    self.basisDer = Numeric.transpose([basis.deriv_all(d).tabulate(points) \
-                                       for d in range(dim)])
-    # We now assume that all simplices of the same dimension are equal
-    self.elementDof = [len(self.element.dual_basis().getNodeIDs(d)[0]) for d in range(dim)]
-    return
-
   def _configure(self):
-    """Set members based using inventory."""
-    self.family = self.inventory.family
-    self.familyOrder = self.inventory.familyOrder
-    self.quadratureOrder = self.inventory.quadratureOrder
-    self.shape = elementShapes[self.inventory.shape]
+    """
+    Set members based using inventory.
+    """
     return
   
 
-# version
-__id__ = "$Id$"
-
 # End of file 

Modified: short/3D/PyLith/trunk/pylith/feassemble/Integrator.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/Integrator.py	2006-11-04 17:52:24 UTC (rev 5178)
+++ short/3D/PyLith/trunk/pylith/feassemble/Integrator.py	2006-11-06 03:42:14 UTC (rev 5179)
@@ -10,46 +10,72 @@
 # ----------------------------------------------------------------------
 #
 
-## @file pyre/feassemble/Integrator.py
-## @brief Python abstract base class for finite-element integration.
+## @file pylith/feassemble/Integrator.py
 
+## @brief Python abstract base class for integration of operator
+## actions with finite-elements.
+
 from pyre.components.Component import Component
 
 # Integrator class
 class Integrator(Component):
-  """Python abstract base class for finite-element integration."""
+  """
+  Python abstract base class for integration of actions with
+  finite-elements.
+  """
 
-  # PUBLIC METHODS /////////////////////////////////////////////////////
+  # INVENTORY //////////////////////////////////////////////////////////
 
-  # INITIALIZE
-  # set mesh, basis, basisDer, quadrature
+  class Inventory(Component.Inventory):
+    """Python object for managing Integrator facilities and properties."""
 
-  def integrateResidual(self, element, state):
-    """Integrate residual for element."""
-    raise NotImplementedError, \
-          "Integrator::integrateResidual() not implemented."
-    return
+    ## @class Inventory
+    ## Python object for managing Integrator facilities and properties.
+    ##
+    ## \b Properties
+    ## @li None
+    ##
+    ## \b Facilities
+    ## @li \b quadrature Quadrature object for integration
 
+    import pyre.inventory
 
-  def integrateJacobian(self, element, state):
-    """Integrate the Jacobian weak form over an element using the given
-    quadrature."""
-    raise NotImplementedError, \
-          "Integrator::integrateJacobian() not implemented."
-    return
+    from Quadrature import Quadrature
+    quadrature = pyre.inventory.facility("quadrature", factory=Quadrature)
+    quadrature.meta['tip'] = "Quadrature object for integration."
 
 
+  # PUBLIC METHODS /////////////////////////////////////////////////////
+
   def __init__(self, name="integrator"):
-    """Constructor."""
+    """
+    Constructor.
+    """
     Component.__init__(self, name, facility="integrator")
-    self.mesh = None
-    self.basis = None
-    self.basisDer = None
+    self.cppHandle = None
     self.quadrature = None
     return
 
 
-# version
-__id__ = "$Id$"
+  def initialize(self):
+    """
+    Initialize C++ integrator object.
+    """
+    q = self.quadrature
+    q.initialize()
+    self.cppHandle.quadrature(q.cppHandle)
+    return
+  
+  
+  # PRIVATE METHODS ////////////////////////////////////////////////////
 
+  def _configure(self):
+    """
+    Set members based using inventory.
+    """
+    Component._configure(self)
+    self.quadrature = self.inventory.quadrature
+    return
+
+
 # End of file 



More information about the cig-commits mailing list