https://mooseframework.inl.gov
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Static Public Attributes | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Types | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
SolidMechanicsPlasticMohrCoulombMulti Class Reference

FiniteStrainMohrCoulombMulti implements rate-independent non-associative mohr-coulomb with hardening/softening in the finite-strain framework, using planar (non-smoothed) surfaces. More...

#include <SolidMechanicsPlasticMohrCoulombMulti.h>

Inheritance diagram for SolidMechanicsPlasticMohrCoulombMulti:
[legend]

Public Types

typedef DataFileName DataFileParameterType
 

Public Member Functions

 SolidMechanicsPlasticMohrCoulombMulti (const InputParameters &parameters)
 
virtual unsigned int numberSurfaces () const override
 The number of yield surfaces for this plasticity model.
 
virtual void yieldFunctionV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &f) const override
 Calculates the yield functions.
 
virtual void dyieldFunction_dstressV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &df_dstress) const override
 The derivative of yield functions with respect to stress.
 
virtual void dyieldFunction_dintnlV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &df_dintnl) const override
 The derivative of yield functions with respect to the internal parameter.
 
virtual void flowPotentialV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &r) const override
 The flow potentials.
 
virtual void dflowPotential_dstressV (const RankTwoTensor &stress, Real intnl, std::vector< RankFourTensor > &dr_dstress) const override
 The derivative of the flow potential with respect to stress.
 
virtual void dflowPotential_dintnlV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &dr_dintnl) const override
 The derivative of the flow potential with respect to the internal parameter.
 
virtual void activeConstraints (const std::vector< Real > &f, const RankTwoTensor &stress, Real intnl, const RankFourTensor &Eijkl, std::vector< bool > &act, RankTwoTensor &returned_stress) const override
 The active yield surfaces, given a vector of yield functions.
 
virtual std::string modelName () const override
 
virtual bool useCustomReturnMap () const override
 Returns false. You will want to override this in your derived class if you write a custom returnMap function.
 
virtual bool returnMap (const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const override
 Performs a custom return-map.
 
void initialize ()
 
void execute ()
 
void finalize ()
 
virtual void hardPotentialV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &h) const
 The hardening potential.
 
virtual void dhardPotential_dstressV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &dh_dstress) const
 The derivative of the hardening potential with respect to stress.
 
virtual void dhardPotential_dintnlV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &dh_dintnl) const
 The derivative of the hardening potential with respect to the internal parameter.
 
virtual bool useCustomCTO () const
 Returns false. You will want to override this in your derived class if you write a custom consistent tangent operator function.
 
virtual RankFourTensor consistentTangentOperator (const RankTwoTensor &trial_stress, Real intnl_old, const RankTwoTensor &stress, Real intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &cumulative_pm) const
 Calculates a custom consistent tangent operator.
 
bool KuhnTuckerSingleSurface (Real yf, Real dpm, Real dpm_tol) const
 Returns true if the Kuhn-Tucker conditions for the single surface are satisfied.
 
virtual Real spatialValue (const Point &) const
 
virtual const std::vector< Point > spatialPoints () const
 
void setPrimaryThreadCopy (UserObject *primary)
 
UserObjectprimaryThreadCopy ()
 
SubProblemgetSubProblem () const
 
bool shouldDuplicateInitialExecution () const
 
void gatherSum (T &value)
 
void gatherMax (T &value)
 
void gatherMin (T &value)
 
void gatherProxyValueMax (T1 &proxy, T2 &value)
 
void gatherProxyValueMin (T1 &proxy, T2 &value)
 
std::set< UserObjectName > getDependObjects () const
 
const std::set< std::string > & getRequestedItems () override
 
const std::set< std::string > & getSuppliedItems () override
 
unsigned int systemNumber () const
 
virtual bool needThreadedCopy () const
 
virtual bool enabled () const
 
std::shared_ptr< MooseObjectgetSharedPtr ()
 
std::shared_ptr< const MooseObjectgetSharedPtr () const
 
bool isKokkosObject () const
 
MooseAppgetMooseApp () const
 
const std::string & type () const
 
const std::string & name () const
 
std::string typeAndName () const
 
MooseObjectParameterName uniqueParameterName (const std::string &parameter_name) const
 
MooseObjectName uniqueName () const
 
const InputParametersparameters () const
 
const hit::Node * getHitNode () const
 
bool hasBase () const
 
const std::string & getBase () const
 
const TgetParam (const std::string &name) const
 
std::vector< std::pair< T1, T2 > > getParam (const std::string &param1, const std::string &param2) const
 
const TqueryParam (const std::string &name) const
 
const TgetRenamedParam (const std::string &old_name, const std::string &new_name) const
 
T getCheckedPointerParam (const std::string &name, const std::string &error_string="") const
 
bool isParamValid (const std::string &name) const
 
bool isParamSetByUser (const std::string &name) const
 
void connectControllableParams (const std::string &parameter, const std::string &object_type, const std::string &object_name, const std::string &object_parameter) const
 
void paramError (const std::string &param, Args... args) const
 
void paramWarning (const std::string &param, Args... args) const
 
void paramWarning (const std::string &param, Args... args) const
 
void paramInfo (const std::string &param, Args... args) const
 
std::string messagePrefix (const bool hit_prefix=true) const
 
std::string errorPrefix (const std::string &) const
 
void mooseError (Args &&... args) const
 
void mooseDocumentedError (const std::string &repo_name, const unsigned int issue_num, Args &&... args) const
 
void mooseErrorNonPrefixed (Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseDeprecatedNoTrace (Args &&... args) const
 
void mooseInfo (Args &&... args) const
 
void callMooseError (std::string msg, const bool with_prefix, const hit::Node *node=nullptr, const bool show_trace=true) const
 
std::string getDataFileName (const std::string &param) const
 
std::string getDataFileNameByName (const std::string &relative_path) const
 
std::string getDataFilePath (const std::string &relative_path) const
 
virtual void initialSetup ()
 
virtual void timestepSetup ()
 
virtual void jacobianSetup ()
 
virtual void residualSetup ()
 
virtual void customSetup (const ExecFlagType &)
 
const ExecFlagEnumgetExecuteOnEnum () const
 
UserObjectName getUserObjectName (const std::string &param_name) const
 
const TgetUserObject (const std::string &param_name, bool is_dependency=true) const
 
const TgetUserObjectByName (const UserObjectName &object_name, bool is_dependency=true) const
 
const UserObjectBasegetUserObjectBase (const std::string &param_name, bool is_dependency=true) const
 
const UserObjectBasegetUserObjectBaseByName (const UserObjectName &object_name, bool is_dependency=true) const
 
const std::vector< MooseVariableScalar * > & getCoupledMooseScalarVars ()
 
const std::set< TagID > & getScalarVariableCoupleableVectorTags () const
 
const std::set< TagID > & getScalarVariableCoupleableMatrixTags () const
 
const GenericMaterialProperty< T, is_ad > & getGenericMaterialProperty (const std::string &name, MaterialData &material_data, const unsigned int state=0)
 
const GenericMaterialProperty< T, is_ad > & getGenericMaterialProperty (const std::string &name, const unsigned int state=0)
 
const GenericMaterialProperty< T, is_ad > & getGenericMaterialProperty (const std::string &name, const unsigned int state=0)
 
const MaterialProperty< T > & getMaterialProperty (const std::string &name, MaterialData &material_data, const unsigned int state=0)
 
const MaterialProperty< T > & getMaterialProperty (const std::string &name, const unsigned int state=0)
 
const MaterialProperty< T > & getMaterialProperty (const std::string &name, const unsigned int state=0)
 
const ADMaterialProperty< T > & getADMaterialProperty (const std::string &name, MaterialData &material_data)
 
const ADMaterialProperty< T > & getADMaterialProperty (const std::string &name)
 
const ADMaterialProperty< T > & getADMaterialProperty (const std::string &name)
 
const MaterialProperty< T > & getMaterialPropertyOld (const std::string &name, MaterialData &material_data)
 
const MaterialProperty< T > & getMaterialPropertyOld (const std::string &name)
 
const MaterialProperty< T > & getMaterialPropertyOld (const std::string &name)
 
const MaterialProperty< T > & getMaterialPropertyOlder (const std::string &name, MaterialData &material_data)
 
const MaterialProperty< T > & getMaterialPropertyOlder (const std::string &name)
 
const MaterialProperty< T > & getMaterialPropertyOlder (const std::string &name)
 
const GenericMaterialProperty< T, is_ad > & getGenericMaterialPropertyByName (const MaterialPropertyName &name, MaterialData &material_data, const unsigned int state)
 
const GenericMaterialProperty< T, is_ad > & getGenericMaterialPropertyByName (const MaterialPropertyName &name, const unsigned int state=0)
 
const GenericMaterialProperty< T, is_ad > & getGenericMaterialPropertyByName (const MaterialPropertyName &name, const unsigned int state=0)
 
const MaterialProperty< T > & getMaterialPropertyByName (const MaterialPropertyName &name, MaterialData &material_data, const unsigned int state=0)
 
const MaterialProperty< T > & getMaterialPropertyByName (const MaterialPropertyName &name, const unsigned int state=0)
 
const MaterialProperty< T > & getMaterialPropertyByName (const MaterialPropertyName &name, const unsigned int state=0)
 
const ADMaterialProperty< T > & getADMaterialPropertyByName (const MaterialPropertyName &name, MaterialData &material_data)
 
const ADMaterialProperty< T > & getADMaterialPropertyByName (const MaterialPropertyName &name)
 
const ADMaterialProperty< T > & getADMaterialPropertyByName (const MaterialPropertyName &name)
 
const MaterialProperty< T > & getMaterialPropertyOldByName (const MaterialPropertyName &name, MaterialData &material_data)
 
const MaterialProperty< T > & getMaterialPropertyOldByName (const MaterialPropertyName &name)
 
const MaterialProperty< T > & getMaterialPropertyOldByName (const MaterialPropertyName &name)
 
const MaterialProperty< T > & getMaterialPropertyOlderByName (const MaterialPropertyName &name, MaterialData &material_data)
 
const MaterialProperty< T > & getMaterialPropertyOlderByName (const MaterialPropertyName &name)
 
const MaterialProperty< T > & getMaterialPropertyOlderByName (const MaterialPropertyName &name)
 
Moose::Kokkos::MaterialProperty< T, dimension > getKokkosMaterialPropertyByName (const std::string &prop_name_in)
 
Moose::Kokkos::MaterialProperty< T, dimension > getKokkosMaterialPropertyOldByName (const std::string &prop_name)
 
Moose::Kokkos::MaterialProperty< T, dimension > getKokkosMaterialPropertyOlderByName (const std::string &prop_name)
 
Moose::Kokkos::MaterialProperty< T, dimension > getKokkosMaterialProperty (const std::string &name)
 
Moose::Kokkos::MaterialProperty< T, dimension > getKokkosMaterialPropertyOld (const std::string &name)
 
Moose::Kokkos::MaterialProperty< T, dimension > getKokkosMaterialPropertyOlder (const std::string &name)
 
std::pair< const MaterialProperty< T > *, std::set< SubdomainID > > getBlockMaterialProperty (const MaterialPropertyName &name)
 
std::pair< Moose::Kokkos::MaterialProperty< T, dimension >, std::set< SubdomainID > > getKokkosBlockMaterialProperty (const MaterialPropertyName &name)
 
const GenericMaterialProperty< T, is_ad > & getGenericZeroMaterialProperty (const std::string &name)
 
const GenericMaterialProperty< T, is_ad > & getGenericZeroMaterialProperty ()
 
const GenericMaterialProperty< T, is_ad > & getGenericZeroMaterialPropertyByName (const std::string &prop_name)
 
const MaterialProperty< T > & getZeroMaterialProperty (Ts... args)
 
std::set< SubdomainIDgetMaterialPropertyBlocks (const std::string &name)
 
std::vector< SubdomainName > getMaterialPropertyBlockNames (const std::string &name)
 
std::set< BoundaryIDgetMaterialPropertyBoundaryIDs (const std::string &name)
 
std::vector< BoundaryName > getMaterialPropertyBoundaryNames (const std::string &name)
 
void checkBlockAndBoundaryCompatibility (std::shared_ptr< MaterialBase > discrete)
 
std::unordered_map< SubdomainID, std::vector< MaterialBase * > > buildRequiredMaterials (bool allow_stateful=true)
 
void statefulPropertiesAllowed (bool)
 
virtual bool getMaterialPropertyCalled () const
 
virtual const std::unordered_set< unsigned int > & getMatPropDependencies () const
 
virtual void resolveOptionalProperties ()
 
const GenericMaterialProperty< T, is_ad > & getPossiblyConstantGenericMaterialPropertyByName (const MaterialPropertyName &prop_name, MaterialData &material_data, const unsigned int state)
 
bool isImplicit ()
 
Moose::StateArg determineState () const
 
virtual void threadJoin (const UserObject &) override
 
virtual void threadJoin (const UserObject &) override
 
virtual void subdomainSetup () override
 
virtual void subdomainSetup () override
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObject (const std::string &param_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
bool hasUserObjectByName (const UserObjectName &object_name) const
 
const GenericOptionalMaterialProperty< T, is_ad > & getGenericOptionalMaterialProperty (const std::string &name, const unsigned int state=0)
 
const GenericOptionalMaterialProperty< T, is_ad > & getGenericOptionalMaterialProperty (const std::string &name, const unsigned int state=0)
 
const OptionalMaterialProperty< T > & getOptionalMaterialProperty (const std::string &name, const unsigned int state=0)
 
const OptionalMaterialProperty< T > & getOptionalMaterialProperty (const std::string &name, const unsigned int state=0)
 
const OptionalADMaterialProperty< T > & getOptionalADMaterialProperty (const std::string &name)
 
const OptionalADMaterialProperty< T > & getOptionalADMaterialProperty (const std::string &name)
 
const OptionalMaterialProperty< T > & getOptionalMaterialPropertyOld (const std::string &name)
 
const OptionalMaterialProperty< T > & getOptionalMaterialPropertyOld (const std::string &name)
 
const OptionalMaterialProperty< T > & getOptionalMaterialPropertyOlder (const std::string &name)
 
const OptionalMaterialProperty< T > & getOptionalMaterialPropertyOlder (const std::string &name)
 
MaterialBasegetMaterial (const std::string &name)
 
MaterialBasegetMaterial (const std::string &name)
 
MaterialBasegetMaterialByName (const std::string &name, bool no_warn=false)
 
MaterialBasegetMaterialByName (const std::string &name, bool no_warn=false)
 
bool hasMaterialProperty (const std::string &name)
 
bool hasMaterialProperty (const std::string &name)
 
bool hasMaterialPropertyByName (const std::string &name)
 
bool hasMaterialPropertyByName (const std::string &name)
 
bool hasADMaterialProperty (const std::string &name)
 
bool hasADMaterialProperty (const std::string &name)
 
bool hasADMaterialPropertyByName (const std::string &name)
 
bool hasADMaterialPropertyByName (const std::string &name)
 
bool hasKokkosMaterialProperty (const std::string &name)
 
bool hasKokkosMaterialProperty (const std::string &name)
 
bool hasKokkosMaterialPropertyByName (const std::string &name)
 
bool hasKokkosMaterialPropertyByName (const std::string &name)
 
bool hasGenericMaterialProperty (const std::string &name)
 
bool hasGenericMaterialProperty (const std::string &name)
 
bool hasGenericMaterialPropertyByName (const std::string &name)
 
bool hasGenericMaterialPropertyByName (const std::string &name)
 
const FunctiongetFunction (const std::string &name) const
 
const FunctiongetFunctionByName (const FunctionName &name) const
 
bool hasFunction (const std::string &param_name) const
 
bool hasFunctionByName (const FunctionName &name) const
 
Moose::Kokkos::Function getKokkosFunction (const std::string &name) const
 
const TgetKokkosFunction (const std::string &name) const
 
Moose::Kokkos::Function getKokkosFunctionByName (const FunctionName &name) const
 
const TgetKokkosFunctionByName (const FunctionName &name) const
 
bool hasKokkosFunction (const std::string &param_name) const
 
bool hasKokkosFunctionByName (const FunctionName &name) const
 
bool isDefaultPostprocessorValue (const std::string &param_name, const unsigned int index=0) const
 
bool hasPostprocessor (const std::string &param_name, const unsigned int index=0) const
 
bool hasPostprocessorByName (const PostprocessorName &name) const
 
std::size_t coupledPostprocessors (const std::string &param_name) const
 
const PostprocessorName & getPostprocessorName (const std::string &param_name, const unsigned int index=0) const
 
const VectorPostprocessorValuegetVectorPostprocessorValue (const std::string &param_name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValue (const std::string &param_name, const std::string &vector_name, bool needs_broadcast) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueByName (const VectorPostprocessorName &name, const std::string &vector_name, bool needs_broadcast) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOld (const std::string &param_name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOld (const std::string &param_name, const std::string &vector_name, bool needs_broadcast) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOldByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
const VectorPostprocessorValuegetVectorPostprocessorValueOldByName (const VectorPostprocessorName &name, const std::string &vector_name, bool needs_broadcast) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValue (const std::string &param_name, const std::string &vector_name) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValueByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValueOld (const std::string &param_name, const std::string &vector_name) const
 
const ScatterVectorPostprocessorValuegetScatterVectorPostprocessorValueOldByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
bool hasVectorPostprocessor (const std::string &param_name, const std::string &vector_name) const
 
bool hasVectorPostprocessor (const std::string &param_name) const
 
bool hasVectorPostprocessorByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
bool hasVectorPostprocessorByName (const VectorPostprocessorName &name) const
 
const VectorPostprocessorName & getVectorPostprocessorName (const std::string &param_name) const
 
TgetSampler (const std::string &name)
 
SamplergetSampler (const std::string &name)
 
TgetSamplerByName (const SamplerName &name)
 
SamplergetSamplerByName (const SamplerName &name)
 
virtual void meshChanged ()
 
virtual void meshDisplaced ()
 
PerfGraphperfGraph ()
 
const PostprocessorValuegetPostprocessorValue (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValue (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOld (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOld (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOlder (const std::string &param_name, const unsigned int index=0) const
 
const PostprocessorValuegetPostprocessorValueOlder (const std::string &param_name, const unsigned int index=0) const
 
virtual const PostprocessorValuegetPostprocessorValueByName (const PostprocessorName &name) const
 
virtual const PostprocessorValuegetPostprocessorValueByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOldByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOldByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOlderByName (const PostprocessorName &name) const
 
const PostprocessorValuegetPostprocessorValueOlderByName (const PostprocessorName &name) const
 
bool isVectorPostprocessorDistributed (const std::string &param_name) const
 
bool isVectorPostprocessorDistributed (const std::string &param_name) const
 
bool isVectorPostprocessorDistributedByName (const VectorPostprocessorName &name) const
 
bool isVectorPostprocessorDistributedByName (const VectorPostprocessorName &name) const
 
const DistributiongetDistribution (const std::string &name) const
 
const TgetDistribution (const std::string &name) const
 
const DistributiongetDistribution (const std::string &name) const
 
const TgetDistribution (const std::string &name) const
 
const DistributiongetDistributionByName (const DistributionName &name) const
 
const TgetDistributionByName (const std::string &name) const
 
const DistributiongetDistributionByName (const DistributionName &name) const
 
const TgetDistributionByName (const std::string &name) const
 
const Parallel::Communicator & comm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static InputParameters validParams ()
 
static void callMooseError (MooseApp *const app, const InputParameters &params, std::string msg, const bool with_prefix, const hit::Node *node, const bool show_trace=true)
 
static void sort (typename std::vector< T > &vector)
 
static void sortDFS (typename std::vector< T > &vector)
 
static void cyclicDependencyError (CyclicDependencyException< T2 > &e, const std::string &header, NameFunc &&name_func)
 
static void cyclicDependencyError (CyclicDependencyException< T2 > &e, const std::string &header)
 

Public Attributes

const Real _f_tol
 Tolerance on yield function.
 
const Real _ic_tol
 Tolerance on internal constraint.
 
 usingCombinedWarningSolutionWarnings
 
const ConsoleStream _console
 

Static Public Attributes

static const std::string type_param
 
static const std::string name_param
 
static const std::string unique_name_param
 
static const std::string app_param
 
static const std::string moose_base_param
 
static const std::string kokkos_object_param
 
static constexpr PropertyValue::id_type default_property_id
 
static constexpr PropertyValue::id_type zero_property_id
 
static constexpr auto SYSTEM
 
static constexpr auto NAME
 

Protected Member Functions

virtual Real cohesion (const Real internal_param) const
 cohesion as a function of residual value, rate, and internal_param
 
virtual Real dcohesion (const Real internal_param) const
 d(cohesion)/d(internal_param) as a function of residual value, rate, and internal_param
 
virtual Real phi (const Real internal_param) const
 phi as a function of residual value, rate, and internal_param
 
virtual Real dphi (const Real internal_param) const
 d(phi)/d(internal_param) as a function of residual value, rate, and internal_param
 
virtual Real psi (const Real internal_param) const
 psi as a function of residual value, rate, and internal_param
 
virtual Real dpsi (const Real internal_param) const
 d(psi)/d(internal_param) as a function of residual value, rate, and internal_param
 
virtual Real yieldFunction (const RankTwoTensor &stress, Real intnl) const
 The following functions are what you should override when building single-plasticity models.
 
virtual RankTwoTensor dyieldFunction_dstress (const RankTwoTensor &stress, Real intnl) const
 The derivative of yield function with respect to stress.
 
virtual Real dyieldFunction_dintnl (const RankTwoTensor &stress, Real intnl) const
 The derivative of yield function with respect to the internal parameter.
 
virtual RankTwoTensor flowPotential (const RankTwoTensor &stress, Real intnl) const
 The flow potential.
 
virtual RankFourTensor dflowPotential_dstress (const RankTwoTensor &stress, Real intnl) const
 The derivative of the flow potential with respect to stress.
 
virtual RankTwoTensor dflowPotential_dintnl (const RankTwoTensor &stress, Real intnl) const
 The derivative of the flow potential with respect to the internal parameter.
 
virtual Real hardPotential (const RankTwoTensor &stress, Real intnl) const
 The hardening potential.
 
virtual RankTwoTensor dhardPotential_dstress (const RankTwoTensor &stress, Real intnl) const
 The derivative of the hardening potential with respect to stress.
 
virtual Real dhardPotential_dintnl (const RankTwoTensor &stress, Real intnl) const
 The derivative of the hardening potential with respect to the internal parameter.
 
virtual void addPostprocessorDependencyHelper (const PostprocessorName &name) const override
 
virtual void addVectorPostprocessorDependencyHelper (const VectorPostprocessorName &name) const override
 
virtual void addUserObjectDependencyHelper (const UserObjectBase &uo) const override
 
void addReporterDependencyHelper (const ReporterName &reporter_name) override
 
void flagInvalidSolutionInternal (const InvalidSolutionID invalid_solution_id) const
 
InvalidSolutionID registerInvalidSolutionInternal (const std::string &message, const bool warning) const
 
const ReporterContextBasegetReporterContextBaseByName (const ReporterName &reporter_name) const
 
const ReporterNamegetReporterName (const std::string &param_name) const
 
TdeclareRestartableData (const std::string &data_name, Args &&... args)
 
ManagedValue< TdeclareManagedRestartableDataWithContext (const std::string &data_name, void *context, Args &&... args)
 
const TgetRestartableData (const std::string &data_name) const
 
TdeclareRestartableDataWithContext (const std::string &data_name, void *context, Args &&... args)
 
TdeclareRecoverableData (const std::string &data_name, Args &&... args)
 
TdeclareRestartableDataWithObjectName (const std::string &data_name, const std::string &object_name, Args &&... args)
 
TdeclareRestartableDataWithObjectNameWithContext (const std::string &data_name, const std::string &object_name, void *context, Args &&... args)
 
std::string restartableName (const std::string &data_name) const
 
const TgetMeshProperty (const std::string &data_name, const std::string &prefix)
 
const TgetMeshProperty (const std::string &data_name)
 
bool hasMeshProperty (const std::string &data_name, const std::string &prefix) const
 
bool hasMeshProperty (const std::string &data_name, const std::string &prefix) const
 
bool hasMeshProperty (const std::string &data_name) const
 
bool hasMeshProperty (const std::string &data_name) const
 
std::string meshPropertyName (const std::string &data_name) const
 
PerfID registerTimedSection (const std::string &section_name, const unsigned int level) const
 
PerfID registerTimedSection (const std::string &section_name, const unsigned int level, const std::string &live_message, const bool print_dots=true) const
 
std::string timedSectionName (const std::string &section_name) const
 
bool isCoupledScalar (const std::string &var_name, unsigned int i=0) const
 
unsigned int coupledScalarComponents (const std::string &var_name) const
 
unsigned int coupledScalar (const std::string &var_name, unsigned int comp=0) const
 
libMesh::Order coupledScalarOrder (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarValue (const std::string &var_name, unsigned int comp=0) const
 
const ADVariableValueadCoupledScalarValue (const std::string &var_name, unsigned int comp=0) const
 
const GenericVariableValue< is_ad > & coupledGenericScalarValue (const std::string &var_name, unsigned int comp=0) const
 
const GenericVariableValue< false > & coupledGenericScalarValue (const std::string &var_name, const unsigned int comp) const
 
const GenericVariableValue< true > & coupledGenericScalarValue (const std::string &var_name, const unsigned int comp) const
 
const VariableValuecoupledVectorTagScalarValue (const std::string &var_name, TagID tag, unsigned int comp=0) const
 
const VariableValuecoupledMatrixTagScalarValue (const std::string &var_name, TagID tag, unsigned int comp=0) const
 
const VariableValuecoupledScalarValueOld (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarValueOlder (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarDot (const std::string &var_name, unsigned int comp=0) const
 
const ADVariableValueadCoupledScalarDot (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarDotDot (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarDotOld (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarDotDotOld (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarDotDu (const std::string &var_name, unsigned int comp=0) const
 
const VariableValuecoupledScalarDotDotDu (const std::string &var_name, unsigned int comp=0) const
 
const MooseVariableScalargetScalarVar (const std::string &var_name, unsigned int comp) const
 
virtual void checkMaterialProperty (const std::string &name, const unsigned int state)
 
virtual void getKokkosMaterialPropertyHook (const std::string &, const unsigned int)
 
void markMatPropRequested (const std::string &)
 
MaterialPropertyName getMaterialPropertyName (const std::string &name) const
 
void checkExecutionStage ()
 
const TgetReporterValue (const std::string &param_name, const std::size_t time_index=0)
 
const TgetReporterValue (const std::string &param_name, ReporterMode mode, const std::size_t time_index=0)
 
const TgetReporterValue (const std::string &param_name, const std::size_t time_index=0)
 
const TgetReporterValue (const std::string &param_name, ReporterMode mode, const std::size_t time_index=0)
 
const TgetReporterValueByName (const ReporterName &reporter_name, const std::size_t time_index=0)
 
const TgetReporterValueByName (const ReporterName &reporter_name, ReporterMode mode, const std::size_t time_index=0)
 
const TgetReporterValueByName (const ReporterName &reporter_name, const std::size_t time_index=0)
 
const TgetReporterValueByName (const ReporterName &reporter_name, ReporterMode mode, const std::size_t time_index=0)
 
bool hasReporterValue (const std::string &param_name) const
 
bool hasReporterValue (const std::string &param_name) const
 
bool hasReporterValue (const std::string &param_name) const
 
bool hasReporterValue (const std::string &param_name) const
 
bool hasReporterValueByName (const ReporterName &reporter_name) const
 
bool hasReporterValueByName (const ReporterName &reporter_name) const
 
bool hasReporterValueByName (const ReporterName &reporter_name) const
 
bool hasReporterValueByName (const ReporterName &reporter_name) const
 
const GenericMaterialProperty< T, is_ad > * defaultGenericMaterialProperty (const std::string &name)
 
const GenericMaterialProperty< T, is_ad > * defaultGenericMaterialProperty (const std::string &name)
 
const MaterialProperty< T > * defaultMaterialProperty (const std::string &name)
 
const MaterialProperty< T > * defaultMaterialProperty (const std::string &name)
 
const ADMaterialProperty< T > * defaultADMaterialProperty (const std::string &name)
 
const ADMaterialProperty< T > * defaultADMaterialProperty (const std::string &name)
 

Static Protected Member Functions

static std::string meshPropertyName (const std::string &data_name, const std::string &prefix)
 

Protected Attributes

const Moose::CoordinateSystemType_coord_sys
 
const THREAD_ID _tid
 
SubProblem_subproblem
 
FEProblemBase_fe_problem
 
SystemBase_sys
 
Assembly_assembly
 
const bool _duplicate_initial_execution
 
std::set< std::string > _depend_uo
 
const bool & _enabled
 
MooseApp_app
 
Factory_factory
 
ActionFactory_action_factory
 
const std::string & _type
 
const std::string & _name
 
const InputParameters_pars
 
const ExecFlagEnum_execute_enum
 
const ExecFlagType_current_execute_flag
 
MooseApp_restartable_app
 
const std::string _restartable_system_name
 
const THREAD_ID _restartable_tid
 
const bool _restartable_read_only
 
FEProblemBase_mci_feproblem
 
FEProblemBase_mdi_feproblem
 
MooseApp_pg_moose_app
 
const std::string _prefix
 
FEProblemBase_sc_fe_problem
 
const THREAD_ID _sc_tid
 
const Real & _real_zero
 
const VariableValue_scalar_zero
 
const Point & _point_zero
 
const InputParameters_mi_params
 
const std::string _mi_name
 
const MooseObjectName _mi_moose_object_name
 
FEProblemBase_mi_feproblem
 
SubProblem_mi_subproblem
 
const THREAD_ID _mi_tid
 
const bool _is_kokkos_object
 
const Moose::MaterialDataType _material_data_type
 
MaterialData_material_data
 
bool _stateful_allowed
 
bool _get_material_property_called
 
std::vector< std::unique_ptr< PropertyValue > > _default_properties
 
std::unordered_set< unsigned int_material_property_dependencies
 
const MaterialPropertyName _get_suffix
 
const bool _use_interpolated_state
 
const InputParameters_ti_params
 
FEProblemBase_ti_feproblem
 
bool _is_implicit
 
Real & _t
 
const Real & _t_old
 
int_t_step
 
Real & _dt
 
Real & _dt_old
 
bool _is_transient
 
const Parallel::Communicator & _communicator
 

Static Protected Attributes

static const std::string _interpolated_old
 
static const std::string _interpolated_older
 

Private Types

enum  return_type {
  tip110100 = 0 , tip010101 = 1 , edge010100 = 2 , edge000101 = 3 ,
  plane000100 = 4
}
 

Private Member Functions

void yieldFunctionEigvals (Real e0, Real e1, Real e2, Real sinphi, Real cohcos, std::vector< Real > &f) const
 Calculates the yield functions given the eigenvalues of stress.
 
void df_dsig (const RankTwoTensor &stress, Real sin_angle, std::vector< RankTwoTensor > &df) const
 this is exactly dyieldFunction_dstress, or flowPotential, depending on whether sin_angle = sin(phi), or sin_angle = sin(psi), respectively
 
void perturbStress (const RankTwoTensor &stress, std::vector< Real > &eigvals, std::vector< RankTwoTensor > &deigvals) const
 perturbs the stress tensor in the case of almost-equal eigenvalues.
 
bool KuhnTuckerOK (const std::vector< Real > &yf, const std::vector< Real > &dpm, Real ep_plastic_tolerance) const
 Returns true if the Kuhn-Tucker conditions are satisfied.
 
bool doReturnMap (const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const
 See doco for returnMap function.
 
bool returnTip (const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
 Tries to return-map to the MC tip using the THREE directions given in n, and THREE dpm values are returned.
 
bool returnPlane (const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
 Tries to return-map to the MC plane using the n[3] direction The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress) will be junk.
 
bool returnEdge000101 (const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, Real mag_E, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
 Tries to return-map to the MC edge using the n[4] and n[6] directions The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress) will be junk.
 
bool returnEdge010100 (const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, Real mag_E, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
 Tries to return-map to the MC edge using the n[1] and n[3] directions The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress) will be junk.
 
const Moose::FunctionBasegetKokkosFunctionByNameHelper (const FunctionName &name) const
 
const UserObjectBasegetUserObjectFromFEProblem (const UserObjectName &object_name, const THREAD_ID tid=0) const
 
const TcastUserObject (const UserObjectBase &uo_base, const std::string &param_name="") const
 
void mooseObjectError (const std::string &param_name, std::stringstream &oss) const
 
const std::string & userObjectType (const UserObjectBase &uo) const
 
const std::string & userObjectName (const UserObjectBase &uo) const
 
const PostprocessorName & getPostprocessorNameInternal (const std::string &param_name, const unsigned int index, const bool allow_default_value=true) const
 
bool isDefaultPostprocessorValueByName (const PostprocessorName &name) const
 
PostprocessorValue getDefaultPostprocessorValueByName (const PostprocessorName &name) const
 
void checkParam (const std::string &param_name, const unsigned int index=std::numeric_limits< unsigned int >::max()) const
 
bool postprocessorsAdded () const
 
const VectorPostprocessorValuegetVectorPostprocessorByNameHelper (const VectorPostprocessorName &name, const std::string &vector_name, bool broadcast, std::size_t t_index) const
 
const VectorPostprocessorContext< VectorPostprocessorValue > & getVectorPostprocessorContextByNameHelper (const VectorPostprocessorName &name, const std::string &vector_name) const
 
bool vectorPostprocessorsAdded () const
 
bool reportersAdded () const
 
void possiblyCheckHasReporter (const ReporterName &reporter_name, const std::string &param_name="") const
 
RestartableDataValueregisterRestartableDataOnApp (std::unique_ptr< RestartableDataValue > data, THREAD_ID tid) const
 
void registerRestartableNameWithFilterOnApp (const std::string &name, Moose::RESTARTABLE_FILTER filter)
 
RestartableData< T > & declareRestartableDataHelper (const std::string &data_name, void *context, Args &&... args) const
 
virtual std::string meshPropertyPrefix (const std::string &data_name) const
 
const RestartableDataValuegetMeshPropertyInternal (const std::string &data_name, const std::string &prefix) const
 
void mooseErrorInternal (Args &&... args) const
 
const VariableValuegetDefaultValue (const std::string &var_name) const
 
const ADVariableValuegetADDefaultValue (const std::string &var_name) const
 
void checkVar (const std::string &var_name) const
 
void validateExecutionerType (const std::string &name, const std::string &fn_name) const
 
Moose::MaterialDataType getMaterialDataType (const std::set< BoundaryID > &boundary_ids) const
 
unsigned int getMaxQps () const
 
void addConsumedPropertyName (const MooseObjectName &obj_name, const std::string &prop_name)
 
const PostprocessorValuegetPostprocessorValueInternal (const std::string &param_name, unsigned int index, std::size_t t_index) const
 
const PostprocessorValuegetPostprocessorValueInternal (const std::string &param_name, unsigned int index, std::size_t t_index) const
 
const PostprocessorValuegetPostprocessorValueByNameInternal (const PostprocessorName &name, std::size_t t_index) const
 
const PostprocessorValuegetPostprocessorValueByNameInternal (const PostprocessorName &name, std::size_t t_index) const
 
void possiblyCheckHasVectorPostprocessor (const std::string &param_name, const std::string &vector_name) const
 
void possiblyCheckHasVectorPostprocessor (const std::string &param_name, const std::string &vector_name) const
 
void possiblyCheckHasVectorPostprocessorByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 
void possiblyCheckHasVectorPostprocessorByName (const VectorPostprocessorName &name, const std::string &vector_name) const
 

Static Private Member Functions

static const hit::Node * getHitNode (const InputParameters &params)
 
static std::string messagePrefix (const InputParameters &params, const bool hit_prefix)
 

Private Attributes

const SolidMechanicsHardeningModel_cohesion
 Hardening model for cohesion.
 
const SolidMechanicsHardeningModel_phi
 Hardening model for phi.
 
const SolidMechanicsHardeningModel_psi
 Hardening model for psi.
 
const unsigned int _max_iters
 Maximum Newton-Raphison iterations in the custom returnMap algorithm.
 
const Real _shift
 yield function is shifted by this amount to avoid problems with stress-derivatives at equal eigenvalues
 
const bool _use_custom_returnMap
 Whether to use the custom return-map algorithm.
 
UserObject_primary_thread_copy
 
std::set< std::string > _supplied_uo
 
const ParallelParamObject_parent
 
const MooseBase_si_moose_base
 
const FEProblemBase_si_problem
 
ExecFlagEnum _empty_execute_enum
 
const MooseObject_fni_object
 
const InputParameters_fni_params
 
FEProblemBase_fni_feproblem
 
const THREAD_ID _fni_tid
 
const MooseObject_uoi_moose_object
 
const FEProblemBase_uoi_feproblem
 
const THREAD_ID _uoi_tid
 
const MooseObject_ppi_moose_object
 
const InputParameters_ppi_params
 
const FEProblemBase_ppi_feproblem
 
std::map< PostprocessorName, std::unique_ptr< PostprocessorValue > > _default_values
 
const bool _broadcast_by_default
 
const MooseObject_vpi_moose_object
 
const FEProblemBase_vpi_feproblem
 
const THREAD_ID _vpi_tid
 
const InputParameters_ri_params
 
FEProblemBase_ri_fe_problem_base
 
const ReporterData_ri_reporter_data
 
const MooseObject_ri_moose_object
 
const InputParameters_dni_params
 
FEProblemBase_dni_feproblem
 
const MooseObject *const _dni_moose_object_ptr
 
const InputParameters_si_params
 
FEProblemBase_si_feproblem
 
THREAD_ID _si_tid
 
const RestartableDataMapName _metaname
 
std::string _restartable_name
 
MooseApp_meta_data_app
 
const MooseObject *const _meta_data_object
 
const InputParameters_sc_parameters
 
const std::string & _sc_name
 
const bool _sc_is_implicit
 
std::unordered_map< std::string, std::vector< MooseVariableScalar * > > _coupled_scalar_vars
 
std::unordered_map< std::string, std::unique_ptr< VariableValue > > _default_value
 
std::unordered_map< std::string, std::unique_ptr< ADVariableValue > > _dual_default_value
 
std::vector< MooseVariableScalar * > _coupled_moose_scalar_vars
 
std::unordered_map< std::string, std::vector< MooseVariableFieldBase * > > _sc_coupled_vars
 
std::set< TagID_sc_coupleable_vector_tags
 
std::set< TagID_sc_coupleable_matrix_tags
 
const MooseObject_mi_moose_object
 
const bool _mi_boundary_restricted
 
const std::set< SubdomainID > & _mi_block_ids
 
const std::set< BoundaryID > & _mi_boundary_ids
 
std::vector< std::unique_ptr< OptionalMaterialPropertyProxyBase< MaterialPropertyInterface > > > _optional_property_proxies
 
const std::string _ti_name
 

Detailed Description

FiniteStrainMohrCoulombMulti implements rate-independent non-associative mohr-coulomb with hardening/softening in the finite-strain framework, using planar (non-smoothed) surfaces.

Definition at line 19 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Member Enumeration Documentation

◆ return_type

Constructor & Destructor Documentation

◆ SolidMechanicsPlasticMohrCoulombMulti()

SolidMechanicsPlasticMohrCoulombMulti::SolidMechanicsPlasticMohrCoulombMulti ( const InputParameters parameters)

Definition at line 61 of file SolidMechanicsPlasticMohrCoulombMulti.C.

64 _cohesion(getUserObject<SolidMechanicsHardeningModel>("cohesion")),
65 _phi(getUserObject<SolidMechanicsHardeningModel>("friction_angle")),
66 _psi(getUserObject<SolidMechanicsHardeningModel>("dilation_angle")),
67 _max_iters(getParam<unsigned int>("max_iterations")),
68 _shift(parameters.isParamValid("shift") ? getParam<Real>("shift") : _f_tol),
69 _use_custom_returnMap(getParam<bool>("use_custom_returnMap"))
70{
71 if (_shift < 0)
72 mooseError("Value of 'shift' in SolidMechanicsPlasticMohrCoulombMulti must not be negative\n");
73 if (_shift > _f_tol)
74 _console << "WARNING: value of 'shift' in SolidMechanicsPlasticMohrCoulombMulti is probably "
75 "set too high"
76 << std::endl;
77 if (LIBMESH_DIM != 3)
78 mooseError("SolidMechanicsPlasticMohrCoulombMulti is only defined for LIBMESH_DIM=3");
80}
const ConsoleStream _console
bool isParamValid(const std::string &name) const
const InputParameters & parameters() const
void mooseError(Args &&... args) const
const T & getParam(const std::string &name) const
void seed(std::size_t i, unsigned int seed)
Plastic Model base class The virtual functions written below must be over-ridden in derived classes t...
const Real _f_tol
Tolerance on yield function.
const Real _shift
yield function is shifted by this amount to avoid problems with stress-derivatives at equal eigenvalu...
const unsigned int _max_iters
Maximum Newton-Raphison iterations in the custom returnMap algorithm.
const SolidMechanicsHardeningModel & _cohesion
Hardening model for cohesion.
const SolidMechanicsHardeningModel & _phi
Hardening model for phi.
const bool _use_custom_returnMap
Whether to use the custom return-map algorithm.
const SolidMechanicsHardeningModel & _psi
Hardening model for psi.

Member Function Documentation

◆ activeConstraints()

void SolidMechanicsPlasticMohrCoulombMulti::activeConstraints ( const std::vector< Real > &  f,
const RankTwoTensor stress,
Real  intnl,
const RankFourTensor Eijkl,
std::vector< bool > &  act,
RankTwoTensor returned_stress 
) const
overridevirtual

The active yield surfaces, given a vector of yield functions.

This is used by FiniteStrainMultiPlasticity to determine the initial set of active constraints at the trial (stress, intnl) configuration. It is up to you (the coder) to determine how accurate you want the returned_stress to be. Currently it is only used by FiniteStrainMultiPlasticity to estimate a good starting value for the Newton-Rahson procedure, so currently it may not need to be super perfect.

Parameters
fvalues of the yield functions
stressstress tensor
intnlinternal parameter
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actact[i] = true if the i_th yield function is active
[out]returned_stressApproximate value of the returned stress

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 248 of file SolidMechanicsPlasticMohrCoulombMulti.C.

254{
255 act.assign(6, false);
256
257 if (f[0] <= _f_tol && f[1] <= _f_tol && f[2] <= _f_tol && f[3] <= _f_tol && f[4] <= _f_tol &&
258 f[5] <= _f_tol)
259 {
260 returned_stress = stress;
261 return;
262 }
263
264 Real returned_intnl;
265 std::vector<Real> dpm(6);
266 RankTwoTensor delta_dp;
267 std::vector<Real> yf(6);
268 bool trial_stress_inadmissible;
269 doReturnMap(stress,
270 intnl,
271 Eijkl,
272 0.0,
273 returned_stress,
274 returned_intnl,
275 dpm,
276 delta_dp,
277 yf,
278 trial_stress_inadmissible);
279
280 for (unsigned i = 0; i < 6; ++i)
281 act[i] = (dpm[i] > 0);
282}
Real f(Real x)
Test function for Brents method.
bool doReturnMap(const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const
See doco for returnMap function.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ cohesion()

Real SolidMechanicsPlasticMohrCoulombMulti::cohesion ( const Real  internal_param) const
protectedvirtual

cohesion as a function of residual value, rate, and internal_param

Definition at line 285 of file SolidMechanicsPlasticMohrCoulombMulti.C.

286{
287 return _cohesion.value(internal_param);
288}
virtual Real value(Real intnl) const

Referenced by doReturnMap(), dyieldFunction_dintnlV(), returnEdge000101(), returnEdge010100(), returnPlane(), returnTip(), and yieldFunctionV().

◆ consistentTangentOperator()

RankFourTensor SolidMechanicsPlasticModel::consistentTangentOperator ( const RankTwoTensor trial_stress,
Real  intnl_old,
const RankTwoTensor stress,
Real  intnl,
const RankFourTensor E_ijkl,
const std::vector< Real > &  cumulative_pm 
) const
virtualinherited

Calculates a custom consistent tangent operator.

You may choose to over-ride this in your derived SolidMechanicsPlasticXXXX class.

(Note, if you over-ride returnMap, you will probably want to override consistentTangentOpertor too, otherwise it will default to E_ijkl.)

Parameters
stress_oldtrial stress before returning
intnl_oldinternal parameter before returning
stresscurrent returned stress state
intnlinternal parameter
E_ijklelasticity tensor
cumulative_pmthe cumulative plastic multipliers
Returns
the consistent tangent operator: E_ijkl if not over-ridden

Reimplemented in SolidMechanicsPlasticDruckerPragerHyperbolic, SolidMechanicsPlasticJ2, SolidMechanicsPlasticMeanCapTC, and SolidMechanicsPlasticTensileMulti.

Definition at line 252 of file SolidMechanicsPlasticModel.C.

259{
260 return E_ijkl;
261}

Referenced by SolidMechanicsPlasticDruckerPragerHyperbolic::consistentTangentOperator(), SolidMechanicsPlasticJ2::consistentTangentOperator(), SolidMechanicsPlasticMeanCapTC::consistentTangentOperator(), and SolidMechanicsPlasticTensileMulti::consistentTangentOperator().

◆ dcohesion()

Real SolidMechanicsPlasticMohrCoulombMulti::dcohesion ( const Real  internal_param) const
protectedvirtual

d(cohesion)/d(internal_param) as a function of residual value, rate, and internal_param

Definition at line 291 of file SolidMechanicsPlasticMohrCoulombMulti.C.

292{
293 return _cohesion.derivative(internal_param);
294}
virtual Real derivative(Real intnl) const

Referenced by dyieldFunction_dintnlV(), returnEdge000101(), returnEdge010100(), returnPlane(), and returnTip().

◆ df_dsig()

void SolidMechanicsPlasticMohrCoulombMulti::df_dsig ( const RankTwoTensor stress,
Real  sin_angle,
std::vector< RankTwoTensor > &  df 
) const
private

this is exactly dyieldFunction_dstress, or flowPotential, depending on whether sin_angle = sin(phi), or sin_angle = sin(psi), respectively

Definition at line 144 of file SolidMechanicsPlasticMohrCoulombMulti.C.

147{
148 std::vector<Real> eigvals;
149 std::vector<RankTwoTensor> deigvals;
150 stress.dsymmetricEigenvalues(eigvals, deigvals);
151
152 if (eigvals[0] > eigvals[1] - 0.1 * _shift || eigvals[1] > eigvals[2] - 0.1 * _shift)
153 perturbStress(stress, eigvals, deigvals);
154
155 df.resize(6);
156 df[0] = 0.5 * (deigvals[0] - deigvals[1]) + 0.5 * (deigvals[0] + deigvals[1]) * sin_angle;
157 df[1] = 0.5 * (deigvals[1] - deigvals[0]) + 0.5 * (deigvals[0] + deigvals[1]) * sin_angle;
158 df[2] = 0.5 * (deigvals[0] - deigvals[2]) + 0.5 * (deigvals[0] + deigvals[2]) * sin_angle;
159 df[3] = 0.5 * (deigvals[2] - deigvals[0]) + 0.5 * (deigvals[0] + deigvals[2]) * sin_angle;
160 df[4] = 0.5 * (deigvals[1] - deigvals[2]) + 0.5 * (deigvals[1] + deigvals[2]) * sin_angle;
161 df[5] = 0.5 * (deigvals[2] - deigvals[1]) + 0.5 * (deigvals[1] + deigvals[2]) * sin_angle;
162}
void dsymmetricEigenvalues(std::vector< T > &eigvals, std::vector< RankTwoTensorTempl< T > > &deigvals) const
void perturbStress(const RankTwoTensor &stress, std::vector< Real > &eigvals, std::vector< RankTwoTensor > &deigvals) const
perturbs the stress tensor in the case of almost-equal eigenvalues.

Referenced by dyieldFunction_dstressV(), and flowPotentialV().

◆ dflowPotential_dintnl()

RankTwoTensor SolidMechanicsPlasticModel::dflowPotential_dintnl ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The derivative of the flow potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
Returns
dr_dintnl(i, j) = dr(i, j)/dintnl

Reimplemented in SolidMechanicsPlasticDruckerPrager, SolidMechanicsPlasticJ2, SolidMechanicsPlasticMeanCap, SolidMechanicsPlasticMeanCapTC, SolidMechanicsPlasticMohrCoulomb, SolidMechanicsPlasticSimpleTester, SolidMechanicsPlasticTensile, SolidMechanicsPlasticWeakPlaneShear, SolidMechanicsPlasticWeakPlaneTensile, and SolidMechanicsPlasticWeakPlaneTensileN.

Definition at line 131 of file SolidMechanicsPlasticModel.C.

133{
134 return RankTwoTensor();
135}

Referenced by SolidMechanicsPlasticModel::dflowPotential_dintnlV().

◆ dflowPotential_dintnlV()

void SolidMechanicsPlasticMohrCoulombMulti::dflowPotential_dintnlV ( const RankTwoTensor stress,
Real  intnl,
std::vector< RankTwoTensor > &  dr_dintnl 
) const
overridevirtual

The derivative of the flow potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
[out]dr_dintnldr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 228 of file SolidMechanicsPlasticMohrCoulombMulti.C.

230{
231 const Real cos_angle = std::cos(psi(intnl));
232 const Real dsin_angle = cos_angle * dpsi(intnl);
233
234 std::vector<Real> eigvals;
235 std::vector<RankTwoTensor> deigvals;
236 stress.dsymmetricEigenvalues(eigvals, deigvals);
237
238 if (eigvals[0] > eigvals[1] - 0.1 * _shift || eigvals[1] > eigvals[2] - 0.1 * _shift)
239 perturbStress(stress, eigvals, deigvals);
240
241 dr_dintnl.resize(6);
242 dr_dintnl[0] = dr_dintnl[1] = 0.5 * (deigvals[0] + deigvals[1]) * dsin_angle;
243 dr_dintnl[2] = dr_dintnl[3] = 0.5 * (deigvals[0] + deigvals[2]) * dsin_angle;
244 dr_dintnl[4] = dr_dintnl[5] = 0.5 * (deigvals[1] + deigvals[2]) * dsin_angle;
245}
virtual Real dpsi(const Real internal_param) const
d(psi)/d(internal_param) as a function of residual value, rate, and internal_param
virtual Real psi(const Real internal_param) const
psi as a function of residual value, rate, and internal_param

◆ dflowPotential_dstress()

RankFourTensor SolidMechanicsPlasticModel::dflowPotential_dstress ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

◆ dflowPotential_dstressV()

void SolidMechanicsPlasticMohrCoulombMulti::dflowPotential_dstressV ( const RankTwoTensor stress,
Real  intnl,
std::vector< RankFourTensor > &  dr_dstress 
) const
overridevirtual

The derivative of the flow potential with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
[out]dr_dstressdr_dstress[alpha](i, j, k, l) = dr[alpha](i, j)/dstress(k, l)

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 204 of file SolidMechanicsPlasticMohrCoulombMulti.C.

206{
207 std::vector<RankFourTensor> d2eigvals;
208 stress.d2symmetricEigenvalues(d2eigvals);
209
210 const Real sinpsi = std::sin(psi(intnl));
211
212 dr_dstress.resize(6);
213 dr_dstress[0] =
214 0.5 * (d2eigvals[0] - d2eigvals[1]) + 0.5 * (d2eigvals[0] + d2eigvals[1]) * sinpsi;
215 dr_dstress[1] =
216 0.5 * (d2eigvals[1] - d2eigvals[0]) + 0.5 * (d2eigvals[0] + d2eigvals[1]) * sinpsi;
217 dr_dstress[2] =
218 0.5 * (d2eigvals[0] - d2eigvals[2]) + 0.5 * (d2eigvals[0] + d2eigvals[2]) * sinpsi;
219 dr_dstress[3] =
220 0.5 * (d2eigvals[2] - d2eigvals[0]) + 0.5 * (d2eigvals[0] + d2eigvals[2]) * sinpsi;
221 dr_dstress[4] =
222 0.5 * (d2eigvals[1] - d2eigvals[2]) + 0.5 * (d2eigvals[1] + d2eigvals[2]) * sinpsi;
223 dr_dstress[5] =
224 0.5 * (d2eigvals[2] - d2eigvals[1]) + 0.5 * (d2eigvals[1] + d2eigvals[2]) * sinpsi;
225}
void d2symmetricEigenvalues(std::vector< RankFourTensorTempl< T > > &deriv) const

◆ dhardPotential_dintnl()

Real SolidMechanicsPlasticModel::dhardPotential_dintnl ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The derivative of the hardening potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
Returns
the derivative

Reimplemented in SolidMechanicsPlasticMeanCapTC.

Definition at line 172 of file SolidMechanicsPlasticModel.C.

174{
175 return 0.0;
176}

Referenced by SolidMechanicsPlasticModel::dhardPotential_dintnlV().

◆ dhardPotential_dintnlV()

void SolidMechanicsPlasticModel::dhardPotential_dintnlV ( const RankTwoTensor stress,
Real  intnl,
std::vector< Real > &  dh_dintnl 
) const
virtualinherited

The derivative of the hardening potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
[out]dh_dintnldh_dintnl[alpha] = dh[alpha]/dintnl

Definition at line 178 of file SolidMechanicsPlasticModel.C.

181{
182 dh_dintnl.resize(numberSurfaces(), dhardPotential_dintnl(stress, intnl));
183}
virtual unsigned int numberSurfaces() const
The number of yield surfaces for this plasticity model.
virtual Real dhardPotential_dintnl(const RankTwoTensor &stress, Real intnl) const
The derivative of the hardening potential with respect to the internal parameter.

◆ dhardPotential_dstress()

RankTwoTensor SolidMechanicsPlasticModel::dhardPotential_dstress ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The derivative of the hardening potential with respect to stress.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
Returns
dh_dstress(i, j) = dh/dstress(i, j)

Reimplemented in SolidMechanicsPlasticMeanCapTC.

Definition at line 158 of file SolidMechanicsPlasticModel.C.

160{
161 return RankTwoTensor();
162}

Referenced by SolidMechanicsPlasticModel::dhardPotential_dstressV().

◆ dhardPotential_dstressV()

void SolidMechanicsPlasticModel::dhardPotential_dstressV ( const RankTwoTensor stress,
Real  intnl,
std::vector< RankTwoTensor > &  dh_dstress 
) const
virtualinherited

The derivative of the hardening potential with respect to stress.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
[out]dh_dstressdh_dstress[alpha](i, j) = dh[alpha]/dstress(i, j)

Definition at line 164 of file SolidMechanicsPlasticModel.C.

167{
168 dh_dstress.assign(numberSurfaces(), dhardPotential_dstress(stress, intnl));
169}
virtual RankTwoTensor dhardPotential_dstress(const RankTwoTensor &stress, Real intnl) const
The derivative of the hardening potential with respect to stress.

◆ doReturnMap()

bool SolidMechanicsPlasticMohrCoulombMulti::doReturnMap ( const RankTwoTensor trial_stress,
Real  intnl_old,
const RankFourTensor E_ijkl,
Real  ep_plastic_tolerance,
RankTwoTensor returned_stress,
Real &  returned_intnl,
std::vector< Real > &  dpm,
RankTwoTensor delta_dp,
std::vector< Real > &  yf,
bool &  trial_stress_inadmissible 
) const
private

See doco for returnMap function.

The interface is identical to this one. This one can be called internally regardless of the value of _use_custom_returnMap

Definition at line 363 of file SolidMechanicsPlasticMohrCoulombMulti.C.

373{
374 mooseAssert(dpm.size() == 6,
375 "SolidMechanicsPlasticMohrCoulombMulti size of dpm should be 6 but it is "
376 << dpm.size());
377
378 std::vector<Real> eigvals;
379 RankTwoTensor eigvecs;
380 trial_stress.symmetricEigenvaluesEigenvectors(eigvals, eigvecs);
381 eigvals[0] += _shift;
382 eigvals[2] -= _shift;
383
384 Real sinphi = std::sin(phi(intnl_old));
385 Real cosphi = std::cos(phi(intnl_old));
386 Real coh = cohesion(intnl_old);
387 Real cohcos = coh * cosphi;
388
389 yieldFunctionEigvals(eigvals[0], eigvals[1], eigvals[2], sinphi, cohcos, yf);
390
391 if (yf[0] <= _f_tol && yf[1] <= _f_tol && yf[2] <= _f_tol && yf[3] <= _f_tol && yf[4] <= _f_tol &&
392 yf[5] <= _f_tol)
393 {
394 // purely elastic (trial_stress, intnl_old)
395 trial_stress_inadmissible = false;
396 return true;
397 }
398
399 trial_stress_inadmissible = true;
400 delta_dp.zero();
401 returned_stress = RankTwoTensor();
402
403 // these are the normals to the 6 yield surfaces, which are const because of the assumption of no
404 // psi hardening
405 std::vector<RealVectorValue> norm(6);
406 const Real sinpsi = std::sin(psi(intnl_old));
407 const Real oneminus = 0.5 * (1 - sinpsi);
408 const Real oneplus = 0.5 * (1 + sinpsi);
409 norm[0](0) = oneplus;
410 norm[0](1) = -oneminus;
411 norm[0](2) = 0;
412 norm[1](0) = -oneminus;
413 norm[1](1) = oneplus;
414 norm[1](2) = 0;
415 norm[2](0) = oneplus;
416 norm[2](1) = 0;
417 norm[2](2) = -oneminus;
418 norm[3](0) = -oneminus;
419 norm[3](1) = 0;
420 norm[3](2) = oneplus;
421 norm[4](0) = 0;
422 norm[4](1) = oneplus;
423 norm[4](2) = -oneminus;
424 norm[5](0) = 0;
425 norm[5](1) = -oneminus;
426 norm[5](2) = oneplus;
427
428 // the flow directions are these norm multiplied by Eijkl.
429 // I call the flow directions "n".
430 // In the following I assume that the Eijkl is
431 // for an isotropic situation. Then I don't have to
432 // rotate to the principal-stress frame, and i don't
433 // have to worry about strange off-diagonal things
434 std::vector<RealVectorValue> n(6);
435 for (unsigned ys = 0; ys < 6; ++ys)
436 for (unsigned i = 0; i < 3; ++i)
437 for (unsigned j = 0; j < 3; ++j)
438 n[ys](i) += E_ijkl(i, i, j, j) * norm[ys](j);
439 const Real mag_E = E_ijkl(0, 0, 0, 0);
440
441 // With non-zero Poisson's ratio and hardening
442 // it is not computationally cheap to know whether
443 // the trial stress will return to the tip, edge,
444 // or plane. The following at least
445 // gives a not-completely-stupid guess
446 // trial_order[0] = type of return to try first
447 // trial_order[1] = type of return to try second
448 // trial_order[2] = type of return to try third
449 // trial_order[3] = type of return to try fourth
450 // trial_order[4] = type of return to try fifth
451 // In the following the "binary" stuff indicates the
452 // deactive (0) and active (1) surfaces, eg
453 // 110100 means that surfaces 0, 1 and 3 are active
454 // and 2, 4 and 5 are deactive
455 const unsigned int number_of_return_paths = 5;
456 std::vector<int> trial_order(number_of_return_paths);
457 if (yf[1] > _f_tol && yf[3] > _f_tol && yf[5] > _f_tol)
458 {
459 trial_order[0] = tip110100;
460 trial_order[1] = edge010100;
461 trial_order[2] = plane000100;
462 trial_order[3] = edge000101;
463 trial_order[4] = tip010101;
464 }
465 else if (yf[1] <= _f_tol && yf[3] > _f_tol && yf[5] > _f_tol)
466 {
467 trial_order[0] = edge000101;
468 trial_order[1] = plane000100;
469 trial_order[2] = tip110100;
470 trial_order[3] = tip010101;
471 trial_order[4] = edge010100;
472 }
473 else if (yf[1] <= _f_tol && yf[3] > _f_tol && yf[5] <= _f_tol)
474 {
475 trial_order[0] = plane000100;
476 trial_order[1] = edge000101;
477 trial_order[2] = edge010100;
478 trial_order[3] = tip110100;
479 trial_order[4] = tip010101;
480 }
481 else
482 {
483 trial_order[0] = edge010100;
484 trial_order[1] = plane000100;
485 trial_order[2] = edge000101;
486 trial_order[3] = tip110100;
487 trial_order[4] = tip010101;
488 }
489
490 unsigned trial;
491 bool nr_converged = false;
492 bool kt_success = false;
493 std::vector<RealVectorValue> ntip(3);
494 std::vector<Real> dpmtip(3);
495
496 for (trial = 0; trial < number_of_return_paths; ++trial)
497 {
498 switch (trial_order[trial])
499 {
500 case tip110100:
501 for (unsigned int i = 0; i < 3; ++i)
502 {
503 ntip[0](i) = n[0](i);
504 ntip[1](i) = n[1](i);
505 ntip[2](i) = n[3](i);
506 }
507 kt_success = returnTip(eigvals,
508 ntip,
509 dpmtip,
510 returned_stress,
511 intnl_old,
512 sinphi,
513 cohcos,
514 0,
515 nr_converged,
516 ep_plastic_tolerance,
517 yf);
518 if (nr_converged && kt_success)
519 {
520 dpm[0] = dpmtip[0];
521 dpm[1] = dpmtip[1];
522 dpm[3] = dpmtip[2];
523 dpm[2] = dpm[4] = dpm[5] = 0;
524 }
525 break;
526
527 case tip010101:
528 for (unsigned int i = 0; i < 3; ++i)
529 {
530 ntip[0](i) = n[1](i);
531 ntip[1](i) = n[3](i);
532 ntip[2](i) = n[5](i);
533 }
534 kt_success = returnTip(eigvals,
535 ntip,
536 dpmtip,
537 returned_stress,
538 intnl_old,
539 sinphi,
540 cohcos,
541 0,
542 nr_converged,
543 ep_plastic_tolerance,
544 yf);
545 if (nr_converged && kt_success)
546 {
547 dpm[1] = dpmtip[0];
548 dpm[3] = dpmtip[1];
549 dpm[5] = dpmtip[2];
550 dpm[0] = dpm[2] = dpm[4] = 0;
551 }
552 break;
553
554 case edge000101:
555 kt_success = returnEdge000101(eigvals,
556 n,
557 dpm,
558 returned_stress,
559 intnl_old,
560 sinphi,
561 cohcos,
562 0,
563 mag_E,
564 nr_converged,
565 ep_plastic_tolerance,
566 yf);
567 break;
568
569 case edge010100:
570 kt_success = returnEdge010100(eigvals,
571 n,
572 dpm,
573 returned_stress,
574 intnl_old,
575 sinphi,
576 cohcos,
577 0,
578 mag_E,
579 nr_converged,
580 ep_plastic_tolerance,
581 yf);
582 break;
583
584 case plane000100:
585 kt_success = returnPlane(eigvals,
586 n,
587 dpm,
588 returned_stress,
589 intnl_old,
590 sinphi,
591 cohcos,
592 0,
593 nr_converged,
594 ep_plastic_tolerance,
595 yf);
596 break;
597 }
598
599 if (nr_converged && kt_success)
600 break;
601 }
602
603 if (trial == number_of_return_paths)
604 {
605 sinphi = std::sin(phi(intnl_old));
606 cosphi = std::cos(phi(intnl_old));
607 coh = cohesion(intnl_old);
608 cohcos = coh * cosphi;
609 yieldFunctionEigvals(eigvals[0], eigvals[1], eigvals[2], sinphi, cohcos, yf);
610 Moose::err << "Trial stress = \n";
611 trial_stress.print(Moose::err);
612 Moose::err << "which has eigenvalues = " << eigvals[0] << " " << eigvals[1] << " " << eigvals[2]
613 << "\n";
614 Moose::err << "and yield functions = " << yf[0] << " " << yf[1] << " " << yf[2] << " " << yf[3]
615 << " " << yf[4] << " " << yf[5] << "\n";
616 Moose::err << "Internal parameter = " << intnl_old << std::endl;
617 mooseError("SolidMechanicsPlasticMohrCoulombMulti: FAILURE! You probably need to implement a "
618 "line search if your hardening is too severe, or you need to tune your tolerances "
619 "(eg, yield_function_tolerance should be a little smaller than (young "
620 "modulus)*ep_plastic_tolerance).\n");
621 return false;
622 }
623
624 // success
625
626 returned_intnl = intnl_old;
627 for (unsigned i = 0; i < 6; ++i)
628 returned_intnl += dpm[i];
629 for (unsigned i = 0; i < 6; ++i)
630 for (unsigned j = 0; j < 3; ++j)
631 delta_dp(j, j) += dpm[i] * norm[i](j);
632 returned_stress = eigvecs * returned_stress * (eigvecs.transpose());
633 delta_dp = eigvecs * delta_dp * (eigvecs.transpose());
634 return true;
635}
RankTwoTensorTempl< T > transpose() const
void print(std::ostream &stm=Moose::out) const
void symmetricEigenvaluesEigenvectors(std::vector< T > &eigvals, RankTwoTensorTempl< T > &eigvecs) const
bool returnEdge010100(const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, Real mag_E, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
Tries to return-map to the MC edge using the n[1] and n[3] directions The return value is true if the...
virtual Real cohesion(const Real internal_param) const
cohesion as a function of residual value, rate, and internal_param
bool returnTip(const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
Tries to return-map to the MC tip using the THREE directions given in n, and THREE dpm values are ret...
bool returnPlane(const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
Tries to return-map to the MC plane using the n[3] direction The return value is true if the internal...
bool returnEdge000101(const std::vector< Real > &eigvals, const std::vector< RealVectorValue > &n, std::vector< Real > &dpm, RankTwoTensor &returned_stress, Real intnl_old, Real &sinphi, Real &cohcos, Real initial_guess, Real mag_E, bool &nr_converged, Real ep_plastic_tolerance, std::vector< Real > &yf) const
Tries to return-map to the MC edge using the n[4] and n[6] directions The return value is true if the...
virtual Real phi(const Real internal_param) const
phi as a function of residual value, rate, and internal_param
void yieldFunctionEigvals(Real e0, Real e1, Real e2, Real sinphi, Real cohcos, std::vector< Real > &f) const
Calculates the yield functions given the eigenvalues of stress.
auto norm(const T &a)

Referenced by activeConstraints(), and returnMap().

◆ dphi()

Real SolidMechanicsPlasticMohrCoulombMulti::dphi ( const Real  internal_param) const
protectedvirtual

d(phi)/d(internal_param) as a function of residual value, rate, and internal_param

Definition at line 303 of file SolidMechanicsPlasticMohrCoulombMulti.C.

304{
305 return _phi.derivative(internal_param);
306}

Referenced by dyieldFunction_dintnlV(), returnEdge000101(), returnEdge010100(), returnPlane(), and returnTip().

◆ dpsi()

Real SolidMechanicsPlasticMohrCoulombMulti::dpsi ( const Real  internal_param) const
protectedvirtual

d(psi)/d(internal_param) as a function of residual value, rate, and internal_param

Definition at line 315 of file SolidMechanicsPlasticMohrCoulombMulti.C.

316{
317 return _psi.derivative(internal_param);
318}

Referenced by dflowPotential_dintnlV().

◆ dyieldFunction_dintnl()

Real SolidMechanicsPlasticModel::dyieldFunction_dintnl ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The derivative of yield function with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
Returns
the derivative

Reimplemented in SolidMechanicsPlasticDruckerPrager, SolidMechanicsPlasticJ2, SolidMechanicsPlasticMeanCap, SolidMechanicsPlasticMeanCapTC, SolidMechanicsPlasticMohrCoulomb, SolidMechanicsPlasticSimpleTester, SolidMechanicsPlasticTensile, SolidMechanicsPlasticWeakPlaneShear, SolidMechanicsPlasticWeakPlaneTensile, and SolidMechanicsPlasticWeakPlaneTensileN.

Definition at line 90 of file SolidMechanicsPlasticModel.C.

92{
93 return 0.0;
94}

Referenced by SolidMechanicsPlasticModel::dyieldFunction_dintnlV().

◆ dyieldFunction_dintnlV()

void SolidMechanicsPlasticMohrCoulombMulti::dyieldFunction_dintnlV ( const RankTwoTensor stress,
Real  intnl,
std::vector< Real > &  df_dintnl 
) const
overridevirtual

The derivative of yield functions with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
[out]df_dintnldf_dintnl[alpha] = df[alpha]/dintnl

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 173 of file SolidMechanicsPlasticMohrCoulombMulti.C.

176{
177 std::vector<Real> eigvals;
178 stress.symmetricEigenvalues(eigvals);
179 eigvals[0] += _shift;
180 eigvals[2] -= _shift;
181
182 const Real sin_angle = std::sin(phi(intnl));
183 const Real cos_angle = std::cos(phi(intnl));
184 const Real dsin_angle = cos_angle * dphi(intnl);
185 const Real dcos_angle = -sin_angle * dphi(intnl);
186 const Real dcohcos = dcohesion(intnl) * cos_angle + cohesion(intnl) * dcos_angle;
187
188 df_dintnl.resize(6);
189 df_dintnl[0] = df_dintnl[1] = 0.5 * (eigvals[0] + eigvals[1]) * dsin_angle - dcohcos;
190 df_dintnl[2] = df_dintnl[3] = 0.5 * (eigvals[0] + eigvals[2]) * dsin_angle - dcohcos;
191 df_dintnl[4] = df_dintnl[5] = 0.5 * (eigvals[1] + eigvals[2]) * dsin_angle - dcohcos;
192}
void symmetricEigenvalues(std::vector< T > &eigvals) const
virtual Real dcohesion(const Real internal_param) const
d(cohesion)/d(internal_param) as a function of residual value, rate, and internal_param
virtual Real dphi(const Real internal_param) const
d(phi)/d(internal_param) as a function of residual value, rate, and internal_param

◆ dyieldFunction_dstress()

RankTwoTensor SolidMechanicsPlasticModel::dyieldFunction_dstress ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The derivative of yield function with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
Returns
df_dstress(i, j) = dyieldFunction/dstress(i, j)

Reimplemented in SolidMechanicsPlasticDruckerPrager, SolidMechanicsPlasticIsotropicSD, SolidMechanicsPlasticJ2, SolidMechanicsPlasticMeanCap, SolidMechanicsPlasticMeanCapTC, SolidMechanicsPlasticMohrCoulomb, SolidMechanicsPlasticOrthotropic, SolidMechanicsPlasticSimpleTester, SolidMechanicsPlasticTensile, SolidMechanicsPlasticWeakPlaneShear, SolidMechanicsPlasticWeakPlaneTensile, and SolidMechanicsPlasticWeakPlaneTensileN.

Definition at line 75 of file SolidMechanicsPlasticModel.C.

77{
78 return RankTwoTensor();
79}

Referenced by SolidMechanicsPlasticModel::dyieldFunction_dstressV().

◆ dyieldFunction_dstressV()

void SolidMechanicsPlasticMohrCoulombMulti::dyieldFunction_dstressV ( const RankTwoTensor stress,
Real  intnl,
std::vector< RankTwoTensor > &  df_dstress 
) const
overridevirtual

The derivative of yield functions with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
[out]df_dstressdf_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 165 of file SolidMechanicsPlasticMohrCoulombMulti.C.

167{
168 const Real sinphi = std::sin(phi(intnl));
169 df_dsig(stress, sinphi, df_dstress);
170}
void df_dsig(const RankTwoTensor &stress, Real sin_angle, std::vector< RankTwoTensor > &df) const
this is exactly dyieldFunction_dstress, or flowPotential, depending on whether sin_angle = sin(phi),...

◆ execute()

void SolidMechanicsPlasticModel::execute ( )
virtualinherited

Implements GeneralUserObject.

Definition at line 45 of file SolidMechanicsPlasticModel.C.

46{
47}

◆ finalize()

void SolidMechanicsPlasticModel::finalize ( )
virtualinherited

Implements GeneralUserObject.

Definition at line 50 of file SolidMechanicsPlasticModel.C.

51{
52}

◆ flowPotential()

RankTwoTensor SolidMechanicsPlasticModel::flowPotential ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

◆ flowPotentialV()

void SolidMechanicsPlasticMohrCoulombMulti::flowPotentialV ( const RankTwoTensor stress,
Real  intnl,
std::vector< RankTwoTensor > &  r 
) const
overridevirtual

The flow potentials.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
[out]rr[alpha] is the flow potential for the "alpha" yield function

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 195 of file SolidMechanicsPlasticMohrCoulombMulti.C.

198{
199 const Real sinpsi = std::sin(psi(intnl));
200 df_dsig(stress, sinpsi, r);
201}

◆ hardPotential()

Real SolidMechanicsPlasticModel::hardPotential ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The hardening potential.

Parameters
stressthe stress at which to calculate the hardening potential
intnlinternal parameter
Returns
the hardening potential

Reimplemented in SolidMechanicsPlasticMeanCapTC.

Definition at line 145 of file SolidMechanicsPlasticModel.C.

146{
147 return -1.0;
148}

Referenced by SolidMechanicsPlasticModel::hardPotentialV().

◆ hardPotentialV()

void SolidMechanicsPlasticModel::hardPotentialV ( const RankTwoTensor stress,
Real  intnl,
std::vector< Real > &  h 
) const
virtualinherited

The hardening potential.

Parameters
stressthe stress at which to calculate the hardening potential
intnlinternal parameter
[out]hh[alpha] is the hardening potential for the "alpha" yield function

Definition at line 150 of file SolidMechanicsPlasticModel.C.

153{
154 h.assign(numberSurfaces(), hardPotential(stress, intnl));
155}
virtual Real hardPotential(const RankTwoTensor &stress, Real intnl) const
The hardening potential.

◆ initialize()

void SolidMechanicsPlasticModel::initialize ( )
virtualinherited

Implements GeneralUserObject.

Definition at line 40 of file SolidMechanicsPlasticModel.C.

41{
42}

◆ KuhnTuckerOK()

bool SolidMechanicsPlasticMohrCoulombMulti::KuhnTuckerOK ( const std::vector< Real > &  yf,
const std::vector< Real > &  dpm,
Real  ep_plastic_tolerance 
) const
private

Returns true if the Kuhn-Tucker conditions are satisfied.

Parameters
yfThe six yield function values
dpmThe six plastic multipliers
ep_plastic_toleranceThe tolerance on the plastic strain (if dpm>-ep_plastic_tolerance then it is grouped as "non-negative" in the Kuhn-Tucker conditions).

Definition at line 1145 of file SolidMechanicsPlasticMohrCoulombMulti.C.

1148{
1149 mooseAssert(
1150 yf.size() == 6,
1151 "SolidMechanicsPlasticMohrCoulomb::KuhnTuckerOK requires yf to be size 6, but your is "
1152 << yf.size());
1153 mooseAssert(
1154 dpm.size() == 6,
1155 "SolidMechanicsPlasticMohrCoulomb::KuhnTuckerOK requires dpm to be size 6, but your is "
1156 << dpm.size());
1157 for (unsigned i = 0; i < 6; ++i)
1158 if (!SolidMechanicsPlasticModel::KuhnTuckerSingleSurface(yf[i], dpm[i], ep_plastic_tolerance))
1159 return false;
1160 return true;
1161}
bool KuhnTuckerSingleSurface(Real yf, Real dpm, Real dpm_tol) const
Returns true if the Kuhn-Tucker conditions for the single surface are satisfied.

◆ KuhnTuckerSingleSurface()

bool SolidMechanicsPlasticModel::KuhnTuckerSingleSurface ( Real  yf,
Real  dpm,
Real  dpm_tol 
) const
inherited

Returns true if the Kuhn-Tucker conditions for the single surface are satisfied.

Parameters
yfYield function value
dpmplastic multiplier
dpm_toltolerance on plastic multiplier: viz dpm>-dpm_tol means "dpm is non-negative"

Definition at line 246 of file SolidMechanicsPlasticModel.C.

247{
248 return (dpm == 0 && yf <= _f_tol) || (dpm > -dpm_tol && yf <= _f_tol && yf >= -_f_tol);
249}

Referenced by SolidMechanicsPlasticTensileMulti::KuhnTuckerOK(), KuhnTuckerOK(), and SolidMechanicsPlasticModel::returnMap().

◆ modelName()

std::string SolidMechanicsPlasticMohrCoulombMulti::modelName ( ) const
overridevirtual

Implements SolidMechanicsPlasticModel.

Definition at line 321 of file SolidMechanicsPlasticMohrCoulombMulti.C.

322{
323 return "MohrCoulombMulti";
324}

◆ numberSurfaces()

unsigned int SolidMechanicsPlasticMohrCoulombMulti::numberSurfaces ( ) const
overridevirtual

The number of yield surfaces for this plasticity model.

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 83 of file SolidMechanicsPlasticMohrCoulombMulti.C.

84{
85 return 6;
86}

◆ perturbStress()

void SolidMechanicsPlasticMohrCoulombMulti::perturbStress ( const RankTwoTensor stress,
std::vector< Real > &  eigvals,
std::vector< RankTwoTensor > &  deigvals 
) const
private

perturbs the stress tensor in the case of almost-equal eigenvalues.

Note that, upon entry, this algorithm assumes that eigvals are the eigvenvalues of stress

Parameters
stressinput stress
[in,out]eigvalseigenvalues after perturbing.
[in,out]deigvalsd(eigenvalues)/dstress_ij after perturbing.

Definition at line 124 of file SolidMechanicsPlasticMohrCoulombMulti.C.

127{
128 Real small_perturbation;
129 RankTwoTensor shifted_stress = stress;
130 while (eigvals[0] > eigvals[1] - 0.1 * _shift || eigvals[1] > eigvals[2] - 0.1 * _shift)
131 {
132 for (unsigned i = 0; i < 3; ++i)
133 for (unsigned j = 0; j <= i; ++j)
134 {
135 small_perturbation = 0.1 * _shift * 2 * (MooseRandom::rand() - 0.5);
136 shifted_stress(i, j) += small_perturbation;
137 shifted_stress(j, i) += small_perturbation;
138 }
139 shifted_stress.dsymmetricEigenvalues(eigvals, deigvals);
140 }
141}
static Real rand()

Referenced by df_dsig(), and dflowPotential_dintnlV().

◆ phi()

Real SolidMechanicsPlasticMohrCoulombMulti::phi ( const Real  internal_param) const
protectedvirtual

phi as a function of residual value, rate, and internal_param

Definition at line 297 of file SolidMechanicsPlasticMohrCoulombMulti.C.

298{
299 return _phi.value(internal_param);
300}

Referenced by doReturnMap(), dyieldFunction_dintnlV(), dyieldFunction_dstressV(), returnEdge000101(), returnEdge010100(), returnPlane(), returnTip(), and yieldFunctionV().

◆ psi()

Real SolidMechanicsPlasticMohrCoulombMulti::psi ( const Real  internal_param) const
protectedvirtual

psi as a function of residual value, rate, and internal_param

Definition at line 309 of file SolidMechanicsPlasticMohrCoulombMulti.C.

310{
311 return _psi.value(internal_param);
312}

Referenced by dflowPotential_dintnlV(), dflowPotential_dstressV(), doReturnMap(), and flowPotentialV().

◆ returnEdge000101()

bool SolidMechanicsPlasticMohrCoulombMulti::returnEdge000101 ( const std::vector< Real > &  eigvals,
const std::vector< RealVectorValue > &  n,
std::vector< Real > &  dpm,
RankTwoTensor returned_stress,
Real  intnl_old,
Real &  sinphi,
Real &  cohcos,
Real  initial_guess,
Real  mag_E,
bool &  nr_converged,
Real  ep_plastic_tolerance,
std::vector< Real > &  yf 
) const
private

Tries to return-map to the MC edge using the n[4] and n[6] directions The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress) will be junk.

Parameters
eigvalsThe three stress eigenvalues, sorted in ascending order
nThe six return directions, n=E_ijkl*r. Note this algorithm assumes isotropic elasticity, so these are vectors in principal stress space
dpm[out]The six plastic multipliers resulting from the return-map to the edge
returned_stress[out]The returned stress. This will be diagonal, with the return-mapped eigenvalues in the diagonal positions, sorted in ascending order
intnl_oldThe internal parameter at stress=eigvals. This algorithm doesn't form the plastic strain, so you will have to use intnl=intnl_old+sum(dpm) if you need the new internal-parameter value at the returned point.
sinphi[out]The new value of sin(friction angle) after returning.
cohcos[out]The new value of cohesion*cos(friction angle) after returning.
mag_EAn approximate value for the magnitude of the Young's modulus. This is used to set appropriate tolerances in the Newton-Raphson procedure
initial_guessA guess for sum(dpm)
nr_converged[out]Whether the internal Newton-Raphson process converged
ep_plastic_toleranceThe user-set tolerance on the plastic strain.
yf[out]The yield functions after return

Definition at line 911 of file SolidMechanicsPlasticMohrCoulombMulti.C.

923{
924 // This returns to the Mohr-Coulomb edge
925 // with 000101 being active. This means that
926 // f3=f5=0. Denoting the returned stress by "a"
927 // (in principal stress space), this means that
928 // a0=a1 = (2Ccosphi + a2(1+sinphi))/(sinphi-1)
929 //
930 // Also, a = eigvals - dpm3*n3 - dpm5*n5,
931 // where the dpm are the plastic multipliers
932 // and the n are the flow directions.
933 //
934 // Hence we have 5 equations and 5 unknowns (a,
935 // dpm3 and dpm5). Eliminating all unknowns
936 // except for x = dpm3+dpm5 gives the residual below
937 // (I used mathematica)
938
939 Real x = initial_guess;
940 sinphi = std::sin(phi(intnl_old + x));
941 Real cosphi = std::cos(phi(intnl_old + x));
942 Real coh = cohesion(intnl_old + x);
943 cohcos = coh * cosphi;
944
945 const Real numer_const =
946 -n[3](2) * eigvals[0] - n[5](1) * eigvals[0] + n[5](2) * eigvals[0] + n[3](2) * eigvals[1] +
947 n[5](0) * eigvals[1] - n[5](2) * eigvals[1] - n[5](0) * eigvals[2] + n[5](1) * eigvals[2] +
948 n[3](0) * (-eigvals[1] + eigvals[2]) - n[3](1) * (-eigvals[0] + eigvals[2]);
949 const Real numer_coeff1 = 2 * (-n[3](0) + n[3](1) + n[5](0) - n[5](1));
950 const Real numer_coeff2 =
951 n[5](2) * (eigvals[0] - eigvals[1]) + n[3](2) * (-eigvals[0] + eigvals[1]) +
952 n[5](1) * (eigvals[0] + eigvals[2]) + (n[3](0) - n[5](0)) * (eigvals[1] + eigvals[2]) -
953 n[3](1) * (eigvals[0] + eigvals[2]);
954 Real numer = numer_const + numer_coeff1 * cohcos + numer_coeff2 * sinphi;
955 const Real denom_const = -n[3](2) * (n[5](0) - n[5](1)) - n[3](1) * (-n[5](0) + n[5](2)) +
956 n[3](0) * (-n[5](1) + n[5](2));
957 const Real denom_coeff = -n[3](2) * (n[5](0) - n[5](1)) - n[3](1) * (n[5](0) + n[5](2)) +
958 n[3](0) * (n[5](1) + n[5](2));
959 Real denom = denom_const + denom_coeff * sinphi;
960 Real residual = -x + numer / denom;
961
962 Real deriv_phi;
963 Real deriv_coh;
964 Real jacobian;
965 const Real tol = Utility::pow<2>(_f_tol / (mag_E * 10.0));
966 unsigned int iter = 0;
967 do
968 {
969 do
970 {
971 deriv_phi = dphi(intnl_old + dpm[3]);
972 deriv_coh = dcohesion(intnl_old + dpm[3]);
973 jacobian = -1 +
974 (numer_coeff1 * deriv_coh * cosphi - numer_coeff1 * coh * sinphi * deriv_phi +
975 numer_coeff2 * cosphi * deriv_phi) /
976 denom -
977 numer * denom_coeff * cosphi * deriv_phi / denom / denom;
978 x -= residual / jacobian;
979 if (iter > _max_iters) // not converging
980 {
981 nr_converged = false;
982 return false;
983 }
984
985 sinphi = std::sin(phi(intnl_old + x));
986 cosphi = std::cos(phi(intnl_old + x));
987 coh = cohesion(intnl_old + x);
988 cohcos = coh * cosphi;
989 numer = numer_const + numer_coeff1 * cohcos + numer_coeff2 * sinphi;
990 denom = denom_const + denom_coeff * sinphi;
991 residual = -x + numer / denom;
992 iter++;
993 } while (residual * residual > tol);
994
995 // now must ensure that yf[3] and yf[5] are both "zero"
996 const Real dpm3minusdpm5 =
997 (2.0 * (eigvals[0] - eigvals[1]) + x * (n[3](1) - n[3](0) + n[5](1) - n[5](0))) /
998 (n[3](0) - n[3](1) + n[5](1) - n[5](0));
999 dpm[3] = (x + dpm3minusdpm5) / 2.0;
1000 dpm[5] = (x - dpm3minusdpm5) / 2.0;
1001
1002 for (unsigned i = 0; i < 3; ++i)
1003 returned_stress(i, i) = eigvals[i] - dpm[3] * n[3](i) - dpm[5] * n[5](i);
1005 returned_stress(0, 0), returned_stress(1, 1), returned_stress(2, 2), sinphi, cohcos, yf);
1006 } while (yf[3] * yf[3] > _f_tol * _f_tol && yf[5] * yf[5] > _f_tol * _f_tol);
1007
1008 // so the NR process converged, but we must
1009 // check Kuhn-Tucker
1010 nr_converged = true;
1011
1012 if (dpm[3] < -ep_plastic_tolerance || dpm[5] < -ep_plastic_tolerance)
1013 // Kuhn-Tucker failure. This is a potential weak-point: if the user has set _f_tol quite
1014 // large, and ep_plastic_tolerance quite small, the above NR process will quickly converge, but
1015 // the solution may be wrong and violate Kuhn-Tucker.
1016 return false;
1017
1018 if (yf[0] > _f_tol || yf[1] > _f_tol || yf[2] > _f_tol || yf[4] > _f_tol)
1019 // Kuhn-Tucker failure
1020 return false;
1021
1022 // success
1023 dpm[0] = dpm[1] = dpm[2] = dpm[4] = 0;
1024 return true;
1025}
const double tol
const std::vector< double > x

Referenced by doReturnMap().

◆ returnEdge010100()

bool SolidMechanicsPlasticMohrCoulombMulti::returnEdge010100 ( const std::vector< Real > &  eigvals,
const std::vector< RealVectorValue > &  n,
std::vector< Real > &  dpm,
RankTwoTensor returned_stress,
Real  intnl_old,
Real &  sinphi,
Real &  cohcos,
Real  initial_guess,
Real  mag_E,
bool &  nr_converged,
Real  ep_plastic_tolerance,
std::vector< Real > &  yf 
) const
private

Tries to return-map to the MC edge using the n[1] and n[3] directions The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress) will be junk.

Parameters
eigvalsThe three stress eigenvalues, sorted in ascending order
nThe six return directions, n=E_ijkl*r. Note this algorithm assumes isotropic elasticity, so these are vectors in principal stress space
dpm[out]The six plastic multipliers resulting from the return-map to the edge
returned_stress[out]The returned stress. This will be diagonal, with the return-mapped eigenvalues in the diagonal positions, sorted in ascending order
intnl_oldThe internal parameter at stress=eigvals. This algorithm doesn't form the plastic strain, so you will have to use intnl=intnl_old+sum(dpm) if you need the new internal-parameter value at the returned point.
sinphi[out]The new value of sin(friction angle) after returning.
cohcos[out]The new value of cohesion*cos(friction angle) after returning.
mag_EAn approximate value for the magnitude of the Young's modulus. This is used to set appropriate tolerances in the Newton-Raphson procedure
initial_guessA guess for sum(dpm)
nr_converged[out]Whether the internal Newton-Raphson process converged
ep_plastic_toleranceThe user-set tolerance on the plastic strain.
yf[out]The yield functions after return

Definition at line 1028 of file SolidMechanicsPlasticMohrCoulombMulti.C.

1040{
1041 // This returns to the Mohr-Coulomb edge
1042 // with 010100 being active. This means that
1043 // f1=f3=0. Denoting the returned stress by "a"
1044 // (in principal stress space), this means that
1045 // a1=a2 = (2Ccosphi + a0(1-sinphi))/(sinphi+1)
1046 //
1047 // Also, a = eigvals - dpm1*n1 - dpm3*n3,
1048 // where the dpm are the plastic multipliers
1049 // and the n are the flow directions.
1050 //
1051 // Hence we have 5 equations and 5 unknowns (a,
1052 // dpm1 and dpm3). Eliminating all unknowns
1053 // except for x = dpm1+dpm3 gives the residual below
1054 // (I used mathematica)
1055
1056 Real x = initial_guess;
1057 sinphi = std::sin(phi(intnl_old + x));
1058 Real cosphi = std::cos(phi(intnl_old + x));
1059 Real coh = cohesion(intnl_old + x);
1060 cohcos = coh * cosphi;
1061
1062 const Real numer_const = -n[1](2) * eigvals[0] - n[3](1) * eigvals[0] + n[3](2) * eigvals[0] -
1063 n[1](0) * eigvals[1] + n[1](2) * eigvals[1] + n[3](0) * eigvals[1] -
1064 n[3](2) * eigvals[1] + n[1](0) * eigvals[2] - n[3](0) * eigvals[2] +
1065 n[3](1) * eigvals[2] - n[1](1) * (-eigvals[0] + eigvals[2]);
1066 const Real numer_coeff1 = 2 * (n[1](1) - n[1](2) - n[3](1) + n[3](2));
1067 const Real numer_coeff2 =
1068 n[3](2) * (-eigvals[0] - eigvals[1]) + n[1](2) * (eigvals[0] + eigvals[1]) +
1069 n[3](1) * (eigvals[0] + eigvals[2]) + (n[1](0) - n[3](0)) * (eigvals[1] - eigvals[2]) -
1070 n[1](1) * (eigvals[0] + eigvals[2]);
1071 Real numer = numer_const + numer_coeff1 * cohcos + numer_coeff2 * sinphi;
1072 const Real denom_const = -n[1](0) * (n[3](1) - n[3](2)) + n[1](2) * (-n[3](0) + n[3](1)) +
1073 n[1](1) * (-n[3](2) + n[3](0));
1074 const Real denom_coeff =
1075 n[1](0) * (n[3](1) - n[3](2)) + n[1](2) * (n[3](0) + n[3](1)) - n[1](1) * (n[3](0) + n[3](2));
1076 Real denom = denom_const + denom_coeff * sinphi;
1077 Real residual = -x + numer / denom;
1078
1079 Real deriv_phi;
1080 Real deriv_coh;
1081 Real jacobian;
1082 const Real tol = Utility::pow<2>(_f_tol / (mag_E * 10.0));
1083 unsigned int iter = 0;
1084 do
1085 {
1086 do
1087 {
1088 deriv_phi = dphi(intnl_old + dpm[3]);
1089 deriv_coh = dcohesion(intnl_old + dpm[3]);
1090 jacobian = -1 +
1091 (numer_coeff1 * deriv_coh * cosphi - numer_coeff1 * coh * sinphi * deriv_phi +
1092 numer_coeff2 * cosphi * deriv_phi) /
1093 denom -
1094 numer * denom_coeff * cosphi * deriv_phi / denom / denom;
1095 x -= residual / jacobian;
1096 if (iter > _max_iters) // not converging
1097 {
1098 nr_converged = false;
1099 return false;
1100 }
1101
1102 sinphi = std::sin(phi(intnl_old + x));
1103 cosphi = std::cos(phi(intnl_old + x));
1104 coh = cohesion(intnl_old + x);
1105 cohcos = coh * cosphi;
1106 numer = numer_const + numer_coeff1 * cohcos + numer_coeff2 * sinphi;
1107 denom = denom_const + denom_coeff * sinphi;
1108 residual = -x + numer / denom;
1109 iter++;
1110 } while (residual * residual > tol);
1111
1112 // now must ensure that yf[1] and yf[3] are both "zero"
1113 Real dpm1minusdpm3 =
1114 (2 * (eigvals[1] - eigvals[2]) + x * (n[1](2) - n[1](1) + n[3](2) - n[3](1))) /
1115 (n[1](1) - n[1](2) + n[3](2) - n[3](1));
1116 dpm[1] = (x + dpm1minusdpm3) / 2.0;
1117 dpm[3] = (x - dpm1minusdpm3) / 2.0;
1118
1119 for (unsigned i = 0; i < 3; ++i)
1120 returned_stress(i, i) = eigvals[i] - dpm[1] * n[1](i) - dpm[3] * n[3](i);
1122 returned_stress(0, 0), returned_stress(1, 1), returned_stress(2, 2), sinphi, cohcos, yf);
1123 } while (yf[1] * yf[1] > _f_tol * _f_tol && yf[3] * yf[3] > _f_tol * _f_tol);
1124
1125 // so the NR process converged, but we must
1126 // check Kuhn-Tucker
1127 nr_converged = true;
1128
1129 if (dpm[1] < -ep_plastic_tolerance || dpm[3] < -ep_plastic_tolerance)
1130 // Kuhn-Tucker failure. This is a potential weak-point: if the user has set _f_tol quite
1131 // large, and ep_plastic_tolerance quite small, the above NR process will quickly converge, but
1132 // the solution may be wrong and violate Kuhn-Tucker
1133 return false;
1134
1135 if (yf[0] > _f_tol || yf[2] > _f_tol || yf[4] > _f_tol || yf[5] > _f_tol)
1136 // Kuhn-Tucker failure
1137 return false;
1138
1139 // success
1140 dpm[0] = dpm[2] = dpm[4] = dpm[5] = 0;
1141 return true;
1142}

Referenced by doReturnMap().

◆ returnMap()

bool SolidMechanicsPlasticMohrCoulombMulti::returnMap ( const RankTwoTensor trial_stress,
Real  intnl_old,
const RankFourTensor E_ijkl,
Real  ep_plastic_tolerance,
RankTwoTensor returned_stress,
Real &  returned_intnl,
std::vector< Real > &  dpm,
RankTwoTensor delta_dp,
std::vector< Real > &  yf,
bool &  trial_stress_inadmissible 
) const
overridevirtual

Performs a custom return-map.

You may choose to over-ride this in your derived SolidMechanicsPlasticXXXX class, and you may implement the return-map algorithm in any way that suits you. Eg, using a Newton-Raphson approach, or a radial-return, etc. This may also be used as a quick way of ascertaining whether (trial_stress, intnl_old) is in fact admissible.

For over-riding this function, please note the following.

(1) Denoting the return value of the function by "successful_return", the only possible output values should be: (A) trial_stress_inadmissible=false, successful_return=true. That is, (trial_stress, intnl_old) is in fact admissible (in the elastic domain). (B) trial_stress_inadmissible=true, successful_return=false. That is (trial_stress, intnl_old) is inadmissible (outside the yield surface), and you didn't return to the yield surface. (C) trial_stress_inadmissible=true, successful_return=true. That is (trial_stress, intnl_old) is inadmissible (outside the yield surface), but you did return to the yield surface. The default implementation only handles case (A) and (B): it does not attempt to do a return-map algorithm.

(2) you must correctly signal "successful_return" using the return value of this function. Don't assume the calling function will do Kuhn-Tucker checking and so forth!

(3) In cases (A) and (B) you needn't set returned_stress, returned_intnl, delta_dp, or dpm. This is for computational efficiency.

(4) In cases (A) and (B), you MUST place the yield function values at (trial_stress, intnl_old) into yf so the calling function can use this information optimally. You will have already calculated these yield function values, which can be quite expensive, and it's not very optimal for the calling function to have to re-calculate them.

(5) In case (C), you need to set: returned_stress (the returned value of stress) returned_intnl (the returned value of the internal variable) delta_dp (the change in plastic strain) dpm (the plastic multipliers needed to bring about the return) yf (yield function values at the returned configuration)

(Note, if you over-ride returnMap, you will probably want to override consistentTangentOpertor too, otherwise it will default to E_ijkl.)

Parameters
trial_stressThe trial stress
intnl_oldValue of the internal parameter
E_ijklElasticity tensor
ep_plastic_toleranceTolerance defined by the user for the plastic strain
[out]returned_stressIn case (C): lies on the yield surface after returning and produces the correct plastic strain (normality condition). Otherwise: not defined
[out]returned_intnlIn case (C): the value of the internal parameter after returning. Otherwise: not defined
[out]dpmIn case (C): the plastic multipliers needed to bring about the return. Otherwise: not defined
[out]delta_dpIn case (C): The change in plastic strain induced by the return process. Otherwise: not defined
[out]yfIn case (C): the yield function at (returned_stress, returned_intnl). Otherwise: the yield function at (trial_stress, intnl_old)
[out]trial_stress_inadmissibleShould be set to false if the trial_stress is admissible, and true if the trial_stress is inadmissible. This can be used by the calling prorgram
Returns
true if a successful return (or a return-map not needed), false if the trial_stress is inadmissible but the return process failed

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 327 of file SolidMechanicsPlasticMohrCoulombMulti.C.

337{
339 return SolidMechanicsPlasticModel::returnMap(trial_stress,
340 intnl_old,
341 E_ijkl,
342 ep_plastic_tolerance,
343 returned_stress,
344 returned_intnl,
345 dpm,
346 delta_dp,
347 yf,
348 trial_stress_inadmissible);
349
350 return doReturnMap(trial_stress,
351 intnl_old,
352 E_ijkl,
353 ep_plastic_tolerance,
354 returned_stress,
355 returned_intnl,
356 dpm,
357 delta_dp,
358 yf,
359 trial_stress_inadmissible);
360}
virtual bool returnMap(const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const
Performs a custom return-map.

◆ returnPlane()

bool SolidMechanicsPlasticMohrCoulombMulti::returnPlane ( const std::vector< Real > &  eigvals,
const std::vector< RealVectorValue > &  n,
std::vector< Real > &  dpm,
RankTwoTensor returned_stress,
Real  intnl_old,
Real &  sinphi,
Real &  cohcos,
Real  initial_guess,
bool &  nr_converged,
Real  ep_plastic_tolerance,
std::vector< Real > &  yf 
) const
private

Tries to return-map to the MC plane using the n[3] direction The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress) will be junk.

Parameters
eigvalsThe three stress eigenvalues, sorted in ascending order
nThe six return directions, n=E_ijkl*r. Note this algorithm assumes isotropic elasticity, so these are vectors in principal stress space
dpm[out]The six plastic multipliers resulting from the return-map to the plane
returned_stress[out]The returned stress. This will be diagonal, with the return-mapped eigenvalues in the diagonal positions, sorted in ascending order
intnl_oldThe internal parameter at stress=eigvals. This algorithm doesn't form the plastic strain, so you will have to use intnl=intnl_old+sum(dpm) if you need the new internal-parameter value at the returned point.
sinphi[out]The new value of sin(friction angle) after returning.
cohcos[out]The new value of cohesion*cos(friction angle) after returning.
initial_guessA guess for sum(dpm)
nr_converged[out]Whether the internal Newton-Raphson process converged
ep_plastic_toleranceThe user-set tolerance on the plastic strain.
yf[out]The yield functions after return

Definition at line 806 of file SolidMechanicsPlasticMohrCoulombMulti.C.

817{
818 // This returns to the Mohr-Coulomb plane using n[3] (ie 000100)
819 //
820 // The returned point is defined by the f[3]=0 and
821 // a = eigvals - dpm[3]*n[3]
822 // where "a" is the returned point and dpm[3] is the plastic multiplier.
823 // This equation is a vector equation in principal stress space.
824 // (Kuhn-Tucker also demands that dpm[3]>=0, but we leave checking
825 // that condition for the end.)
826 // Since f[3]=0, we must have
827 // a[2]*(1+sinphi) + a[0]*(-1+sinphi) - 2*coh*cosphi = 0
828 // which gives dpm[3] as the solution of
829 // alpha*dpm[3] + eigvals[2] - eigvals[0] + beta*sinphi - 2*coh*cosphi = 0
830 // with alpha = n[3](0) - n[3](2) - (n[3](2) + n[3](0))*sinphi
831 // beta = eigvals[2] + eigvals[0]
832
833 mooseAssert(n.size() == 6,
834 "SolidMechanicsPlasticMohrCoulombMulti: Custom plane-return algorithm must be "
835 "supplied with n of size 6, whereas yours is "
836 << n.size());
837 mooseAssert(dpm.size() == 6,
838 "SolidMechanicsPlasticMohrCoulombMulti: Custom plane-return algorithm must be "
839 "supplied with dpm of size 6, whereas yours is "
840 << dpm.size());
841 mooseAssert(yf.size() == 6,
842 "SolidMechanicsPlasticMohrCoulombMulti: Custom tip-return algorithm must be "
843 "supplied with yf of size 6, whereas yours is "
844 << yf.size());
845
846 dpm[3] = initial_guess;
847 sinphi = std::sin(phi(intnl_old + dpm[3]));
848 Real cosphi = std::cos(phi(intnl_old + dpm[3]));
849 Real coh = cohesion(intnl_old + dpm[3]);
850 cohcos = coh * cosphi;
851
852 Real alpha = n[3](0) - n[3](2) - (n[3](2) + n[3](0)) * sinphi;
853 Real deriv_phi;
854 Real dalpha;
855 const Real beta = eigvals[2] + eigvals[0];
856 Real deriv_coh;
857
858 Real residual =
859 alpha * dpm[3] + eigvals[2] - eigvals[0] + beta * sinphi - 2.0 * cohcos; // this is 2*yf[3]
860 Real jacobian;
861
862 const Real f_tol2 = Utility::pow<2>(_f_tol);
863 unsigned int iter = 0;
864 do
865 {
866 deriv_phi = dphi(intnl_old + dpm[3]);
867 dalpha = -(n[3](2) + n[3](0)) * cosphi * deriv_phi;
868 deriv_coh = dcohesion(intnl_old + dpm[3]);
869 jacobian = alpha + dalpha * dpm[3] + beta * cosphi * deriv_phi - 2.0 * deriv_coh * cosphi +
870 2.0 * coh * sinphi * deriv_phi;
871
872 dpm[3] -= residual / jacobian;
873 if (iter > _max_iters) // not converging
874 {
875 nr_converged = false;
876 return false;
877 }
878
879 sinphi = std::sin(phi(intnl_old + dpm[3]));
880 cosphi = std::cos(phi(intnl_old + dpm[3]));
881 coh = cohesion(intnl_old + dpm[3]);
882 cohcos = coh * cosphi;
883 alpha = n[3](0) - n[3](2) - (n[3](2) + n[3](0)) * sinphi;
884 residual = alpha * dpm[3] + eigvals[2] - eigvals[0] + beta * sinphi - 2.0 * cohcos;
885 iter++;
886 } while (residual * residual > f_tol2);
887
888 // so the NR process converged, but we must
889 // check Kuhn-Tucker
890 nr_converged = true;
891 if (dpm[3] < -ep_plastic_tolerance)
892 // Kuhn-Tucker failure
893 return false;
894
895 for (unsigned i = 0; i < 3; ++i)
896 returned_stress(i, i) = eigvals[i] - dpm[3] * n[3](i);
898 returned_stress(0, 0), returned_stress(1, 1), returned_stress(2, 2), sinphi, cohcos, yf);
899
900 // by construction abs(yf[3]) = abs(residual/2) < _f_tol/2
901 if (yf[0] > _f_tol || yf[1] > _f_tol || yf[2] > _f_tol || yf[4] > _f_tol || yf[5] > _f_tol)
902 // Kuhn-Tucker failure
903 return false;
904
905 // success!
906 dpm[0] = dpm[1] = dpm[2] = dpm[4] = dpm[5] = 0;
907 return true;
908}

Referenced by doReturnMap().

◆ returnTip()

bool SolidMechanicsPlasticMohrCoulombMulti::returnTip ( const std::vector< Real > &  eigvals,
const std::vector< RealVectorValue > &  n,
std::vector< Real > &  dpm,
RankTwoTensor returned_stress,
Real  intnl_old,
Real &  sinphi,
Real &  cohcos,
Real  initial_guess,
bool &  nr_converged,
Real  ep_plastic_tolerance,
std::vector< Real > &  yf 
) const
private

Tries to return-map to the MC tip using the THREE directions given in n, and THREE dpm values are returned.

Note that you must supply THREE suitale n vectors out of the total of SIX flow directions, and then interpret the THREE dpm values appropriately. The return value is true if the internal Newton-Raphson process has converged and Kuhn-Tucker is satisfied, otherwise it is false If the return value is false and/or nr_converged=false then the "out" parameters (sinphi, cohcos, yf, returned_stress, dpm) will be junk.

Parameters
eigvalsThe three stress eigenvalues, sorted in ascending order
nThe three return directions, n=E_ijkl*r. Note this algorithm assumes isotropic elasticity, so these are 3 vectors in principal stress space
dpm[out]The three plastic multipliers resulting from the return-map to the tip.
returned_stress[out]The returned stress. This will be diagonal, with the return-mapped eigenvalues in the diagonal positions, sorted in ascending order
intnl_oldThe internal parameter at stress=eigvals. This algorithm doesn't form the plastic strain, so you will have to use intnl=intnl_old+sum(dpm) if you need the new internal-parameter value at the returned point.
sinphi[out]The new value of sin(friction angle) after returning.
cohcos[out]The new value of cohesion*cos(friction angle) after returning.
initial_guessA guess for sum(dpm)
nr_converged[out]Whether the internal Newton-Raphson process converged
ep_plastic_toleranceThe user-set tolerance on the plastic strain.
yf[out]The yield functions after return

Definition at line 638 of file SolidMechanicsPlasticMohrCoulombMulti.C.

649{
650 // This returns to the Mohr-Coulomb tip using the THREE directions
651 // given in n, and yields the THREE dpm values. Note that you
652 // must supply THREE suitable n vectors out of the total of SIX
653 // flow directions, and then interpret the THREE dpm values appropriately.
654 //
655 // Eg1. You supply the flow directions n[0], n[1] and n[3] as
656 // the "n" vectors. This is return-to-the-tip via 110100.
657 // Then the three returned dpm values will be dpm[0], dpm[1] and dpm[3].
658
659 // Eg2. You supply the flow directions n[1], n[3] and n[5] as
660 // the "n" vectors. This is return-to-the-tip via 010101.
661 // Then the three returned dpm values will be dpm[1], dpm[3] and dpm[5].
662
663 // The returned point is defined by the three yield functions (corresonding
664 // to the three supplied flow directions) all being zero.
665 // that is, returned_stress = diag(cohcot, cohcot, cohcot), where
666 // cohcot = cohesion*cosphi/sinphi
667 // where intnl = intnl_old + dpm[0] + dpm[1] + dpm[2]
668 // The 3 plastic multipliers, dpm, are defiend by the normality condition
669 // eigvals - cohcot = dpm[0]*n[0] + dpm[1]*n[1] + dpm[2]*n[2]
670 // (Kuhn-Tucker demands that all dpm are non-negative, but we leave
671 // that checking for the end.)
672 // (Remember that these "n" vectors and "dpm" values must be interpreted
673 // differently depending on what you pass into this function.)
674 // This is a vector equation with solution (A):
675 // dpm[0] = triple(eigvals - cohcot, n[1], n[2])/trip;
676 // dpm[1] = triple(eigvals - cohcot, n[2], n[0])/trip;
677 // dpm[2] = triple(eigvals - cohcot, n[0], n[1])/trip;
678 // where trip = triple(n[0], n[1], n[2]).
679 // By adding the three components of that solution together
680 // we can get an equation for x = dpm[0] + dpm[1] + dpm[2],
681 // and then our Newton-Raphson only involves one variable (x).
682 // In the following, i specialise to the isotropic situation.
683
684 mooseAssert(n.size() == 3,
685 "SolidMechanicsPlasticMohrCoulombMulti: Custom tip-return algorithm must be "
686 "supplied with n of size 3, whereas yours is "
687 << n.size());
688 mooseAssert(dpm.size() == 3,
689 "SolidMechanicsPlasticMohrCoulombMulti: Custom tip-return algorithm must be "
690 "supplied with dpm of size 3, whereas yours is "
691 << dpm.size());
692 mooseAssert(yf.size() == 6,
693 "SolidMechanicsPlasticMohrCoulombMulti: Custom tip-return algorithm must be "
694 "supplied with yf of size 6, whereas yours is "
695 << yf.size());
696
697 Real x = initial_guess;
698 const Real trip = triple_product(n[0], n[1], n[2]);
699 sinphi = std::sin(phi(intnl_old + x));
700 Real cosphi = std::cos(phi(intnl_old + x));
701 Real coh = cohesion(intnl_old + x);
702 cohcos = coh * cosphi;
703 Real cohcot = cohcos / sinphi;
704
705 if (_cohesion.modelName().compare("Constant") != 0 || _phi.modelName().compare("Constant") != 0)
706 {
707 // Finding x is expensive. Therefore
708 // although x!=0 for Constant Hardening, solution (A)
709 // demonstrates that we don't
710 // actually need to know x to find the dpm for
711 // Constant Hardening.
712 //
713 // However, for nontrivial Hardening, the following
714 // is necessary
715 // cohcot_coeff = [1,1,1].(Cross[n[1], n[2]] + Cross[n[2], n[0]] + Cross[n[0], n[1]])/trip
716 Real cohcot_coeff =
717 (n[0](0) * (n[1](1) - n[1](2) - n[2](1)) + (n[1](2) - n[1](1)) * n[2](0) +
718 (n[1](0) - n[1](2)) * n[2](1) + n[0](2) * (n[1](0) - n[1](1) - n[2](0) + n[2](1)) +
719 n[0](1) * (n[1](2) - n[1](0) + n[2](0) - n[2](2)) +
720 (n[0](0) - n[1](0) + n[1](1)) * n[2](2)) /
721 trip;
722 // eig_term = eigvals.(Cross[n[1], n[2]] + Cross[n[2], n[0]] + Cross[n[0], n[1]])/trip
723 Real eig_term = eigvals[0] *
724 (-n[0](2) * n[1](1) + n[0](1) * n[1](2) + n[0](2) * n[2](1) -
725 n[1](2) * n[2](1) - n[0](1) * n[2](2) + n[1](1) * n[2](2)) /
726 trip;
727 eig_term += eigvals[1] *
728 (n[0](2) * n[1](0) - n[0](0) * n[1](2) - n[0](2) * n[2](0) + n[1](2) * n[2](0) +
729 n[0](0) * n[2](2) - n[1](0) * n[2](2)) /
730 trip;
731 eig_term += eigvals[2] *
732 (n[0](0) * n[1](1) - n[1](1) * n[2](0) + n[0](1) * n[2](0) - n[0](1) * n[1](0) -
733 n[0](0) * n[2](1) + n[1](0) * n[2](1)) /
734 trip;
735 // and finally, the equation we want to solve is:
736 // x - eig_term + cohcot*cohcot_coeff = 0
737 // but i divide by cohcot_coeff so the result has the units of
738 // stress, so using _f_tol as a convergence check is reasonable
739 eig_term /= cohcot_coeff;
740 Real residual = x / cohcot_coeff - eig_term + cohcot;
741 Real jacobian;
742 Real deriv_phi;
743 Real deriv_coh;
744 unsigned int iter = 0;
745 do
746 {
747 deriv_phi = dphi(intnl_old + x);
748 deriv_coh = dcohesion(intnl_old + x);
749 jacobian = 1.0 / cohcot_coeff + deriv_coh * cosphi / sinphi -
750 coh * deriv_phi / Utility::pow<2>(sinphi);
751 x += -residual / jacobian;
752
753 if (iter > _max_iters) // not converging
754 {
755 nr_converged = false;
756 return false;
757 }
758
759 sinphi = std::sin(phi(intnl_old + x));
760 cosphi = std::cos(phi(intnl_old + x));
761 coh = cohesion(intnl_old + x);
762 cohcos = coh * cosphi;
763 cohcot = cohcos / sinphi;
764 residual = x / cohcot_coeff - eig_term + cohcot;
765 iter++;
766 } while (residual * residual > _f_tol * _f_tol / 100);
767 }
768
769 // so the NR process converged, but we must
770 // calculate the individual dpm values and
771 // check Kuhn-Tucker
772 nr_converged = true;
773 if (x < -3 * ep_plastic_tolerance)
774 // obviously at least one of the dpm are < -ep_plastic_tolerance. No point in proceeding. This
775 // is a potential weak-point: if the user has set _f_tol quite large, and ep_plastic_tolerance
776 // quite small, the above NR process will quickly converge, but the solution may be wrong and
777 // violate Kuhn-Tucker.
778 return false;
779
780 // The following is the solution (A) written above
781 // (dpm[0] = triple(eigvals - cohcot, n[1], n[2])/trip, etc)
782 // in the isotropic situation
784 v(0) = eigvals[0] - cohcot;
785 v(1) = eigvals[1] - cohcot;
786 v(2) = eigvals[2] - cohcot;
787 dpm[0] = triple_product(v, n[1], n[2]) / trip;
788 dpm[1] = triple_product(v, n[2], n[0]) / trip;
789 dpm[2] = triple_product(v, n[0], n[1]) / trip;
790
791 if (dpm[0] < -ep_plastic_tolerance || dpm[1] < -ep_plastic_tolerance ||
792 dpm[2] < -ep_plastic_tolerance)
793 // Kuhn-Tucker failure. No point in proceeding
794 return false;
795
796 // Kuhn-Tucker has succeeded: just need returned_stress and yf values
797 // I do not use the dpm to calculate returned_stress, because that
798 // might add error (and computational effort), simply:
799 returned_stress(0, 0) = returned_stress(1, 1) = returned_stress(2, 2) = cohcot;
800 // So by construction the yield functions are all zero
801 yf[0] = yf[1] = yf[2] = yf[3] = yf[4] = yf[5] = 0;
802 return true;
803}
const double v
virtual std::string modelName() const =0
T triple_product(const TypeVector< T > &a, const TypeVector< T > &b, const TypeVector< T > &c)
VectorValue< Real > RealVectorValue

Referenced by doReturnMap().

◆ useCustomCTO()

bool SolidMechanicsPlasticModel::useCustomCTO ( ) const
virtualinherited

Returns false. You will want to override this in your derived class if you write a custom consistent tangent operator function.

Reimplemented in SolidMechanicsPlasticDruckerPragerHyperbolic, SolidMechanicsPlasticJ2, SolidMechanicsPlasticMeanCapTC, and SolidMechanicsPlasticTensileMulti.

Definition at line 213 of file SolidMechanicsPlasticModel.C.

214{
215 return false;
216}

◆ useCustomReturnMap()

bool SolidMechanicsPlasticMohrCoulombMulti::useCustomReturnMap ( ) const
overridevirtual

Returns false. You will want to override this in your derived class if you write a custom returnMap function.

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 1164 of file SolidMechanicsPlasticMohrCoulombMulti.C.

1165{
1166 return _use_custom_returnMap;
1167}

◆ validParams()

InputParameters SolidMechanicsPlasticMohrCoulombMulti::validParams ( )
static

Definition at line 24 of file SolidMechanicsPlasticMohrCoulombMulti.C.

25{
27 params.addClassDescription("Non-associative Mohr-Coulomb plasticity with hardening/softening");
28 params.addRequiredParam<UserObjectName>(
29 "cohesion", "A SolidMechanicsHardening UserObject that defines hardening of the cohesion");
30 params.addRequiredParam<UserObjectName>("friction_angle",
31 "A SolidMechanicsHardening UserObject "
32 "that defines hardening of the "
33 "friction angle (in radians)");
34 params.addRequiredParam<UserObjectName>("dilation_angle",
35 "A SolidMechanicsHardening UserObject "
36 "that defines hardening of the "
37 "dilation angle (in radians)");
38 params.addParam<unsigned int>("max_iterations",
39 10,
40 "Maximum number of Newton-Raphson iterations "
41 "allowed in the custom return-map algorithm. "
42 " For highly nonlinear hardening this may "
43 "need to be higher than 10.");
44 params.addParam<Real>("shift",
45 "Yield surface is shifted by this amount to avoid problems with "
46 "defining derivatives when eigenvalues are equal. If this is "
47 "larger than f_tol, a warning will be issued. This may be set "
48 "very small when using the custom returnMap. Default = f_tol.");
49 params.addParam<bool>("use_custom_returnMap",
50 true,
51 "Use a custom return-map algorithm for this "
52 "plasticity model, which may speed up "
53 "computations considerably. Set to true "
54 "only for isotropic elasticity with no "
55 "hardening of the dilation angle. In this "
56 "case you may set 'shift' very small.");
57
58 return params;
59}
void addRequiredParam(const std::string &name, const std::string &doc_string)
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
void addClassDescription(const std::string &doc_string)
static InputParameters validParams()

◆ yieldFunction()

Real SolidMechanicsPlasticModel::yieldFunction ( const RankTwoTensor stress,
Real  intnl 
) const
protectedvirtualinherited

The following functions are what you should override when building single-plasticity models.

The yield function

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
Returns
the yield function

Reimplemented in SolidMechanicsPlasticDruckerPrager, SolidMechanicsPlasticDruckerPragerHyperbolic, SolidMechanicsPlasticIsotropicSD, SolidMechanicsPlasticJ2, SolidMechanicsPlasticMeanCap, SolidMechanicsPlasticMeanCapTC, SolidMechanicsPlasticMohrCoulomb, SolidMechanicsPlasticOrthotropic, SolidMechanicsPlasticSimpleTester, SolidMechanicsPlasticTensile, SolidMechanicsPlasticWeakPlaneShear, SolidMechanicsPlasticWeakPlaneTensile, and SolidMechanicsPlasticWeakPlaneTensileN.

Definition at line 61 of file SolidMechanicsPlasticModel.C.

62{
63 return 0.0;
64}

Referenced by SolidMechanicsPlasticModel::yieldFunctionV().

◆ yieldFunctionEigvals()

void SolidMechanicsPlasticMohrCoulombMulti::yieldFunctionEigvals ( Real  e0,
Real  e1,
Real  e2,
Real  sinphi,
Real  cohcos,
std::vector< Real > &  f 
) const
private

Calculates the yield functions given the eigenvalues of stress.

Parameters
e0Smallest eigenvalue
e1Middle eigenvalue
e2Largest eigenvalue
sinphisin(friction angle)
cohcoscohesion*cos(friction angle)
[out]fthe yield functions

Definition at line 106 of file SolidMechanicsPlasticMohrCoulombMulti.C.

108{
109 // Naively it seems a shame to have 6 yield functions active instead of just
110 // 3. But 3 won't do. Eg, think of a loading with eigvals[0]=eigvals[1]=eigvals[2]
111 // Then to return to the yield surface would require 2 positive plastic multipliers
112 // and one negative one. Boo hoo.
113
114 f.resize(6);
115 f[0] = 0.5 * (e0 - e1) + 0.5 * (e0 + e1) * sinphi - cohcos;
116 f[1] = 0.5 * (e1 - e0) + 0.5 * (e0 + e1) * sinphi - cohcos;
117 f[2] = 0.5 * (e0 - e2) + 0.5 * (e0 + e2) * sinphi - cohcos;
118 f[3] = 0.5 * (e2 - e0) + 0.5 * (e0 + e2) * sinphi - cohcos;
119 f[4] = 0.5 * (e1 - e2) + 0.5 * (e1 + e2) * sinphi - cohcos;
120 f[5] = 0.5 * (e2 - e1) + 0.5 * (e1 + e2) * sinphi - cohcos;
121}

Referenced by doReturnMap(), returnEdge000101(), returnEdge010100(), returnPlane(), and yieldFunctionV().

◆ yieldFunctionV()

void SolidMechanicsPlasticMohrCoulombMulti::yieldFunctionV ( const RankTwoTensor stress,
Real  intnl,
std::vector< Real > &  f 
) const
overridevirtual

Calculates the yield functions.

Note that for single-surface plasticity you don't want to override this - override the private yieldFunction below

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
[out]fthe yield functions

Reimplemented from SolidMechanicsPlasticModel.

Definition at line 89 of file SolidMechanicsPlasticMohrCoulombMulti.C.

92{
93 std::vector<Real> eigvals;
94 stress.symmetricEigenvalues(eigvals);
95 eigvals[0] += _shift;
96 eigvals[2] -= _shift;
97
98 const Real sinphi = std::sin(phi(intnl));
99 const Real cosphi = std::cos(phi(intnl));
100 const Real cohcos = cohesion(intnl) * cosphi;
101
102 yieldFunctionEigvals(eigvals[0], eigvals[1], eigvals[2], sinphi, cohcos, f);
103}

Member Data Documentation

◆ _cohesion

const SolidMechanicsHardeningModel& SolidMechanicsPlasticMohrCoulombMulti::_cohesion
private

Hardening model for cohesion.

Definition at line 95 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Referenced by cohesion(), dcohesion(), and returnTip().

◆ _f_tol

const Real SolidMechanicsPlasticModel::_f_tol
inherited

◆ _ic_tol

const Real SolidMechanicsPlasticModel::_ic_tol
inherited

Tolerance on internal constraint.

Definition at line 173 of file SolidMechanicsPlasticModel.h.

◆ _max_iters

const unsigned int SolidMechanicsPlasticMohrCoulombMulti::_max_iters
private

Maximum Newton-Raphison iterations in the custom returnMap algorithm.

Definition at line 104 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Referenced by returnEdge000101(), returnEdge010100(), returnPlane(), and returnTip().

◆ _phi

const SolidMechanicsHardeningModel& SolidMechanicsPlasticMohrCoulombMulti::_phi
private

Hardening model for phi.

Definition at line 98 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Referenced by dphi(), phi(), and returnTip().

◆ _psi

const SolidMechanicsHardeningModel& SolidMechanicsPlasticMohrCoulombMulti::_psi
private

Hardening model for psi.

Definition at line 101 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Referenced by dpsi(), and psi().

◆ _shift

const Real SolidMechanicsPlasticMohrCoulombMulti::_shift
private

yield function is shifted by this amount to avoid problems with stress-derivatives at equal eigenvalues

Definition at line 107 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Referenced by df_dsig(), dflowPotential_dintnlV(), doReturnMap(), dyieldFunction_dintnlV(), perturbStress(), SolidMechanicsPlasticMohrCoulombMulti(), and yieldFunctionV().

◆ _use_custom_returnMap

const bool SolidMechanicsPlasticMohrCoulombMulti::_use_custom_returnMap
private

Whether to use the custom return-map algorithm.

Definition at line 110 of file SolidMechanicsPlasticMohrCoulombMulti.h.

Referenced by returnMap(), and useCustomReturnMap().


The documentation for this class was generated from the following files: