[cig-commits] r19135 - short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults
brad at geodynamics.org
brad at geodynamics.org
Mon Oct 31 10:08:18 PDT 2011
Author: brad
Date: 2011-10-31 10:08:17 -0700 (Mon, 31 Oct 2011)
New Revision: 19135
Modified:
short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.cc
short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.hh
Log:
Merge from trunk.
Modified: short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.cc
===================================================================
--- short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.cc 2011-10-31 16:57:52 UTC (rev 19134)
+++ short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.cc 2011-10-31 17:08:17 UTC (rev 19135)
@@ -454,7 +454,8 @@
(scalar_array*,
const scalar_array&,
const scalar_array&,
- const scalar_array&);
+ const scalar_array&,
+ const bool);
assert(0 != fields);
assert(0 != _quadrature);
@@ -580,10 +581,12 @@
// Use fault constitutive model to compute traction associated with
// friction.
dLagrangeTpdtVertex = 0.0;
+ const bool iterating = true; // Iterating to get friction
CALL_MEMBER_FN(*this,
constrainSolnSpaceFn)(&dLagrangeTpdtVertex,
slipVertex, slipRateVertex,
- tractionTpdtVertex);
+ tractionTpdtVertex,
+ iterating);
// Rotate traction back to global coordinate system.
dLagrangeTpdtVertexGlobal = 0.0;
@@ -820,7 +823,8 @@
(scalar_array*,
const scalar_array&,
const scalar_array&,
- const scalar_array&);
+ const scalar_array&,
+ const bool);
assert(0 != fields);
assert(0 != _quadrature);
@@ -1040,10 +1044,12 @@
// Use fault constitutive model to compute traction associated with
// friction.
dLagrangeTpdtVertex = 0.0;
+ const bool iterating = false; // No iteration for friction in lumped soln
CALL_MEMBER_FN(*this,
constrainSolnSpaceFn)(&dLagrangeTpdtVertex,
slipVertex, slipRateVertex,
- tractionTpdtVertex);
+ tractionTpdtVertex,
+ iterating);
// Rotate traction back to global coordinate system.
dLagrangeTpdtVertexGlobal = 0.0;
@@ -1918,103 +1924,13 @@
} // _sensitivityUpdateSoln
// ----------------------------------------------------------------------
-// Solve slip/Lagrange multiplier sensitivity problem for case of lumped Jacobian in 1-D.
-void
-pylith::faults::FaultCohesiveDyn::_sensitivitySolveLumped1D(
- scalar_array* slip,
- const scalar_array& dLagrangeTpdt,
- const scalar_array& jacobianN,
- const scalar_array& jacobianP)
-{ // _sensitivitySolveLumped1D
- assert(0 != slip);
-
- // Sensitivity of slip to changes in the Lagrange multipliers
- const PylithScalar Spp = 1.0 / jacobianN[0] + 1.0
- / jacobianP[0];
-
- const PylithScalar dlp = dLagrangeTpdt[0];
- (*slip)[0] -= Spp * dlp;
-
- PetscLogFlops(2);
-} // _sensitivitySolveLumped1D
-
-// ----------------------------------------------------------------------
-// Solve slip/Lagrange multiplier sensitivity problem for case of lumped Jacobian in 2-D.
-void
-pylith::faults::FaultCohesiveDyn::_sensitivitySolveLumped2D(
- scalar_array* slip,
- const scalar_array& dLagrangeTpdt,
- const scalar_array& jacobianN,
- const scalar_array& jacobianP)
-{ // _sensitivitySolveLumped2D
- assert(0 != slip);
-
- // Sensitivity of slip to changes in the Lagrange multipliers
- // Spp = Sqq = 1.0/Aix + 1.0/Ajx
- assert(jacobianN[0] > 0.0);
- assert(jacobianP[0] > 0.0);
-
- // Check to make sure Jacobian is same at all DOF for
- // vertices i and j (means S is diagonal with equal enties).
- assert(jacobianN[0] == jacobianN[1]);
- assert(jacobianP[0] == jacobianP[1]);
-
- const PylithScalar Spp = 1.0 / jacobianN[0] + 1.0 / jacobianP[0];
- const PylithScalar Sqq = Spp;
-
- const PylithScalar dlp = dLagrangeTpdt[0];
- const PylithScalar dlq = dLagrangeTpdt[1];
- (*slip)[0] -= Spp * dlp;
- (*slip)[1] -= Sqq * dlq;
-
- PetscLogFlops(7);
-} // _sensitivitySolveLumped2D
-
-// ----------------------------------------------------------------------
-// Solve slip/Lagrange multiplier sensitivity problem for case of lumped Jacobian in 3-D.
-void
-pylith::faults::FaultCohesiveDyn::_sensitivitySolveLumped3D(
- scalar_array* slip,
- const scalar_array& dLagrangeTpdt,
- const scalar_array& jacobianN,
- const scalar_array& jacobianP)
-{ // _sensitivitySolveLumped3D
- assert(0 != slip);
-
- // Sensitivity of slip to changes in the Lagrange multipliers
- // Aixjx = 1.0/Aix + 1.0/Ajx
- assert(jacobianN[0] > 0.0);
- assert(jacobianP[0] > 0.0);
-
- // Check to make sure Jacobian is same at all DOF for
- // vertices i and j (means S is diagonal with equal enties).
- assert(jacobianN[0] == jacobianN[1] &&
- jacobianN[0] == jacobianN[2]);
- assert(jacobianP[0] == jacobianP[1] &&
- jacobianP[0] == jacobianP[2]);
-
-
- const PylithScalar Spp = 1.0 / jacobianN[0] + 1.0 / jacobianP[0];
- const PylithScalar Sqq = Spp;
- const PylithScalar Srr = Spp;
-
- const PylithScalar dlp = dLagrangeTpdt[0];
- const PylithScalar dlq = dLagrangeTpdt[1];
- const PylithScalar dlr = dLagrangeTpdt[2];
- (*slip)[0] -= Spp * dlp;
- (*slip)[1] -= Sqq * dlq;
- (*slip)[2] -= Srr * dlr;
-
- PetscLogFlops(9);
-} // _sensitivitySolveLumped3D
-
-// ----------------------------------------------------------------------
// Constrain solution space in 1-D.
void
pylith::faults::FaultCohesiveDyn::_constrainSolnSpace1D(scalar_array* dLagrangeTpdt,
const scalar_array& slip,
const scalar_array& sliprate,
- const scalar_array& tractionTpdt)
+ const scalar_array& tractionTpdt,
+ const bool iterating)
{ // _constrainSolnSpace1D
assert(0 != dLagrangeTpdt);
@@ -2036,7 +1952,8 @@
pylith::faults::FaultCohesiveDyn::_constrainSolnSpace2D(scalar_array* dLagrangeTpdt,
const scalar_array& slip,
const scalar_array& slipRate,
- const scalar_array& tractionTpdt)
+ const scalar_array& tractionTpdt,
+ const bool iterating)
{ // _constrainSolnSpace2D
assert(0 != dLagrangeTpdt);
@@ -2052,11 +1969,7 @@
// if in compression and no opening
const PylithScalar frictionStress = _friction->calcFriction(slipMag, slipRateMag,
tractionNormal);
-#if 0
- if (tractionShearMag > frictionStress || slipRateMag > 0.0) {
-#else
- if (tractionShearMag > frictionStress) {
-#endif
+ if (tractionShearMag > frictionStress || (iterating && slipRateMag > 0.0)) {
// traction is limited by friction, so have sliding OR
// friction exceeds traction due to overshoot in slip
@@ -2093,7 +2006,8 @@
pylith::faults::FaultCohesiveDyn::_constrainSolnSpace3D(scalar_array* dLagrangeTpdt,
const scalar_array& slip,
const scalar_array& slipRate,
- const scalar_array& tractionTpdt)
+ const scalar_array& tractionTpdt,
+ const bool iterating)
{ // _constrainSolnSpace3D
assert(0 != dLagrangeTpdt);
@@ -2113,11 +2027,7 @@
// if in compression and no opening
const PylithScalar frictionStress =
_friction->calcFriction(slipShearMag, slipRateMag, tractionNormal);
-#if 0
- if (tractionShearMag > frictionStress || slipRateMag > 0.0) {
-#else
- if (tractionShearMag > frictionStress) {
-#endif
+ if (tractionShearMag > frictionStress || (iterating && slipRateMag > 0.0)) {
// traction is limited by friction, so have sliding OR
// friction exceeds traction due to overshoot in slip
Modified: short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.hh
===================================================================
--- short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.hh 2011-10-31 16:57:52 UTC (rev 19134)
+++ short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/faults/FaultCohesiveDyn.hh 2011-10-31 17:08:17 UTC (rev 19135)
@@ -199,42 +199,6 @@
*/
void _sensitivityUpdateSoln(const bool negativeSide);
- /** Solve slip/Lagrange multiplier sensitivity problem for case of lumped Jacobian in 1-D.
- *
- * @param slip Adjustment to slip assoc. w/Lagrange multiplier vertex.
- * @param dLagrangeTpdt Change in Lagrange multiplier.
- * @param jacobianN Jacobian for vertex on - side of the fault.
- * @param jacobianP Jacobian for vertex on + side of the fault.
- */
- void _sensitivitySolveLumped1D(scalar_array* slip,
- const scalar_array& dLagrangeTpdt,
- const scalar_array& jacobianN,
- const scalar_array& jacobianP);
-
- /** Solve slip/Lagrange multiplier sensitivity problem for case of lumped Jacobian in 2-D.
- *
- * @param slip Adjustment to slip assoc. w/Lagrange multiplier vertex.
- * @param dLagrangeTpdt Change in Lagrange multiplier.
- * @param jacobianN Jacobian for vertex on - side of the fault.
- * @param jacobianP Jacobian for vertex on + side of the fault.
- */
- void _sensitivitySolveLumped2D(scalar_array* slip,
- const scalar_array& dLagrangeTpdt,
- const scalar_array& jacobianN,
- const scalar_array& jacobianP);
-
- /** Solve slip/Lagrange multiplier sensitivity problem for case of lumped Jacobian in 3-D.
- *
- * @param slip Adjustment to slip assoc. w/Lagrange multiplier vertex.
- * @param dLagrangeTpdt Change in Lagrange multiplier.
- * @param jacobianN Jacobian for vertex on - side of the fault.
- * @param jacobianP Jacobian for vertex on + side of the fault.
- */
- void _sensitivitySolveLumped3D(scalar_array* slip,
- const scalar_array& dLagrangeTpdt,
- const scalar_array& jacobianN,
- const scalar_array& jacobianP);
-
/** Constrain solution space with lumped Jacobian in 1-D.
*
* @param dLagrangeTpdt Adjustment to Lagrange multiplier.
@@ -243,9 +207,10 @@
* @param tractionTpdt Fault traction assoc. w/Lagrange multiplier vertex.
*/
void _constrainSolnSpace1D(scalar_array* dLagrangeTpdt,
- const scalar_array& slip,
- const scalar_array& slipRate,
- const scalar_array& tractionTpdt);
+ const scalar_array& slip,
+ const scalar_array& slipRate,
+ const scalar_array& tractionTpdt,
+ const bool iterating =true);
/** Constrain solution space with lumped Jacobian in 2-D.
*
@@ -255,9 +220,10 @@
* @param tractionTpdt Fault traction assoc. w/Lagrange multiplier vertex.
*/
void _constrainSolnSpace2D(scalar_array* dLagrangeTpdt,
- const scalar_array& slip,
- const scalar_array& slipRate,
- const scalar_array& tractionTpdt);
+ const scalar_array& slip,
+ const scalar_array& slipRate,
+ const scalar_array& tractionTpdt,
+ const bool iterating =true);
/** Constrain solution space with lumped Jacobian in 3-D.
*
@@ -267,9 +233,10 @@
* @param tractionTpdt Fault traction assoc. w/Lagrange multiplier vertex.
*/
void _constrainSolnSpace3D(scalar_array* dLagrangeTpdt,
- const scalar_array& slip,
- const scalar_array& slipRate,
- const scalar_array& tractionTpdt);
+ const scalar_array& slip,
+ const scalar_array& slipRate,
+ const scalar_array& tractionTpdt,
+ const bool iterating =true);
// PRIVATE MEMBERS ////////////////////////////////////////////////////
private :
More information about the CIG-COMMITS
mailing list