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

brad at geodynamics.org brad at geodynamics.org
Wed May 23 15:36:15 PDT 2007


Author: brad
Date: 2007-05-23 15:36:14 -0700 (Wed, 23 May 2007)
New Revision: 6948

Added:
   short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc
Modified:
   short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
   short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh
   short/3D/PyLith/trunk/libsrc/materials/Makefile.am
   short/3D/PyLith/trunk/libsrc/materials/Material.cc
   short/3D/PyLith/trunk/libsrc/materials/Material.hh
   short/3D/PyLith/trunk/libsrc/materials/Material.icc
   short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStrain.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStress.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStrain1D.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStress1D.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.cc
   short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.hh
Log:
Added timestep() and needNewJacobian() to Material and useElasticBehavior to ElasticMaterial. Implemented corresponding C++ unit tests.

Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -124,13 +124,14 @@
 // ----------------------------------------------------------------------
 // Update state variables (for next time step).
 void
-pylith::materials::ElasticMaterial::updateState(void)
+pylith::materials::ElasticMaterial::updateState(
+				const std::vector<double_array>& totalStrain)
 { // updateState
   const int numQuadPts = _numQuadPts;
   assert(_paramsCell.size() == numQuadPts);
 
   for (int iQuad=0; iQuad < numQuadPts; ++iQuad)
-    _updateState(&_paramsCell[iQuad]);
+    _updateState(&_paramsCell[iQuad], totalStrain[iQuad]);
 } // updateState
 
 // ----------------------------------------------------------------------
@@ -168,8 +169,11 @@
 // ----------------------------------------------------------------------
 // Update parameters (for next time step).
 void
-pylith::materials::ElasticMaterial::_updateState(std::vector<double_array>* const parameters)
+pylith::materials::ElasticMaterial::_updateState(
+				 std::vector<double_array>* const parameters,
+				 const double_array& totalStrain)
 { // _updateState
 } // _updateState
 
+
 // End of file 

Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh	2007-05-23 22:36:14 UTC (rev 6948)
@@ -118,9 +118,26 @@
   const std::vector<double_array>&
   calcDerivElastic(const std::vector<double_array>& totalStrain);
 
-  /// Update state variables (for next time step).
-  void updateState(void);
+  /** Update state variables (for next time step).
+   *
+   * @param totalStrain Total strain tensor at quadrature points
+   *    [numQuadPts][tensorSize]
+   */
+  void updateState(const std::vector<double_array>& totalStrain);
 
+  /** Set whether elastic or inelastic constitutive relations are used.
+   *
+   * @param flag True to use elastic, false to use inelastic.
+   */
+  void useElasticBehavior(const bool flag);
+
+  /** Get flag indicating whether elastic or inelastic constitutive
+   * relations are used.
+   *
+   * @returns True if using elastic and false if using inelastic.
+   */
+  bool useElasticBehavior(void) const;
+
   // PROTECTED METHODS //////////////////////////////////////////////////
 protected :
 
@@ -178,9 +195,11 @@
   /** Update parameters (for next time step).
    *
    * @param parameters Parameters at location.
+   * @param totalStrain Total strain at location.
    */
   virtual
-  void _updateState(std::vector<double_array>* const parameters);
+  void _updateState(std::vector<double_array>* const parameters,
+		    const double_array& totalStrain);
 
   // NOT IMPLEMENTED ////////////////////////////////////////////////////
 private :
@@ -230,8 +249,16 @@
    */
   std::vector<double_array> _elasticConsts;
 
+  /** Flag indicating whether elastic (true) or inelastic (false)
+   * constitutive model should be used. Always true for purely elastic
+   * materials.
+   */
+  bool _useElasticBehavior;
+
 }; // class ElasticMaterial
 
+#include "ElasticMaterial.icc" // inline methods
+
 #endif // pylith_materials_elasticmaterial_hh
 
 

Added: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -0,0 +1,32 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+//                           Brad T. Aagaard
+//                        U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ----------------------------------------------------------------------
+//
+
+#if !defined(pylith_materials_elasticmaterial_hh)
+#error "ElasticMaterial.icc can only be included from ElasticMaterial.hh"
+#endif
+
+// Set whether elastic or inelastic constitutive relations are used.
+inline
+void
+pylith::materials::ElasticMaterial::useElasticBehavior(const bool flag) {
+  _useElasticBehavior = flag;
+} // useElasticBehavior
+
+// Get flag indicating whether elastic or inelastic constitutive
+inline
+bool
+pylith::materials::ElasticMaterial::useElasticBehavior(void) const {
+  return _useElasticBehavior;
+} // useElasticBehavior
+
+
+// End of file 

Modified: short/3D/PyLith/trunk/libsrc/materials/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Makefile.am	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/Makefile.am	2007-05-23 22:36:14 UTC (rev 6948)
@@ -21,6 +21,7 @@
 	ElasticIsotropic3D.hh \
 	ElasticIsotropic3D.icc \
 	ElasticMaterial.hh \
+	ElasticMaterial.icc \
 	ElasticPlaneStrain.hh \
 	ElasticPlaneStrain.icc \
 	ElasticPlaneStress.hh \

Modified: short/3D/PyLith/trunk/libsrc/materials/Material.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -29,8 +29,10 @@
 // ----------------------------------------------------------------------
 // Default constructor.
 pylith::materials::Material::Material(void) :
+  _dt(0.0),
   _parameters(0),
   _dimension(0),
+  _needNewJacobian(false),
   _db(0),
   _id(0),
   _label("")
@@ -50,8 +52,10 @@
 // ----------------------------------------------------------------------
 // Copy constructor.
 pylith::materials::Material::Material(const Material& m) :
+  _dt(m._dt),
   _parameters(m._parameters),
   _dimension(m._dimension),
+  _needNewJacobian(m._needNewJacobian),
   _db(m._db),
   _id(m._id),
   _label(m._label)

Modified: short/3D/PyLith/trunk/libsrc/materials/Material.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.hh	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.hh	2007-05-23 22:36:14 UTC (rev 6948)
@@ -108,6 +108,13 @@
    */
   const std::string& label(void) const;
 
+  /** Set current time step.
+   *
+   * @param dt Current time step.
+   */
+  virtual
+  void timestep(const double dt);
+
   /** Initialize material by getting physical property parameters from
    * database.
    *
@@ -119,6 +126,13 @@
 		  const spatialdata::geocoords::CoordSys* cs,
 		  pylith::feassemble::Quadrature* quadrature);
   
+  /** Get flag indicating whether Jacobian matrix must be reformed for
+   * current state.
+   *
+   * @returns True if Jacobian matrix must be reformed, false otherwise.
+   */
+  bool needNewJacobian(void) const;
+
   // PROTECTED METHODS //////////////////////////////////////////////////
 protected :
 
@@ -178,10 +192,13 @@
   // PROTECTED MEMBERS //////////////////////////////////////////////////
 protected :
 
+  double _dt; ///< Current time step
+
   ///< Manager of parameters for physical properties of material
   topology::FieldsManager* _parameters;
 
   int _dimension; ///< Spatial dimension associated with material
+  bool _needNewJacobian; ///< True if need to reform Jacobian, false otherwise
 
   // PRIVATE MEMBERS ////////////////////////////////////////////////////
 private :

Modified: short/3D/PyLith/trunk/libsrc/materials/Material.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.icc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.icc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -56,5 +56,21 @@
   return _label;
 }
 
+// Set current time step.
+inline
+void
+pylith::materials::Material::timestep(const double dt) {
+  _dt = dt;
+}
 
+/* Get flag indicating whether Jacobian matrix must be reformed for
+ * current state.
+ */
+inline
+bool
+pylith::materials::Material::needNewJacobian(void) const {
+  return _needNewJacobian;
+}
+
+
 // End of file 

Modified: short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -202,12 +202,6 @@
   const double lambda2mu = lambda + 2.0 * mu;
   const double bulkmodulus = lambda + 2.0 * mu/3.0;
 
-  // The following three definitions are dummies until we figure out how to get the
-  // information into the function.
-  const bool elastic = true;
-  const double dt = 1.0;
-  const bool updatestate = true;
-  
   const double e11 = totalStrain[0];
   const double e22 = totalStrain[1];
   const double e33 = totalStrain[2];
@@ -219,7 +213,7 @@
   const double emeanTpdt = etraceTpdt/3.0;
   const double s123 = lambda * etraceTpdt;
 
-  if (elastic) {
+  if (useElasticBehavior()) {
     (*stress)[0] = s123 + 2.0*mu*e11;
     (*stress)[1] = s123 + 2.0*mu*e22;
     (*stress)[2] = s123 + 2.0*mu*e33;
@@ -227,33 +221,27 @@
     (*stress)[4] = 2.0 * mu * e23;
     (*stress)[5] = 2.0 * mu * e13;
   } else {
-    const double diag[6] = {1.0, 1.0, 1.0, 0.0, 0.0, 0.0};
+    const double diag[] = { 1.0, 1.0, 1.0, 0.0, 0.0, 0.0 };
     const double smeanTpdt = bulkmodulus * etraceTpdt;
     // The code below should probably be in a separate function since it
     // is used more than once.  I should also probably cover the possibility
     // that Maxwell time is zero (although this should never happen).
     const double timeFrac = 1.0e-5;
-    const double numTerms = 5.0;
+    const int numTerms = 5;
     double dq = 0.0;
-    if(maxwelltime < timeFrac*dt) {
+    if(maxwelltime < timeFrac*_dt) {
       double fSign = 1.0;
       double factorial = 1.0;
       double fraction = 1.0;
       dq = 1.0;
-      double term = 2.0;
-      while (term <= numTerms) {
-	factorial = factorial*term;
-	fSign = -1.0 * fSign;
-	fraction = fraction * (dt/maxwelltime);
-	dq = dq + fSign*fraction/factorial;
-	term++;
-      }
-    } else {
+      for (int iTerm=0; iTerm < numTerms; ++iTerm) {
+	factorial *= iTerm;
+	fSign *= -1.0;
+	fraction *= _dt/maxwelltime;
+	dq += fSign*fraction/factorial;
+      } // for
+    } else
       dq = maxwelltime*(1.0-exp(-dt/maxwelltime))/dt;
-    }
-
-
-
 } // _calcStress
 
 // ----------------------------------------------------------------------

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -87,7 +87,9 @@
 pylith::materials::TestElasticIsotropic3D::testUpdateState(void)
 { // testUpdateState
   ElasticIsotropic3D material;
-  material.updateState();
+
+  std::vector<double_array> totalStrain;
+  material.updateState(totalStrain);
 } // testUpdateState
 
 

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -44,6 +44,22 @@
 } // testClone
 
 // ----------------------------------------------------------------------
+// Test useElasticBehavior()
+void
+pylith::materials::TestElasticMaterial::testUseElasticBehavior(void)
+{ // testUseElasticBehavior
+  ElasticIsotropic3D material;
+
+  bool flag = true;
+  material.useElasticBehavior(flag);
+  CPPUNIT_ASSERT_EQUAL(flag, material._useElasticBehavior);
+
+  flag = false;
+  material.useElasticBehavior(flag);
+  CPPUNIT_ASSERT_EQUAL(flag, material._useElasticBehavior);
+} // testUseElasticBehavior
+
+// ----------------------------------------------------------------------
 // Test calcDensity()
 void
 pylith::materials::TestElasticMaterial::testCalcDensity(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh	2007-05-23 22:36:14 UTC (rev 6948)
@@ -39,6 +39,7 @@
   // CPPUNIT TEST SUITE /////////////////////////////////////////////////
   CPPUNIT_TEST_SUITE( TestElasticMaterial );
   CPPUNIT_TEST( testClone );
+  CPPUNIT_TEST( testUseElasticBehavior );
   CPPUNIT_TEST( testCalcDensity );
   CPPUNIT_TEST( testCalcStress );
   CPPUNIT_TEST( testCalcDerivElastic );
@@ -50,6 +51,9 @@
   /// Test clone()
   void testClone(void);
 
+  /// Test useElasticBehavior()
+  void testUseElasticBehavior(void);
+
   /// Test calcDensity()
   void testCalcDensity(void);
 

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStrain.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStrain.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStrain.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -87,7 +87,9 @@
 pylith::materials::TestElasticPlaneStrain::testUpdateState(void)
 { // testUpdateState
   ElasticPlaneStrain material;
-  material.updateState();
+
+  std::vector<double_array> totalStrain;
+  material.updateState(totalStrain);
 } // testUpdateState
 
 

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStress.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStress.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticPlaneStress.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -87,7 +87,9 @@
 pylith::materials::TestElasticPlaneStress::testUpdateState(void)
 { // testUpdateState
   ElasticPlaneStress material;
-  material.updateState();
+
+  std::vector<double_array> totalStrain;
+  material.updateState(totalStrain);
 } // testUpdateState
 
 

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStrain1D.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStrain1D.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStrain1D.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -87,7 +87,9 @@
 pylith::materials::TestElasticStrain1D::testUpdateState(void)
 { // testUpdateState
   ElasticStrain1D material;
-  material.updateState();
+
+  std::vector<double_array> totalStrain;
+  material.updateState(totalStrain);
 } // testUpdateState
 
 

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStress1D.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStress1D.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticStress1D.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -87,7 +87,9 @@
 pylith::materials::TestElasticStress1D::testUpdateState(void)
 { // testUpdateState
   ElasticStress1D material;
-  material.updateState();
+
+  std::vector<double_array> totalStrain;
+  material.updateState(totalStrain);
 } // testUpdateState
 
 

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.cc	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.cc	2007-05-23 22:36:14 UTC (rev 6948)
@@ -73,6 +73,34 @@
 } // testLabel
     
 // ----------------------------------------------------------------------
+// Test timestep()
+void
+pylith::materials::TestMaterial::testTimestep(void) 
+{ // testTimestep
+  const double dt = 2.0;
+  ElasticIsotropic3D material;
+  material.timestep(dt);
+  
+  CPPUNIT_ASSERT_EQUAL(dt, material._dt);
+} // testTimestep
+
+// ----------------------------------------------------------------------
+// Test needNewJacobian()
+void
+pylith::materials::TestMaterial::testNeedNewJacobian(void)
+{ // testNeedNewJacobian
+  ElasticIsotropic3D material;
+
+  bool flag = false;
+  material._needNewJacobian = flag;
+  CPPUNIT_ASSERT_EQUAL(flag, material.needNewJacobian());
+
+  flag = true;
+  material._needNewJacobian = flag;
+  CPPUNIT_ASSERT_EQUAL(flag, material.needNewJacobian());
+} // testNeedNewJacobian
+
+// ----------------------------------------------------------------------
 // Test initialize()
 void
 pylith::materials::TestMaterial::testInitialize(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.hh	2007-05-23 20:53:23 UTC (rev 6947)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestMaterial.hh	2007-05-23 22:36:14 UTC (rev 6948)
@@ -41,6 +41,8 @@
   CPPUNIT_TEST( testDB );
   CPPUNIT_TEST( testID );
   CPPUNIT_TEST( testLabel );
+  CPPUNIT_TEST( testTimestep );
+  CPPUNIT_TEST( testNeedNewJacobian );
   CPPUNIT_TEST( testInitialize );
   CPPUNIT_TEST_SUITE_END();
 
@@ -56,6 +58,12 @@
   /// Test label()
   void testLabel(void);
 
+  /// Test timestep()
+  void testTimestep(void);
+
+  /// Test needNewJacobian()
+  void testNeedNewJacobian(void);
+
   /// Test initialize()
   void testInitialize(void);
 



More information about the cig-commits mailing list