[cig-commits] r17181 - in short/3D/PyLith/trunk: libsrc/topology unittests/libtests/topology

brad at geodynamics.org brad at geodynamics.org
Wed Sep 8 16:20:47 PDT 2010


Author: brad
Date: 2010-09-08 16:20:47 -0700 (Wed, 08 Sep 2010)
New Revision: 17181

Modified:
   short/3D/PyLith/trunk/libsrc/topology/FieldsNew.cc
   short/3D/PyLith/trunk/libsrc/topology/FieldsNew.hh
   short/3D/PyLith/trunk/libsrc/topology/UniformSectionDS.cc
   short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsNewMesh.cc
Log:
Fixed some bugs in FieldsNew and UniformSectionDS.

Modified: short/3D/PyLith/trunk/libsrc/topology/FieldsNew.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/topology/FieldsNew.cc	2010-09-08 20:38:25 UTC (rev 17180)
+++ short/3D/PyLith/trunk/libsrc/topology/FieldsNew.cc	2010-09-08 23:20:47 UTC (rev 17181)
@@ -94,7 +94,11 @@
   // Set fiber dimension
   const int fiberDim = _fiberDim();
   assert(fiberDim > 0);
+#if defined(USE_UNIFORMSECTION)
+  _section = new section_type(_mesh.comm(), fiberDim, _mesh.debug());
+#else
   _section = new section_type(_mesh.comm(), _mesh.debug());
+#endif
   assert(!_section.isNull());
 
   // Set spaces
@@ -104,6 +108,7 @@
        ++f_iter)
     _section->addSpace();
 
+  assert(!points.isNull());
   if (points->size() > 0) {
     const point_type pointMin = 
       *std::min_element(points->begin(), points->end());
@@ -142,7 +147,11 @@
   // Set fiber dimension
   const int fiberDim = _fiberDim();
   assert(fiberDim > 0);
+#if defined(USE_UNIFORMSECTION)
+  _section = new section_type(_mesh.comm(), fiberDim, _mesh.debug());
+#else
   _section = new section_type(_mesh.comm(), _mesh.debug());
+#endif
   assert(!_section.isNull());
 
   // Set spaces
@@ -205,7 +214,7 @@
 // ----------------------------------------------------------------------
 // Get field.
 template<typename mesh_type>
-const pylith::topology::Field<mesh_type>&
+const pylith::topology::Field<mesh_type, typename pylith::topology::FieldsNew<mesh_type>::section_type>&
 pylith::topology::FieldsNew<mesh_type>::get(const char* name) const
 { // get
   typename map_type::const_iterator f_iter = _fields.find(name);
@@ -219,14 +228,12 @@
   const int fibration = f_iter->second.fibration;
   assert(fibration >= 0 && fibration < _fields.size());
 
-  if (0 == f_iter->second.field) {
-    assert(!_section.isNull());
-    assert(!_section->getFibration(fibration).isNull());
-    f_iter->second.field = 
-      new Field<mesh_type>(_mesh, _section->getFibration(fibration), 
-			   f_iter->second.metadata);
-    assert(0 != f_iter->second.field);
-  } // if
+  delete f_iter->second.field; f_iter->second.field = 0;
+  assert(!_section.isNull());
+  f_iter->second.field = 
+    new Field<mesh_type>(_mesh, _section->getFibration(fibration), 
+			 f_iter->second.metadata);
+  assert(0 != f_iter->second.field);
 
   return *f_iter->second.field;
 } // get
@@ -234,7 +241,7 @@
 // ----------------------------------------------------------------------
 // Get field.
 template<typename mesh_type>
-pylith::topology::Field<mesh_type>&
+pylith::topology::Field<mesh_type, typename pylith::topology::FieldsNew<mesh_type>::section_type>&
 pylith::topology::FieldsNew<mesh_type>::get(const char* name)
 { // get
   typename map_type::iterator f_iter = _fields.find(name);
@@ -247,13 +254,11 @@
   const int fibration = f_iter->second.fibration;
   assert(fibration >= 0 && fibration < _fields.size());
 
-  if (0 == f_iter->second.field) {
-    assert(!_section.isNull());
-    assert(!_section->getFibration(fibration).isNull());
-    f_iter->second.field = 
-      new Field<mesh_type>(_mesh, _section->getFibration(fibration), 
-			   f_iter->second.metadata);
-  } // if
+  delete f_iter->second.field; f_iter->second.field = 0;
+  assert(!_section.isNull());
+  f_iter->second.field = 
+    new Field<mesh_type, section_type>(_mesh, _section->getFibration(fibration), 
+				       f_iter->second.metadata);
   assert(0 != f_iter->second.field);
 
   return *f_iter->second.field;

Modified: short/3D/PyLith/trunk/libsrc/topology/FieldsNew.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/topology/FieldsNew.hh	2010-09-08 20:38:25 UTC (rev 17180)
+++ short/3D/PyLith/trunk/libsrc/topology/FieldsNew.hh	2010-09-08 23:20:47 UTC (rev 17181)
@@ -27,6 +27,8 @@
 #if !defined(pylith_topology_fields_hh)
 #define pylith_topology_fields_hh
 
+#define USE_UNIFORMSECTION
+
 // Include directives ---------------------------------------------------
 #include "topologyfwd.hh" // forward declarations
 
@@ -42,7 +44,11 @@
   friend class TestFieldsNewMesh; // unit testing
   friend class TestFieldsNewSubMesh; // unit testing
 
+#if defined(USE_UNIFORMSECTION)
+  typedef typename mesh_type::RealUniformSection section_type;
+#else
   typedef typename mesh_type::RealSection section_type;
+#endif
 
 // PUBLIC MEMBERS ///////////////////////////////////////////////////////
 public :
@@ -105,13 +111,13 @@
    *
    * @param name Name of field.
    */
-  const Field<mesh_type>& get(const char* name) const;
+  const Field<mesh_type, section_type>& get(const char* name) const;
 	   
   /** Get field.
    *
    * @param name Name of field.
    */
-  Field<mesh_type>& get(const char* name);
+  Field<mesh_type, section_type>& get(const char* name);
 	   
   /** Get mesh associated with fields.
    *
@@ -140,7 +146,7 @@
     FieldBase::Metadata metadata;
     int fiberDim;
     int fibration; 
-    Field<mesh_type>* field;
+    Field<mesh_type, section_type>* field;
   }; // FieldInfo
 
 // PROTECTED TYPEDEFS ///////////////////////////////////////////////////

Modified: short/3D/PyLith/trunk/libsrc/topology/UniformSectionDS.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/topology/UniformSectionDS.cc	2010-09-08 20:38:25 UTC (rev 17180)
+++ short/3D/PyLith/trunk/libsrc/topology/UniformSectionDS.cc	2010-09-08 23:20:47 UTC (rev 17181)
@@ -16,6 +16,8 @@
 // ======================================================================
 //
 
+#include <cassert> // USES assert()
+
 // ----------------------------------------------------------------------
 template<typename point_type, 
 	 typename value_type, 
@@ -183,6 +185,10 @@
   this->_atlas->setChart(chart);
   int dim = _fiberDim;
   this->_atlas->updatePoint(*this->getChart().begin(), &dim);
+
+  const int numSpaces = this->_spaces.size();
+  for(int i=0; i < numSpaces; ++i)
+    this->_spaces[i]->setChart(chart);
 } // setChart
 
 // ----------------------------------------------------------------------
@@ -392,6 +398,7 @@
 { // allocatePoint
   this->_array = this->_allocator.allocate(this->sizeWithBC());
   this->_array -= this->getChart().min()*_fiberDim;
+  assert(this->_array);
   const index_type chartEnd = this->getChart().max()*_fiberDim;
   for(index_type i=this->getChart().min()*_fiberDim;
       i < chartEnd;
@@ -467,6 +474,7 @@
 void
 ALE::IUniformSectionDS<point_type, value_type, alloc_type>::zero(void)
 { // zero
+  assert(this->_array);
   memset(this->_array+(this->getChart().min()*_fiberDim), 0,
 	 this->sizeWithBC()*sizeof(value_type));
 } // zero
@@ -479,6 +487,7 @@
 const typename ALE::IUniformSectionDS<point_type, value_type, alloc_type>::values_type& 
 ALE::IUniformSectionDS<point_type, value_type, alloc_type>::restrictSpace(void) const
 { // restrictSpace
+  assert(this->_array);
   return this->_array;
 } // restrictSpace
 
@@ -490,6 +499,7 @@
 const typename ALE::IUniformSectionDS<point_type, value_type, alloc_type>::value_type*
 ALE::IUniformSectionDS<point_type, value_type, alloc_type>::restrictPoint(const point_type& p) const
 { // restrictPoint
+  assert(this->_array);
   if (!this->hasPoint(p))
     return this->_emptyValue.v;
   return &this->_array[p*_fiberDim];
@@ -504,6 +514,7 @@
 ALE::IUniformSectionDS<point_type, value_type, alloc_type>::updatePoint(const point_type& p,
 				    const value_type v[])
 { // updatePoint
+  assert(this->_array);
   for(int i = 0, idx = p*_fiberDim; i < _fiberDim; ++i, ++idx)
     this->_array[idx] = v[i];
 } // updatePoint
@@ -517,6 +528,7 @@
 ALE::IUniformSectionDS<point_type, value_type, alloc_type>::updateAddPoint(const point_type& p,
 				       const value_type v[])
 { // updateAddPoint
+  assert(this->_array);
   for(int i = 0, idx = p*_fiberDim; i < _fiberDim; ++i, ++idx)
     this->_array[idx] += v[i];
 } // updateAddPoint
@@ -607,6 +619,10 @@
 ALE::IUniformSectionDS<point_type, value_type, alloc_type>::addSpace(void)
 { // addSpace
   Obj<atlas_type> space = new atlas_type(this->comm(), this->debug());
+  assert(!space.isNull());
+
+  assert(!this->_atlas.isNull());
+  space->setChart(this->_atlas->getChart());
   this->_spaces.push_back(space);
 } // addSpace
 
@@ -619,6 +635,9 @@
 						    const point_type& p,
 						    const int space) const
 { // getFiberDimension
+  assert(space < this->_spaces.size());
+  assert(!this->_spaces[space].isNull());
+
   return *this->_spaces[space]->restrictPoint(p);
 } // getFiberDimension
 
@@ -632,9 +651,11 @@
 					    int dim,
 					    const int space)
 { // setFiberDimension
-  const index_type idx = -1;
+  assert(space < this->_spaces.size());
+  assert(!this->_spaces[space].isNull());
+
   this->_spaces[space]->addPoint(p);
-  this->_spaces[space]->updatePoint(p, &idx);
+  this->_spaces[space]->updatePoint(p, &dim);
 } // setFiberDimension
 
 // ----------------------------------------------------------------------
@@ -708,10 +729,9 @@
   assert(fiberDim > 0);
   Obj<IUniformSectionDS> field = new IUniformSectionDS(this->comm(),
 						       fiberDim,
-						       *chart.begin(),
-						       *chart.end(),
 						       this->debug());
-  
+  field->_atlas->setChart(chart);
+
   // Copy sizes
   const typename chart_type::const_iterator chartEnd = chart.end();
   for(typename chart_type::const_iterator c_iter=chart.begin();
@@ -721,14 +741,23 @@
     if (fDim)
       field->setFiberDimension(*c_iter, fDim);
   } // for
-  
+  field->allocatePoint();
+
   // Copy values
+  value_type* spaceValues = (fiberDim > 0) ? new value_type[fiberDim] : NULL;
   const chart_type& newChart = field->getChart();
   const typename chart_type::const_iterator newChartEnd = newChart.end();
   for(typename chart_type::const_iterator c_iter = newChart.begin();
       c_iter != newChartEnd;
       ++c_iter)
-    field->updatePoint(*c_iter, this->restrictPoint(*c_iter));
+    if (field->getFiberDimension(*c_iter) > 0) {
+      assert(fiberDim <= field->getFiberDimension(*c_iter));
+      const value_type* allValues = this->restrictPoint(*c_iter);
+      for (int i=0; i < fiberDim; ++i)
+	spaceValues[i] = allValues[i];
+      field->updatePoint(*c_iter, spaceValues);
+    } // if
+  delete[] spaceValues; spaceValues = NULL;
 
   return field;
 } // getFibration

Modified: short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsNewMesh.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsNewMesh.cc	2010-09-08 20:38:25 UTC (rev 17180)
+++ short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsNewMesh.cc	2010-09-08 23:20:47 UTC (rev 17181)
@@ -31,9 +31,14 @@
 
 // ----------------------------------------------------------------------
 typedef pylith::topology::FieldsNew<pylith::topology::Mesh> FieldsNewMesh;
-typedef pylith::topology::Mesh::RealSection RealSection;
 typedef pylith::topology::Mesh::SieveMesh SieveMesh;
 
+#if defined(USE_UNIFORMSECTION)
+typedef pylith::topology::Mesh::RealUniformSection section_type;
+#else
+typedef pylith::topology::Mesh::RealSection section_type;
+#endif
+
 // ----------------------------------------------------------------------
 void
 pylith::topology::TestFieldsNewMesh::setUp(void)
@@ -121,7 +126,7 @@
   const size_t size = 2;
   CPPUNIT_ASSERT_EQUAL(size, fields._fields.size());
 
-  const ALE::Obj<RealSection>& section = fields.section();
+  const ALE::Obj<section_type>& section = fields.section();
   CPPUNIT_ASSERT(!section.isNull());
   for (SieveMesh::label_sequence::iterator v_iter=vertices->begin();
        v_iter != vertices->end();
@@ -163,7 +168,7 @@
   const size_t size = 2;
   CPPUNIT_ASSERT_EQUAL(size, fields._fields.size());
 
-  const ALE::Obj<RealSection>& section = fields.section();
+  const ALE::Obj<section_type>& section = fields.section();
   CPPUNIT_ASSERT(!section.isNull());
   for (int i=0; i < nptsIn; ++i)
     CPPUNIT_ASSERT_EQUAL(fiberDim, section->getFiberDimension(verticesIn[i]));
@@ -188,7 +193,7 @@
   const size_t size = 2;
   CPPUNIT_ASSERT_EQUAL(size, fields._fields.size());
 
-  const ALE::Obj<RealSection>& section = fields.section();
+  const ALE::Obj<section_type>& section = fields.section();
   CPPUNIT_ASSERT(!section.isNull());
   const ALE::Obj<SieveMesh>& sieveMesh = _mesh->sieveMesh();
   CPPUNIT_ASSERT(!sieveMesh.isNull());
@@ -213,14 +218,14 @@
   fields.add("field B", "displacement", 4, FieldBase::OTHER, 2.0, true);
   fields.allocate(FieldBase::VERTICES_FIELD);
 
-  const Field<Mesh>& fieldA = fields.get("field A");
+  Field<Mesh, section_type>& fieldA = fields.get("field A");
   CPPUNIT_ASSERT_EQUAL(std::string("velocity"), std::string(fieldA.label()));
   CPPUNIT_ASSERT_EQUAL(FieldBase::VECTOR,
 		       fieldA.vectorFieldType());
   CPPUNIT_ASSERT_EQUAL(1.0, fieldA.scale());
   CPPUNIT_ASSERT_EQUAL(false, fieldA.addDimensionOkay());
 
-  const Field<Mesh>& fieldB = fields.get("field B");
+  Field<Mesh, section_type>& fieldB = fields.get("field B");
   CPPUNIT_ASSERT_EQUAL(std::string("displacement"), 
 		       std::string(fieldB.label()));
   CPPUNIT_ASSERT_EQUAL(FieldBase::OTHER,
@@ -241,14 +246,14 @@
   fields.add("field B", "displacement", 4, FieldBase::OTHER, 2.0, true);
   fields.allocate(FieldBase::VERTICES_FIELD);
 
-  const Field<Mesh>& fieldA = fields.get("field A");
+  const Field<Mesh, section_type>& fieldA = fields.get("field A");
   CPPUNIT_ASSERT_EQUAL(std::string("velocity"), std::string(fieldA.label()));
   CPPUNIT_ASSERT_EQUAL(FieldBase::VECTOR,
 		       fieldA.vectorFieldType());
   CPPUNIT_ASSERT_EQUAL(1.0, fieldA.scale());
   CPPUNIT_ASSERT_EQUAL(false, fieldA.addDimensionOkay());
 
-  const Field<Mesh>& fieldB = fields.get("field B");
+  const Field<Mesh, section_type>& fieldB = fields.get("field B");
   CPPUNIT_ASSERT_EQUAL(std::string("displacement"), 
 		       std::string(fieldB.label()));
   CPPUNIT_ASSERT_EQUAL(FieldBase::OTHER,



More information about the CIG-COMMITS mailing list