Skip to content

Commit

Permalink
Add option to check scale of material properties when getting values …
Browse files Browse the repository at this point in the history 
…from spatial database.
  • Loading branch information
@baagaard-usgs
baagaard-usgs committed Nov 5, 2024
1 parent 77bfcc8 commit b8c03dc
Show file tree
Hide file tree
Showing 17 changed files with 263 additions and 103 deletions.
1 change: 1 addition & 0 deletions CHANGES.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@ The version numbers are in the form `MAJOR.MINOR.PATCH`, where major releases in
* **Added**
* Default filenames for progress monitor and parameters file are set from the simulation name like
the other output files.
* Consistency check for material properties and scales used in nondimensionalization (currently just the shear modulus).
* Add section to User Guide on troubleshooting solver issues.
* **Changed**
* Switch CI from Azure Pipelines to GitHub Actions.
Expand Down
1 change: 1 addition & 0 deletions libsrc/pylith/materials/AuxiliaryFactoryElastic.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -55,6 +55,7 @@ pylith::materials::AuxiliaryFactoryElastic::addShearModulus(void) {
description.componentNames[0] = subfieldName;
description.scale = pressureScale;
description.validator = pylith::topology::FieldQuery::validatorNonnegative;
description.validatorTolerance = 25.0;

_field->subfieldAdd(description, getSubfieldDiscretization(subfieldName));
pylith::materials::Query::shearModulusFromVM(subfieldName, this);
Expand Down
6 changes: 5 additions & 1 deletion libsrc/pylith/materials/AuxiliaryFactoryElasticity.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,11 @@ pylith::materials::AuxiliaryFactoryElasticity::addBodyForce(void) {
PYLITH_JOURNAL_DEBUG("addBodyForce(void)");

const char* subfieldName = "body_force";
const char* componentNames[3] = { "body_force_x", "body_force_y", "body_force_z" };
const char* componentNames[3] = {
"body_force_x",
"body_force_y",
"body_force_z",
};

const PylithReal forceScale = _normalizer->getPressureScale() / _normalizer->getLengthScale();

Expand Down
23 changes: 14 additions & 9 deletions libsrc/pylith/materials/AuxiliaryFactoryPoroelastic.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::~AuxiliaryFactoryPoroelastic(voi
// ----------------------------------------------------------------------------
// Add isotropic permeability subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addIsotropicPermeability(void) { // isotropicPermeablity
pylith::materials::AuxiliaryFactoryPoroelastic::addIsotropicPermeability(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addIsotropicPermeability(void)");

Expand All @@ -56,7 +56,8 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addIsotropicPermeability(void) {
description.componentNames.resize(1);
description.componentNames[0] = subfieldName;
description.scale = permeabilityScale;
description.validator = NULL;
// description.validator = pylith::topology::FieldQuery::validatorScale;
// description.validatorTolerance = 10.0;

_field->subfieldAdd(description, getSubfieldDiscretization(subfieldName));
this->setSubfieldQuery(subfieldName);
Expand All @@ -68,7 +69,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addIsotropicPermeability(void) {
// ----------------------------------------------------------------------------
// Add isotropic permeability subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addTensorPermeability(void) { // tensorPermeablity
pylith::materials::AuxiliaryFactoryPoroelastic::addTensorPermeability(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addTensorPermeability(void)");

Expand All @@ -95,7 +96,8 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addTensorPermeability(void) { //
description.componentNames[i] = componentNames[i];
} // for
description.scale = permeabilityScale;
description.validator = NULL;
// description.validator = pylith::topology::FieldQuery::validatorScale;
// description.validatorTolerance = 10.0;

_field->subfieldAdd(description, getSubfieldDiscretization(subfieldName));
this->setSubfieldQuery(subfieldName);
Expand All @@ -107,7 +109,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addTensorPermeability(void) { //
// --------------------------------------------------------------------
// Add drained bulk modulus subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addDrainedBulkModulus(void) { // DrainedBulkModulus
pylith::materials::AuxiliaryFactoryPoroelastic::addDrainedBulkModulus(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addDrainedBulkModulus(void)");

Expand All @@ -123,6 +125,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addDrainedBulkModulus(void) { //
description.componentNames[0] = subfieldName;
description.scale = pressureScale;
description.validator = pylith::topology::FieldQuery::validatorPositive;
description.validatorTolerance = 25.0;

_field->subfieldAdd(description, getSubfieldDiscretization(subfieldName));
this->setSubfieldQuery(subfieldName);
Expand All @@ -134,7 +137,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addDrainedBulkModulus(void) { //
// ---------------------------------------------------------------------
// Add biot coefficient subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addBiotCoefficient(void) { // biotCoefficient
pylith::materials::AuxiliaryFactoryPoroelastic::addBiotCoefficient(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addBiotCoefficient(void)");

Expand All @@ -149,6 +152,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addBiotCoefficient(void) { // bi
description.componentNames[0] = subfieldName;
description.scale = 1.0;
description.validator = pylith::topology::FieldQuery::validatorPositive;
description.validatorTolerance = 10.0;

_field->subfieldAdd(description, getSubfieldDiscretization(subfieldName));
this->setSubfieldQuery(subfieldName);
Expand All @@ -160,7 +164,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addBiotCoefficient(void) { // bi
// ---------------------------------------------------------------------
// Add biot modulus subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addBiotModulus(void) { // biotCoefficient
pylith::materials::AuxiliaryFactoryPoroelastic::addBiotModulus(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addBiotModulus(void)");

Expand All @@ -176,6 +180,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addBiotModulus(void) { // biotCo
description.componentNames[0] = subfieldName;
description.scale = pressureScale;
description.validator = pylith::topology::FieldQuery::validatorPositive;
description.validatorTolerance = 25.0;

_field->subfieldAdd(description, getSubfieldDiscretization(subfieldName));
pylith::materials::Query::biotModulusFromInput(subfieldName, this);
Expand All @@ -187,7 +192,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addBiotModulus(void) { // biotCo
// ----------------------------------------------------------------------
// Add reference stress subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addReferenceStress(void) { // referenceStress
pylith::materials::AuxiliaryFactoryPoroelastic::addReferenceStress(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addReferenceStress(void)");

Expand Down Expand Up @@ -225,7 +230,7 @@ pylith::materials::AuxiliaryFactoryPoroelastic::addReferenceStress(void) { // re
// ----------------------------------------------------------------------
// Add reference strain subfield to auxiliary fields.
void
pylith::materials::AuxiliaryFactoryPoroelastic::addReferenceStrain(void) { // addReferenceStrain
pylith::materials::AuxiliaryFactoryPoroelastic::addReferenceStrain(void) {
PYLITH_METHOD_BEGIN;
PYLITH_JOURNAL_DEBUG("addRefrenceStrain(void)");

Expand Down
2 changes: 1 addition & 1 deletion libsrc/pylith/materials/AuxiliaryFactoryPoroelasticity.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -116,7 +116,7 @@ pylith::materials::AuxiliaryFactoryPoroelasticity::addPorosity(void) { // porosi
PYLITH_JOURNAL_DEBUG("addPorosity(void)");

const char* subfieldName = "porosity";
const PylithReal noScale = 1;
const PylithReal noScale = 1.0;

pylith::topology::Field::Description description;
description.label = subfieldName;
Expand Down
4 changes: 2 additions & 2 deletions libsrc/pylith/materials/AuxiliaryFactoryPoroelasticity.hh
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,11 +10,11 @@
#pragma once

#include "pylith/materials/materialsfwd.hh" // forward declarations
#include "pylith/feassemble/AuxiliaryFactory.hh" // ISA AuxiliaryFactory
#include "pylith/materials/AuxiliaryFactoryElasticity.hh" // ISA AuxiliaryFactory

#include "spatialdata/spatialdb/spatialdbfwd.hh" // USES GravityField

class pylith::materials::AuxiliaryFactoryPoroelasticity : public pylith::feassemble::AuxiliaryFactory {
class pylith::materials::AuxiliaryFactoryPoroelasticity : public pylith::materials::AuxiliaryFactoryElasticity {
friend class TestAuxiliaryFactoryPoroelasticity; // unit testing

// PUBLIC METHODS //////////////////////////////////////////////////////////////////////////////////////////////////
Expand Down
7 changes: 5 additions & 2 deletions libsrc/pylith/materials/AuxiliaryFactoryViscoelastic.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -114,8 +114,11 @@ pylith::materials::AuxiliaryFactoryViscoelastic::addShearModulusRatioGeneralized
description.vectorFieldType = pylith::topology::Field::OTHER;
description.numComponents = 3;
description.componentNames.resize(3);
const char* componentNames[3] = { "shear_modulus_ratio_1", "shear_modulus_ratio_2",
"shear_modulus_ratio_3" };
const char* componentNames[3] = {
"shear_modulus_ratio_1",
"shear_modulus_ratio_2",
"shear_modulus_ratio_3",
};
for (int i = 0; i < 3; ++i) {
description.componentNames[i] = componentNames[i];
} // for
Expand Down
38 changes: 19 additions & 19 deletions libsrc/pylith/materials/Query.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -222,9 +222,9 @@ pylith::materials::_Query::vmToShearModulus(PylithScalar valueSubfield[],

std::ostringstream msg;
if (density <= 0) {
msg << "Found negative density (" << density << ").";
msg << "Found non-positive density (" << density << ").";
} // if
if (vs <= 0) {
if (vs < 0) {
msg << "Found negative shear wave speed (" << vs << ").";
} // if

Expand Down Expand Up @@ -257,13 +257,13 @@ pylith::materials::_Query::vmToBulkModulus(PylithScalar valueSubfield[],

std::ostringstream msg;
if (density <= 0) {
msg << "Found nonpositive density (" << density << ").";
msg << "Found non-positive density (" << density << ").";
} // if
if (vs < 0) {
msg << "Found negative shear wave speed (" << vs << ").";
} // if
if (vp <= 0) {
msg << "Found nonpositive dilatational wave speed (" << vp << ").";
msg << "Found non-positive dilatational wave speed (" << vp << ").";
} // if

PYLITH_METHOD_RETURN(msg.str());
Expand Down Expand Up @@ -296,13 +296,13 @@ pylith::materials::_Query::vmToMaxwellTime(PylithScalar valueSubfield[],

std::ostringstream msg;
if (density <= 0) {
msg << "Found negative density (" << density << ").";
msg << "Found non-positive density (" << density << ").";
} // if
if (vs <= 0) {
msg << "Found negative shear wave speed (" << vs << ").";
msg << "Found non-positive shear wave speed (" << vs << ").";
} // if
if (viscosity <= 0) {
msg << "Found nonpositive viscosity (" << viscosity << ").";
msg << "Found non-positive viscosity (" << viscosity << ").";
} // if

PYLITH_METHOD_RETURN(msg.str());
Expand Down Expand Up @@ -351,28 +351,28 @@ pylith::materials::_Query::vmToGeneralizedMaxwellTimes(PylithScalar valueSubfiel

std::ostringstream msg;
if (density <= 0) {
msg << "Found negative density (" << density << ").";
msg << "Found non-positive density (" << density << ").";
} // if
if (vs <= 0) {
msg << "Found negative shear wave speed (" << vs << ").";
msg << "Found non-positive shear wave speed (" << vs << ").";
} // if
if (viscosity1 <= 0) {
msg << "Found nonpositive viscosity 1 (" << viscosity1 << ").";
msg << "Found non-positive viscosity 1 (" << viscosity1 << ").";
} // if
if (viscosity2 <= 0) {
msg << "Found nonpositive viscosity 2 (" << viscosity2 << ").";
msg << "Found non-positive viscosity 2 (" << viscosity2 << ").";
} // if
if (viscosity3 <= 0) {
msg << "Found nonpositive viscosity 3 (" << viscosity3 << ").";
msg << "Found non-positive viscosity 3 (" << viscosity3 << ").";
} // if
if (shearModulusRatio1 <= 0) {
msg << "Found negative shear modulus ratio 1 (" << shearModulusRatio1 << ").";
msg << "Found non-positive shear modulus ratio 1 (" << shearModulusRatio1 << ").";
} // if
if (shearModulusRatio2 <= 0) {
msg << "Found negative shear modulus ratio 2 (" << shearModulusRatio2 << ").";
msg << "Found non-positive shear modulus ratio 2 (" << shearModulusRatio2 << ").";
} // if
if (shearModulusRatio3 <= 0) {
msg << "Found negative shear modulus ratio 3 (" << shearModulusRatio3 << ").";
msg << "Found non-positive shear modulus ratio 3 (" << shearModulusRatio3 << ").";
} // if

const double ratioSum = shearModulusRatio1 + shearModulusRatio2 + shearModulusRatio3;
Expand Down Expand Up @@ -409,13 +409,13 @@ pylith::materials::_Query::vmToGeneralizedMaxwellShearModulusRatios(PylithScalar

std::ostringstream msg;
if (shearModulusRatio1 <= 0) {
msg << "Found negative shear modulus ratio 1 (" << shearModulusRatio1 << ").";
msg << "Found non-positive shear modulus ratio 1 (" << shearModulusRatio1 << ").";
} // if
if (shearModulusRatio2 <= 0) {
msg << "Found negative shear modulus ratio 2 (" << shearModulusRatio2 << ").";
msg << "Found non-positive shear modulus ratio 2 (" << shearModulusRatio2 << ").";
} // if
if (shearModulusRatio3 <= 0) {
msg << "Found negative shear modulus ratio 3 (" << shearModulusRatio3 << ").";
msg << "Found non-positive shear modulus ratio 3 (" << shearModulusRatio3 << ").";
} // if
const double ratioSum = shearModulusRatio1 + shearModulusRatio2 + shearModulusRatio3;
if (ratioSum > 1) {
Expand Down Expand Up @@ -498,7 +498,7 @@ pylith::materials::_Query::inputToBiotModulus(PylithScalar valueSubfield[],
msg << "Found negative porosity (" << porosity << ").";
} // if
if (biot_coefficient <= 0) {
msg << "Found negative biot coefficient (" << biot_coefficient << ").";
msg << "Found non-positive biot coefficient (" << biot_coefficient << ").";
} // if

// Debug
Expand Down
8 changes: 7 additions & 1 deletion libsrc/pylith/topology/FieldBase.hh
View file
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,9 @@ public:
public:

/// Function prototype for validator functions.
typedef const char* (*validatorfn_type)(const PylithReal);
typedef std::string (*validatorfn_type)(const PylithReal /*value */,
const PylithReal /* scale */,
const PylithReal /* tolerance */);

// PUBLIC STRUCTS ///////////////////////////////////////////////////////
public:
Expand All @@ -58,6 +60,7 @@ public:
pylith::string_vector componentNames; ///< Names of components.
size_t numComponents; ///< Number of components.
PylithReal scale; ///< Dimension scale associated with values.
PylithReal validatorTolerance; ///< Tolerance relative to scale for validation.
validatorfn_type validator; ///< Validator for values in field;
bool hasHistory; ///< Has subfields with history, i.e., state variables.
size_t historySize; ///< Number of points in time history (currently only).
Expand All @@ -68,6 +71,7 @@ public:
vectorFieldType(OTHER),
numComponents(0),
scale(1.0),
validatorTolerance(0.0),
validator(NULL),
hasHistory(false),
historySize(0) {}
Expand All @@ -79,6 +83,7 @@ public:
const size_t numComponentsValue=0,
const VectorFieldEnum vectorFieldTypeValue=SCALAR,
const PylithReal scaleValue=1.0,
const PylithReal validatorToleranceValue=0.0,
const validatorfn_type validatorValue=NULL,
bool isFaultOnlyValue=false,
bool hasHistoryValue=false,
Expand All @@ -89,6 +94,7 @@ public:
componentNames(componentNamesValue),
numComponents(numComponentsValue),
scale(scaleValue),
validatorTolerance(validatorToleranceValue),
validator(validatorValue),
hasHistory(hasHistoryValue),
historySize(historySizeValue) {}
Expand Down
Loading

0 comments on commit b8c03dc

Please sign in to comment.