KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
|
This class is the base class which define the Plastic Damage model developed by Luccioni B. and Oller S. More...
#include <generic_small_strain_plastic_damage_model.h>
Classes | |
struct | PlasticDamageParameters |
Public Member Functions | |
Life Cycle | |
GenericSmallStrainPlasticDamageModel () | |
ConstitutiveLaw::Pointer | Clone () const override |
GenericSmallStrainPlasticDamageModel (const GenericSmallStrainPlasticDamageModel &rOther) | |
~GenericSmallStrainPlasticDamageModel () override | |
Operations | |
void | CalculateMaterialResponsePK1 (ConstitutiveLaw::Parameters &rValues) override |
Computes the material response in terms of 1st Piola-Kirchhoff stresses and constitutive tensor. More... | |
void | CalculateMaterialResponsePK2 (ConstitutiveLaw::Parameters &rValues) override |
Computes the material response in terms of 2nd Piola-Kirchhoff stresses and constitutive tensor. More... | |
void | CalculateMaterialResponseKirchhoff (ConstitutiveLaw::Parameters &rValues) override |
Computes the material response in terms of Kirchhoff stresses and constitutive tensor. More... | |
void | CalculateMaterialResponseCauchy (ConstitutiveLaw::Parameters &rValues) override |
Computes the material response in terms of Cauchy stresses and constitutive tensor. More... | |
void | InitializeMaterial (const Properties &rMaterialProperties, const GeometryType &rElementGeometry, const Vector &rShapeFunctionsValues) override |
This is to be called at the very beginning of the calculation (e.g. from InitializeElement) in order to initialize all relevant attributes of the constitutive law. More... | |
void | InitializeMaterialResponsePK1 (ConstitutiveLaw::Parameters &rValues) override |
Initialize the material response in terms of 1st Piola-Kirchhoff stresses. More... | |
void | InitializeMaterialResponsePK2 (ConstitutiveLaw::Parameters &rValues) override |
Initialize the material response in terms of 2nd Piola-Kirchhoff stresses. More... | |
void | InitializeMaterialResponseKirchhoff (ConstitutiveLaw::Parameters &rValues) override |
Initialize the material response in terms of Kirchhoff stresses. More... | |
void | InitializeMaterialResponseCauchy (ConstitutiveLaw::Parameters &rValues) override |
Initialize the material response in terms of Cauchy stresses. More... | |
void | FinalizeMaterialResponsePK1 (ConstitutiveLaw::Parameters &rValues) override |
Finalize the material response in terms of 1st Piola-Kirchhoff stresses. More... | |
void | FinalizeMaterialResponsePK2 (ConstitutiveLaw::Parameters &rValues) override |
Finalize the material response in terms of 2nd Piola-Kirchhoff stresses. More... | |
void | FinalizeMaterialResponseKirchhoff (ConstitutiveLaw::Parameters &rValues) override |
Finalize the material response in terms of Kirchhoff stresses. More... | |
void | FinalizeMaterialResponseCauchy (ConstitutiveLaw::Parameters &rValues) override |
bool | Has (const Variable< double > &rThisVariable) override |
Returns whether this constitutive Law has specified variable (double) More... | |
bool | Has (const Variable< Vector > &rThisVariable) override |
Returns whether this constitutive Law has specified variable (Vector) More... | |
bool | Has (const Variable< Matrix > &rThisVariable) override |
Returns whether this constitutive Law has specified variable (Matrix) More... | |
void | SetValue (const Variable< double > &rThisVariable, const double &rValue, const ProcessInfo &rCurrentProcessInfo) override |
Sets the value of a specified variable (double) More... | |
void | SetValue (const Variable< Vector > &rThisVariable, const Vector &rValue, const ProcessInfo &rCurrentProcessInfo) override |
Sets the value of a specified variable (Vector) More... | |
double & | GetValue (const Variable< double > &rThisVariable, double &rValue) override |
Returns the value of a specified variable (double) More... | |
Vector & | GetValue (const Variable< Vector > &rThisVariable, Vector &rValue) override |
Returns the value of a specified variable (Vector) More... | |
Matrix & | GetValue (const Variable< Matrix > &rThisVariable, Matrix &rValue) override |
Returns the value of a specified variable (matrix) More... | |
bool | RequiresInitializeMaterialResponse () override |
If the CL requires to initialize the material response, called by the element in InitializeSolutionStep. More... | |
bool | RequiresFinalizeMaterialResponse () override |
If the CL requires to initialize the material response, called by the element in InitializeSolutionStep. More... | |
double & | CalculateValue (ConstitutiveLaw::Parameters &rParameterValues, const Variable< double > &rThisVariable, double &rValue) override |
Returns the value of a specified variable (double) More... | |
Vector & | CalculateValue (ConstitutiveLaw::Parameters &rParameterValues, const Variable< Vector > &rThisVariable, Vector &rValue) override |
Returns the value of a specified variable (vector) More... | |
Matrix & | CalculateValue (ConstitutiveLaw::Parameters &rParameterValues, const Variable< Matrix > &rThisVariable, Matrix &rValue) override |
Returns the value of a specified variable (matrix) More... | |
int | Check (const Properties &rMaterialProperties, const GeometryType &rElementGeometry, const ProcessInfo &rCurrentProcessInfo) const override |
This function provides the place to perform checks on the completeness of the input. More... | |
double | CalculateDamageParameters (PlasticDamageParameters &rParameters, const Matrix &rElasticMatrix, ConstitutiveLaw::Parameters &rValues) |
This method works as the damage integrator im the isotropic damage CL and computes the associated parameters. More... | |
void | CalculateIndicatorsFactors (const array_1d< double, 6 > &rPredictiveStressVector, double &rTensileIndicatorFactor, double &rCompressionIndicatorFactor, double &rSumPrincipalStresses, array_1d< double, 3 > &rPrincipalStresses) |
This method computes the tensile/compressive indicators. More... | |
void | CheckInternalVariable (double &rInternalVariable) |
This method checks the value of some of the variables used in the CL and guarantees that their value is between 0 and 1. More... | |
void | CalculateIncrementsPlasticDamageCase (PlasticDamageParameters &rParameters, const Matrix &rElasticMatrix) |
This method computes the increments for the damage internal variable and for the plastic consistency parameter. More... | |
double | CalculatePlasticParameters (PlasticDamageParameters &rParameters, const Matrix &rConstitutiveMatrix, ConstitutiveLaw::Parameters &rValues) |
This method works as the damage integrator im the isotropic damage CL and computes the associated parameters. More... | |
void | CalculatePlasticDenominator (const array_1d< double, VoigtSize > &rFFlux, const array_1d< double, VoigtSize > &rGFlux, const Matrix &rConstitutiveMatrix, double &rHardeningParameter, const double Damage, double &rPlasticDenominator) |
This method computes the plastic denominator needed to obtain the plastic consistency factor to compute the plastic consistency factor. More... | |
Protected Member Functions | |
Protected Operators | |
double & | GetThresholdPlasticity () |
double & | GetPlasticDissipation () |
Vector & | GetPlasticStrain () |
void | SetThresholdPlasticity (const double ThresholdPlasticity) |
void | SetPlasticDissipation (const double PlasticDissipation) |
void | SetPlasticStrain (const array_1d< double, VoigtSize > &rPlasticStrain) |
Type Definitions | |
typedef std::size_t | IndexType |
typedef std::conditional< VoigtSize==6, ElasticIsotropic3D, LinearPlaneStrain >::type | BaseType |
Definition of the base class. More... | |
typedef array_1d< double, VoigtSize > | BoundedArrayType |
The definition of the Voigt array type. More... | |
typedef BoundedMatrix< double, Dimension, Dimension > | BoundedMatrixType |
The definition of the bounded matrix type. More... | |
typedef Node | NodeType |
The node definition. More... | |
typedef Geometry< NodeType > | GeometryType |
The geometry definition. More... | |
static constexpr SizeType | Dimension = TPlasticityIntegratorType::Dimension |
The define the working dimension size, already defined in the integrator. More... | |
static constexpr SizeType | VoigtSize = TPlasticityIntegratorType::VoigtSize |
The define the Voigt size, already defined in the integrator. More... | |
static constexpr double | tolerance = std::numeric_limits<double>::epsilon() |
The machine precision tolerance. More... | |
KRATOS_CLASS_POINTER_DEFINITION (GenericSmallStrainPlasticDamageModel) | |
Counted pointer of GenericSmallStrainPlasticDamageModel. More... | |
Un accessible methods | |
class | Serializer |
This class is the base class which define the Plastic Damage model developed by Luccioni B. and Oller S.
This class considers a constitutive law integrator for the damage and another one for the plasticity process
TPlasticityIntegratorType | and TDamageIntegratorType The constitutive law integrators considered |
typedef std::conditional<VoigtSize == 6, ElasticIsotropic3D, LinearPlaneStrain >::type Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::BaseType |
Definition of the base class.
typedef array_1d<double, VoigtSize> Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::BoundedArrayType |
The definition of the Voigt array type.
typedef BoundedMatrix<double, Dimension, Dimension> Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::BoundedMatrixType |
The definition of the bounded matrix type.
typedef Geometry<NodeType> Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::GeometryType |
The geometry definition.
typedef std::size_t Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::IndexType |
typedef Node Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::NodeType |
The node definition.
|
inline |
Default constructor.
|
inline |
Copy constructor.
|
inlineoverride |
Destructor.
double Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::CalculateDamageParameters | ( | PlasticDamageParameters & | rParameters, |
const Matrix & | rElasticMatrix, | ||
ConstitutiveLaw::Parameters & | rValues | ||
) |
This method works as the damage integrator im the isotropic damage CL and computes the associated parameters.
It is used to call the damage integrator and apply the particularities of this CL.
rPredictiveStressVector | The predictive stress vector S = C:(E-Ep) |
rStrainVector | The strain vector |
rUniaxialStress | The uniaxial stress of the damage model |
rDamageThreshold | The maximum uniaxial stress achieved previously |
rDamageDissipation | The internal variable of energy dissipation due to damage |
rConstitutiveMatrix | The elastic constitutive matrix |
rValues | Parameters of the constitutive law |
CharacteristicLength | The equivalent length of the FE |
rDamageFlux | The derivative of the yield surface used for damage |
rPlasticStrain | The plastic component of the strain |
Damage | The internal variable of the damage model |
DamageIncrement | The increment of the internal variable of the damage model at this time step |
UndamagedFreeEnergy | The undamaged free energy |
rHardd | Hardening parameter for the damage model |
rDamageDissipationIncrement | Increment of the internal variable of energy dissipation due to damage |
void Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::CalculateIncrementsPlasticDamageCase | ( | PlasticDamageParameters & | rParameters, |
const Matrix & | rElasticMatrix | ||
) |
This method computes the increments for the damage internal variable and for the plastic consistency parameter.
rFluxDamageYield | The derivative of the yield surface used for damage |
rStrainVector | The strain vector |
Damage | The internal variable of the damage model |
rPlasticityFlux | The derivative of the yield surface used for plasticity |
rPlasticityGFlux | The derivative of the potential used for plasticity |
rElasticMatrix | The elastic constitutive matrix |
DamageIndicator | The difference between the uniaxial stress and the damage threshold |
PlasticityIndicator | The difference between the uniaxial stress and the plasticity threshold |
rPlasticStrain | The plastic component of the strain |
rDamageIncrement | The increment of the internal variable of the damage model at this time step |
rPlasticConsistencyIncrement | The increment of the internal plastic consistency variable of the plasticity model at this time step |
PlasticDenominator | The plasticity numerical value to obtain the pastic consistency factor |
UniaxialStressPlast | The uniaxial stress of the plasticity model |
Hardd | Hardening parameter for the damage model |
rDamageDissipationIncrement | Increment of the internal variable of energy dissipation due to damage |
void Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::CalculateIndicatorsFactors | ( | const array_1d< double, 6 > & | rPredictiveStressVector, |
double & | rTensileIndicatorFactor, | ||
double & | rCompressionIndicatorFactor, | ||
double & | rSumPrincipalStresses, | ||
array_1d< double, 3 > & | rPrincipalStresses | ||
) |
This method computes the tensile/compressive indicators.
rPredictiveStressVector | The predictive stress vector S = C:(E-Ep) |
rTensileIndicatorFactor | The tensile indicator |
rCompressionIndicatorFactor | The compressive indicator |
rSumPrincipalStresses | The sum of the principal stresses |
rPrincipalStresses | The principal stresses |
|
override |
Computes the material response in terms of Cauchy stresses and constitutive tensor.
|
override |
Computes the material response in terms of Kirchhoff stresses and constitutive tensor.
|
override |
Computes the material response in terms of 1st Piola-Kirchhoff stresses and constitutive tensor.
|
override |
Computes the material response in terms of 2nd Piola-Kirchhoff stresses and constitutive tensor.
void Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::CalculatePlasticDenominator | ( | const array_1d< double, VoigtSize > & | rFFlux, |
const array_1d< double, VoigtSize > & | rGFlux, | ||
const Matrix & | rConstitutiveMatrix, | ||
double & | rHardeningParameter, | ||
const double | Damage, | ||
double & | rPlasticDenominator | ||
) |
This method computes the plastic denominator needed to obtain the plastic consistency factor to compute the plastic consistency factor.
rFflux | The derivative of the yield surface |
rGflux | The derivative of the plastic potential |
rConstitutiveMatrix | The elastic constitutive matrix |
rHardeningParameter | The hardening parameter needed for the plasticity algorithm |
rPlasticDenominator | The plasticity numerical value to obtain the plastic consistency factor |
double Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::CalculatePlasticParameters | ( | PlasticDamageParameters & | rParameters, |
const Matrix & | rConstitutiveMatrix, | ||
ConstitutiveLaw::Parameters & | rValues | ||
) |
This method works as the damage integrator im the isotropic damage CL and computes the associated parameters.
rPredictiveStressVector | The predictive stress vector S = C:(E-Ep) |
rStrainVector | The strain vector |
rUniaxialStress | The uniaxial stress of the plasticity model |
rThreshold | The maximum uniaxial stress achieved previously |
rPlasticDenominator | The plasticity numerical value to obtain the plastic consistency factor |
rFflux | The derivative of the yield surface used for plasticity |
rGflux | The derivative of the potential used for plasticity |
rPlasticDissipation | The internal variable of energy dissipation due to plasticity |
rPlasticStrainIncrement | The increment of the plastic strain |
rConstitutiveMatrix | The elastic constitutive matrix |
rValues | Parameters of the constitutive law |
CharacteristicLength | The equivalent length of the FE |
rPlasticStrain | The plastic component of the strain |
Damage | The internal variable of the damage model |
|
override |
Returns the value of a specified variable (double)
rParameterValues | the needed parameters for the CL calculation |
rThisVariable | the variable to be returned |
rValue | a reference to the returned value |
rValue | output: the value of the specified variable |
|
override |
Returns the value of a specified variable (matrix)
rParameterValues | the needed parameters for the CL calculation |
rThisVariable | the variable to be returned |
rValue | a reference to the returned value |
rValue | output: the value of the specified variable |
|
override |
Returns the value of a specified variable (vector)
rParameterValues | the needed parameters for the CL calculation |
rThisVariable | the variable to be returned |
rValue | a reference to the returned value |
rValue | output: the value of the specified variable |
|
override |
This function provides the place to perform checks on the completeness of the input.
It is designed to be called only once (or anyway, not often) typically at the beginning of the calculations, so to verify that nothing is missing from the input or that no common error is found.
rMaterialProperties | The properties of the material |
rElementGeometry | The geometry of the element |
rCurrentProcessInfo | The current process info instance |
void Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::CheckInternalVariable | ( | double & | rInternalVariable | ) |
This method checks the value of some of the variables used in the CL and guarantees that their value is between 0 and 1.
rInternalVariable | Internal variable of the CL (damage, rDamageDissipationIncrement and rDamageDissipation) |
|
inlineoverride |
Clone.
|
override |
Finalize the material response in terms of Cauchy stresses
|
override |
Finalize the material response in terms of Kirchhoff stresses.
|
override |
Finalize the material response in terms of 1st Piola-Kirchhoff stresses.
|
override |
Finalize the material response in terms of 2nd Piola-Kirchhoff stresses.
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
override |
Returns the value of a specified variable (double)
rThisVariable | the variable to be returned |
rValue | a reference to the returned value |
|
override |
Returns the value of a specified variable (matrix)
rThisVariable | the variable to be returned |
rValue | a reference to the returned value |
|
override |
Returns the value of a specified variable (Vector)
rThisVariable | the variable to be returned |
rValue | a reference to the returned value |
|
override |
Returns whether this constitutive Law has specified variable (double)
rThisVariable | the variable to be checked for |
|
override |
Returns whether this constitutive Law has specified variable (Matrix)
rThisVariable | the variable to be checked for |
|
override |
Returns whether this constitutive Law has specified variable (Vector)
rThisVariable | the variable to be checked for |
|
override |
This is to be called at the very beginning of the calculation (e.g. from InitializeElement) in order to initialize all relevant attributes of the constitutive law.
rMaterialProperties | the Properties instance of the current element |
rElementGeometry | the geometry of the current element |
rShapeFunctionsValues | the shape functions values in the current integration point |
|
override |
Initialize the material response in terms of Cauchy stresses.
|
override |
Initialize the material response in terms of Kirchhoff stresses.
|
override |
Initialize the material response in terms of 1st Piola-Kirchhoff stresses.
|
override |
Initialize the material response in terms of 2nd Piola-Kirchhoff stresses.
Kratos::GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType >::KRATOS_CLASS_POINTER_DEFINITION | ( | GenericSmallStrainPlasticDamageModel< TPlasticityIntegratorType, TDamageIntegratorType > | ) |
Counted pointer of GenericSmallStrainPlasticDamageModel.
|
inlineoverride |
If the CL requires to initialize the material response, called by the element in InitializeSolutionStep.
|
inlineoverride |
If the CL requires to initialize the material response, called by the element in InitializeSolutionStep.
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
override |
Sets the value of a specified variable (double)
rVariable | the variable to be returned |
rValue | new value of the specified variable |
rCurrentProcessInfo | the process info |
|
override |
Sets the value of a specified variable (Vector)
rThisVariable | the variable to be returned |
rValue | new value of the specified variable |
rCurrentProcessInfo | the process info |
|
friend |
|
staticconstexpr |
The define the working dimension size, already defined in the integrator.
|
staticconstexpr |
The machine precision tolerance.
|
staticconstexpr |
The define the Voigt size, already defined in the integrator.