[cig-commits] r6480 - in short/3D/PyLith/trunk: . libsrc/feassemble
libsrc/materials pylith/problems unittests/libtests/materials
unittests/libtests/materials/data
brad at geodynamics.org
brad at geodynamics.org
Fri Mar 30 18:41:17 PDT 2007
Author: brad
Date: 2007-03-30 18:41:15 -0700 (Fri, 30 Mar 2007)
New Revision: 6480
Removed:
short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc
Modified:
short/3D/PyLith/trunk/TODO
short/3D/PyLith/trunk/libsrc/feassemble/ExplicitElasticity.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.hh
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/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/pylith/problems/Explicit.py
short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc
short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.hh
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/data/ElasticIsotropic3D.py
short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.cc
short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.hh
short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialApp.py
short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.cc
short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.hh
short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.cc
short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.hh
Log:
Reworked material interface for compatibility with linear and nonlinear materials.
Modified: short/3D/PyLith/trunk/TODO
===================================================================
--- short/3D/PyLith/trunk/TODO 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/TODO 2007-03-31 01:41:15 UTC (rev 6480)
@@ -5,6 +5,8 @@
Error checking
add isNull() assertions before using ALE::Obj.
+ add check to material::initialize material dimension must match cell dimension
+
0. Add unit tests for ElasticIsotropic1D, ElasticPlaneStrain,
ElasticPlaneStress.
Modified: short/3D/PyLith/trunk/libsrc/feassemble/ExplicitElasticity.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/ExplicitElasticity.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/feassemble/ExplicitElasticity.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -107,15 +107,10 @@
const double* basisDeriv = _quadrature->basisDeriv();
const double* jacobianDet = _quadrature->jacobianDet();
- // Get material physical properties at quadrature points for this cell
- _material->calcProperties(*cellIter, numQuadPts);
- const double* density = _material->density();
- const double* elasticConsts = _material->elasticConsts();
- const int numElasticConsts = _material->numElasticConsts();
-
// Compute action for cell
// Compute action for inertial terms
+ const double* density = _material->calcDensity(*cellIter, numQuadPts);
for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
const double wt =
quadWts[iQuad] * jacobianDet[iQuad] * density[iQuad] / dt2;
@@ -133,7 +128,7 @@
} // for
} // for
PetscErrorCode err =
- PetscLogFlops(numQuadPts*(3+numBasis*(1+numBasis*(1+4*spaceDim))));
+ PetscLogFlops(numQuadPts*(3+numBasis*(1+numBasis*(5*spaceDim))));
if (err)
throw std::runtime_error("Logging PETSc flops failed.");
@@ -149,141 +144,124 @@
// Compute action for elastic terms
if (1 == cellDim) {
- assert(1 == numElasticConsts);
+ // Compute total strains
+ const int stressSize = _material->stressSize();
+ assert(numQuadPts == stressSize);
+ const int strainSize = stressSize;;
+ double* totalStrain = (strainSize > 0) ? new double[strainSize] : 0;
+ memset(totalStrain, 0, strainSize*sizeof(double));
for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis)
+ totalStrain[iQuad] += basisDeriv[iQ+iBasis] * dispTCell[iBasis];
+ } // for
+ const double* stress =
+ _material->calcStress(*cellIter, totalStrain, numQuadPts,
+ spaceDim);
+
+ // Compute elastic action
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
const double wt = quadWts[iQuad] * jacobianDet[iQuad];
- const double C1111 = elasticConsts[iQuad];
+ const int iStress = iQuad*stressSize;
+ const double s11 = stress[iStress ];
for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
const int iBlock = iBasis * spaceDim;
- const double valI = wt*basisDeriv[iQ+iBasis]*C1111;
- for (int jBasis=0; jBasis < numBasis; ++jBasis) {
- const int jBlock = jBasis * spaceDim;
- const double valIJ = valI * basisDeriv[iQ+jBasis];
- _cellVector[iBlock] -= valIJ * dispTCell[jBlock];
- } // for
+ const double N1 = wt*basisDeriv[iQ+iBasis*cellDim ];
+ _cellVector[iBlock ] -= N1*s11;
} // for
- } // for
- PetscErrorCode err =
- PetscLogFlops(numQuadPts*(1+numBasis*(2+numBasis*3)));
+ } // for
+ delete[] totalStrain; totalStrain = 0;
+ PetscErrorCode err = PetscLogFlops(numQuadPts*(1+numBasis*5));
if (err)
throw std::runtime_error("Logging PETSc flops failed.");
} else if (2 == cellDim) {
- assert(6 == numElasticConsts);
+ // Compute total strains
+ const int stressSize = _material->stressSize();
+ assert(3*numQuadPts == stressSize);
+ const int strainSize = stressSize;;
+ double* totalStrain = (strainSize > 0) ? new double[strainSize] : 0;
+ memset(totalStrain, 0, strainSize*sizeof(double));
for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
+ totalStrain[iQuad ] +=
+ basisDeriv[iQ+iBasis ] * dispTCell[iBasis ];
+ totalStrain[iQuad+1] +=
+ basisDeriv[iQ+iBasis+1] * dispTCell[iBasis+1];
+ totalStrain[iQuad+2] +=
+ 0.5 * (basisDeriv[iQ+iBasis+1] * dispTCell[iBasis ] +
+ basisDeriv[iQ+iBasis ] * dispTCell[iBasis+1]);
+ } // for
+ } // for
+ const double* stress =
+ _material->calcStress(*cellIter, totalStrain, numQuadPts,
+ spaceDim);
+ // Compute elastic action
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
const double wt = quadWts[iQuad] * jacobianDet[iQuad];
- const int iConst = iQuad*numElasticConsts;
- const double C1111 = elasticConsts[iConst+0];
- const double C1122 = elasticConsts[iConst+1];
- const double C1112 = elasticConsts[iConst+2];
- const double C2222 = elasticConsts[iConst+3];
- const double C2212 = elasticConsts[iConst+4];
- const double C1212 = elasticConsts[iConst+5];
+ const int iStress = iQuad*stressSize;
+ const double s11 = stress[iStress ];
+ const double s22 = stress[iStress+1];
+ const double s12 = stress[iStress+2];
for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
const int iBlock = iBasis * spaceDim;
- const double Nip = wt*basisDeriv[iQ+iBasis*cellDim+0];
- const double Niq = wt*basisDeriv[iQ+iBasis*cellDim+1];
- for (int jBasis=0; jBasis < numBasis; ++jBasis) {
- const int jBlock = jBasis * spaceDim;
- const double Njp = basisDeriv[iQ+jBasis*cellDim+0];
- const double Njq = basisDeriv[iQ+jBasis*cellDim+1];
- const double ki0j0 =
- C1111 * Nip * Njp + C1112 * Niq * Njp +
- C1112 * Nip * Njq + C1212 * Niq * Njq;
- const double ki0j1 =
- C1122 * Nip * Njq + C2212 * Niq * Njq +
- C1112 * Nip * Njp + C1212 * Niq * Njp;
- const double ki1j1 =
- C2222 * Niq * Njq + C2212 * Nip * Njq +
- C2212 * Niq * Njp + C1212 * Nip * Njp;
- _cellVector[iBlock ] -=
- ki0j0 * dispTCell[jBlock ] + ki0j1 * dispTCell[jBlock+1];
- _cellVector[iBlock+1] -=
- ki0j1 * dispTCell[jBlock ] + ki1j1 * dispTCell[jBlock+1];
- } // for
+ const double N1 = wt*basisDeriv[iQ+iBasis*cellDim ];
+ const double N2 = wt*basisDeriv[iQ+iBasis*cellDim+1];
+ _cellVector[iBlock ] -= N1*s11 + N2*s12;
+ _cellVector[iBlock+1] -= N1*s12 + N2*s22;
} // for
} // for
- PetscErrorCode err =
- PetscLogFlops(numQuadPts*(1+numBasis*(2+numBasis*(2+3*11+2*4))));
+ err = PetscLogFlops(numQuadPts*(1+numBasis*(8+2+9)));
if (err)
throw std::runtime_error("Logging PETSc flops failed.");
} else if (3 == cellDim) {
- assert(21 == numElasticConsts);
+ // Compute total strains
+ const int stressSize = _material->stressSize();
+ assert(6*numQuadPts == stressSize);
+ const int strainSize = stressSize;;
+ double* totalStrain = (strainSize > 0) ? new double[strainSize] : 0;
+ memset(totalStrain, 0, strainSize*sizeof(double));
for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
+ totalStrain[iQuad ] +=
+ basisDeriv[iQ+iBasis ] * dispTCell[iBasis ];
+ totalStrain[iQuad+1] +=
+ basisDeriv[iQ+iBasis+1] * dispTCell[iBasis+1];
+ totalStrain[iQuad+2] +=
+ basisDeriv[iQ+iBasis+2] * dispTCell[iBasis+2];
+ totalStrain[iQuad+3] +=
+ 0.5 * (basisDeriv[iQ+iBasis+1] * dispTCell[iBasis ] +
+ basisDeriv[iQ+iBasis ] * dispTCell[iBasis+1]);
+ totalStrain[iQuad+4] +=
+ 0.5 * (basisDeriv[iQ+iBasis+2] * dispTCell[iBasis+1] +
+ basisDeriv[iQ+iBasis+1] * dispTCell[iBasis+2]);
+ totalStrain[iQuad+5] +=
+ 0.5 * (basisDeriv[iQ+iBasis+2] * dispTCell[iBasis ] +
+ basisDeriv[iQ+iBasis ] * dispTCell[iBasis+2]);
+ } // for
+ } // for
+ const double* stress =
+ _material->calcStress(*cellIter, totalStrain, numQuadPts, spaceDim);
+ // Compute elastic action
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
const double wt = quadWts[iQuad] * jacobianDet[iQuad];
- const int iConst = iQuad*numElasticConsts;
- const double C1111 = elasticConsts[iConst+ 0];
- const double C1122 = elasticConsts[iConst+ 1];
- const double C1133 = elasticConsts[iConst+ 2];
- const double C1112 = elasticConsts[iConst+ 3];
- const double C1123 = elasticConsts[iConst+ 4];
- const double C1113 = elasticConsts[iConst+ 5];
- const double C2222 = elasticConsts[iConst+ 6];
- const double C2233 = elasticConsts[iConst+ 7];
- const double C2212 = elasticConsts[iConst+ 8];
- const double C2223 = elasticConsts[iConst+ 9];
- const double C2213 = elasticConsts[iConst+10];
- const double C3333 = elasticConsts[iConst+11];
- const double C3312 = elasticConsts[iConst+12];
- const double C3323 = elasticConsts[iConst+13];
- const double C3313 = elasticConsts[iConst+14];
- const double C1212 = elasticConsts[iConst+15];
- const double C1223 = elasticConsts[iConst+16];
- const double C1213 = elasticConsts[iConst+17];
- const double C2323 = elasticConsts[iConst+18];
- const double C2313 = elasticConsts[iConst+19];
- const double C1313 = elasticConsts[iConst+20];
+ const int iStress = iQuad*stressSize;
+ const double s11 = stress[iStress ];
+ const double s22 = stress[iStress+1];
+ const double s33 = stress[iStress+2];
+ const double s12 = stress[iStress+3];
+ const double s23 = stress[iStress+4];
+ const double s13 = stress[iStress+5];
+
for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
const int iBlock = iBasis * spaceDim;
- const double Nip = wt*basisDeriv[iQ+iBasis*cellDim+0];
- const double Niq = wt*basisDeriv[iQ+iBasis*cellDim+1];
- const double Nir = wt*basisDeriv[iQ+iBasis*cellDim+2];
- for (int jBasis=0; jBasis < numBasis; ++jBasis) {
- const int jBlock = jBasis * spaceDim;
- const double Njp = basisDeriv[iQ+jBasis*cellDim+0];
- const double Njq = basisDeriv[iQ+jBasis*cellDim+1];
- const double Njr = basisDeriv[iQ+jBasis*cellDim+2];
- const double ki0j0 =
- C1111 * Nip * Njp + C1112 * Niq * Njp + C1113 * Nir * Njp +
- C1112 * Nip * Njq + C1212 * Niq * Njq + C1213 * Nir * Njq +
- C1113 * Nip * Njr + C1213 * Niq * Njr + C1313 * Nir * Njr;
- const double ki0j1 =
- C1122 * Nip * Njq + C2212 * Niq * Njq + C2213 * Nir * Njq +
- C1112 * Nip * Njp + C1212 * Niq * Njp + C1213 * Nir * Njp +
- C1123 * Nip * Njr + C1223 * Niq * Njr + C2313 * Nir * Njr;
- const double ki0j2 =
- C1133 * Nip * Njr + C3312 * Niq * Njr + C3313 * Nir * Njr +
- C1123 * Nip * Njq + C1223 * Niq * Njq + C2313 * Nir * Njq +
- C1113 * Nip * Njp + C1213 * Niq * Njp + C1313 * Nir * Njp;
- const double ki1j1 =
- C2222 * Niq * Njq + C2212 * Nip * Njq + C2223 * Nir * Njq +
- C2212 * Niq * Njp + C1212 * Nip * Njp + C1223 * Nir * Njp +
- C2223 * Niq * Njr + C1223 * Nip * Njr + C2323 * Nir * Njr;
- const double ki1j2 =
- C2233 * Niq * Njr + C3312 * Nip * Njr + C3323 * Nir * Njr +
- C2223 * Niq * Njq + C1223 * Nip * Njq + C2323 * Nir * Njq +
- C2213 * Niq * Njp + C1213 * Nip * Njp + C2313 * Nir * Njp;
- const double ki2j2 =
- C3333 * Nir * Njr + C3323 * Niq * Njr + C3313 * Nip * Njr +
- C3323 * Nir * Njq + C2323 * Niq * Njq + C2313 * Nip * Njq +
- C3313 * Nir * Njp + C2313 * Niq * Njp + C1313 * Nip * Njp;
-
- _cellVector[iBlock ] -=
- ki0j0 * dispTCell[jBlock ] +
- ki0j1 * dispTCell[jBlock+1] +
- ki0j2 * dispTCell[jBlock+2];
- _cellVector[iBlock+1] -=
- ki0j1 * dispTCell[jBlock ] +
- ki1j1 * dispTCell[jBlock+1] +
- ki1j2 * dispTCell[jBlock+2];
- _cellVector[iBlock+2] -=
- ki0j2 * dispTCell[jBlock ] +
- ki1j2 * dispTCell[jBlock+1] +
- ki2j2 * dispTCell[jBlock+2];
- } // for
+ const double N1 = wt*basisDeriv[iQ+iBasis*cellDim+0];
+ const double N2 = wt*basisDeriv[iQ+iBasis*cellDim+1];
+ const double N3 = wt*basisDeriv[iQ+iBasis*cellDim+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;
} // for
} // for
- PetscErrorCode err =
- PetscLogFlops(numQuadPts*(1+numBasis*(3+numBasis*(3+6*26+3*6))));
+ err = PetscLogFlops(numQuadPts*(1+numBasis*(3+12)));
if (err)
throw std::runtime_error("Logging PETSc flops failed.");
} // if/else
@@ -334,8 +312,7 @@
const double* jacobianDet = _quadrature->jacobianDet();
// Get material physical properties at quadrature points for this cell
- _material->calcProperties(*cellIter, numQuadPts);
- const double* density = _material->density();
+ const double* density = _material->calcDensity(*cellIter, numQuadPts);
// Compute Jacobian for cell
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -28,7 +28,10 @@
class pylith::materials::_ElasticIsotropic1D {
public:
- // Number of elastic constants (for general 3-D elastic material)
+ // Number of entries in stress tensor.
+ static const int stressSize;
+
+ // Number of entries in derivative of elasticity matrix.
static const int numElasticConsts;
// Values expected in spatial database
@@ -48,6 +51,7 @@
static const int pidLambda2Mu;
}; // _ElasticIsotropic1D
+const int pylith::materials::_ElasticIsotropic1D::stressSize = 1;
const int pylith::materials::_ElasticIsotropic1D::numElasticConsts = 1;
const int pylith::materials::_ElasticIsotropic1D::numDBValues = 2;
const char* pylith::materials::_ElasticIsotropic1D::namesDBValues[] =
@@ -64,6 +68,7 @@
// Default constructor.
pylith::materials::ElasticIsotropic1D::ElasticIsotropic1D(void)
{ // constructor
+ _dimension = 1;
} // constructor
// ----------------------------------------------------------------------
@@ -81,8 +86,16 @@
} // copy constructor
// ----------------------------------------------------------------------
+// Get number of entries in stress tensor.
+int
+pylith::materials::ElasticIsotropic1D::stressSize(void) const
+{ // stressSize
+ return _ElasticIsotropic1D::stressSize;
+} // stressSize
+
+// ----------------------------------------------------------------------
// Get number of elastic constants for material.
-const int
+int
pylith::materials::ElasticIsotropic1D::numElasticConsts(void) const
{ // numElasticConsts
return _ElasticIsotropic1D::numElasticConsts;
@@ -148,38 +161,67 @@
pylith::materials::ElasticIsotropic1D::_calcDensity(const double* parameters,
const int numParameters,
const int numLocs)
-{ // calcDensity
+{ // _calcDensity
assert(0 != _density);
assert(0 != parameters);
for (int iLoc=0, index=0; iLoc < numLocs; ++iLoc, index+=numParameters)
_density[iLoc] =
parameters[index+_ElasticIsotropic1D::pidDensity];
-} // calcDensity
+} // _calcDensity
// ----------------------------------------------------------------------
-// Compute density at location from parameters.
+// Compute stress tensor at location from parameters.
void
+pylith::materials::ElasticIsotropic1D::_calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
+{ // _calcStress
+ assert(0 != _stress);
+ assert(0 != parameters);
+ assert(_ElasticIsotropic1D::numParameters == numParameters);
+ assert(spaceDim == _dimension);
+
+ for (int iLoc=0, indexP=0, indexS=0;
+ iLoc < numLocs;
+ ++iLoc,
+ indexP+=_ElasticIsotropic1D::numParameters,
+ indexS+=_ElasticIsotropic1D::stressSize) {
+ const double density = parameters[indexP+_ElasticIsotropic1D::pidDensity];
+ const double lambda2mu = parameters[indexP+_ElasticIsotropic1D::pidLambda2Mu];
+
+ const double e11 = totalStrain[indexS ];
+ _stress[indexS ] = lambda2mu * e11;
+ } // for
+} // _calcStress
+
+// ----------------------------------------------------------------------
+// Compute derivative of elasticity matrix at location from parameters.
+void
pylith::materials::ElasticIsotropic1D::_calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs)
-{ // calcElasticConsts
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
+{ // _calcElasticConsts
assert(0 != _elasticConsts);
assert(0 != parameters);
assert(_ElasticIsotropic1D::numParameters == numParameters);
-
+ assert(spaceDim == _dimension);
+
for (int iLoc=0, indexP=0, indexC=0;
iLoc < numLocs;
++iLoc,
- indexP+=_ElasticIsotropic1D::numParameters,
- indexC+=_ElasticIsotropic1D::numElasticConsts) {
+ indexP+=_ElasticIsotropic1D::numParameters,
+ indexC+=_ElasticIsotropic1D::numElasticConsts) {
const double density = parameters[indexP+_ElasticIsotropic1D::pidDensity];
- const double lambda2mu =
- parameters[indexP+_ElasticIsotropic1D::pidLambda2Mu];
-
+ const double lambda2mu = parameters[indexP+_ElasticIsotropic1D::pidLambda2Mu];
+
_elasticConsts[indexC+ 0] = lambda2mu; // C1111
} // for
-} // calcElasticConsts
+} // _calcElasticConsts
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic1D.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -54,6 +54,16 @@
*/
ElasticMaterial* clone(void) const;
+ /** Get number of entries in stress tensor.
+ *
+ * 1-D = 1
+ * 2-D = 3
+ * 3-D = 6
+ *
+ * @returns Number of entries in stress tensor.
+ */
+ int stressSize(void) const;
+
/** Get number of elastic constants for material.
*
* 1-D = 1
@@ -62,7 +72,7 @@
*
* @returns Number of elastic constants
*/
- const int numElasticConsts(void) const;
+ int numElasticConsts(void) const;
// PROTECTED METHODS //////////////////////////////////////////////////
protected :
@@ -128,19 +138,39 @@
const int numParameters,
const int numLocs);
- /** Compute density at locations from parameters.
+ /** Compute stress tensor at locations from parameters.
*
- * Results are stored in _elasticConsts.
+ * Results are stored in _stress.
*
* Index into parameters = iLoc*numParameters + iParam
*
* @param parameters Parameters at location
* @param numParameters Number of parameters
+ * @param totalStrain Total strain at locations
* @param numLocs Number of locations
+ * @param spaceDim Spatial dimension for locations.
*/
+ void _calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
+
+ /** Compute derivatives of elasticity matrix at locations from parameters.
+ *
+ * Results are stored in _elasticConsts.
+ *
+ * @param parameters Parameters at locations.
+ * @param numParameters Number of parameters.
+ * @param totalStrain Total strain at locations.
+ * @param numLocs Number of locations.
+ * @param spaceDim Spatial dimension for locations.
+ */
void _calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs);
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -28,7 +28,10 @@
class pylith::materials::_ElasticIsotropic3D {
public:
- // Number of elastic constants (for general 3-D elastic material)
+ // Number of entries in stress tensor.
+ static const int stressSize;
+
+ // Number of entries in derivative of elasticity matrix.
static const int numElasticConsts;
// Values expected in spatial database
@@ -50,6 +53,7 @@
static const int pidLambda;
}; // _ElasticIsotropic3D
+const int pylith::materials::_ElasticIsotropic3D::stressSize = 6;
const int pylith::materials::_ElasticIsotropic3D::numElasticConsts = 21;
const int pylith::materials::_ElasticIsotropic3D::numDBValues = 3;
const char* pylith::materials::_ElasticIsotropic3D::namesDBValues[] =
@@ -68,6 +72,7 @@
// Default constructor.
pylith::materials::ElasticIsotropic3D::ElasticIsotropic3D(void)
{ // constructor
+ _dimension = 3;
} // constructor
// ----------------------------------------------------------------------
@@ -85,8 +90,16 @@
} // copy constructor
// ----------------------------------------------------------------------
-// Get number of elastic constants for material.
-const int
+// Get number of entries in stress tensor.
+int
+pylith::materials::ElasticIsotropic3D::stressSize(void) const
+{ // stressSize
+ return _ElasticIsotropic3D::stressSize;
+} // stressSize
+
+// ----------------------------------------------------------------------
+// Get number of entries in elasticity matrix for material.
+int
pylith::materials::ElasticIsotropic3D::numElasticConsts(void) const
{ // numElasticConsts
return _ElasticIsotropic3D::numElasticConsts;
@@ -155,58 +168,106 @@
pylith::materials::ElasticIsotropic3D::_calcDensity(const double* parameters,
const int numParameters,
const int numLocs)
-{ // calcDensity
+{ // _calcDensity
assert(0 != _density);
assert(0 != parameters);
for (int iLoc=0, index=0; iLoc < numLocs; ++iLoc, index+=numParameters)
_density[iLoc] =
parameters[index+_ElasticIsotropic3D::pidDensity];
-} // calcDensity
+} // _calcDensity
// ----------------------------------------------------------------------
-// Compute density at location from parameters.
+// Compute stress tensor at location from parameters.
void
+pylith::materials::ElasticIsotropic3D::_calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
+{ // _calcStress
+ assert(0 != _stress);
+ assert(0 != parameters);
+ assert(_ElasticIsotropic3D::numParameters == numParameters);
+ assert(spaceDim == _dimension);
+
+ for (int iLoc=0, indexP=0, indexS=0;
+ iLoc < numLocs;
+ ++iLoc,
+ indexP+=_ElasticIsotropic3D::numParameters,
+ indexS+=_ElasticIsotropic3D::stressSize) {
+ const double density = parameters[indexP+_ElasticIsotropic3D::pidDensity];
+ const double mu = parameters[indexP+_ElasticIsotropic3D::pidMu];
+ const double lambda = parameters[indexP+_ElasticIsotropic3D::pidLambda];
+
+ const double lambda2mu = lambda + 2.0 * mu;
+
+ const double e11 = totalStrain[indexS ];
+ const double e22 = totalStrain[indexS+1];
+ const double e33 = totalStrain[indexS+2];
+ const double e12 = totalStrain[indexS+3];
+ const double e23 = totalStrain[indexS+4];
+ const double e13 = totalStrain[indexS+5];
+
+ const double s123 = lambda * (e11 + e22 + e33);
+
+ _stress[indexS ] = s123 + 2.0*mu*e11;
+ _stress[indexS+1] = s123 + 2.0*mu*e22;
+ _stress[indexS+2] = s123 + 2.0*mu*e33;
+ _stress[indexS+3] = 2.0 * mu * e12;
+ _stress[indexS+4] = 2.0 * mu * e23;
+ _stress[indexS+5] = 2.0 * mu * e13;
+ } // for
+} // _calcStress
+
+// ----------------------------------------------------------------------
+// Compute derivative of elasticity matrix at location from parameters.
+void
pylith::materials::ElasticIsotropic3D::_calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs)
-{ // calcElasticConsts
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
+{ // _calcElasticConsts
assert(0 != _elasticConsts);
assert(0 != parameters);
assert(_ElasticIsotropic3D::numParameters == numParameters);
-
+ assert(spaceDim == _dimension);
+
for (int iLoc=0, indexP=0, indexC=0;
iLoc < numLocs;
++iLoc,
- indexP+=_ElasticIsotropic3D::numParameters,
- indexC+=_ElasticIsotropic3D::numElasticConsts) {
+ indexP+=_ElasticIsotropic3D::numParameters,
+ indexC+=_ElasticIsotropic3D::numElasticConsts) {
const double density = parameters[indexP+_ElasticIsotropic3D::pidDensity];
const double mu = parameters[indexP+_ElasticIsotropic3D::pidMu];
const double lambda = parameters[indexP+_ElasticIsotropic3D::pidLambda];
-
- _elasticConsts[indexC+ 0] = lambda + 2.0*mu; // C1111
+
+ const double lambda2mu = lambda + 2.0 * mu;
+
+ _elasticConsts[indexC+ 0] = lambda2mu; // C1111
_elasticConsts[indexC+ 1] = lambda; // C1122
_elasticConsts[indexC+ 2] = lambda; // C1133
_elasticConsts[indexC+ 3] = 0; // C1112
_elasticConsts[indexC+ 4] = 0; // C1123
_elasticConsts[indexC+ 5] = 0; // C1113
- _elasticConsts[indexC+ 6] = lambda + 2.0*mu; // C2222
+ _elasticConsts[indexC+ 6] = lambda2mu; // C2222
_elasticConsts[indexC+ 7] = lambda; // C2233
_elasticConsts[indexC+ 8] = 0; // C2212
_elasticConsts[indexC+ 9] = 0; // C2223
_elasticConsts[indexC+10] = 0; // C2213
- _elasticConsts[indexC+11] = lambda + 2.0*mu; // C3333
+ _elasticConsts[indexC+11] = lambda2mu; // C3333
_elasticConsts[indexC+12] = 0; // C3312
_elasticConsts[indexC+13] = 0; // C3323
_elasticConsts[indexC+14] = 0; // C3313
- _elasticConsts[indexC+15] = mu; // C1212
+ _elasticConsts[indexC+15] = 2.0 * mu; // C1212
_elasticConsts[indexC+16] = 0; // C1223
_elasticConsts[indexC+17] = 0; // C1213
- _elasticConsts[indexC+18] = mu; // C2323
+ _elasticConsts[indexC+18] = 2.0 * mu; // C2323
_elasticConsts[indexC+19] = 0; // C2313
- _elasticConsts[indexC+20] = mu; // C1313
+ _elasticConsts[indexC+20] = 2.0 * mu; // C1313
} // for
-} // calcElasticConsts
+} // _calcElasticConsts
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -55,15 +55,25 @@
*/
ElasticMaterial* clone(void) const;
- /** Get number of elastic constants for material.
+ /** Get number of entries in stress tensor.
*
* 1-D = 1
+ * 2-D = 3
+ * 3-D = 6
+ *
+ * @returns Number of entries in stress tensor.
+ */
+ int stressSize(void) const;
+
+ /** Get number of entries in derivative of elasticity matrix.
+ *
+ * 1-D = 1
* 2-D = 6
* 3-D = 21
*
- * @returns Number of elastic constants
+ * @returns Number of entries in derivative of elasticity matrix.
*/
- const int numElasticConsts(void) const;
+ int numElasticConsts(void) const;
// PROTECTED METHODS //////////////////////////////////////////////////
protected :
@@ -129,19 +139,39 @@
const int numParameters,
const int numLocs);
- /** Compute density at locations from parameters.
+ /** Compute stress tensor at locations from parameters.
*
- * Results are stored in _elasticConsts.
+ * Results are stored in _stress.
*
* Index into parameters = iLoc*numParameters + iParam
*
* @param parameters Parameters at location
* @param numParameters Number of parameters
+ * @param totalStrain Total strain at locations
* @param numLocs Number of locations
+ * @param spaceDim Spatial dimension for locations.
*/
+ void _calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
+
+ /** Compute derivatives of elasticity matrix at locations from parameters.
+ *
+ * Results are stored in _elasticConsts.
+ *
+ * @param parameters Parameters at locations.
+ * @param numParameters Number of parameters.
+ * @param totalStrain Total strain at locations.
+ * @param numLocs Number of locations.
+ * @param spaceDim Spatial dimension for locations.
+ */
void _calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs);
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -24,6 +24,7 @@
// Default constructor.
pylith::materials::ElasticMaterial::ElasticMaterial(void) :
_density(0),
+ _stress(0),
_elasticConsts(0)
{ // constructor
} // constructor
@@ -33,6 +34,7 @@
pylith::materials::ElasticMaterial::~ElasticMaterial(void)
{ // destructor
delete[] _density; _density = 0;
+ delete[] _stress; _stress = 0;
delete[] _elasticConsts; _elasticConsts = 0;
} // destructor
@@ -41,49 +43,121 @@
pylith::materials::ElasticMaterial::ElasticMaterial(const ElasticMaterial& m) :
Material(m),
_density(0),
+ _stress(0),
_elasticConsts(0)
{ // copy constructor
} // copy constructor
// ----------------------------------------------------------------------
-// Compute physical properties of cell at quadrature points.
-void
-pylith::materials::ElasticMaterial::calcProperties(
+// Compute density for cell at quadrature points.
+const double*
+pylith::materials::ElasticMaterial::calcDensity(
const Mesh::point_type& cell,
const int numQuadPts)
-{ // calcProperties
- _initCellData(numQuadPts);
+{ // calcDensity
+ assert(0 != _density);
+ int size = numQuadPts;
+ memset(_density, 0, size*sizeof(double));
+ double* paramsCell = 0;
+ _getParameters(¶msCell, cell, patch, numQuadPts);
const int numParams = _numParameters();
- const char** paramNames = _parameterNames();
+ _calcDensity(paramsCell, numParams, numQuadPts);
+ delete[] paramsCell; paramsCell = 0;
- int size = numParams * numQuadPts;
- double* paramsCell = (size > 0) ? new double[size] : 0;
- for (int iParam=0; iParam < numParams; ++iParam) {
- const ALE::Obj<real_section_type> parameter =
- _parameters->getReal(paramNames[iParam]);
+ return _density;
+} // calcDensity
- assert(numQuadPts == parameter->getFiberDimension(cell));
- const real_section_type::value_type* parameterCell =
- parameter->restrictPoint(cell);
- for (int iQuadPt=0; iQuadPt < numQuadPts; ++iQuadPt)
- paramsCell[iQuadPt*numParams+iParam] = parameterCell[iQuadPt];
- } // for
- _calcDensity(paramsCell, numParams, numQuadPts);
- _calcElasticConsts(paramsCell, numParams, numQuadPts);
+// ----------------------------------------------------------------------
+// Compute stress tensor for cell at quadrature points.
+const double*
+pylith::materials::ElasticMaterial::calcStress(
+ const topology_type::point_type& cell,
+ const topology_type::patch_type& patch,
+ const double* totalStrain,
+ const int numQuadPts,
+ const int spaceDim)
+{ // calcStress
+ assert(0 != totalStrain);
+ assert(0 != _stress);
+
+ int size = numQuadPts * stressSize();
+ memset(_stress, 0, size*sizeof(double));
+
+ double* paramsCell = 0;
+ _getParameters(¶msCell, cell, patch, numQuadPts);
+ const int numParams = _numParameters();
+ _calcStress(paramsCell, numParams, totalStrain, numQuadPts, spaceDim);
delete[] paramsCell; paramsCell = 0;
-} // calcProperties
+ return _stress;
+} // calcStress
+
// ----------------------------------------------------------------------
+// Compute derivative of elasticity matrix for cell at quadrature points.
+const double*
+pylith::materials::ElasticMaterial::calcDerivElastic(
+ const topology_type::point_type& cell,
+ const topology_type::patch_type& patch,
+ const double* totalStrain,
+ const int numQuadPts,
+ const int spaceDim)
+{ // calcDerivElastic
+ assert(0 != totalStrain);
+ assert(0 != _elasticConsts);
+
+ int size = numQuadPts * numElasticConsts();
+ memset(_elasticConsts, 0, size*sizeof(double));
+
+ double* paramsCell = 0;
+ _getParameters(¶msCell, cell, patch, numQuadPts);
+ const int numParams = _numParameters();
+ _calcElasticConsts(paramsCell, numParams, totalStrain, numQuadPts, spaceDim);
+ delete[] paramsCell; paramsCell = 0;
+
+ return _elasticConsts;
+} // calcDerivElastic
+
+// ----------------------------------------------------------------------
// Initialize arrays holding cell data.
void
pylith::materials::ElasticMaterial::_initCellData(const int numQuadPts)
{ // _initCellData
int size = numQuadPts;
delete[] _density; _density = (size > 0) ? new double[size] : 0;
+
+ size = numQuadPts * stressSize();
+ delete[] _stress; _stress = (size > 0) ? new double[size] : 0;
+
size = numQuadPts * numElasticConsts();
delete[] _elasticConsts; _elasticConsts = (size > 0) ? new double[size] : 0;
} // _initCellData
+// ----------------------------------------------------------------------
+// Get parameters for cell.
+void
+pylith::materials::ElasticMaterial::_getParameters(double** paramsCell,
+ const topology_type::point_type& cell,
+ const topology_type::patch_type& patch,
+ const int numQuadPts)
+{ // _getParameters
+ const int numParams = _numParameters();
+ const char** paramNames = _parameterNames();
+
+ const int size = numParams * numQuadPts;
+ delete[] *paramsCell; *paramsCell = (size > 0) ? new double[size] : 0;
+ for (int iParam=0; iParam < numParams; ++iParam) {
+ const ALE::Obj<real_section_type> parameter =
+ _parameters->getReal(paramNames[iParam]);
+
+ assert(numQuadPts == parameter->getFiberDimension(patch, cell));
+ const real_section_type::value_type* parameterCell =
+ parameter->restrict(patch, cell);
+ for (int iQuadPt=0; iQuadPt < numQuadPts; ++iQuadPt)
+ (*paramsCell)[iQuadPt*numParams+iParam] = parameterCell[iQuadPt];
+ } // for
+} // _getParameters
+
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -14,7 +14,7 @@
*
* @brief C++ ElasticMaterial object
*
- * Interface definition for 3-D linear and nonlinear elastic materials.
+ * Interface definition for linear and nonlinear elastic materials.
*/
#if !defined(pylith_materials_elasticmaterial_hh)
@@ -61,14 +61,59 @@
virtual
ElasticMaterial* clone(void) const = 0;
- /** Get density values for cell at quadrature points.
+ /** Compute density for cell at quadrature points.
*
+ * @param cell Finite-element cell
+ * @param patch Finite-element patch
+ * @param numQuadPts Number of quadrature points (consistency check)
+ *
* @returns Array of density values at cell's quadrature points.
*/
- const double* density(void) const;
+ const double* calcDensity(const topology_type::point_type& cell,
+ const topology_type::patch_type& patch,
+ const int numQuadPts);
+
+ /** Get stress tensor at quadrature points.
+ *
+ * Index into array of elasticity constants:
+ * index = iQuadPt*STRESSSIZE + iStress
+ *
+ * Order of stresses for 3-D:
+ * 0: S11, 1: S22, 2: S33, 3: S12, 4: S23, 5: S13,
+ *
+ * Order of stresses for 2-D:
+ * 0: S11, 1: S22, 2: S12,
+ *
+ * Order of elasticity constants for 1-D:
+ * 0: S11
+ *
+ * @param cell Finite-element cell
+ * @param patch Finite-element patch
+ * @param totalStrain Total strain tensor at quadrature points
+ * @param numQuadPts Number of quadrature points (consistency check)
+ * @param spaceDim Spatial dimension (consistency check)
+ *
+ * @returns Array of stresses at cell's quadrature points.
+ */
+ const double* calcStress(const topology_type::point_type& cell,
+ const topology_type::patch_type& patch,
+ const double* totalStrain,
+ const int numQuadPts,
+ const int spaceDim);
- /** Get elasticity constants for cell at quadrature points.
+ /** Get number of entries in stress tensor for material.
*
+ * 1-D = 1
+ * 2-D = 3
+ * 3-D = 6
+ *
+ * @returns Number of entries in stress tensor
+ */
+ virtual
+ int stressSize(void) const = 0;
+
+ /** Compute derivative of elasticity matrix for cell at quadrature points.
+ *
* Index into array of elasticity constants:
* index = iQuadPt*NUMELASTCONSTS + iConstant
*
@@ -88,25 +133,32 @@
* Order of elasticity constants for 1-D:
* 0: C1111
*
- * @returns Array of elasticity constants at cell's quadrature points.
+ * @param cell Finite-element cell
+ * @param patch Finite-element patch
+ * @param totalStrain Total strain tensor at quadrature points
+ * @param numQuadPts Number of quadrature points (consistency check)
+ * @param spaceDim Spatial dimension (consistency check)
*/
- const double* elasticConsts(void) const;
+ const double* calcDerivElastic(const topology_type::point_type& cell,
+ const topology_type::patch_type& patch,
+ const double* totalStrain,
+ const int numQuadPts,
+ const int spaceDim);
- /** Get number of elastic constants for material.
+ /** Get number of entries in derivatives of elasticity matrix for material.
*
* 1-D = 1
* 2-D = 6
* 3-D = 21
*
- * @returns Number of elastic constants
+ * @returns Number of entries in derivative of elasticity matrix
*/
virtual
- const int numElasticConsts(void) const = 0;
+ int numElasticConsts(void) const = 0;
/** Compute physical properties of cell at quadrature points.
*
* @param cell Finite-element cell
- * @param patch Finite-element patch
* @param numQuadPts Number of quadrature points (consistency check)
*/
void calcProperties(const Mesh::point_type& cell,
@@ -140,18 +192,39 @@
const int numParameters,
const int numLocs) = 0;
- /** Compute density at locations from parameters.
+ /** Compute stress density at locations from parameters.
*
+ * Results are stored in _stress.
+ *
+ * @param parameters Parameters at locations.
+ * @param numParameters Number of parameters.
+ * @param totalStrain Total strain at locations.
+ * @param numLocs Number of locations.
+ * @param spaceDim Spatial dimension for locations.
+ */
+ virtual
+ void _calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim) = 0;
+
+ /** Compute derivatives of elasticity matrix at locations from parameters.
+ *
* Results are stored in _elasticConsts.
*
- * @param parameters Parameters at location
- * @param numParameters Number of parameters
- * @param numLocs Number of locations
+ * @param parameters Parameters at locations.
+ * @param numParameters Number of parameters.
+ * @param totalStrain Total strain at locations.
+ * @param numLocs Number of locations.
+ * @param spaceDim Spatial dimension for locations.
*/
virtual
void _calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs) = 0;
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim) = 0;
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
@@ -169,17 +242,35 @@
*/
double* _density;
+ /** Array of stress tensor at quadrature points for current cell.
+ *
+ * size = stressSize*numQuadPts
+ * index = iQuadPt*stressSize+iStress
+ */
+ double* _stress;
+
/** Array of elasticity constants at quadrature points for current cell.
*
* size = numElasticConsts*numQuadPts
- * index = iQuadPt*numElasticConsts+iElasticConst
+ * index = iQuadPt*numElasticConsts+iConstant
*/
double* _elasticConsts;
+ // PRIVATE METHODS ////////////////////////////////////////////////////
+private :
+
+ /** Get parameters for cell.
+ *
+ * @param paramsCell Array of parameters for cell
+ * @param cell Finite-element cell
+ * @param numQuadPts Number of quadrature points (consistency check)
+ */
+ void _getParameters(double** paramsCells,
+ const Mesh::point_type& cell,
+ const int numQuadPts);
+
}; // class ElasticMaterial
-#include "ElasticMaterial.icc" // inline methods
-
#endif // pylith_materials_elasticmaterial_hh
Deleted: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.icc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -1,32 +0,0 @@
-// -*- 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
-
-// Get elasticity constants for cell at quadrature points.
-inline
-const double*
-pylith::materials::ElasticMaterial::elasticConsts(void) const {
- return _elasticConsts;
-}
-
-// Get density values for cell at quadrature points.
-inline
-const double*
-pylith::materials::ElasticMaterial::density(void) const {
- return _density;
-}
-
-
-// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -28,6 +28,9 @@
class pylith::materials::_ElasticPlaneStrain {
public:
+ // Number of entries in stress tensor.
+ static const int stressSize;
+
// Number of elastic constants (for general 3-D elastic material)
static const int numElasticConsts;
@@ -50,6 +53,7 @@
static const int pidLambda;
}; // _ElasticPlaneStrain
+const int pylith::materials::_ElasticPlaneStrain::stressSize = 3;
const int pylith::materials::_ElasticPlaneStrain::numElasticConsts = 6;
const int pylith::materials::_ElasticPlaneStrain::numDBValues = 3;
const char* pylith::materials::_ElasticPlaneStrain::namesDBValues[] =
@@ -68,6 +72,7 @@
// Default constructor.
pylith::materials::ElasticPlaneStrain::ElasticPlaneStrain(void)
{ // constructor
+ _dimension = 2;
} // constructor
// ----------------------------------------------------------------------
@@ -85,8 +90,16 @@
} // copy constructor
// ----------------------------------------------------------------------
+// Get number of entries in stress tensor.
+int
+pylith::materials::ElasticPlaneStrain::stressSize(void) const
+{ // stressSize
+ return _ElasticPlaneStrain::stressSize;
+} // stressSize
+
+// ----------------------------------------------------------------------
// Get number of elastic constants for material.
-const int
+int
pylith::materials::ElasticPlaneStrain::numElasticConsts(void) const
{ // numElasticConsts
return _ElasticPlaneStrain::numElasticConsts;
@@ -165,11 +178,50 @@
} // calcDensity
// ----------------------------------------------------------------------
+// Compute stress tensor at location from parameters.
+void
+pylith::materials::ElasticPlaneStrain::_calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
+{ // _calcStress
+ assert(0 != _stress);
+ assert(0 != parameters);
+ assert(_ElasticPlaneStrain::numParameters == numParameters);
+ assert(spaceDim == _dimension);
+
+ for (int iLoc=0, indexP=0, indexS=0;
+ iLoc < numLocs;
+ ++iLoc,
+ indexP+=_ElasticPlaneStrain::numParameters,
+ indexS+=_ElasticPlaneStrain::stressSize) {
+ const double density = parameters[indexP+_ElasticPlaneStrain::pidDensity];
+ const double mu = parameters[indexP+_ElasticPlaneStrain::pidMu];
+ const double lambda = parameters[indexP+_ElasticPlaneStrain::pidLambda];
+
+ const double lambda2mu = lambda + 2.0 * mu;
+
+ const double e11 = totalStrain[indexS ];
+ const double e22 = totalStrain[indexS+1];
+ const double e12 = totalStrain[indexS+3];
+
+ const double s12 = lambda * (e11 + e22);
+
+ _stress[indexS ] = s12 + 2.0*mu*e11;
+ _stress[indexS+1] = s12 + 2.0*mu*e22;
+ _stress[indexS+2] = 2.0 * mu * e12;
+ } // for
+} // _calcStress
+
+// ----------------------------------------------------------------------
// Compute density at location from parameters.
void
pylith::materials::ElasticPlaneStrain::_calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs)
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
{ // calcElasticConsts
assert(0 != _elasticConsts);
assert(0 != parameters);
@@ -183,13 +235,15 @@
const double density = parameters[indexP+_ElasticPlaneStrain::pidDensity];
const double mu = parameters[indexP+_ElasticPlaneStrain::pidMu];
const double lambda = parameters[indexP+_ElasticPlaneStrain::pidLambda];
+
+ const double lambda2mu = lambda + 2.0*mu;
- _elasticConsts[indexC+ 0] = lambda + 2.0*mu; // C1111
+ _elasticConsts[indexC+ 0] = lambda2mu; // C1111
_elasticConsts[indexC+ 1] = lambda; // C1122
_elasticConsts[indexC+ 3] = 0; // C1112
- _elasticConsts[indexC+ 6] = lambda + 2.0*mu; // C2222
+ _elasticConsts[indexC+ 6] = lambda2mu; // C2222
_elasticConsts[indexC+ 8] = 0; // C2212
- _elasticConsts[indexC+15] = mu; // C1212
+ _elasticConsts[indexC+15] = 2.0*mu; // C1212
} // for
} // calcElasticConsts
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -55,6 +55,16 @@
*/
ElasticMaterial* clone(void) const;
+ /** Get number of entries in stress tensor.
+ *
+ * 1-D = 1
+ * 2-D = 3
+ * 3-D = 6
+ *
+ * @returns Number of entries in stress tensor.
+ */
+ int stressSize(void) const;
+
/** Get number of elastic constants for material.
*
* 1-D = 1
@@ -63,7 +73,7 @@
*
* @returns Number of elastic constants
*/
- const int numElasticConsts(void) const;
+ int numElasticConsts(void) const;
// PROTECTED METHODS //////////////////////////////////////////////////
protected :
@@ -129,19 +139,39 @@
const int numParameters,
const int numLocs);
- /** Compute density at locations from parameters.
+ /** Compute stress tensor at locations from parameters.
*
- * Results are stored in _elasticConsts.
+ * Results are stored in _stress.
*
* Index into parameters = iLoc*numParameters + iParam
*
* @param parameters Parameters at location
* @param numParameters Number of parameters
+ * @param totalStrain Total strain at locations
* @param numLocs Number of locations
+ * @param spaceDim Spatial dimension for locations.
*/
+ void _calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
+
+ /** Compute derivatives of elasticity matrix at locations from parameters.
+ *
+ * Results are stored in _elasticConsts.
+ *
+ * @param parameters Parameters at locations.
+ * @param numParameters Number of parameters.
+ * @param totalStrain Total strain at locations.
+ * @param numLocs Number of locations.
+ * @param spaceDim Spatial dimension for locations.
+ */
void _calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs);
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -28,6 +28,9 @@
class pylith::materials::_ElasticPlaneStress {
public:
+ // Number of entries in stress tensor.
+ static const int stressSize;
+
// Number of elastic constants (for general 3-D elastic material)
static const int numElasticConsts;
@@ -50,6 +53,7 @@
static const int pidLambda;
}; // _ElasticPlaneStress
+const int pylith::materials::_ElasticPlaneStress::stressSize = 3;
const int pylith::materials::_ElasticPlaneStress::numElasticConsts = 6;
const int pylith::materials::_ElasticPlaneStress::numDBValues = 3;
const char* pylith::materials::_ElasticPlaneStress::namesDBValues[] =
@@ -68,6 +72,7 @@
// Default constructor.
pylith::materials::ElasticPlaneStress::ElasticPlaneStress(void)
{ // constructor
+ _dimension = 2;
} // constructor
// ----------------------------------------------------------------------
@@ -85,8 +90,16 @@
} // copy constructor
// ----------------------------------------------------------------------
+// Get number of entries in stress tensor.
+int
+pylith::materials::ElasticPlaneStress::stressSize(void) const
+{ // stressSize
+ return _ElasticPlaneStress::stressSize;
+} // stressSize
+
+// ----------------------------------------------------------------------
// Get number of elastic constants for material.
-const int
+int
pylith::materials::ElasticPlaneStress::numElasticConsts(void) const
{ // numElasticConsts
return _ElasticPlaneStress::numElasticConsts;
@@ -165,11 +178,51 @@
} // calcDensity
// ----------------------------------------------------------------------
+// Compute stress tensor at location from parameters.
+void
+pylith::materials::ElasticPlaneStress::_calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
+{ // _calcStress
+ assert(0 != _stress);
+ assert(0 != parameters);
+ assert(_ElasticPlaneStress::numParameters == numParameters);
+ assert(spaceDim == _dimension);
+
+ for (int iLoc=0, indexP=0, indexS=0;
+ iLoc < numLocs;
+ ++iLoc,
+ indexP+=_ElasticPlaneStress::numParameters,
+ indexS+=_ElasticPlaneStress::stressSize) {
+ const double density = parameters[indexP+_ElasticPlaneStress::pidDensity];
+ const double mu = parameters[indexP+_ElasticPlaneStress::pidMu];
+ const double lambda = parameters[indexP+_ElasticPlaneStress::pidLambda];
+
+ const double lambda2mu = lambda + 2.0 * mu;
+ const double lambdamu = lambda + mu;
+
+ const double e11 = totalStrain[indexS ];
+ const double e22 = totalStrain[indexS+1];
+ const double e12 = totalStrain[indexS+3];
+
+ _stress[indexS ] =
+ (4.0*mu*lambdamu * e11 + 2.0*mu*lambda * e22)/lambda2mu;
+ _stress[indexS+1] =
+ (2.0*mu*lambda * e11 + 4.0*mu*lambdamu * e22)/lambda2mu;
+ _stress[indexS+2] = 2.0 * mu * e12;
+ } // for
+} // _calcStress
+
+// ----------------------------------------------------------------------
// Compute density at location from parameters.
void
pylith::materials::ElasticPlaneStress::_calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs)
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim)
{ // calcElasticConsts
assert(0 != _elasticConsts);
assert(0 != parameters);
@@ -184,12 +237,15 @@
const double mu = parameters[indexP+_ElasticPlaneStress::pidMu];
const double lambda = parameters[indexP+_ElasticPlaneStress::pidLambda];
- _elasticConsts[indexC+ 0] = lambda + 2.0*mu; // C1111
- _elasticConsts[indexC+ 1] = lambda; // C1122
- _elasticConsts[indexC+ 2] = 0; // C1112
- _elasticConsts[indexC+ 3] = lambda + 2.0*mu; // C2222
- _elasticConsts[indexC+ 4] = 0; // C2212
- _elasticConsts[indexC+15] = mu; // C1212
+ const double lambda2mu = lambda + 2.0 * mu;
+ const double c11 = 4.0 * mu * (lambda + mu) / lambda2mu;
+
+ _elasticConsts[indexC+0] = c11; // C1111
+ _elasticConsts[indexC+1] = 2.0 * mu * lambda / lambda2mu; // C1122
+ _elasticConsts[indexC+2] = 0; // C1112
+ _elasticConsts[indexC+3] = c11; // C2222
+ _elasticConsts[indexC+4] = 0; // C2212
+ _elasticConsts[indexC+5] = 2.0 * mu; // C1212
} // for
} // calcElasticConsts
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -55,6 +55,16 @@
*/
ElasticMaterial* clone(void) const;
+ /** Get number of entries in stress tensor.
+ *
+ * 1-D = 1
+ * 2-D = 3
+ * 3-D = 6
+ *
+ * @returns Number of entries in stress tensor.
+ */
+ int stressSize(void) const;
+
/** Get number of elastic constants for material.
*
* 1-D = 1
@@ -63,7 +73,7 @@
*
* @returns Number of elastic constants
*/
- const int numElasticConsts(void) const;
+ int numElasticConsts(void) const;
// PROTECTED METHODS //////////////////////////////////////////////////
protected :
@@ -129,19 +139,39 @@
const int numParameters,
const int numLocs);
- /** Compute density at locations from parameters.
+ /** Compute stress tensor at locations from parameters.
*
- * Results are stored in _elasticConsts.
+ * Results are stored in _stress.
*
* Index into parameters = iLoc*numParameters + iParam
*
* @param parameters Parameters at location
* @param numParameters Number of parameters
+ * @param totalStrain Total strain at locations
* @param numLocs Number of locations
+ * @param spaceDim Spatial dimension for locations.
*/
+ void _calcStress(const double* parameters,
+ const int numParameters,
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
+
+ /** Compute derivatives of elasticity matrix at locations from parameters.
+ *
+ * Results are stored in _elasticConsts.
+ *
+ * @param parameters Parameters at locations.
+ * @param numParameters Number of parameters.
+ * @param totalStrain Total strain at locations.
+ * @param numLocs Number of locations.
+ * @param spaceDim Spatial dimension for locations.
+ */
void _calcElasticConsts(const double* parameters,
const int numParameters,
- const int numLocs);
+ const double* totalStrain,
+ const int numLocs,
+ const int spaceDim);
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Makefile.am 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/Makefile.am 2007-03-31 01:41:15 UTC (rev 6480)
@@ -19,7 +19,6 @@
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-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -29,6 +29,7 @@
// Default constructor.
pylith::materials::Material::Material(void) :
_parameters(0),
+ _dimension(0),
_db(0),
_id(0),
_label("")
@@ -49,6 +50,7 @@
// Copy constructor.
pylith::materials::Material::Material(const Material& m) :
_parameters(m._parameters),
+ _dimension(m._dimension),
_db(m._db),
_id(m._id),
_label(m._label)
@@ -163,7 +165,16 @@
// Close database
_db->close();
+
+ // Initialize cell data
+ _initCellData(numQuadPts);
} // initialize
+// ----------------------------------------------------------------------
+// Initialize arrays holding cell data.
+void
+pylith::materials::Material::_initCellData(const int numQuadPts)
+{}
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/materials/Material.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -68,6 +68,12 @@
virtual
~Material(void);
+ /** Get spatial dimension of material.
+ *
+ * @returns Spatial dimension.
+ */
+ int dimension(void) const;
+
/** Set database for physical property parameters.
*
* @param value Pointer to database.
@@ -161,6 +167,13 @@
const double* dbValues,
const int numValues) const = 0;
+ /** Initialize arrays holding cell data.
+ *
+ * @param numQuadPts Number of quadrature points
+ */
+ virtual
+ void _initCellData(const int numQuadPts);
+
// NOT IMPLEMENTED ////////////////////////////////////////////////////
private :
@@ -173,6 +186,8 @@
///< Manager of parameters for physical properties of material
pylith::feassemble::ParameterManager* _parameters;
+ int _dimension; ///< Spatial dimension associated with material
+
// PRIVATE MEMBERS ////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/materials/Material.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.icc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.icc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -14,6 +14,13 @@
#error "Material.icc can only be included from Material.hh"
#endif
+// Get spatial dimension of material.
+inline
+int
+pylith::materials::Material::dimension(void) const {
+ return _dimension;
+}
+
// Set database for material parameters.
inline
void
Modified: short/3D/PyLith/trunk/pylith/problems/Explicit.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/Explicit.py 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/pylith/problems/Explicit.py 2007-03-31 01:41:15 UTC (rev 6480)
@@ -138,7 +138,10 @@
"""
Hook for doing stuff after advancing time step.
"""
- self._info.log("WARNING: Explicit::poststep() not implemented.")
+ tmp = self.dispTmdt
+ self.dispTmdt = self.dispT
+ self.dispT = dispTpdt
+ self.dispTpdt = tmp
return
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -61,6 +61,16 @@
} // testCalcDensity
// ----------------------------------------------------------------------
+// Test calcStress()
+void
+pylith::materials::TestElasticIsotropic3D::testCalcStress(void)
+{ // testCalcStress
+ ElasticIsotropic3D material;
+ ElasticIsotropic3DData data;
+ _testCalcStress(&material, data);
+} // testCalcStress
+
+// ----------------------------------------------------------------------
// Test calcElasticConsts()
void
pylith::materials::TestElasticIsotropic3D::testCalcElasticConsts(void)
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticIsotropic3D.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -42,6 +42,7 @@
CPPUNIT_TEST( testDBValues );
CPPUNIT_TEST( testParameters );
CPPUNIT_TEST( testCalcDensity );
+ CPPUNIT_TEST( testCalcStress );
CPPUNIT_TEST( testCalcElasticConsts );
CPPUNIT_TEST_SUITE_END();
@@ -60,6 +61,9 @@
/// Test calcDensity()
void testCalcDensity(void);
+ /// Test calcStress()
+ void testCalcStress(void);
+
/// Test calcElasticConsts()
void testCalcElasticConsts(void);
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -43,48 +43,201 @@
} // testClone
// ----------------------------------------------------------------------
-// Test density()
+// Test calcDensity()
void
-pylith::materials::TestElasticMaterial::testDensity(void)
-{ // testDensity
- const double densityE[] = { 1.0, 8.0 };
- const int numQuadPts = 2;
+pylith::materials::TestElasticMaterial::testCalcDensity(void)
+{ // testCalcDensity
+ typedef ALE::Mesh::topology_type topology_type;
+ typedef topology_type::sieve_type sieve_type;
+ typedef ALE::Mesh::real_section_type real_section_type;
+ ALE::Obj<ALE::Mesh> mesh;
+ { // create mesh
+ const int cellDim = 1;
+ const int numCorners = 2;
+ const int spaceDim = 1;
+ const int numVertices = 2;
+ const int numCells = 1;
+ const double vertCoords[] = { -1.0, 1.0};
+ const int cells[] = { 0, 1};
+ CPPUNIT_ASSERT(0 != vertCoords);
+ CPPUNIT_ASSERT(0 != cells);
+
+ mesh = new ALE::Mesh(PETSC_COMM_WORLD, cellDim);
+ ALE::Obj<sieve_type> sieve = new sieve_type(mesh->comm());
+ ALE::Obj<topology_type> topology = new topology_type(mesh->comm());
+
+ const bool interpolate = false;
+ ALE::New::SieveBuilder<sieve_type>::buildTopology(sieve, cellDim, numCells,
+ const_cast<int*>(cells), numVertices, interpolate, numCorners);
+ sieve->stratify();
+ topology->setPatch(0, sieve);
+ topology->stratify();
+ mesh->setTopology(topology);
+ ALE::New::SieveBuilder<sieve_type>::buildCoordinates(
+ mesh->getRealSection("coordinates"), spaceDim, vertCoords);
+ } // create mesh
+
+ // Get cells associated with material
+ const ALE::Mesh::int_section_type::patch_type patch = 0;
+ const ALE::Obj<real_section_type>& coordinates =
+ mesh->getRealSection("coordinates");
+ const ALE::Obj<topology_type>& topology = coordinates->getTopology();
+ const ALE::Obj<topology_type::label_sequence>& cells =
+ topology->heightStratum(patch, 0);
+
ElasticIsotropic3D material;
- material._density = (numQuadPts > 0) ? new double[numQuadPts] : 0;
- CPPUNIT_ASSERT(0 != material._density);
- memcpy(material._density, densityE, numQuadPts*sizeof(double));
-
- const double* density = material.density();
- for (int i=0; i < numQuadPts; ++i)
- CPPUNIT_ASSERT_EQUAL(densityE[i], density[i]);
-} // testDensity
+ ElasticIsotropic3DData data;
+ delete material._parameters;
+ material._parameters = new feassemble::ParameterManager(mesh);
+ const int numQuadPts = 2;
+ const int fiberDim = numQuadPts; // number of values in field per cell
+ topology_type::label_sequence::iterator cellIter=cells->begin();
+
+ material._parameters->addReal("density");
+ const ALE::Obj<real_section_type>& parameterDensity =
+ material._parameters->getReal("density");
+ parameterDensity->setFiberDimension(patch, cells, fiberDim);
+ parameterDensity->allocate();
+ double cellData[numQuadPts];
+ cellData[0] = data.parameterData[0];
+ cellData[1] = data.parameterData[3];
+ parameterDensity->updateAdd(patch, *cellIter, cellData);
+
+ material._parameters->addReal("mu");
+ const ALE::Obj<real_section_type>& parameterMu =
+ material._parameters->getReal("mu");
+ parameterMu->setFiberDimension(patch, cells, fiberDim);
+ parameterMu->allocate();
+ cellData[0] = data.parameterData[1];
+ cellData[1] = data.parameterData[4];
+ parameterMu->updateAdd(patch, *cellIter, cellData);
+
+ material._parameters->addReal("lambda");
+ const ALE::Obj<real_section_type>& parameterLambda =
+ material._parameters->getReal("lambda");
+ parameterLambda->setFiberDimension(patch, cells, fiberDim);
+ parameterLambda->allocate();
+ cellData[0] = data.parameterData[2];
+ cellData[1] = data.parameterData[5];
+ parameterLambda->updateAdd(patch, *cellIter, cellData);
+
+ material._initCellData(numQuadPts);
+ const double* density = material.calcDensity(*cellIter, patch, numQuadPts);
+
+ const double tolerance = 1.0e-06;
+ const double* densityE = data.density;
+ CPPUNIT_ASSERT(0 != density);
+ CPPUNIT_ASSERT(0 != densityE);
+ const double size = numQuadPts;
+ for (int i=0; i < size; ++i)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, density[i]/densityE[i], tolerance);
+} // testCalcProperties
+
// ----------------------------------------------------------------------
-// Test numElasticConsts() and elasticConsts()
+// Test calcStress()
void
-pylith::materials::TestElasticMaterial::testElasticConsts(void)
-{ // testElasticConsts
- const double elasticConstsE[] = { 4.0, 1.0, 2.0 };
- const int numQuadPts = 3;
+pylith::materials::TestElasticMaterial::testCalcStress(void)
+{ // testCalcProperties
+ typedef ALE::Mesh::topology_type topology_type;
+ typedef topology_type::sieve_type sieve_type;
+ typedef ALE::Mesh::real_section_type real_section_type;
+ ALE::Obj<ALE::Mesh> mesh;
+ { // create mesh
+ const int cellDim = 1;
+ const int numCorners = 2;
+ const int spaceDim = 1;
+ const int numVertices = 2;
+ const int numCells = 1;
+ const double vertCoords[] = { -1.0, 1.0};
+ const int cells[] = { 0, 1};
+ CPPUNIT_ASSERT(0 != vertCoords);
+ CPPUNIT_ASSERT(0 != cells);
+
+ mesh = new ALE::Mesh(PETSC_COMM_WORLD, cellDim);
+ ALE::Obj<sieve_type> sieve = new sieve_type(mesh->comm());
+ ALE::Obj<topology_type> topology = new topology_type(mesh->comm());
+
+ const bool interpolate = false;
+ ALE::New::SieveBuilder<sieve_type>::buildTopology(sieve, cellDim, numCells,
+ const_cast<int*>(cells), numVertices, interpolate, numCorners);
+ sieve->stratify();
+ topology->setPatch(0, sieve);
+ topology->stratify();
+ mesh->setTopology(topology);
+ ALE::New::SieveBuilder<sieve_type>::buildCoordinates(
+ mesh->getRealSection("coordinates"), spaceDim, vertCoords);
+ } // create mesh
+
+ // Get cells associated with material
+ const ALE::Mesh::int_section_type::patch_type patch = 0;
+ const ALE::Obj<real_section_type>& coordinates =
+ mesh->getRealSection("coordinates");
+ const ALE::Obj<topology_type>& topology = coordinates->getTopology();
+ const ALE::Obj<topology_type::label_sequence>& cells =
+ topology->heightStratum(patch, 0);
+
ElasticIsotropic3D material;
- material._elasticConsts = (numQuadPts > 0) ? new double[numQuadPts] : 0;
- CPPUNIT_ASSERT(0 != material._elasticConsts);
- memcpy(material._elasticConsts, elasticConstsE, numQuadPts*sizeof(double));
-
- const double* elasticConsts = material.elasticConsts();
- for (int i=0; i < numQuadPts; ++i)
- CPPUNIT_ASSERT_EQUAL(elasticConstsE[i], elasticConsts[i]);
+ ElasticIsotropic3DData data;
+ delete material._parameters;
+ material._parameters = new feassemble::ParameterManager(mesh);
+ const int numQuadPts = 2;
+ const int fiberDim = numQuadPts; // number of values in field per cell
- CPPUNIT_ASSERT_EQUAL(21, material.numElasticConsts());
-} // testElasticConsts
+ topology_type::label_sequence::iterator cellIter=cells->begin();
+ material._parameters->addReal("density");
+ const ALE::Obj<real_section_type>& parameterDensity =
+ material._parameters->getReal("density");
+ parameterDensity->setFiberDimension(patch, cells, fiberDim);
+ parameterDensity->allocate();
+ double cellData[numQuadPts];
+ cellData[0] = data.parameterData[0];
+ cellData[1] = data.parameterData[3];
+ parameterDensity->updateAdd(patch, *cellIter, cellData);
+
+ material._parameters->addReal("mu");
+ const ALE::Obj<real_section_type>& parameterMu =
+ material._parameters->getReal("mu");
+ parameterMu->setFiberDimension(patch, cells, fiberDim);
+ parameterMu->allocate();
+ cellData[0] = data.parameterData[1];
+ cellData[1] = data.parameterData[4];
+ parameterMu->updateAdd(patch, *cellIter, cellData);
+
+ material._parameters->addReal("lambda");
+ const ALE::Obj<real_section_type>& parameterLambda =
+ material._parameters->getReal("lambda");
+ parameterLambda->setFiberDimension(patch, cells, fiberDim);
+ parameterLambda->allocate();
+ cellData[0] = data.parameterData[2];
+ cellData[1] = data.parameterData[5];
+ parameterLambda->updateAdd(patch, *cellIter, cellData);
+
+ material._initCellData(numQuadPts);
+ const double* stress = material.calcStress(*cellIter, patch, data.strain,
+ numQuadPts, data.dimension);
+
+ const double tolerance = 1.0e-06;
+ const double* stressE = data.stress;
+ CPPUNIT_ASSERT(0 != stress);
+ CPPUNIT_ASSERT(0 != stressE);
+ const double size = numQuadPts * material.stressSize();
+ for (int i=0; i < size; ++i)
+ if (fabs(stressE[i]) > tolerance)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, stress[i]/stressE[i],
+ tolerance);
+ else
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(stressE[i], stress[i], tolerance);
+} // testCalcStress
+
// ----------------------------------------------------------------------
-// Test calcProperties()
+// Test calcDerivElastic()
void
-pylith::materials::TestElasticMaterial::testCalcProperties(void)
-{ // testCalcProperties
+pylith::materials::TestElasticMaterial::testCalcDerivElastic(void)
+{ // testCalcDerivElastic
typedef ALE::Field::Mesh Mesh;
typedef Mesh::sieve_type sieve_type;
typedef Mesh::real_section_type real_section_type;
@@ -154,22 +307,12 @@
cellData[1] = data.parameterData[5];
parameterLambda->updateAddPoint(*cellIter, cellData);
- material.calcProperties(*cellIter, numQuadPts);
+ material._initCellData(numQuadPts);
+ const double* elasticConsts =
+ material.calcDerivElastic(*cellIter, patch, data.strain,
+ numQuadPts, data.dimension);
const double tolerance = 1.0e-06;
-
- { // check density
- const double* density = material.density();
- const double* densityE = data.density;
- CPPUNIT_ASSERT(0 != density);
- CPPUNIT_ASSERT(0 != densityE);
- const double size = numQuadPts;
- for (int i=0; i < size; ++i)
- CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, density[i]/densityE[i], tolerance);
- } // check density
-
- { // check elasticConsts
- const double* elasticConsts = material.elasticConsts();
const double* elasticConstsE = data.elasticConsts;
CPPUNIT_ASSERT(0 != elasticConsts);
CPPUNIT_ASSERT(0 != elasticConstsE);
@@ -181,8 +324,7 @@
else
CPPUNIT_ASSERT_DOUBLES_EQUAL(elasticConstsE[i], elasticConsts[i],
tolerance);
- } // check elasticConsts
-} // testCalcProperties
+} // testCalcDerivElastic
// ----------------------------------------------------------------------
// Test _initCellData()
@@ -196,6 +338,7 @@
const int numQuadPts = 4;
material._initCellData(numQuadPts);
CPPUNIT_ASSERT(0 != material._density);
+ CPPUNIT_ASSERT(0 != material._stress);
CPPUNIT_ASSERT(0 != material._elasticConsts);
} // testInitCellData
@@ -209,8 +352,8 @@
const int numQuadPts = data.numLocs;
material->_initCellData(numQuadPts);
material->_calcDensity(data.parameterData, data.numParameters, data.numLocs);
- const double* density = material->density();
const double* densityE = data.density;
+ const double* density = material->_density;
CPPUNIT_ASSERT(0 != density);
CPPUNIT_ASSERT(0 != densityE);
@@ -222,6 +365,34 @@
} // _testCalcDensity
// ----------------------------------------------------------------------
+// Test _calcStress()
+void
+pylith::materials::TestElasticMaterial::_testCalcStress(
+ ElasticMaterial* material,
+ const ElasticMaterialData& data) const
+{ // _testCalcElasticConsts
+ const int numQuadPts = data.numLocs;
+ material->_initCellData(numQuadPts);
+ material->_calcStress(data.parameterData, data.numParameters, data.strain,
+ data.numLocs, data.dimension);
+ const double* stressE = data.stress;
+ const double* stress = material->_stress;
+
+ CPPUNIT_ASSERT(0 != stress);
+ CPPUNIT_ASSERT(0 != stressE);
+
+ const double tolerance = 1.0e-06;
+ const double size = numQuadPts * material->stressSize();
+ for (int i=0; i < size; ++i)
+ if (fabs(stressE[i]) > tolerance)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, stress[i]/stressE[i],
+ tolerance);
+ else
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(stressE[i], stress[i],
+ tolerance);
+} // _testCalcStress
+
+// ----------------------------------------------------------------------
// Test _calcElasticConsts()
void
pylith::materials::TestElasticMaterial::_testCalcElasticConsts(
@@ -231,9 +402,9 @@
const int numQuadPts = data.numLocs;
material->_initCellData(numQuadPts);
material->_calcElasticConsts(data.parameterData, data.numParameters,
- data.numLocs);
- const double* elasticConsts = material->elasticConsts();
+ data.strain, data.numLocs, data.dimension);
const double* elasticConstsE = data.elasticConsts;
+ const double* elasticConsts = material->_elasticConsts;
CPPUNIT_ASSERT(0 != elasticConsts);
CPPUNIT_ASSERT(0 != elasticConstsE);
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/TestElasticMaterial.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -39,9 +39,9 @@
// CPPUNIT TEST SUITE /////////////////////////////////////////////////
CPPUNIT_TEST_SUITE( TestElasticMaterial );
CPPUNIT_TEST( testClone );
- CPPUNIT_TEST( testDensity );
- CPPUNIT_TEST( testElasticConsts );
- CPPUNIT_TEST( testCalcProperties );
+ CPPUNIT_TEST( testCalcDensity );
+ CPPUNIT_TEST( testCalcStress );
+ CPPUNIT_TEST( testCalcDerivElastic );
CPPUNIT_TEST( testInitCellData );
CPPUNIT_TEST_SUITE_END();
@@ -51,14 +51,14 @@
/// Test clone()
void testClone(void);
- /// Test density()
- void testDensity(void);
+ /// Test calcDensity()
+ void testCalcDensity(void);
- /// Test elaticConsts()
- void testElasticConsts(void);
+ /// Test calcStress()
+ void testCalcStress(void);
- /// Test calcProperties()
- void testCalcProperties(void);
+ /// Test calcDerivElastic()
+ void testCalcDerivElastic(void);
/// Test _initCellData()
void testInitCellData(void);
@@ -74,6 +74,14 @@
void _testCalcDensity(ElasticMaterial* material,
const ElasticMaterialData& data) const;
+ /** Test _calcStress()
+ *
+ * @param material Pointer to material
+ * @param data Data for testing material
+ */
+ void _testCalcStress(ElasticMaterial* material,
+ const ElasticMaterialData& data) const;
+
/** Test _calcElasticConsts()
*
* @param material Pointer to material
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3D.py
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3D.py 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3D.py 2007-03-31 01:41:15 UTC (rev 6480)
@@ -35,6 +35,8 @@
"""
ElasticMaterialApp.__init__(self, name)
+ self.dimension = 3
+
self.numDBValues = 3
self.dbValues = ["density", "vs", "vp"]
self.numParameters = 3
@@ -43,9 +45,13 @@
densityA = 2500.0
vsA = 3000.0
vpA = vsA*3**0.5
+ strainA = [1.1e-4, 2.2e-4, 3.3e-4, 4.4e-4, 5.5e-4, 6.6e-4]
+
densityB = 2000.0
vsB = 1200.0
vpB = vsB*3**0.5
+ strainB = [1.2e-4, 2.3e-4, 3.4e-4, 4.5e-4, 5.6e-4, 6.7e-4]
+
self.dbData = numpy.array([ [densityA, vsA, vpA],
[densityB, vsB, vpB] ],
dtype=numpy.float64)
@@ -61,16 +67,23 @@
numElasticConsts = 21
self.density = numpy.array([densityA, densityB],
dtype=numpy.float64)
+
+ self.strain = numpy.array([strainA, strainB],
+ dtype=numpy.float64)
+ self.stress = numpy.zeros( (self.numLocs, 6), dtype=numpy.float64)
self.elasticConsts = numpy.zeros( (self.numLocs, numElasticConsts),
dtype=numpy.float64)
- self.elasticConsts[0,:] = self._calcElasticConsts(densityA, muA, lambdaA)
- self.elasticConsts[1,:] = self._calcElasticConsts(densityB, muB, lambdaB)
+
+ (self.elasticConsts[0,:], self.stress[0,:]) = \
+ self._calcStress(strainA, densityA, muA, lambdaA)
+ (self.elasticConsts[1,:], self.stress[1,:]) = \
+ self._calcStress(strainB, densityB, muB, lambdaB)
return
- def _calcElasticConsts(self, densityV, muV, lambdaV):
+ def _calcStress(self, strainV, densityV, muV, lambdaV):
"""
- Compute elastic constants.
+ Compute stress and derivative of elasticity matrix.
"""
C1111 = lambdaV + 2.0*muV
C1122 = lambdaV
@@ -87,19 +100,29 @@
C3312 = 0.0
C3323 = 0.0
C3313 = 0.0
- C1212 = muV
+ C1212 = 2.0*muV
C1223 = 0.0
C1213 = 0.0
- C2323 = muV
+ C2323 = 2.0*muV
C2313 = 0.0
- C1313 = muV
+ C1313 = 2.0*muV
elasticConsts = numpy.array([C1111, C1122, C1133, C1112, C1123, C1113,
C2222, C2233, C2212, C2223, C2213,
C3333, C3312, C3323, C3313,
C1212, C1223, C1213,
C2323, C2313,
C1313], dtype=numpy.float64)
- return elasticConsts
+
+ strain = numpy.reshape(strainV, (6,1))
+ elastic = numpy.array([ [C1111, C1122, C1133, C1112, C1123, C1113],
+ [C1122, C2222, C2233, C2212, C2223, C2213],
+ [C1133, C2233, C3333, C3312, C3323, C3313],
+ [C1112, C2212, C3312, C1212, C1223, C1213],
+ [C1123, C2223, C3323, C1223, C2323, C2313],
+ [C1113, C2213, C3313, C1213, C2313, C1313] ],
+ dtype=numpy.float64)
+ stress = numpy.dot(elastic, strain)
+ return (elasticConsts, numpy.ravel(stress))
# MAIN /////////////////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -15,6 +15,8 @@
#include "ElasticIsotropic3DData.hh"
+const int pylith::materials::ElasticIsotropic3DData::_dimension = 3;
+
const int pylith::materials::ElasticIsotropic3DData::_numDBValues = 3;
const int pylith::materials::ElasticIsotropic3DData::_numParameters = 3;
@@ -56,6 +58,36 @@
2.00000000e+03,
};
+const double pylith::materials::ElasticIsotropic3DData::_strain[] = {
+ 1.10000000e-04,
+ 2.20000000e-04,
+ 3.30000000e-04,
+ 4.40000000e-04,
+ 5.50000000e-04,
+ 6.60000000e-04,
+ 1.20000000e-04,
+ 2.30000000e-04,
+ 3.40000000e-04,
+ 4.50000000e-04,
+ 5.60000000e-04,
+ 6.70000000e-04,
+};
+
+const double pylith::materials::ElasticIsotropic3DData::_stress[] = {
+ 1.98000000e+07,
+ 2.47500000e+07,
+ 2.97000000e+07,
+ 1.98000000e+07,
+ 2.47500000e+07,
+ 2.97000000e+07,
+ 2.67840000e+06,
+ 3.31200000e+06,
+ 3.94560000e+06,
+ 2.59200000e+06,
+ 3.22560000e+06,
+ 3.85920000e+06,
+};
+
const double pylith::materials::ElasticIsotropic3DData::_elasticConsts[] = {
6.75000000e+10,
2.25000000e+10,
@@ -72,12 +104,12 @@
0.00000000e+00,
0.00000000e+00,
0.00000000e+00,
- 2.25000000e+10,
+ 4.50000000e+10,
0.00000000e+00,
0.00000000e+00,
- 2.25000000e+10,
+ 4.50000000e+10,
0.00000000e+00,
- 2.25000000e+10,
+ 4.50000000e+10,
8.64000000e+09,
2.88000000e+09,
2.88000000e+09,
@@ -93,16 +125,17 @@
0.00000000e+00,
0.00000000e+00,
0.00000000e+00,
- 2.88000000e+09,
+ 5.76000000e+09,
0.00000000e+00,
0.00000000e+00,
- 2.88000000e+09,
+ 5.76000000e+09,
0.00000000e+00,
- 2.88000000e+09,
+ 5.76000000e+09,
};
pylith::materials::ElasticIsotropic3DData::ElasticIsotropic3DData(void)
{ // constructor
+ dimension = _dimension;
numDBValues = _numDBValues;
numParameters = _numParameters;
numLocs = _numLocs;
@@ -111,6 +144,8 @@
dbData = const_cast<double*>(_dbData);
parameterData = const_cast<double*>(_parameterData);
density = const_cast<double*>(_density);
+ strain = const_cast<double*>(_strain);
+ stress = const_cast<double*>(_stress);
elasticConsts = const_cast<double*>(_elasticConsts);
} // constructor
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticIsotropic3DData.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -37,6 +37,8 @@
private:
+ static const int _dimension;
+
static const int _numDBValues;
static const int _numParameters;
@@ -53,6 +55,10 @@
static const double _density[];
+ static const double _strain[];
+
+ static const double _stress[];
+
static const double _elasticConsts[];
};
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialApp.py
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialApp.py 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialApp.py 2007-03-31 01:41:15 UTC (rev 6480)
@@ -58,6 +58,7 @@
Script.__init__(self, name)
# Material information
+ self.dimension = None
self.numDBValues = None
self.numParameters = None
self.dbValues = None
@@ -68,7 +69,9 @@
# Elastic material information
self.numLocs = None
self.density = None
- self.parameters = None
+ self.strain = None
+ self.stress = None
+ self.elasticConsts = None
return
@@ -93,6 +96,9 @@
def _initData(self):
+ self.data.addScalar(vtype="int", name="_dimension",
+ value=self.dimension,
+ format="%d")
self.data.addScalar(vtype="int", name="_numDBValues",
value=self.numDBValues,
format="%d")
@@ -115,6 +121,12 @@
self.data.addArray(vtype="double", name="_density",
values=self.density,
format="%16.8e", ncols=1)
+ self.data.addArray(vtype="double", name="_strain",
+ values=self.strain,
+ format="%16.8e", ncols=1)
+ self.data.addArray(vtype="double", name="_stress",
+ values=self.stress,
+ format="%16.8e", ncols=1)
self.data.addArray(vtype="double", name="_elasticConsts",
values=self.elasticConsts,
format="%16.8e", ncols=1)
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -15,6 +15,8 @@
pylith::materials::ElasticMaterialData::ElasticMaterialData(void) :
numLocs(0),
density(0),
+ strain(0),
+ stress(0),
elasticConsts(0)
{ // constructor
} // constructor
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/ElasticMaterialData.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -39,6 +39,8 @@
int numLocs; ///< Number of locations
double* density; ///< Density at locations
+ double* strain; ///< Strain at locations
+ double* stress; ///< Stress at locations
double* elasticConsts; ///< Elastic constants at locations
};
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.cc 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.cc 2007-03-31 01:41:15 UTC (rev 6480)
@@ -15,6 +15,7 @@
// ----------------------------------------------------------------------
// Constructor
pylith::materials::MaterialData::MaterialData(void) :
+ dimension(0),
numDBValues(0),
numParameters(0),
dbValues(0),
Modified: short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.hh 2007-03-31 01:24:06 UTC (rev 6479)
+++ short/3D/PyLith/trunk/unittests/libtests/materials/data/MaterialData.hh 2007-03-31 01:41:15 UTC (rev 6480)
@@ -34,6 +34,8 @@
// PUBLIC MEMBERS ///////////////////////////////////////////////////////
public:
+ int dimension; ///< Number of dimensions
+
int numDBValues; ///< Number of database values
int numParameters; ///< Number of parameters
More information about the cig-commits
mailing list