[cig-commits] r15056 - in short/3D/PyLith/trunk: . libsrc/feassemble libsrc/materials libsrc/meshio pylith/feassemble
brad at geodynamics.org
brad at geodynamics.org
Tue May 26 11:58:34 PDT 2009
Author: brad
Date: 2009-05-26 11:58:33 -0700 (Tue, 26 May 2009)
New Revision: 15056
Modified:
short/3D/PyLith/trunk/TODO
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc
short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
short/3D/PyLith/trunk/libsrc/materials/Material.cc
short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc
short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py
Log:
Fixed nondimensionalization bugs. Setting normalizer for materials was missing. Fixed setting scales for material properties and state variables.
Modified: short/3D/PyLith/trunk/TODO
===================================================================
--- short/3D/PyLith/trunk/TODO 2009-05-26 16:27:18 UTC (rev 15055)
+++ short/3D/PyLith/trunk/TODO 2009-05-26 18:58:33 UTC (rev 15056)
@@ -9,7 +9,6 @@
Brad
Test precompute geometry [DONE]
Nondimensionalization
- Fault coordinates - need coordinates_dimensioned
Material state vars - use dimstatevars and dimproperties to get scales
cleanup
full-scale testing
@@ -44,22 +43,12 @@
1. SNES
- Cleanup SolutionFields. [full scale test]
-
Power-law rheology
2. Nondimensionalization
- Add Quadrature::computeGeometry(coordinatesCell)
- Use #define PRECOMPUTE_GEOMETRY [TEST]
+ const slip rate in output (slip rate not slip)
- DataWriteVTK needs dimensionalizer
-
- Done coordinates in output (nondimensional instead of dimensioned)
- create dimensionalized copy (NAME_dimensioned) solution in output
- (nondimensional instead of dimensioned) const slip rate in output
- (slip rate not slip)
-
Ask constraints if a block matrix is okay. If okay and matrix type
is "unknown" (not set by user), then use "baij" rather than
"aij". Do this in Python.
@@ -200,7 +189,6 @@
Cleanup SlipTimeFn tests (refactor test/initialize stuff)
6. Memory model [Matt]
- playpen/memcheck/*
7. 2-D Plane strain Maxwell viscoelastic rheology [Charles]
Modified: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc 2009-05-26 16:27:18 UTC (rev 15055)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc 2009-05-26 18:58:33 UTC (rev 15056)
@@ -328,10 +328,10 @@
_allocateTensorField(mesh);
topology::Field<topology::Mesh>& buffer =
_outputFields->get("buffer (tensor)");
- _calcStrainStressField(&buffer, name, fields);
buffer.label("total_strain");
buffer.scale(1.0);
buffer.addDimensionOkay(true);
+ _calcStrainStressField(&buffer, name, fields);
return buffer;
} else if (!_material->hasStateVars() && 0 == strcasecmp(name, "stress")) {
@@ -339,10 +339,10 @@
_allocateTensorField(mesh);
topology::Field<topology::Mesh>& buffer =
_outputFields->get("buffer (tensor)");
- _calcStrainStressField(&buffer, name, fields);
buffer.label("stress");
buffer.scale(_normalizer->pressureScale());
buffer.addDimensionOkay(true);
+ _calcStrainStressField(&buffer, name, fields);
return buffer;
} else if (0 == strcasecmp(name, "stress")) {
@@ -351,10 +351,10 @@
topology::Field<topology::Mesh>& buffer =
_outputFields->get("buffer (tensor)");
_material->getField(&buffer, "total_strain");
- _calcStressFromStrain(&buffer);
buffer.label(name);
buffer.scale(_normalizer->pressureScale());
buffer.addDimensionOkay(true);
+ _calcStressFromStrain(&buffer);
return buffer;
} else {
@@ -446,10 +446,6 @@
double_array stressCell(numQuadPts*tensorSize);
stressCell = 0.0;
- // Get normalizer
- assert(0 != _normalizer);
- const double pressureScale = _normalizer->pressureScale();
-
// Get cell information
const ALE::Obj<SieveMesh>& sieveMesh = field->mesh().sieveMesh();
assert(!sieveMesh.isNull());
@@ -509,8 +505,6 @@
else {
_material->retrievePropsAndVars(*c_iter);
stressCell = _material->calcStress(strainCell);
- _normalizer->dimensionalize(&stressCell[0], stressCell.size(),
- pressureScale);
fieldSection->updatePoint(*c_iter, &stressCell[0]);
} // else
} // for
@@ -537,10 +531,6 @@
double_array stressCell(numQuadPts*tensorSize);
stressCell = 0.0;
- // Get normalizer
- assert(0 != _normalizer);
- const double pressureScale = _normalizer->pressureScale();
-
// Get cell information
const ALE::Obj<SieveMesh>& sieveMesh = field->mesh().sieveMesh();
assert(!sieveMesh.isNull());
@@ -562,8 +552,6 @@
fieldSection->restrictPoint(*c_iter, &strainCell[0], strainCell.size());
_material->retrievePropsAndVars(*c_iter);
stressCell = _material->calcStress(strainCell);
- _normalizer->dimensionalize(&stressCell[0], stressCell.size(),
- pressureScale);
fieldSection->updatePoint(*c_iter, &stressCell[0]);
} // for
} // _calcStressFromStrain
Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc 2009-05-26 16:27:18 UTC (rev 15055)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticMaterial.cc 2009-05-26 18:58:33 UTC (rev 15056)
@@ -608,10 +608,6 @@
} // if
} // for
- // Nondimensionalize strain
- _normalizer->nondimensionalize(&strainCell[0], strainCell.size(),
- pressureScale);
-
initialStrainSection->updatePoint(*c_iter, &strainCell[0]);
} // for
Modified: short/3D/PyLith/trunk/libsrc/materials/Material.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/Material.cc 2009-05-26 16:27:18 UTC (rev 15055)
+++ short/3D/PyLith/trunk/libsrc/materials/Material.cc 2009-05-26 18:58:33 UTC (rev 15056)
@@ -314,8 +314,6 @@
_metadata.fiberDim(properties[i].c_str(), Metadata::PROPERTY);
const int fiberDim = _metadata.fiberDim(name, Metadata::PROPERTY);
- // :TODO: Get scale information
-
// Get properties section
const ALE::Obj<RealSection>& propertiesSection = _properties->section();
assert(!propertiesSection.isNull());
@@ -331,6 +329,11 @@
assert(totalPropsFiberDim == numQuadPts * numPropsQuadPt);
const int totalFiberDim = numQuadPts * fiberDim;
+ // Get dimension scale information for properties.
+ double_array propertyScales(numPropsQuadPt);
+ propertyScales = 1.0;
+ _dimProperties(&propertyScales[0], propertyScales.size());
+
ALE::MemoryLogger& logger = ALE::MemoryLogger::singleton();
logger.stagePush("Materials");
// Allocate buffer for property field if necessary.
@@ -353,6 +356,7 @@
} // if
assert(!fieldSection.isNull());
field->label(name);
+ field->scale(propertyScales[propOffset]);
fieldType = _metadata.fieldType(name, Metadata::PROPERTY);
logger.stagePop();
@@ -367,13 +371,10 @@
propertiesSection->restrictPoint(*c_iter,
&propertiesCell[0], propertiesCell.size());
- for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
- _dimProperties(&propertiesCell[iQuad*numPropsQuadPt],
- numPropsQuadPt);
- memcpy(&fieldCell[iQuad*fiberDim],
- &propertiesCell[iQuad*numPropsQuadPt+propOffset],
- fiberDim*sizeof(double));
- } // for
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad)
+ for (int i=0; i < fiberDim; ++i)
+ fieldCell[iQuad*fiberDim+i] =
+ propertiesCell[iQuad*numPropsQuadPt+propOffset+i];
fieldSection->updatePoint(*c_iter, &fieldCell[0]);
} // for
@@ -388,8 +389,6 @@
_metadata.fiberDim(stateVars[i].c_str(), Metadata::STATEVAR);
const int fiberDim = _metadata.fiberDim(name, Metadata::STATEVAR);
- // :TODO: Get scale information
-
// Get state variables section
const ALE::Obj<RealSection>& stateVarsSection = _stateVars->section();
assert(!stateVarsSection.isNull());
@@ -405,6 +404,11 @@
assert(totalVarsFiberDim == numQuadPts * numVarsQuadPt);
const int totalFiberDim = numQuadPts * fiberDim;
+ // Get dimension scale information for state variables.
+ double_array stateVarScales(numVarsQuadPt);
+ stateVarScales = 1.0;
+ _dimStateVars(&stateVarScales[0], stateVarScales.size());
+
ALE::MemoryLogger& logger = ALE::MemoryLogger::singleton();
logger.stagePush("Materials");
// Allocate buffer for state variable field if necessary.
@@ -428,6 +432,7 @@
assert(!fieldSection.isNull());
fieldType = _metadata.fieldType(name, Metadata::STATEVAR);
field->label(name);
+ field->scale(stateVarScales[varOffset]);
logger.stagePop();
// Buffer for state variable at cell's quadrature points
@@ -441,14 +446,11 @@
stateVarsSection->restrictPoint(*c_iter,
&stateVarsCell[0], stateVarsCell.size());
- for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
- _dimStateVars(&stateVarsCell[iQuad*numVarsQuadPt],
- numVarsQuadPt);
- memcpy(&fieldCell[iQuad*fiberDim],
- &stateVarsCell[iQuad*numVarsQuadPt+varOffset],
- fiberDim*sizeof(double));
- } // for
-
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad)
+ for (int i=0; i < fiberDim; ++i)
+ fieldCell[iQuad*fiberDim+i] =
+ stateVarsCell[iQuad*numVarsQuadPt+varOffset+i];
+
fieldSection->updatePoint(*c_iter, &fieldCell[0]);
} // for
} // if/else
Modified: short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc 2009-05-26 16:27:18 UTC (rev 15055)
+++ short/3D/PyLith/trunk/libsrc/meshio/CellFilterAvg.cc 2009-05-26 18:58:33 UTC (rev 15056)
@@ -101,11 +101,6 @@
_fieldAvg->allocate();
} // else
- std::cout << "FIELD " << fieldIn.label()
- << ", fiberDim: " << totalFiberDim << "\n"
- << "AVG FIELD: " << _fieldAvg->label()
- << ", fiberDim: " << fiberDim << std::endl;
-
assert(0 != _fieldAvg);
switch (fieldIn.vectorFieldType())
{ // switch
Modified: short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py 2009-05-26 16:27:18 UTC (rev 15055)
+++ short/3D/PyLith/trunk/pylith/feassemble/IntegratorElasticity.py 2009-05-26 18:58:33 UTC (rev 15056)
@@ -86,6 +86,7 @@
self._info.log("Initializing integrator for material '%s'." % \
self.materialObj.label)
Integrator.initialize(self, totalTime, numTimeSteps, normalizer)
+ self.materialObj.normalizer(normalizer)
self._logger.eventEnd(logEvent)
return
More information about the CIG-COMMITS
mailing list