Skip to content

Commit

Permalink
[GeoMechanicsApplication] Removed the strain calculation from geo con…
Browse files Browse the repository at this point in the history
…stitutive laws (#12719)

* Removed the strain calculation from the linear elastic laws in geomechanics
* Removed the strain calculation from the geo UDSM laws
* Removed strain calculation of the UMAT constitutive laws in geo
  • Loading branch information
rfaasse authored Oct 4, 2024
1 parent bd503d9 commit 3e10ca1
Show file tree
Hide file tree
Showing 26 changed files with 26 additions and 288 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -117,22 +117,4 @@ void ElasticIsotropicK03DLaw::CalculatePK2Stress(const Vector& rS
}
}

void ElasticIsotropicK03DLaw::CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector)
{
const SizeType space_dimension = this->WorkingSpaceDimension();

// 1.-Compute total deformation gradient
const Matrix& F = rValues.GetDeformationGradientF();
KRATOS_DEBUG_ERROR_IF(F.size1() != space_dimension || F.size2() != space_dimension)
<< "expected size of F " << space_dimension << "x" << space_dimension << ", got "
<< F.size1() << "x" << F.size2() << std::endl;

Matrix E_tensor = prod(trans(F), F);
for (unsigned int i = 0; i < space_dimension; ++i)
E_tensor(i, i) -= 1.0;
E_tensor *= 0.5;

noalias(rStrainVector) = MathUtils<double>::StrainTensorToVector(E_tensor);
}

} // Namespace Kratos
Original file line number Diff line number Diff line change
Expand Up @@ -138,13 +138,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) ElasticIsotropicK03DLaw : public Geo
Vector& rStressVector,
ConstitutiveLaw::Parameters& rValues) override;

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}

private:
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -52,21 +52,19 @@ void GeoLinearElasticLaw::CalculateMaterialResponsePK2(ConstitutiveLaw::Paramete

const Flags& r_options = rValues.GetOptions();

Vector& r_strain_vector = rValues.GetStrainVector();
KRATOS_DEBUG_ERROR_IF(r_options.IsNot(ConstitutiveLaw::USE_ELEMENT_PROVIDED_STRAIN))
<< "The GeoLinearElasticLaw needs an element provided strain" << std::endl;

// NOTE: SINCE THE ELEMENT IS IN SMALL STRAINS WE CAN USE ANY STRAIN MEASURE. HERE EMPLOYING THE CAUCHY_GREEN
if (r_options.IsNot(ConstitutiveLaw::USE_ELEMENT_PROVIDED_STRAIN)) {
CalculateCauchyGreenStrain(rValues, r_strain_vector);
}
KRATOS_ERROR_IF(!rValues.IsSetStrainVector() || rValues.GetStrainVector().size() != GetStrainSize())
<< "Constitutive laws in the geomechanics application need a valid provided strain"
<< std::endl;

if (r_options.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
Matrix& r_constitutive_matrix = rValues.GetConstitutiveMatrix();
CalculateElasticMatrix(r_constitutive_matrix, rValues);
CalculateElasticMatrix(rValues.GetConstitutiveMatrix(), rValues);
}

if (r_options.Is(ConstitutiveLaw::COMPUTE_STRESS)) {
Vector& r_stress_vector = rValues.GetStressVector();
CalculatePK2Stress(r_strain_vector, r_stress_vector, rValues);
CalculatePK2Stress(rValues.GetStrainVector(), rValues.GetStressVector(), rValues);
}

KRATOS_CATCH("")
Expand All @@ -89,7 +87,6 @@ double& GeoLinearElasticLaw::CalculateValue(ConstitutiveLaw::Parameters& rParame
if (rThisVariable == STRAIN_ENERGY) {
Vector& r_strain_vector = rParameterValues.GetStrainVector();
Vector& r_stress_vector = rParameterValues.GetStressVector();
this->CalculateCauchyGreenStrain(rParameterValues, r_strain_vector);
this->CalculatePK2Stress(r_strain_vector, r_stress_vector, rParameterValues);

rValue = 0.5 * inner_prod(r_strain_vector, r_stress_vector); // Strain energy = 0.5*E:C:E
Expand All @@ -102,11 +99,8 @@ Vector& GeoLinearElasticLaw::CalculateValue(ConstitutiveLaw::Parameters& rParame
const Variable<Vector>& rThisVariable,
Vector& rValue)
{
if (rThisVariable == STRAIN || rThisVariable == GREEN_LAGRANGE_STRAIN_VECTOR ||
rThisVariable == ALMANSI_STRAIN_VECTOR) {
this->CalculateCauchyGreenStrain(rParameterValues, rValue);
} else if (rThisVariable == STRESSES || rThisVariable == CAUCHY_STRESS_VECTOR ||
rThisVariable == KIRCHHOFF_STRESS_VECTOR || rThisVariable == PK2_STRESS_VECTOR) {
if (rThisVariable == STRESSES || rThisVariable == CAUCHY_STRESS_VECTOR ||
rThisVariable == KIRCHHOFF_STRESS_VECTOR || rThisVariable == PK2_STRESS_VECTOR) {
// Get Values to compute the constitutive law:
Flags& r_flags = rParameterValues.GetOptions();

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoLinearElasticLaw : public Constit

protected:
virtual void CalculateElasticMatrix(Matrix& rConstitutiveMatrix, ConstitutiveLaw::Parameters& rValues) = 0;
virtual void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) = 0;
virtual void CalculatePK2Stress(const Vector& rStrainVector,
Vector& rStressVector,
ConstitutiveLaw::Parameters& rValues) = 0;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -132,24 +132,4 @@ void LinearPlaneStrainK0Law::CalculatePK2Stress(const Vector& rSt
KRATOS_CATCH("")
}

void LinearPlaneStrainK0Law::CalculateCauchyGreenStrain(Parameters& rValues, Vector& rStrainVector)
{
// 1.-Compute total deformation gradient
const Matrix& F = rValues.GetDeformationGradientF();

// for shells/membranes in case the DeformationGradient is of size 3x3
BoundedMatrix<double, 2, 2> F2x2;
for (unsigned int i = 0; i < 2; ++i)
for (unsigned int j = 0; j < 2; ++j)
F2x2(i, j) = F(i, j);

Matrix E_tensor = prod(trans(F2x2), F2x2);

for (unsigned int i = 0; i < 2; ++i)
E_tensor(i, i) -= 1.0;

E_tensor *= 0.5;
noalias(rStrainVector) = MathUtils<double>::StrainTensorToVector(E_tensor);
}

} // Namespace Kratos
Original file line number Diff line number Diff line change
Expand Up @@ -159,13 +159,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) LinearPlaneStrainK0Law : public GeoL
Vector& rStressVector,
ConstitutiveLaw::Parameters& rValues) override;

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}

private:
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -76,25 +76,4 @@ void GeoLinearElasticPlaneStress2DLaw::CalculatePK2Stress(const Vector& rStrainV
noalias(rStressVector) = prod(C, rStrainVector);
}

void GeoLinearElasticPlaneStress2DLaw::CalculateCauchyGreenStrain(Parameters& rValues, Vector& rStrainVector)
{
// 1.-Compute total deformation gradient
const Matrix& F = rValues.GetDeformationGradientF();

// for shells/membranes in case the DeformationGradient is of size 3x3
BoundedMatrix<double, 2, 2> F2x2;
for (unsigned int i = 0; i < 2; ++i)
for (unsigned int j = 0; j < 2; ++j)
F2x2(i, j) = F(i, j);

Matrix E_tensor = prod(trans(F2x2), F2x2);

for (unsigned int i = 0; i < 2; ++i)
E_tensor(i, i) -= 1.0;

E_tensor *= 0.5;

noalias(rStrainVector) = MathUtils<double>::StrainTensorToVector(E_tensor);
}

} // namespace Kratos
Original file line number Diff line number Diff line change
Expand Up @@ -152,13 +152,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoLinearElasticPlaneStress2DLaw : p
Vector& rStressVector,
ConstitutiveLaw::Parameters& rValues) override;

/**
* It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}

private:
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -90,13 +90,6 @@ indexStress3D SmallStrainUDSM2DInterfaceLaw::getIndex3D(const indexStress2DInter
}
}

void SmallStrainUDSM2DInterfaceLaw::CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues,
Vector& rStrainVector)
{
KRATOS_ERROR << "CalculateCauchyGreenStrain is not implemented in SmallStrainUDSM2DInterfaceLaw"
<< std::endl;
}

Vector& SmallStrainUDSM2DInterfaceLaw::GetValue(const Variable<Vector>& rThisVariable, Vector& rValue)
{
if (rThisVariable == STATE_VARIABLES) {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -100,13 +100,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SmallStrainUDSM2DInterfaceLaw : publ
*/
StressMeasure GetStressMeasure() override { return StressMeasure_Cauchy; }

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}
///@name Inquiry
///@{
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -84,27 +84,6 @@ void SmallStrainUDSM2DPlaneStrainLaw::CopyConstitutiveMatrix(ConstitutiveLaw::Pa
KRATOS_CATCH("")
}

void SmallStrainUDSM2DPlaneStrainLaw::CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues,
Vector& rStrainVector)
{
// 1.-Compute total deformation gradient
const Matrix& F = rValues.GetDeformationGradientF();

// for shells/membranes in case the DeformationGradient is of size 3x3
BoundedMatrix<double, 2, 2> F2x2;
for (unsigned int i = 0; i < 2; ++i)
for (unsigned int j = 0; j < 2; ++j)
F2x2(i, j) = F(i, j);

Matrix E_tensor = prod(trans(F2x2), F2x2);

for (unsigned int i = 0; i < 2; ++i)
E_tensor(i, i) -= 1.0;

E_tensor *= 0.5;
noalias(rStrainVector) = MathUtils<double>::StrainTensorToVector(E_tensor);
}

Vector& SmallStrainUDSM2DPlaneStrainLaw::GetValue(const Variable<Vector>& rThisVariable, Vector& rValue)
{
if (rThisVariable == STATE_VARIABLES) {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -100,13 +100,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SmallStrainUDSM2DPlaneStrainLaw : pu
*/
StressMeasure GetStressMeasure() override { return StressMeasure_Cauchy; }

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}
///@name Inquiry
///@{
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -101,13 +101,6 @@ indexStress3D SmallStrainUDSM3DInterfaceLaw::getIndex3D(const indexStress3DInter
}
}

void SmallStrainUDSM3DInterfaceLaw::CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues,
Vector& rStrainVector)
{
KRATOS_ERROR << "CalculateCauchyGreenStrain is not implemented in SmallStrainUDSM3DInterfaceLaw"
<< std::endl;
}

Vector& SmallStrainUDSM3DInterfaceLaw::GetValue(const Variable<Vector>& rThisVariable, Vector& rValue)
{
if (rThisVariable == STATE_VARIABLES) {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -99,13 +99,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SmallStrainUDSM3DInterfaceLaw : publ
*/
StressMeasure GetStressMeasure() override { return StressMeasure_Cauchy; }

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}
///@name Inquiry
///@{
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -588,11 +588,12 @@ void SmallStrainUDSM3DLaw::CalculateMaterialResponseCauchy(ConstitutiveLaw::Para
// Get Values to compute the constitutive law:
const Flags& rOptions = rValues.GetOptions();

// NOTE: SINCE THE ELEMENT IS IN SMALL STRAINS WE CAN USE ANY STRAIN MEASURE. HERE EMPLOYING THE CAUCHY_GREEN
if (rOptions.IsNot(ConstitutiveLaw::USE_ELEMENT_PROVIDED_STRAIN)) {
Vector& rStrainVector = rValues.GetStrainVector();
CalculateCauchyGreenStrain(rValues, rStrainVector);
}
KRATOS_DEBUG_ERROR_IF(rOptions.IsNot(ConstitutiveLaw::USE_ELEMENT_PROVIDED_STRAIN))
<< "The GeoLinearElasticLaw needs an element provided strain" << std::endl;

KRATOS_ERROR_IF(!rValues.IsSetStrainVector() || rValues.GetStrainVector().size() != GetStrainSize())
<< "Constitutive laws in the geomechanics application need a valid provided strain"
<< std::endl;

if (rOptions.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
// Constitutive matrix (D matrix)
Expand Down Expand Up @@ -815,31 +816,12 @@ void SmallStrainUDSM3DLaw::UpdateInternalStrainVectorFinalized(ConstitutiveLaw::
this->SetInternalStrainVector(rStrainVector);
}

void SmallStrainUDSM3DLaw::CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector)
{
const SizeType space_dimension = this->WorkingSpaceDimension();

//-Compute total deformation gradient
const Matrix& F = rValues.GetDeformationGradientF();
KRATOS_DEBUG_ERROR_IF(F.size1() != space_dimension || F.size2() != space_dimension)
<< "expected size of F " << space_dimension << "x" << space_dimension << ", got "
<< F.size1() << "x" << F.size2() << std::endl;

Matrix E_tensor = prod(trans(F), F);
for (unsigned int i = 0; i < space_dimension; ++i)
E_tensor(i, i) -= 1.0;
E_tensor *= 0.5;

noalias(rStrainVector) = MathUtils<double>::StrainTensorToVector(E_tensor);
}

double& SmallStrainUDSM3DLaw::CalculateValue(ConstitutiveLaw::Parameters& rParameterValues,
const Variable<double>& rThisVariable,
double& rValue)
{
if (rThisVariable == STRAIN_ENERGY) {
Vector& rStrainVector = rParameterValues.GetStrainVector();
this->CalculateCauchyGreenStrain(rParameterValues, rStrainVector);
Vector& rStressVector = rParameterValues.GetStressVector();
this->CalculateStress(rParameterValues, rStressVector);

Expand All @@ -852,11 +834,8 @@ Vector& SmallStrainUDSM3DLaw::CalculateValue(ConstitutiveLaw::Parameters& rParam
const Variable<Vector>& rThisVariable,
Vector& rValue)
{
if (rThisVariable == STRAIN || rThisVariable == GREEN_LAGRANGE_STRAIN_VECTOR ||
rThisVariable == ALMANSI_STRAIN_VECTOR) {
this->CalculateCauchyGreenStrain(rParameterValues, rValue);
} else if (rThisVariable == STRESSES || rThisVariable == CAUCHY_STRESS_VECTOR ||
rThisVariable == KIRCHHOFF_STRESS_VECTOR || rThisVariable == PK2_STRESS_VECTOR) {
if (rThisVariable == STRESSES || rThisVariable == CAUCHY_STRESS_VECTOR ||
rThisVariable == KIRCHHOFF_STRESS_VECTOR || rThisVariable == PK2_STRESS_VECTOR) {
// Get Values to compute the constitutive law:
Flags& rFlags = rParameterValues.GetOptions();

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -284,13 +284,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SmallStrainUDSM3DLaw : public Consti
void InitializeMaterialResponsePK2(ConstitutiveLaw::Parameters& rValues) override;
void InitializeMaterialResponseKirchhoff(ConstitutiveLaw::Parameters& rValues) override;

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
virtual void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector);

/**
* This can be used in order to reset all internal variables of the
* constitutive law (e.g. if a model should be reset to its reference state)
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -86,13 +86,6 @@ indexStress3D SmallStrainUMAT2DInterfaceLaw::getIndex3D(indexStress2DInterface i
}
}

void SmallStrainUMAT2DInterfaceLaw::CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues,
Vector& rStrainVector)
{
KRATOS_ERROR << "CalculateCauchyGreenStrain is not implemented in SmallStrainUMAT2DInterfaceLaw"
<< std::endl;
}

Vector& SmallStrainUMAT2DInterfaceLaw::GetValue(const Variable<Vector>& rThisVariable, Vector& rValue)
{
if (rThisVariable == STATE_VARIABLES) {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -98,13 +98,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SmallStrainUMAT2DInterfaceLaw : publ
*/
StressMeasure GetStressMeasure() override { return StressMeasure_Cauchy; }

/**
* @brief It calculates the strain vector
* @param rValues The internal values of the law
* @param rStrainVector The strain vector in Voigt notation
*/
void CalculateCauchyGreenStrain(ConstitutiveLaw::Parameters& rValues, Vector& rStrainVector) override;

///@}
///@name Inquiry
///@{
Expand Down
Loading

0 comments on commit 3e10ca1

Please sign in to comment.