[cig-commits] r12385 - in short/3D/PyLith/trunk: libsrc/feassemble libsrc/materials pylith/problems unittests/libtests/feassemble unittests/pytests/bc unittests/pytests/faults unittests/pytests/feassemble unittests/pytests/problems
brad at geodynamics.org
brad at geodynamics.org
Mon Jul 7 16:14:51 PDT 2008
Author: brad
Date: 2008-07-07 16:14:50 -0700 (Mon, 07 Jul 2008)
New Revision: 12385
Modified:
short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.cc
short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.hh
short/3D/PyLith/trunk/libsrc/feassemble/Integrator.cc
short/3D/PyLith/trunk/libsrc/feassemble/Integrator.hh
short/3D/PyLith/trunk/libsrc/feassemble/Integrator.icc
short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh
short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh
short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh
short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh
short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.hh
short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.hh
short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.cc
short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.hh
short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc
short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.hh
short/3D/PyLith/trunk/pylith/problems/Formulation.py
short/3D/PyLith/trunk/pylith/problems/TimeStep.py
short/3D/PyLith/trunk/pylith/problems/TimeStepAdapt.py
short/3D/PyLith/trunk/pylith/problems/TimeStepUniform.py
short/3D/PyLith/trunk/pylith/problems/TimeStepUser.py
short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityExplicit.cc
short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityImplicit.cc
short/3D/PyLith/trunk/unittests/pytests/bc/TestAbsorbingDampers.py
short/3D/PyLith/trunk/unittests/pytests/bc/TestNeumann.py
short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py
short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityExplicit.py
short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityImplicit.py
short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStep.py
short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStepAdapt.py
Log:
Added stability_factor to adaptive time stepping to permit flexibility in adpative time stepping. Added calculation of stable time step for integrators and materials. Default is an extremely large time step. Still need to add unit tests for ElasticMaterial::stableTimeStepImplicit().
Modified: short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -57,6 +57,15 @@
} // timeStep
// ----------------------------------------------------------------------
+// Get stable time step for advancing from time t to time t+dt.
+double
+pylith::feassemble::ElasticityImplicit::stableTimeStep(void) const
+{ // stableTimeStep
+ assert(0 != _material);
+ return _material->stableTimeStepImplicit();
+} // stableTimeStep
+
+// ----------------------------------------------------------------------
// Set flag for setting constraints for total field solution or
// incremental field solution.
void
Modified: short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/feassemble/ElasticityImplicit.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -85,6 +85,14 @@
*/
void timeStep(const double dt);
+ /** Get stable time step for advancing from time t to time t+dt.
+ *
+ * Default is current time step.
+ *
+ * @returns Time step
+ */
+ double stableTimeStep(void) const;
+
/** Set flag for setting constraints for total field solution or
* incremental field solution.
*
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Integrator.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Integrator.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Integrator.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -20,6 +20,11 @@
#include <assert.h> // USES assert()
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
// Constructor
pylith::feassemble::Integrator::Integrator(void) :
@@ -61,10 +66,18 @@
pylith::feassemble::Integrator::gravityField(spatialdata::spatialdb::GravityField* const gravityField)
{ // gravityField
_gravityField = gravityField;
-
} // gravityField
// ----------------------------------------------------------------------
+// Get stable time step for advancing from time t to time t+dt.
+double
+pylith::feassemble::Integrator::stableTimeStep(void) const
+{ // stableTimeStep
+ // Assume any time step will work.
+ return MAXFLOAT;
+} // stableTimeStep
+
+// ----------------------------------------------------------------------
// Initialize vector containing result of integration action for cell.
void
pylith::feassemble::Integrator::_initCellVector(void)
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Integrator.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Integrator.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Integrator.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -85,7 +85,7 @@
/** Get stable time step for advancing from time t to time t+dt.
*
- * Default is current time step.
+ * Default is MAXFLOAT (or 1.0e+30 if MAXFLOAT is not defined in math.h).
*
* @returns Time step
*/
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Integrator.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Integrator.icc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Integrator.icc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -21,14 +21,6 @@
_dt = dt;
} // timeStep
-// Get stable time step for advancing from time t to time t+dt.
-inline
-double
-pylith::feassemble::Integrator::stableTimeStep(void) const {
- // Default is current time step
- return _dt;
-} // stableTimeStep
-
// Check whether Jacobian needs to be recomputed.
inline
bool
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -22,6 +22,11 @@
#include <sstream> // USES std::ostringstream
#include <stdexcept> // USES std::runtime_error
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
namespace pylith {
namespace materials {
@@ -223,5 +228,14 @@
PetscLogFlops(2);
} // _calcElasticConsts
+// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::ElasticIsotropic3D::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ return MAXFLOAT;
+} // _stableTimeStepImplicit
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -110,6 +110,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -20,6 +20,11 @@
#include <string.h> // USES memcpy()
#include <assert.h> // USES assert()
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
// Default constructor.
pylith::materials::ElasticMaterial::ElasticMaterial(const int tensorSize,
@@ -102,6 +107,27 @@
} // calcDerivElastic
// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::ElasticMaterial::stableTimeStepImplicit(void) const
+{ // stableTimeStepImplicit
+ const int numQuadPts = _numQuadPts;
+ const int totalPropsQuadPt = _totalPropsQuadPt;
+ assert(_propertiesCell.size() == numQuadPts*totalPropsQuadPt);
+
+ double dtStable = MAXFLOAT;
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ const double dt =
+ _stableTimeStepImplicit(&_propertiesCell[iQuad*totalPropsQuadPt],
+ totalPropsQuadPt);
+ if (dt < dtStable)
+ dtStable = dt;
+ } // for
+
+ return dtStable;
+} // stableTimeStepImplicit
+
+// ----------------------------------------------------------------------
// Get cell's property information from material's sections.
void
pylith::materials::ElasticMaterial::getPropertiesCell(const Mesh::point_type& cell,
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -125,6 +125,15 @@
const double_array&
calcDerivElastic(const double_array& totalStrain);
+ /** Get stable time step for implicit time integration.
+ *
+ * Default is MAXFLOAT (or 1.0e+30 if MAXFLOAT is not defined in math.h).
+ *
+ * @returns Time step
+ */
+ virtual
+ double stableTimeStepImplicit(void) const;
+
/** Update properties (for next time step).
*
* @param totalStrain Total strain tensor at quadrature points
@@ -202,6 +211,14 @@
const double* totalStrain,
const int strainSize) = 0;
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ virtual
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const = 0;
+
/** Update properties (for next time step).
*
* @param properties Properties at location.
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -20,6 +20,11 @@
#include <sstream> // USES std::ostringstream
#include <stdexcept> // USES std::runtime_error
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
namespace pylith {
namespace materials {
@@ -200,5 +205,14 @@
PetscLogFlops(2);
} // calcElasticConsts
+// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::ElasticPlaneStrain::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ return MAXFLOAT;
+} // _stableTimeStepImplicit
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -110,6 +110,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -20,6 +20,11 @@
#include <sstream> // USES std::ostringstream
#include <stdexcept> // USES std::runtime_error
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
namespace pylith {
namespace materials {
@@ -198,5 +203,14 @@
PetscLogFlops(8);
} // calcElasticConsts
+// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::ElasticPlaneStress::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ return MAXFLOAT;
+} // _stableTimeStepImplicit
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -110,6 +110,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -18,6 +18,11 @@
#include <assert.h> // USES assert()
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
namespace pylith {
namespace materials {
@@ -180,5 +185,14 @@
PetscLogFlops(2);
} // _calcElasticConsts
+// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::ElasticStrain1D::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ return MAXFLOAT;
+} // _stableTimeStepImplicit
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -109,6 +109,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -18,6 +18,11 @@
#include <assert.h> // USES assert()
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
namespace pylith {
namespace materials {
@@ -177,5 +182,14 @@
PetscLogFlops(6);
} // _calcElasticConsts
+// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::ElasticStress1D::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ return MAXFLOAT;
+} // _stableTimeStepImplicit
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -110,6 +110,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -25,6 +25,11 @@
#include <sstream> // USES std::ostringstream
#include <stdexcept> // USES std::runtime_error
+#include <math.h> // USES MAXFLOAT
+#if !defined(MAXFLOAT)
+#define MAXFLOAT 1.0e+30
+#endif
+
// ----------------------------------------------------------------------
namespace pylith {
namespace materials {
@@ -536,6 +541,29 @@
} // _calcElasticConstsViscoelastic
// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::GenMaxwellIsotropic3D::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ assert(0 != properties);
+ assert(_totalPropsQuadPt == numProperties);
+
+ double dtStable = MAXFLOAT;
+
+ const int numMaxwellModels = _GenMaxwellIsotropic3D::numMaxwellModels;
+ for (int i=0; i < numMaxwellModels; ++i) {
+ const double maxwellTime =
+ properties[_GenMaxwellIsotropic3D::pidMaxwellTime+i];
+ const double dt = 0.1*maxwellTime;
+ if (dt < dtStable)
+ dtStable = dt;
+ } // for
+
+ return dtStable;
+} // _stableTimeStepImplicit
+
+// ----------------------------------------------------------------------
// Update state variables.
void
pylith::materials::GenMaxwellIsotropic3D::_updatePropertiesElastic(
Modified: short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/GenMaxwellIsotropic3D.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -139,6 +139,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
/** Update properties (for next time step).
*
* @param properties Properties at location.
Modified: short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -411,6 +411,22 @@
} // _calcElasticConstsViscoelastic
// ----------------------------------------------------------------------
+// Get stable time step for implicit time integration.
+double
+pylith::materials::MaxwellIsotropic3D::_stableTimeStepImplicit(const double* properties,
+ const int numProperties) const
+{ // _stableTimeStepImplicit
+ assert(0 != properties);
+ assert(_totalPropsQuadPt == numProperties);
+
+ const double maxwellTime =
+ properties[_MaxwellIsotropic3D::pidMaxwellTime];
+ const double dtStable = 0.1*maxwellTime;
+
+ return dtStable;
+} // _stableTimeStepImplicit
+
+// ----------------------------------------------------------------------
// Update state variables.
void
pylith::materials::MaxwellIsotropic3D::_updatePropertiesElastic(
Modified: short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.hh 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.hh 2008-07-07 23:14:50 UTC (rev 12385)
@@ -128,6 +128,13 @@
const double* totalStrain,
const int strainSize);
+ /** Get stable time step for implicit time integration.
+ *
+ * @returns Time step
+ */
+ double _stableTimeStepImplicit(const double* properties,
+ const int numProperties) const;
+
/** Update properties (for next time step).
*
* @param properties Properties at location.
Modified: short/3D/PyLith/trunk/pylith/problems/Formulation.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/Formulation.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/pylith/problems/Formulation.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -216,14 +216,8 @@
logEvent = "%stimestep" % self._loggingPrefix
self._logger.eventBegin(logEvent)
- dt = 1.0e+30*second
- for integrator in self.integrators:
- stableDt = integrator.stableTimeStep()
- if dt < stableDt:
- dt = stableDt
+ dt = self.timeStep.timeStep(self.integrators)
- dt = self.timeStep.timeStep(dt)
-
self._logger.eventEnd(logEvent)
return dt
Modified: short/3D/PyLith/trunk/pylith/problems/TimeStep.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/TimeStep.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/pylith/problems/TimeStep.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -83,11 +83,12 @@
return 0
- def timeStep(self, dtStable):
+ def timeStep(self, integrators):
"""
- Adjust stable time step for advancing forward in time.
+ Get stable time step for advancing forward in time.
"""
- return dtStable
+ # Default is to use the current time step.
+ return self.dt
def currentStep(self):
Modified: short/3D/PyLith/trunk/pylith/problems/TimeStepAdapt.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/TimeStepAdapt.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/pylith/problems/TimeStepAdapt.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -44,6 +44,7 @@
## @li \b total_time Time duration for simulation.
## @li \b max_dt Maximum time step.
## @li \b adapt_skip Number of time steps to skip between adjusting value.
+ ## @li \b stability_factor "Safety factor" for stable time step.
##
## \b Facilities
## @li None
@@ -64,7 +65,11 @@
adaptSkip.meta['tip'] = "Number of time steps to skip between " \
"adjusting value."
+ stabilityFactor = pyre.inventory.float("stability_factor", default=1.2,
+ validator=pyre.inventory.greater(0.0))
+ stabilityFactor.meta['tip'] = "'Safety factor' for stable time step."
+
# PUBLIC METHODS /////////////////////////////////////////////////////
def __init__(self, name="timestepadapt"):
@@ -85,17 +90,23 @@
return nsteps
- def timeStep(self, dtStable):
+ def timeStep(self, integrators):
"""
Adjust stable time step for advancing forward in time.
"""
from pyre.units.time import second
+ dtStable = 1.0e+30*second
+ for integrator in integrators:
+ dt = integrator.stableTimeStep()
+ if dt < dtStable:
+ dtStable = dt
+
if self.skipped < self.adaptSkip and \
self.dt != 0.0*second and \
self.dt < dtStable:
self.skipped += 1
else:
- self.dt = min(dtStable, self.maxDt)
+ self.dt = min(dtStable/self.stabilityFactor, self.maxDt)
self.skipped = 0
return self.dt
@@ -110,6 +121,7 @@
self.totalTime = self.inventory.totalTime
self.maxDt = self.inventory.maxDt
self.adaptSkip = self.inventory.adaptSkip
+ self.stabilityFactor = self.inventory.stabilityFactor
self.dt = self.maxDt
return
Modified: short/3D/PyLith/trunk/pylith/problems/TimeStepUniform.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/TimeStepUniform.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/pylith/problems/TimeStepUniform.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -74,7 +74,7 @@
return nsteps
- def timeStep(self, dtStable):
+ def timeStep(self, integrators):
"""
Adjust stable time step for advancing forward in time.
"""
Modified: short/3D/PyLith/trunk/pylith/problems/TimeStepUser.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/TimeStepUser.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/pylith/problems/TimeStepUser.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -137,17 +137,15 @@
return nsteps
- def timeStep(self, dtStable):
+ def timeStep(self, integrators):
"""
- Adjust stable time step for advancing forward in time.
+ Get time step for advancing forward in time.
"""
self.dt = self.steps[self.index]
- self.index += 1
- if self.index >= len(self.steps):
- if self.loopSteps:
- self.index = 0
- else:
- self.index -= 1
+ if self.index+1 < len(self.steps):
+ self.index += 1
+ elif self.loopSteps:
+ self.index = 0
return self.dt
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityExplicit.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityExplicit.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityExplicit.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -85,10 +85,8 @@
{ // testStableTimeStep
ElasticityExplicit integrator;
- const double dt1 = 2.0;
- integrator.timeStep(dt1);
const double stableTimeStep = integrator.stableTimeStep();
- CPPUNIT_ASSERT_EQUAL(dt1, stableTimeStep);
+ CPPUNIT_ASSERT_EQUAL(1.0e+30, stableTimeStep);
} // testStableTimeStep
// ----------------------------------------------------------------------
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityImplicit.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityImplicit.cc 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/TestElasticityImplicit.cc 2008-07-07 23:14:50 UTC (rev 12385)
@@ -83,10 +83,15 @@
{ // testStableTimeStep
ElasticityImplicit integrator;
- const double dt1 = 2.0;
- integrator.timeStep(dt1);
+ materials::ElasticIsotropic3D material;
+ const int id = 3;
+ const std::string label("my material");
+ material.id(id);
+ material.label(label.c_str());
+ integrator.material(&material);
+
const double stableTimeStep = integrator.stableTimeStep();
- CPPUNIT_ASSERT_EQUAL(dt1, stableTimeStep);
+ CPPUNIT_ASSERT_EQUAL(1.0e+30, stableTimeStep);
} // testStableTimeStep
// ----------------------------------------------------------------------
Modified: short/3D/PyLith/trunk/unittests/pytests/bc/TestAbsorbingDampers.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/bc/TestAbsorbingDampers.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/bc/TestAbsorbingDampers.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -79,9 +79,7 @@
"""
(mesh, bc, fields) = self._initialize()
- dt = 0.24*second
- bc.timeStep(dt)
- self.assertEqual(dt, bc.stableTimeStep())
+ self.assertEqual(1.0e+30*second, bc.stableTimeStep())
return
Modified: short/3D/PyLith/trunk/unittests/pytests/bc/TestNeumann.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/bc/TestNeumann.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/bc/TestNeumann.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -79,9 +79,7 @@
"""
(mesh, bc, fields) = self._initialize()
- dt = 0.24*second
- bc.timeStep(dt)
- self.assertEqual(dt, bc.stableTimeStep())
+ self.assertEqual(1.0e+30*second, bc.stableTimeStep())
return
Modified: short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -107,11 +107,10 @@
"""
Test stableTimeStep().
"""
- from pyre.units.time import second
- dt = 2.4*second
(mesh, fault, fields) = self._initialize()
- fault.timeStep(dt)
- self.assertEqual(dt, fault.stableTimeStep())
+
+ from pyre.units.time import second
+ self.assertEqual(1.0e+30*second, fault.stableTimeStep())
return
Modified: short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityExplicit.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityExplicit.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityExplicit.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -102,11 +102,10 @@
"""
Test stableTimeStep().
"""
- from pyre.units.time import second
- dt = 2.3*second
(mesh, integrator, fields) = self._initialize()
- integrator.timeStep(dt)
- self.assertEqual(dt, integrator.stableTimeStep())
+
+ from pyre.units.time import second
+ self.assertEqual(1.0e+30*second, integrator.stableTimeStep())
return
Modified: short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityImplicit.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityImplicit.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/feassemble/TestElasticityImplicit.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -102,11 +102,10 @@
"""
Test stableTimeStep().
"""
- from pyre.units.time import second
- dt = 2.3*second
(mesh, integrator, fields) = self._initialize()
- integrator.timeStep(dt)
- self.assertEqual(dt, integrator.stableTimeStep())
+
+ from pyre.units.time import second
+ self.assertEqual(1.0e+30*second, integrator.stableTimeStep())
return
Modified: short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStep.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStep.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStep.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -89,7 +89,8 @@
tstep = TimeStep()
tstep._configure()
- self.assertEqual(0.5*second, tstep.timeStep(0.5*second))
+ integrators = None
+ self.assertEqual(0.0*second, tstep.timeStep(integrators))
return
Modified: short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStepAdapt.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStepAdapt.py 2008-07-07 17:07:57 UTC (rev 12384)
+++ short/3D/PyLith/trunk/unittests/pytests/problems/TestTimeStepAdapt.py 2008-07-07 23:14:50 UTC (rev 12385)
@@ -20,6 +20,17 @@
from pyre.units.time import second
# ----------------------------------------------------------------------
+class Integrator:
+
+ def __init__(self, dt):
+ self.dt = dt
+
+
+ def stableTimeStep(self):
+ return self.dt
+
+
+# ----------------------------------------------------------------------
class TestTimeStepAdapt(unittest.TestCase):
"""
Unit testing of TimeStepAdapt object.
@@ -47,39 +58,46 @@
"""
tstep = TimeStepAdapt()
tstep._configure()
-
+
tstep.adaptSkip = 2
+ integrators = [Integrator(2.0*second),
+ Integrator(0.5*second)]
# Set time step
- dt = 0.5*second
- self.assertEqual(dt, tstep.timeStep(dt))
+ dt = 0.5*second / 1.2
+ self.assertEqual(dt, tstep.timeStep(integrators))
- # Time step should not change
- self.assertEqual(dt, tstep.timeStep(2.0*second))
+ # Increase stable time step, but time step should not change (skipped)
+ integrators[1].dt = 0.8*second
+ self.assertEqual(dt, tstep.timeStep(integrators))
- # Reduce time step even if should have skipped
- dt = 0.2*second
- self.assertEqual(dt, tstep.timeStep(dt))
+ # Reduce time step even if though should have skipped
+ integrators[1].dt = 0.2*second
+ dt = 0.2*second / 1.2
+ self.assertEqual(dt, tstep.timeStep(integrators))
# Skip adjusting time step
- self.assertEqual(dt, tstep.timeStep(0.5*second))
+ integrators[1].dt = 0.8*second
+ self.assertEqual(dt, tstep.timeStep(integrators))
# Skip adjusting time step
- self.assertEqual(dt, tstep.timeStep(0.5*second))
+ self.assertEqual(dt, tstep.timeStep(integrators))
- # Adjust time step
- dt = 0.8*second
- self.assertEqual(dt, tstep.timeStep(dt))
+ # Adjust time step and stability factor
+ tstep.stabilityFactor = 2.0
+ dt = 0.8*second / 2.0
+ self.assertEqual(dt, tstep.timeStep(integrators))
# Skip adjusting time step
- self.assertEqual(dt, tstep.timeStep(2.0*second))
+ integrators[1].dt = 2.0*second
+ self.assertEqual(dt, tstep.timeStep(integrators))
# Skip adjusting time step
- self.assertEqual(dt, tstep.timeStep(2.0*second))
+ self.assertEqual(dt, tstep.timeStep(integrators))
# Adjust time step with value bigger than max
dt = 1.0*second
- self.assertEqual(dt, tstep.timeStep(2.0*second))
+ self.assertEqual(dt, tstep.timeStep(integrators))
return
@@ -95,8 +113,11 @@
self.assertEqual(10.0*second, tstep.currentStep())
- dt = 4.0*second
- tstep.timeStep(dt)
+
+ integrators = [Integrator(3.0*second),
+ Integrator(2.4*second)]
+ dt = 2.4*second / 1.2
+ tstep.timeStep(integrators)
self.assertEqual(dt, tstep.currentStep())
return
More information about the cig-commits
mailing list