[cig-commits] r15952 - in short/3D/PyLith/branches/pylith-friction/libsrc: . feassemble
brad at geodynamics.org
brad at geodynamics.org
Tue Nov 10 19:29:33 PST 2009
Author: brad
Date: 2009-11-10 19:29:33 -0800 (Tue, 10 Nov 2009)
New Revision: 15952
Modified:
short/3D/PyLith/branches/pylith-friction/libsrc/Makefile.am
short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.cc
short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.hh
short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.cc
short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.hh
short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Makefile.am
Log:
Updated integrating residual for large deformations (not tested).
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/Makefile.am 2009-11-11 02:14:10 UTC (rev 15951)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/Makefile.am 2009-11-11 03:29:33 UTC (rev 15952)
@@ -71,9 +71,10 @@
feassemble/Quadrature2Din3D.cc \
feassemble/Quadrature3D.cc \
feassemble/IntegratorElasticity.cc \
- feassemble/IntegratorElasticityLgDeform.cc \
feassemble/ElasticityImplicit.cc \
feassemble/ElasticityExplicit.cc \
+ feassemble/IntegratorElasticityLgDeform.cc \
+ feassemble/ElasticityImplicitLgDeform.cc \
materials/Metadata.cc \
materials/Material.cc \
materials/ElasticMaterial.cc \
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.cc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.cc 2009-11-11 02:14:10 UTC (rev 15951)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.cc 2009-11-11 03:29:33 UTC (rev 15952)
@@ -12,7 +12,7 @@
#include <portinfo>
-#include "ElasticityImplicit.hh" // implementation of class methods
+#include "ElasticityImplicitLgDeform.hh" // implementation of class methods
#include "Quadrature.hh" // USES Quadrature
#include "CellGeometry.hh" // USES CellGeometry
@@ -44,14 +44,14 @@
// ----------------------------------------------------------------------
// Constructor
-pylith::feassemble::ElasticityImplicit::ElasticityImplicit(void) :
+pylith::feassemble::ElasticityImplicitLgDeform::ElasticityImplicitLgDeform(void) :
_dtm1(-1.0)
{ // constructor
} // constructor
// ----------------------------------------------------------------------
// Destructor
-pylith::feassemble::ElasticityImplicit::~ElasticityImplicit(void)
+pylith::feassemble::ElasticityImplicitLgDeform::~ElasticityImplicitLgDeform(void)
{ // destructor
deallocate();
} // destructor
@@ -59,15 +59,15 @@
// ----------------------------------------------------------------------
// Deallocate PETSc and local data structures.
void
-pylith::feassemble::ElasticityImplicit::deallocate(void)
+pylith::feassemble::ElasticityImplicitLgDeform::deallocate(void)
{ // deallocate
- IntegratorElasticity::deallocate();
+ IntegratorElasticityLgDeform::deallocate();
} // deallocate
// ----------------------------------------------------------------------
// Set time step for advancing from time t to time t+dt.
void
-pylith::feassemble::ElasticityImplicit::timeStep(const double dt)
+pylith::feassemble::ElasticityImplicitLgDeform::timeStep(const double dt)
{ // timeStep
if (_dt != -1.0)
_dtm1 = _dt;
@@ -81,7 +81,7 @@
// ----------------------------------------------------------------------
// Get stable time step for advancing from time t to time t+dt.
double
-pylith::feassemble::ElasticityImplicit::stableTimeStep(const topology::Mesh& mesh) const
+pylith::feassemble::ElasticityImplicitLgDeform::stableTimeStep(const topology::Mesh& mesh) const
{ // stableTimeStep
assert(0 != _material);
return _material->stableTimeStepImplicit(mesh);
@@ -91,7 +91,7 @@
// Set flag for setting constraints for total field solution or
// incremental field solution.
void
-pylith::feassemble::ElasticityImplicit::useSolnIncr(const bool flag)
+pylith::feassemble::ElasticityImplicitLgDeform::useSolnIncr(const bool flag)
{ // useSolnIncr
assert(0 != _material);
_useSolnIncr = flag;
@@ -101,14 +101,14 @@
// ----------------------------------------------------------------------
// Integrate constributions to residual term (r) for operator.
void
-pylith::feassemble::ElasticityImplicit::integrateResidual(
+pylith::feassemble::ElasticityImplicitLgDeform::integrateResidual(
const topology::Field<topology::Mesh>& residual,
const double t,
topology::SolutionFields* const fields)
{ // integrateResidual
/// Member prototype for _elasticityResidualXD()
- typedef void (pylith::feassemble::ElasticityImplicit::*elasticityResidual_fn_type)
- (const double_array&);
+ typedef void (pylith::feassemble::ElasticityImplicitLgDeform::*elasticityResidual_fn_type)
+ (const double_array&, const double_array&);
assert(0 != _quadrature);
assert(0 != _material);
@@ -127,6 +127,7 @@
// Get cell geometry information that doesn't depend on cell
const int numQuadPts = _quadrature->numQuadPts();
+ const double_array& vertices = _quadrature->refGeometry().vertices();
const double_array& quadWts = _quadrature->quadWts();
assert(quadWts.size() == numQuadPts);
const int numBasis = _quadrature->numBasis();
@@ -143,19 +144,19 @@
elasticityResidual_fn_type elasticityResidualFn;
if (1 == cellDim) {
elasticityResidualFn =
- &pylith::feassemble::ElasticityImplicit::_elasticityResidual1D;
+ &pylith::feassemble::ElasticityImplicitLgDeform::_elasticityResidual1D;
calcTotalStrainFn =
- &pylith::feassemble::IntegratorElasticity::_calcTotalStrain1D;
+ &pylith::feassemble::IntegratorElasticityLgDeform::_calcTotalStrain1D;
} else if (2 == cellDim) {
elasticityResidualFn =
- &pylith::feassemble::ElasticityImplicit::_elasticityResidual2D;
+ &pylith::feassemble::ElasticityImplicitLgDeform::_elasticityResidual2D;
calcTotalStrainFn =
- &pylith::feassemble::IntegratorElasticity::_calcTotalStrain2D;
+ &pylith::feassemble::IntegratorElasticityLgDeform::_calcTotalStrain2D;
} else if (3 == cellDim) {
elasticityResidualFn =
- &pylith::feassemble::ElasticityImplicit::_elasticityResidual3D;
+ &pylith::feassemble::ElasticityImplicitLgDeform::_elasticityResidual3D;
calcTotalStrainFn =
- &pylith::feassemble::IntegratorElasticity::_calcTotalStrain3D;
+ &pylith::feassemble::IntegratorElasticityLgDeform::_calcTotalStrain3D;
} else
assert(0);
@@ -163,8 +164,8 @@
double_array dispTCell(numBasis*spaceDim);
double_array dispTIncrCell(numBasis*spaceDim);
double_array dispTpdtCell(numBasis*spaceDim);
+ double_array deformCell(numQuadPts*spaceDim*spaceDim);
double_array strainCell(numQuadPts*tensorSize);
- strainCell = 0.0;
double_array gravVec(spaceDim);
double_array quadPtsGlobal(numQuadPts*spaceDim);
@@ -287,15 +288,17 @@
_logger->eventEnd(computeEvent);
} // if
- // Compute B(transpose) * sigma, first computing strains
+ // Compute B(transpose) * sigma, first computing deformation
+ // tensor and strains
_logger->eventBegin(stressEvent);
- calcTotalStrainFn(&strainCell, basisDeriv, dispTpdtCell,
- numBasis, numQuadPts);
+ _calcDeformation(&deformCell, basisDeriv, vertices, dispTpdtCell,
+ numBasis, numQuadPts, spaceDim);
+ calcTotalStrainFn(&strainCell, deformCell, numQuadPts);
const double_array& stressCell = _material->calcStress(strainCell, true);
_logger->eventEnd(stressEvent);
_logger->eventBegin(computeEvent);
- CALL_MEMBER_FN(*this, elasticityResidualFn)(stressCell);
+ CALL_MEMBER_FN(*this, elasticityResidualFn)(stressCell, dispTpdtCell);
_logger->eventEnd(computeEvent);
#if 0 // DEBUGGING
@@ -315,14 +318,14 @@
// ----------------------------------------------------------------------
// Compute stiffness matrix.
void
-pylith::feassemble::ElasticityImplicit::integrateJacobian(
+pylith::feassemble::ElasticityImplicitLgDeform::integrateJacobian(
topology::Jacobian* jacobian,
const double t,
topology::SolutionFields* fields)
{ // integrateJacobian
/// Member prototype for _elasticityJacobianXD()
- typedef void (pylith::feassemble::ElasticityImplicit::*elasticityJacobian_fn_type)
- (const double_array&);
+ typedef void (pylith::feassemble::ElasticityImplicitLgDeform::*elasticityJacobian_fn_type)
+ (const double_array&, const double_array&, const double_array&);
assert(0 != _quadrature);
assert(0 != _material);
@@ -341,6 +344,7 @@
// Get cell geometry information that doesn't depend on cell
const int numQuadPts = _quadrature->numQuadPts();
+ const double_array& vertices = _quadrature->refGeometry().vertices();
const double_array& quadWts = _quadrature->quadWts();
assert(quadWts.size() == numQuadPts);
const int numBasis = _quadrature->numBasis();
@@ -357,19 +361,19 @@
elasticityJacobian_fn_type elasticityJacobianFn;
if (1 == cellDim) {
elasticityJacobianFn =
- &pylith::feassemble::ElasticityImplicit::_elasticityJacobian1D;
+ &pylith::feassemble::ElasticityImplicitLgDeform::_elasticityJacobian1D;
calcTotalStrainFn =
- &pylith::feassemble::IntegratorElasticity::_calcTotalStrain1D;
+ &pylith::feassemble::IntegratorElasticityLgDeform::_calcTotalStrain1D;
} else if (2 == cellDim) {
elasticityJacobianFn =
- &pylith::feassemble::ElasticityImplicit::_elasticityJacobian2D;
+ &pylith::feassemble::ElasticityImplicitLgDeform::_elasticityJacobian2D;
calcTotalStrainFn =
- &pylith::feassemble::IntegratorElasticity::_calcTotalStrain2D;
+ &pylith::feassemble::IntegratorElasticityLgDeform::_calcTotalStrain2D;
} else if (3 == cellDim) {
elasticityJacobianFn =
- &pylith::feassemble::ElasticityImplicit::_elasticityJacobian3D;
+ &pylith::feassemble::ElasticityImplicitLgDeform::_elasticityJacobian3D;
calcTotalStrainFn =
- &pylith::feassemble::IntegratorElasticity::_calcTotalStrain3D;
+ &pylith::feassemble::IntegratorElasticityLgDeform::_calcTotalStrain3D;
} else
assert(0);
@@ -377,8 +381,9 @@
double_array dispTCell(numBasis*spaceDim);
double_array dispTIncrCell(numBasis*spaceDim);
double_array dispTpdtCell(numBasis*spaceDim);
+ double_array deformCell(numQuadPts*spaceDim*spaceDim);
double_array strainCell(numQuadPts*tensorSize);
- strainCell = 0.0;
+ double_array stressCell(numQuadPts*tensorSize);
// Get cell information
const ALE::Obj<SieveMesh>& sieveMesh = fields->mesh().sieveMesh();
@@ -473,15 +478,20 @@
dispTpdtCell = dispTCell + dispTIncrCell;
_logger->eventBegin(computeEvent);
- // Compute strains
- calcTotalStrainFn(&strainCell, basisDeriv, dispTpdtCell,
- numBasis, numQuadPts);
-
+ // Compute deformation tensor, strains, and stresses
+ _calcDeformation(&deformCell, basisDeriv, vertices, dispTpdtCell,
+ numBasis, numQuadPts, spaceDim);
+ calcTotalStrainFn(&strainCell, deformCell, numQuadPts);
+
// Get "elasticity" matrix at quadrature points for this cell
const double_array& elasticConsts =
_material->calcDerivElastic(strainCell);
- CALL_MEMBER_FN(*this, elasticityJacobianFn)(elasticConsts);
+ // Get Second Priola-Kirchoff stress tensor
+ const double_array& stressCell = _material->calcStress(strainCell, true);
+
+ CALL_MEMBER_FN(*this, elasticityJacobianFn)(elasticConsts, stressCell,
+ dispTpdtCell);
_logger->eventEnd(computeEvent);
if (_quadrature->checkConditioning()) {
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.hh 2009-11-11 02:14:10 UTC (rev 15951)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/ElasticityImplicitLgDeform.hh 2009-11-11 03:29:33 UTC (rev 15952)
@@ -11,10 +11,10 @@
//
/**
- * @file libsrc/feassemble/ElasticityImplicit.hh
+ * @file libsrc/feassemble/ElasticityImplicitLgDeform.hh
*
* @brief Implicit time integration of quasi-static elasticity equation
- * using finite-elements.
+ * with large rigid body motion and small strains using finite-elements.
*
* Note: This object operates on a single finite-element family, which
* is defined by the quadrature and a database of material property
@@ -44,25 +44,26 @@
* information.
*/
-#if !defined(pylith_feassemble_elasticityimplicit_hh)
-#define pylith_feassemble_elasticityimplicit_hh
+#if !defined(pylith_feassemble_elasticityimplicitlgdeform_hh)
+#define pylith_feassemble_elasticityimplicitlgdeform_hh
// Include directives ---------------------------------------------------
-#include "IntegratorElasticity.hh" // ISA IntegratorElasticity
+#include "IntegratorElasticityLgDeform.hh" // ISA IntegratorElasticityLgDeform
-// ElasticityImplicit ---------------------------------------------------
-class pylith::feassemble::ElasticityImplicit : public IntegratorElasticity
-{ // ElasticityImplicit
- friend class TestElasticityImplicit; // unit testing
+// ElasticityImplicitLgDeform -------------------------------------------
+class pylith::feassemble::ElasticityImplicitLgDeform :
+ public IntegratorElasticityLgDeform
+{ // ElasticityImplicitLgDeform
+ friend class TestElasticityImplicitLgDeform; // unit testing
// PUBLIC MEMBERS ///////////////////////////////////////////////////////
public :
/// Constructor
- ElasticityImplicit(void);
+ ElasticityImplicitLgDeform(void);
/// Destructor
- ~ElasticityImplicit(void);
+ ~ElasticityImplicitLgDeform(void);
/// Deallocate PETSc and local data structures.
void deallocate(void);
@@ -122,19 +123,19 @@
private :
/// Not implemented.
- ElasticityImplicit(const ElasticityImplicit&);
+ ElasticityImplicitLgDeform(const ElasticityImplicitLgDeform&);
/// Not implemented
- const ElasticityImplicit& operator=(const ElasticityImplicit&);
+ const ElasticityImplicitLgDeform& operator=(const ElasticityImplicitLgDeform&);
// PRIVATE MEMBERS //////////////////////////////////////////////////////
private :
double _dtm1; ///< Time step for t-dt1 -> t
-}; // ElasticityImplicit
+}; // ElasticityImplicitLgDeform
-#endif // pylith_feassemble_elasticityimplicit_hh
+#endif // pylith_feassemble_elasticityimplicitlgdeform_hh
// End of file
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.cc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.cc 2009-11-11 02:14:10 UTC (rev 15951)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.cc 2009-11-11 03:29:33 UTC (rev 15952)
@@ -199,9 +199,7 @@
double_array dispCell(numBasis*spaceDim);
double_array deformCell(numQuadPts*spaceDim*spaceDim);
double_array strainCell(numQuadPts*tensorSize);
- strainCell = 0.0;
double_array stressCell(numQuadPts*tensorSize);
- stressCell = 0.0;
// Get cell information
const ALE::Obj<SieveMesh>& sieveMesh = field->mesh().sieveMesh();
@@ -318,7 +316,8 @@
// Integrate elasticity term in residual for 1-D cells.
void
pylith::feassemble::IntegratorElasticityLgDeform::_elasticityResidual1D(
- const double_array& stress)
+ const double_array& stress,
+ const double_array& disp)
{ // _elasticityResidual1D
const int numQuadPts = _quadrature->numQuadPts();
const int numBasis = _quadrature->numBasis();
@@ -335,18 +334,22 @@
const double wt = quadWts[iQuad] * jacobianDet[iQuad];
const double s11 = stress[iQuad];
for (int iBasis=0; iBasis < numBasis; ++iBasis) {
- const double N1 = wt*basisDeriv[iQuad*numBasis+iBasis ];
+ double l11 = 0.0;
+ for (int kBasis=0; kBasis < numBasis; ++kBasis)
+ l11 += basisDeriv[iQuad*numBasis+kBasis ] * disp[kBasis ];
+ const double N1 = wt * (1.0 + l11) * basisDeriv[iQuad*numBasis+iBasis ];
_cellVector[iBasis*spaceDim ] -= N1*s11;
} // for
} // for
- PetscLogFlops(numQuadPts*(1+numBasis*5));
+ PetscLogFlops(numQuadPts*(1+numBasis*5+numBasis*numBasis*2));
} // _elasticityResidual1D
// ----------------------------------------------------------------------
// Integrate elasticity term in residual for 2-D cells.
void
pylith::feassemble::IntegratorElasticityLgDeform::_elasticityResidual2D(
- const double_array& stress)
+ const double_array& stress,
+ const double_array& disp)
{ // _elasticityResidual2D
const int numQuadPts = _quadrature->numQuadPts();
const int numBasis = _quadrature->numBasis();
@@ -368,20 +371,35 @@
for (int iBasis=0, iQ=iQuad*numBasis*spaceDim;
iBasis < numBasis;
++iBasis) {
- const double N1 = wt*basisDeriv[iQ+iBasis*spaceDim ];
- const double N2 = wt*basisDeriv[iQ+iBasis*spaceDim+1];
- _cellVector[iBasis*spaceDim ] -= N1*s11 + N2*s12;
- _cellVector[iBasis*spaceDim+1] -= N1*s12 + N2*s22;
+ double l11 = 0.0;
+ double l12 = 0.0;
+ double l21 = 0.0;
+ double l22 = 0.0;
+ for (int kBasis=0; kBasis < numBasis; ++kBasis) {
+ const int kB = kBasis*spaceDim;
+ l11 += basisDeriv[iQ+kB ] * disp[kB ];
+ l12 += basisDeriv[iQ+kB+1] * disp[kB ];
+ l21 += basisDeriv[iQ+kB ] * disp[kB+1];
+ l22 += basisDeriv[iQ+kB+1] * disp[kB+1];
+ } // for
+ const int iB = iBasis*spaceDim;
+ const double N1 = basisDeriv[iQ+iB ];
+ const double N2 = basisDeriv[iQ+iB+1];
+ _cellVector[iB ] -=
+ wt*((1.0+l11)*N1*s11 + N2*l12*s22 + ((1.0+l11)*N2 + l12*N1)*s12);
+ _cellVector[iB+1] -=
+ wt*(l21*N1*s11 + (1.0+l22)*N2*s22 + ((1.0+l22)*N1 + l21*N2)*s12);
} // for
} // for
- PetscLogFlops(numQuadPts*(1+numBasis*(8+2+9)));
+ PetscLogFlops(numQuadPts*(1+numBasis*(numBasis*8+14)));
} // _elasticityResidual2D
// ----------------------------------------------------------------------
// Integrate elasticity term in residual for 3-D cells.
void
pylith::feassemble::IntegratorElasticityLgDeform::_elasticityResidual3D(
- const double_array& stress)
+ const double_array& stress,
+ const double_array& disp)
{ // _elasticityResidual3D
const int numQuadPts = _quadrature->numQuadPts();
const int numBasis = _quadrature->numBasis();
@@ -407,24 +425,59 @@
for (int iBasis=0, iQ=iQuad*numBasis*spaceDim;
iBasis < numBasis;
++iBasis) {
- const int iBlock = iBasis*spaceDim;
- const double N1 = wt*basisDeriv[iQ+iBasis*spaceDim+0];
- const double N2 = wt*basisDeriv[iQ+iBasis*spaceDim+1];
- const double N3 = wt*basisDeriv[iQ+iBasis*spaceDim+2];
+ double l11 = 0.0;
+ double l12 = 0.0;
+ double l13 = 0.0;
+ double l21 = 0.0;
+ double l22 = 0.0;
+ double l23 = 0.0;
+ double l31 = 0.0;
+ double l32 = 0.0;
+ double l33 = 0.0;
+ for (int kBasis=0; kBasis < numBasis; ++kBasis) {
+ const int kB = kBasis*spaceDim;
+ l11 += basisDeriv[iQ+kB ] * disp[kB ];
+ l12 += basisDeriv[iQ+kB+1] * disp[kB ];
+ l13 += basisDeriv[iQ+kB+2] * disp[kB ];
+ l21 += basisDeriv[iQ+kB ] * disp[kB+1];
+ l22 += basisDeriv[iQ+kB+1] * disp[kB+1];
+ l23 += basisDeriv[iQ+kB+2] * disp[kB+1];
+ l31 += basisDeriv[iQ+kB ] * disp[kB+2];
+ l32 += basisDeriv[iQ+kB+1] * disp[kB+2];
+ l33 += basisDeriv[iQ+kB+2] * disp[kB+2];
+ } // for
+ const int iB = iBasis*spaceDim;
+ const double N1 = basisDeriv[iQ+iB ];
+ const double N2 = basisDeriv[iQ+iB+1];
+ const double N3 = basisDeriv[iQ+iB+2];
- _cellVector[iBlock ] -= N1*s11 + N2*s12 + N3*s13;
- _cellVector[iBlock+1] -= N1*s12 + N2*s22 + N3*s23;
- _cellVector[iBlock+2] -= N1*s13 + N2*s23 + N3*s33;
+ _cellVector[iB ] -=
+ wt*((1.0+l11)*N1*s11 + l12*N2*s22 + l13*N3*s33 +
+ ((1.0+l11)*N2 + l12*N1)*s12 +
+ (l12*N3 + l13*N2)*s23 +
+ ((1.0+l11)*N3 + l13*N1)*s13);
+ _cellVector[iB+1] -=
+ wt*(l21*N1*s11 + (1.0+l22)*N2*s22 + l23*N3*s33 +
+ ((1.0+l22)*N1 + l21*N2)*s12 +
+ ((1.0+l22)*N3 + l23*N2)*s23 +
+ (l21*N3 + l23*N1)*s13);
+ _cellVector[iB+2] -=
+ wt*(l31*N1*s11 + l32*N2*s22 + (1.0+l33)*N3*s33 +
+ (l31*N2 + l32*N1)*s12 +
+ ((1.0+l33)*N2 + l32*N3)*s23 +
+ ((1.0+l33)*N1 + l31*N3)*s13);
} // for
} // for
- PetscLogFlops(numQuadPts*(1+numBasis*(3+12)));
+ PetscLogFlops(numQuadPts*(1+numBasis*(numBasis*18+3*27)));
} // _elasticityResidual3D
// ----------------------------------------------------------------------
// Integrate elasticity term in Jacobian for 1-D cells.
void
pylith::feassemble::IntegratorElasticityLgDeform::_elasticityJacobian1D(
- const double_array& elasticConsts)
+ const double_array& elasticConsts,
+ const double_array& stress,
+ const double_array& disp)
{ // _elasticityJacobian1D
const int numQuadPts = _quadrature->numQuadPts();
const int numBasis = _quadrature->numBasis();
@@ -457,7 +510,9 @@
// Integrate elasticity term in Jacobian for 2-D cells.
void
pylith::feassemble::IntegratorElasticityLgDeform::_elasticityJacobian2D(
- const double_array& elasticConsts)
+ const double_array& elasticConsts,
+ const double_array& stress,
+ const double_array& disp)
{ // _elasticityJacobian2D
const int numQuadPts = _quadrature->numQuadPts();
const int numBasis = _quadrature->numBasis();
@@ -521,7 +576,9 @@
// Integrate elasticity term in Jacobian for 3-D cells.
void
pylith::feassemble::IntegratorElasticityLgDeform::_elasticityJacobian3D(
- const double_array& elasticConsts)
+ const double_array& elasticConsts,
+ const double_array& stress,
+ const double_array& disp)
{ // _elasticityJacobian3D
const int numQuadPts = _quadrature->numQuadPts();
const int numBasis = _quadrature->numBasis();
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.hh 2009-11-11 02:14:10 UTC (rev 15951)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/IntegratorElasticityLgDeform.hh 2009-11-11 03:29:33 UTC (rev 15952)
@@ -91,38 +91,56 @@
/** Integrate elasticity term in residual for 1-D cells.
*
* @param stress Stress tensor for cell at quadrature points.
+ * @param disp Displacement field at cell's DOF.
*/
- void _elasticityResidual1D(const double_array& stress);
+ void _elasticityResidual1D(const double_array& stress,
+ const double_array& disp);
/** Integrate elasticity term in residual for 2-D cells.
*
* @param stress Stress tensor for cell at quadrature points.
+ * @param disp Displacement field at cell's DOF.
*/
- void _elasticityResidual2D(const double_array& stress);
+ void _elasticityResidual2D(const double_array& stress,
+ const double_array& disp);
/** Integrate elasticity term in residual for 3-D cells.
*
* @param stress Stress tensor for cell at quadrature points.
+ * @param disp Displacement field at cell's DOF.
*/
- void _elasticityResidual3D(const double_array& stress);
+ void _elasticityResidual3D(const double_array& stress,
+ const double_array& disp);
/** Integrate elasticity term in Jacobian for 1-D cells.
*
* @param elasticConsts Matrix of elasticity constants at quadrature points.
+ * @param stress Stress tensor for cell at quadrature points.
+ * @param disp Displacement field at cell's DOF.
*/
- void _elasticityJacobian1D(const double_array& elasticConsts);
+ void _elasticityJacobian1D(const double_array& elasticConsts,
+ const double_array& stress,
+ const double_array& disp);
/** Integrate elasticity term in Jacobian for 2-D cells.
*
* @param elasticConsts Matrix of elasticity constants at quadrature points.
+ * @param stress Stress tensor for cell at quadrature points.
+ * @param disp Displacement field at cell's DOF.
*/
- void _elasticityJacobian2D(const double_array& elasticConsts);
+ void _elasticityJacobian2D(const double_array& elasticConsts,
+ const double_array& stress,
+ const double_array& disp);
/** Integrate elasticity term in Jacobian for 3-D cells.
*
* @param elasticConsts Matrix of elasticity constants at quadrature points.
+ * @param stress Stress tensor for cell at quadrature points.
+ * @param disp Displacement field at cell's DOF.
*/
- void _elasticityJacobian3D(const double_array& elasticConsts);
+ void _elasticityJacobian3D(const double_array& elasticConsts,
+ const double_array& stress,
+ const double_array& disp);
/** Calculate Green-Lagrange strain tensor at quadrature points of a
* 1-D cell.
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Makefile.am 2009-11-11 02:14:10 UTC (rev 15951)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Makefile.am 2009-11-11 03:29:33 UTC (rev 15952)
@@ -20,6 +20,7 @@
Constraint.icc \
ElasticityExplicit.hh \
ElasticityImplicit.hh \
+ ElasticityImplicitLgDeform.hh \
Integrator.hh \
Integrator.icc \
Integrator.cc \
More information about the CIG-COMMITS
mailing list