www.mooseframework.org
Public Member Functions | Static Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
ComputeMultiPlasticityStress Class Reference

ComputeMultiPlasticityStress performs the return-map algorithm and associated stress updates for plastic models defined by a General User Objects. More...

#include <ComputeMultiPlasticityStress.h>

Inheritance diagram for ComputeMultiPlasticityStress:
[legend]

Public Member Functions

 ComputeMultiPlasticityStress (const InputParameters &parameters)
 
void outputAndCheckDebugParameters ()
 Outputs the debug parameters: _fspb_debug_stress, _fspd_debug_pm, etc and checks that they are sized correctly. More...
 
void checkDerivatives ()
 Checks the derivatives, eg dyieldFunction_dstress by using finite difference approximations. More...
 
void checkJacobian (const RankFourTensor &E_inv, const std::vector< Real > &intnl_old)
 Checks the full Jacobian, which is just certain linear combinations of the dyieldFunction_dstress, etc, by using finite difference approximations. More...
 
void checkSolution (const RankFourTensor &E_inv)
 Checks that Ax does equal b in the NR procedure. More...
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Types

enum  TangentOperatorEnum { elastic, linear, nonlinear }
 The type of tangent operator to return. tangent operator = d(stress_rate)/d(strain_rate). More...
 
enum  DeactivationSchemeEnum {
  optimized, safe, dumb, optimized_to_safe,
  safe_to_dumb, optimized_to_safe_to_dumb, optimized_to_dumb
}
 
enum  quickStep_called_from_t { computeQpStress_function, returnMap_function }
 The functions from which quickStep can be called. More...
 

Protected Member Functions

virtual void computeQpStress ()
 Compute the stress and store it in the _stress material property for the current quadrature point. More...
 
virtual void initQpStatefulProperties ()
 
virtual bool reinstateLinearDependentConstraints (std::vector< bool > &deactivated_due_to_ld)
 makes all deactivated_due_to_ld false, and if >0 of them were initially true, returns true More...
 
virtual unsigned int numberActive (const std::vector< bool > &active)
 counts the number of active constraints More...
 
virtual Real residual2 (const std::vector< Real > &pm, const std::vector< Real > &f, const RankTwoTensor &epp, const std::vector< Real > &ic, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld)
 The residual-squared. More...
 
virtual bool returnMap (const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &f, unsigned int &iter, bool can_revert_to_dumb, bool &linesearch_needed, bool &ld_encountered, bool &constraints_added, bool final_step, RankFourTensor &consistent_tangent_operator, std::vector< Real > &cumulative_pm)
 Implements the return map. More...
 
virtual bool lineSearch (Real &nr_res2, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, const RankFourTensor &E_inv, RankTwoTensor &delta_dp, const RankTwoTensor &dstress, const std::vector< Real > &dpm, const std::vector< Real > &dintnl, std::vector< Real > &f, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, bool &linesearch_needed)
 Performs a line search. More...
 
virtual bool singleStep (Real &nr_res2, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, RankTwoTensor &delta_dp, const RankFourTensor &E_inv, std::vector< Real > &f, RankTwoTensor &epp, std::vector< Real > &ic, std::vector< bool > &active, DeactivationSchemeEnum deactivation_scheme, bool &linesearch_needed, bool &ld_encountered)
 Performs a single Newton-Raphson + linesearch step Constraints are deactivated and the step is re-done if deactivation_scheme is set appropriately. More...
 
virtual bool checkAdmissible (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &all_f)
 Checks whether the yield functions are in the admissible region. More...
 
void buildDumbOrder (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< unsigned int > &dumb_order)
 Builds the order which "dumb" activation will take. More...
 
virtual void incrementDumb (int &dumb_iteration, const std::vector< unsigned int > &dumb_order, std::vector< bool > &act)
 Increments "dumb_iteration" by 1, and sets "act" appropriately (act[alpha] = true iff alpha_th bit of dumb_iteration == 1) More...
 
virtual bool checkKuhnTucker (const std::vector< Real > &f, const std::vector< Real > &pm, const std::vector< bool > &active)
 Checks Kuhn-Tucker conditions, and alters "active" if appropriate. More...
 
virtual void applyKuhnTucker (const std::vector< Real > &f, const std::vector< Real > &pm, std::vector< bool > &active)
 Checks Kuhn-Tucker conditions, and alters "active" if appropriate. More...
 
virtual void preReturnMap ()
 
virtual void postReturnMap ()
 
virtual bool quickStep (const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, std::vector< Real > &cumulative_pm, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &yf, unsigned int &iterations, RankFourTensor &consistent_tangent_operator, const quickStep_called_from_t called_from, bool final_step)
 Attempts to find an admissible (stress, intnl) by using the customized return-map algorithms defined through the TensorMechanicsPlasticXXXX.returnMap functions. More...
 
virtual bool plasticStep (const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &yf, unsigned int &iterations, bool &linesearch_needed, bool &ld_encountered, bool &constraints_added, RankFourTensor &consistent_tangent_operator)
 performs a plastic step More...
 
bool canChangeScheme (DeactivationSchemeEnum current_deactivation_scheme, bool can_revert_to_dumb)
 
bool canIncrementDumb (int dumb_iteration)
 
void changeScheme (const std::vector< bool > &initial_act, bool can_revert_to_dumb, const RankTwoTensor &initial_stress, const std::vector< Real > &intnl_old, DeactivationSchemeEnum &current_deactivation_scheme, std::vector< bool > &act, int &dumb_iteration, std::vector< unsigned int > &dumb_order)
 
bool canAddConstraints (const std::vector< bool > &act, const std::vector< Real > &all_f)
 
unsigned int activeCombinationNumber (const std::vector< bool > &act)
 
RankFourTensor consistentTangentOperator (const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &pm_this_step, const std::vector< Real > &cumulative_pm)
 Computes the consistent tangent operator (another name for the jacobian = d(stress_rate)/d(strain_rate) More...
 
virtual void computeQpProperties () override
 
virtual void calculateConstraints (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &f, std::vector< RankTwoTensor > &r, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active)
 The constraints. More...
 
virtual void calculateRHS (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &rhs, const std::vector< bool > &active, bool eliminate_ld, std::vector< bool > &deactivated_due_to_ld)
 Calculate the RHS which is rhs = -(epp(0,0), epp(1,0), epp(1,1), epp(2,0), epp(2,1), epp(2,2), f[0], f[1], ..., f[num_f], ic[0], ic[1], ..., ic[num_ic]) More...
 
virtual void calculateJacobian (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, std::vector< std::vector< Real >> &jac)
 d(rhs)/d(dof) More...
 
virtual void nrStep (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const RankTwoTensor &delta_dp, RankTwoTensor &dstress, std::vector< Real > &dpm, std::vector< Real > &dintnl, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
 Performs one Newton-Raphson step. More...
 
virtual void yieldFunction (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &f)
 The active yield function(s) More...
 
virtual void dyieldFunction_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &df_dstress)
 The derivative of the active yield function(s) with respect to stress. More...
 
virtual void dyieldFunction_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &df_dintnl)
 The derivative of active yield function(s) with respect to their internal parameters (the user objects assume there is exactly one internal param per yield function) More...
 
virtual void flowPotential (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &r)
 The active flow potential(s) - one for each yield function. More...
 
virtual void dflowPotential_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankFourTensor > &dr_dstress)
 The derivative of the active flow potential(s) with respect to stress. More...
 
virtual void dflowPotential_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dr_dintnl)
 The derivative of the active flow potentials with respect to the active internal parameters The UserObjects explicitly assume that r[alpha] is only dependent on intnl[alpha]. More...
 
virtual void hardPotential (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &h)
 The active hardening potentials (one for each internal parameter and for each yield function) by assumption in the Userobjects, the h[a][alpha] is nonzero only if the surface alpha is part of model a, so we only calculate those here. More...
 
virtual void dhardPotential_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dh_dstress)
 The derivative of the active hardening potentials with respect to stress By assumption in the Userobjects, the h[a][alpha] is nonzero only for a = alpha, so we only calculate those here. More...
 
virtual void dhardPotential_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &dh_dintnl)
 The derivative of the active hardening potentials with respect to the active internal parameters. More...
 
virtual void buildActiveConstraints (const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
 Constructs a set of active constraints, given the yield functions, f. More...
 
unsigned int modelNumber (unsigned int surface)
 returns the model number, given the surface number More...
 
bool anyActiveSurfaces (int model, const std::vector< bool > &active)
 returns true if any internal surfaces of the given model are active according to 'active' More...
 
void activeModelSurfaces (int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces_of_model)
 Returns the internal surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model. More...
 
void activeSurfaces (int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces)
 Returns the external surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model. More...
 
bool returnMapAll (const RankTwoTensor &trial_stress, const std::vector< Real > &intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &stress, std::vector< Real > &intnl, std::vector< Real > &pm, std::vector< Real > &cumulative_pm, RankTwoTensor &delta_dp, std::vector< Real > &yf, unsigned &num_successful_plastic_returns, unsigned &custom_model)
 Performs a returnMap for each plastic model using their inbuilt returnMap functions. More...
 

Protected Attributes

unsigned int _max_iter
 Maximum number of Newton-Raphson iterations allowed. More...
 
Real _min_stepsize
 Minimum fraction of applied strain that may be applied during adaptive stepsizing. More...
 
Real _max_stepsize_for_dumb
 "dumb" deactivation will only be used if the stepsize falls below this quantity More...
 
bool _ignore_failures
 Even if the returnMap fails, return the best values found for stress and internal parameters. More...
 
enum ComputeMultiPlasticityStress::TangentOperatorEnum _tangent_operator_type
 
Real _epp_tol
 Tolerance on the plastic strain increment ("direction") constraint. More...
 
std::vector< Real > _dummy_pm
 dummy "consistency parameters" (plastic multipliers) used in quickStep when called from computeQpStress More...
 
std::vector< Real > _cumulative_pm
 the sum of the plastic multipliers over all the sub-steps. More...
 
enum ComputeMultiPlasticityStress::DeactivationSchemeEnum _deactivation_scheme
 
bool _n_supplied
 User supplied the transverse direction vector. More...
 
RealVectorValue _n_input
 the supplied transverse direction vector More...
 
RealTensorValue _rot
 rotation matrix that takes _n to (0, 0, 1) More...
 
bool _perform_finite_strain_rotations
 whether to perform the rotations necessary in finite-strain simulations More...
 
const std::string _elasticity_tensor_name
 Name of the elasticity tensor material property. More...
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 Elasticity tensor material property. More...
 
MaterialProperty< RankTwoTensor > & _plastic_strain
 plastic strain More...
 
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
 Old value of plastic strain. More...
 
MaterialProperty< std::vector< Real > > & _intnl
 internal parameters More...
 
const MaterialProperty< std::vector< Real > > & _intnl_old
 old values of internal parameters More...
 
MaterialProperty< std::vector< Real > > & _yf
 yield functions More...
 
MaterialProperty< Real > & _iter
 Number of Newton-Raphson iterations used in the return-map. More...
 
MaterialProperty< Real > & _linesearch_needed
 Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise) More...
 
MaterialProperty< Real > & _ld_encountered
 Whether linear-dependence was encountered in the latest Newton-Raphson process (1 if true, 0 otherwise) More...
 
MaterialProperty< Real > & _constraints_added
 Whether constraints were added in during the latest Newton-Raphson process (1 if true, 0 otherwise) More...
 
MaterialProperty< RealVectorValue > & _n
 current value of transverse direction More...
 
const MaterialProperty< RealVectorValue > & _n_old
 old value of transverse direction More...
 
const MaterialProperty< RankTwoTensor > & _strain_increment
 strain increment (coming from ComputeIncrementalSmallStrain, for example) More...
 
const MaterialProperty< RankTwoTensor > & _total_strain_old
 Old value of total strain (coming from ComputeIncrementalSmallStrain, for example) More...
 
const MaterialProperty< RankTwoTensor > & _rotation_increment
 Rotation increment (coming from ComputeIncrementalSmallStrain, for example) More...
 
const MaterialProperty< RankTwoTensor > & _stress_old
 Old value of stress. More...
 
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
 Old value of elastic strain. More...
 
bool _cosserat
 whether Cosserat mechanics should be used More...
 
const MaterialProperty< RankTwoTensor > * _curvature
 The Cosserat curvature strain. More...
 
const MaterialProperty< RankFourTensor > * _elastic_flexural_rigidity_tensor
 The Cosserat elastic flexural rigidity tensor. More...
 
MaterialProperty< RankTwoTensor > * _couple_stress
 the Cosserat couple-stress More...
 
const MaterialProperty< RankTwoTensor > * _couple_stress_old
 the old value of Cosserat couple-stress More...
 
MaterialProperty< RankFourTensor > * _Jacobian_mult_couple
 derivative of couple-stress w.r.t. curvature More...
 
RankFourTensor _my_elasticity_tensor
 Elasticity tensor that can be rotated by this class (ie, its not const) More...
 
RankTwoTensor _my_strain_increment
 Strain increment that can be rotated by this class, and split into multiple increments (ie, its not const) More...
 
RankFourTensor _my_flexural_rigidity_tensor
 Flexual rigidity tensor that can be rotated by this class (ie, its not const) More...
 
RankTwoTensor _my_curvature
 Curvature that can be rotated by this class, and split into multiple increments (ie, its not const) More...
 
const std::string _base_name
 Base name prepended to all material property names to allow for multi-material systems. More...
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 Mechanical strain material property. More...
 
MaterialProperty< RankTwoTensor > & _stress
 Stress material property. More...
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 Elastic strain material property. More...
 
const MaterialProperty< RankTwoTensor > & _extra_stress
 Extra stress tensor. More...
 
std::vector< const Function * > _initial_stress_fcn
 initial stress components More...
 
MaterialProperty< RankFourTensor > & _Jacobian_mult
 derivative of stress w.r.t. strain (_dstress_dstrain) More...
 
MooseEnum _fspb_debug
 none - don't do any debugging crash - currently inactive jacobian - check the jacobian entries jacobian_and_linear_system - check entire jacobian and check that Ax=b More...
 
RankTwoTensor _fspb_debug_stress
 Debug the Jacobian entries at this stress. More...
 
std::vector< Real > _fspb_debug_pm
 Debug the Jacobian entires at these plastic multipliers. More...
 
std::vector< Real > _fspb_debug_intnl
 Debug the Jacobian entires at these internal parameters. More...
 
Real _fspb_debug_stress_change
 Debug finite-differencing parameter for the stress. More...
 
std::vector< Real > _fspb_debug_pm_change
 Debug finite-differencing parameters for the plastic multipliers. More...
 
std::vector< Real > _fspb_debug_intnl_change
 Debug finite-differencing parameters for the internal parameters. More...
 
Real _svd_tol
 Tolerance on the minimum ratio of singular values before flow-directions are deemed linearly dependent. More...
 
Real _min_f_tol
 Minimum value of the _f_tol parameters for the Yield Function User Objects. More...
 
const InputParameters & _params
 
unsigned int _num_models
 Number of plastic models for this material. More...
 
unsigned int _num_surfaces
 Number of surfaces within the plastic models. More...
 
std::vector< std::vector< unsigned int > > _surfaces_given_model
 _surfaces_given_model[model_number] = vector of surface numbers for this model More...
 
MooseEnum _specialIC
 Allows initial set of active constraints to be chosen optimally. More...
 
std::vector< const TensorMechanicsPlasticModel * > _f
 User objects that define the yield functions, flow potentials, etc. More...
 

Private Member Functions

RankTwoTensor rot (const RankTwoTensor &tens)
 
void fddyieldFunction_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &df_dstress)
 The finite-difference derivative of yield function(s) with respect to stress. More...
 
void fddyieldFunction_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &df_dintnl)
 The finite-difference derivative of yield function(s) with respect to internal parameter(s) More...
 
virtual void fddflowPotential_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankFourTensor > &dr_dstress)
 The finite-difference derivative of the flow potential(s) with respect to stress. More...
 
virtual void fddflowPotential_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &dr_dintnl)
 The finite-difference derivative of the flow potentials with respect to internal parameters. More...
 
virtual void fdJacobian (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, const RankFourTensor &E_inv, bool eliminate_ld, std::vector< std::vector< Real >> &jac)
 The Jacobian calculated using finite differences. More...
 
bool dof_included (unsigned int dof, const std::vector< bool > &deactivated_due_to_ld)
 
virtual int singularValuesOfR (const std::vector< RankTwoTensor > &r, std::vector< Real > &s)
 Performs a singular-value decomposition of r and returns the singular values. More...
 
virtual void eliminateLinearDependence (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &f, const std::vector< RankTwoTensor > &r, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
 Performs a number of singular-value decompositions to check for linear-dependence of the active directions "r" If linear dependence is found, then deactivated_due_to_ld will contain 'true' entries where surfaces need to be deactivated_due_to_ld. More...
 
void buildActiveConstraintsRock (const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
 "Rock" version Constructs a set of active constraints, given the yield functions, f. More...
 
void buildActiveConstraintsJoint (const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
 "Joint" version Constructs a set of active constraints, given the yield functions, f. More...
 

Private Attributes

std::vector< unsigned int > _model_given_surface
 given a surface number, this returns the model number More...
 
std::vector< unsigned int > _model_surface_given_surface
 given a surface number, this returns the corresponding-model's internal surface number More...
 

Detailed Description

ComputeMultiPlasticityStress performs the return-map algorithm and associated stress updates for plastic models defined by a General User Objects.

Note that if run in debug mode you might have to use the –no-trap-fpe flag because PETSc-LAPACK-BLAS explicitly compute 0/0 and 1/0, and this causes Libmesh to trap the floating-point exceptions

Definition at line 30 of file ComputeMultiPlasticityStress.h.

Member Enumeration Documentation

◆ DeactivationSchemeEnum

enum ComputeMultiPlasticityStress::DeactivationSchemeEnum
protected
Enumerator
optimized 
safe 
dumb 
optimized_to_safe 
safe_to_dumb 
optimized_to_safe_to_dumb 
optimized_to_dumb 

Definition at line 79 of file ComputeMultiPlasticityStress.h.

80  {
81  optimized,
82  safe,
83  dumb,
84  optimized_to_safe,
85  safe_to_dumb,
86  optimized_to_safe_to_dumb,
87  optimized_to_dumb
88  } _deactivation_scheme;

◆ quickStep_called_from_t

enum ComputeMultiPlasticityStress::quickStep_called_from_t
protected

The functions from which quickStep can be called.

Enumerator
computeQpStress_function 
returnMap_function 

Definition at line 445 of file ComputeMultiPlasticityStress.h.

446  {
447  computeQpStress_function,
448  returnMap_function
449  };

◆ TangentOperatorEnum

enum ComputeMultiPlasticityStress::TangentOperatorEnum
protected

The type of tangent operator to return. tangent operator = d(stress_rate)/d(strain_rate).

Enumerator
elastic 
linear 
nonlinear 

Definition at line 54 of file ComputeMultiPlasticityStress.h.

55  {
56  elastic,
57  linear,
58  nonlinear
59  } _tangent_operator_type;

Constructor & Destructor Documentation

◆ ComputeMultiPlasticityStress()

ComputeMultiPlasticityStress::ComputeMultiPlasticityStress ( const InputParameters &  parameters)

Definition at line 99 of file ComputeMultiPlasticityStress.C.

100  : ComputeStressBase(parameters),
101  MultiPlasticityDebugger(this),
102  _max_iter(getParam<unsigned int>("max_NR_iterations")),
103  _min_stepsize(getParam<Real>("min_stepsize")),
104  _max_stepsize_for_dumb(getParam<Real>("max_stepsize_for_dumb")),
105  _ignore_failures(getParam<bool>("ignore_failures")),
106 
107  _tangent_operator_type((TangentOperatorEnum)(int)getParam<MooseEnum>("tangent_operator")),
108 
109  _epp_tol(getParam<Real>("ep_plastic_tolerance")),
110 
111  _dummy_pm(0),
112 
113  _cumulative_pm(0),
114 
115  _deactivation_scheme((DeactivationSchemeEnum)(int)getParam<MooseEnum>("deactivation_scheme")),
116 
117  _n_supplied(parameters.isParamValid("transverse_direction")),
118  _n_input(_n_supplied ? getParam<RealVectorValue>("transverse_direction") : RealVectorValue()),
119  _rot(RealTensorValue()),
120 
121  _perform_finite_strain_rotations(getParam<bool>("perform_finite_strain_rotations")),
122 
123  _elasticity_tensor_name(_base_name + "elasticity_tensor"),
124  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
125  _plastic_strain(declareProperty<RankTwoTensor>("plastic_strain")),
126  _plastic_strain_old(getMaterialPropertyOld<RankTwoTensor>("plastic_strain")),
127  _intnl(declareProperty<std::vector<Real>>("plastic_internal_parameter")),
128  _intnl_old(getMaterialPropertyOld<std::vector<Real>>("plastic_internal_parameter")),
129  _yf(declareProperty<std::vector<Real>>("plastic_yield_function")),
130  _iter(declareProperty<Real>("plastic_NR_iterations")), // this is really an unsigned int, but
131  // for visualisation i convert it to Real
132  _linesearch_needed(declareProperty<Real>("plastic_linesearch_needed")), // this is really a
133  // boolean, but for
134  // visualisation i
135  // convert it to Real
136  _ld_encountered(declareProperty<Real>(
137  "plastic_linear_dependence_encountered")), // this is really a boolean, but for
138  // visualisation i convert it to Real
139  _constraints_added(declareProperty<Real>("plastic_constraints_added")), // this is really a
140  // boolean, but for
141  // visualisation i
142  // convert it to Real
143  _n(declareProperty<RealVectorValue>("plastic_transverse_direction")),
144  _n_old(getMaterialPropertyOld<RealVectorValue>("plastic_transverse_direction")),
145 
146  _strain_increment(getMaterialPropertyByName<RankTwoTensor>(_base_name + "strain_increment")),
147  _total_strain_old(getMaterialPropertyOldByName<RankTwoTensor>(_base_name + "total_strain")),
148  _rotation_increment(
149  getMaterialPropertyByName<RankTwoTensor>(_base_name + "rotation_increment")),
150 
151  _stress_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "stress")),
152  _elastic_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "elastic_strain")),
153 
154  // TODO: This design does NOT work. It makes these materials construction order dependent and it
155  // disregards block restrictions.
156  _cosserat(hasMaterialProperty<RankTwoTensor>("curvature") &&
157  hasMaterialProperty<RankFourTensor>("elastic_flexural_rigidity_tensor")),
158  _curvature(_cosserat ? &getMaterialPropertyByName<RankTwoTensor>("curvature") : nullptr),
159  _elastic_flexural_rigidity_tensor(
160  _cosserat ? &getMaterialPropertyByName<RankFourTensor>("elastic_flexural_rigidity_tensor")
161  : nullptr),
162  _couple_stress(_cosserat ? &declareProperty<RankTwoTensor>("couple_stress") : nullptr),
163  _couple_stress_old(_cosserat ? &getMaterialPropertyOld<RankTwoTensor>("couple_stress")
164  : nullptr),
165  _Jacobian_mult_couple(_cosserat ? &declareProperty<RankFourTensor>("couple_Jacobian_mult")
166  : nullptr),
167 
168  _my_elasticity_tensor(RankFourTensor()),
169  _my_strain_increment(RankTwoTensor()),
170  _my_flexural_rigidity_tensor(RankFourTensor()),
171  _my_curvature(RankTwoTensor())
172 {
173  if (_epp_tol <= 0)
174  mooseError("ComputeMultiPlasticityStress: ep_plastic_tolerance must be positive");
175 
176  if (_n_supplied)
177  {
178  // normalise the inputted transverse_direction
179  if (_n_input.norm() == 0)
180  mooseError(
181  "ComputeMultiPlasticityStress: transverse_direction vector must not have zero length");
182  else
183  _n_input /= _n_input.norm();
184  }
185 
186  if (_num_surfaces == 1)
187  _deactivation_scheme = safe;
188 }

Member Function Documentation

◆ activeCombinationNumber()

unsigned int ComputeMultiPlasticityStress::activeCombinationNumber ( const std::vector< bool > &  act)
protected

Definition at line 1574 of file ComputeMultiPlasticityStress.C.

1575 {
1576  unsigned num = 0;
1577  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1578  if (act[surface])
1579  num += (1 << surface); // (1 << x) = 2^x
1580 
1581  return num;
1582 }

Referenced by returnMap().

◆ activeModelSurfaces()

void MultiPlasticityRawComponentAssembler::activeModelSurfaces ( int  model,
const std::vector< bool > &  active,
std::vector< unsigned int > &  active_surfaces_of_model 
)
protectedinherited

Returns the internal surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model.

Parameters
modelthe model number
activearray with entries being 'true' if the surface is active
[out]active_surfaces_of_modelthe output

Definition at line 811 of file MultiPlasticityRawComponentAssembler.C.

815 {
816  active_surfaces_of_model.resize(0);
817  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
818  if (active[_surfaces_given_model[model][model_surface]])
819  active_surfaces_of_model.push_back(model_surface);
820 }

Referenced by MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ activeSurfaces()

void MultiPlasticityRawComponentAssembler::activeSurfaces ( int  model,
const std::vector< bool > &  active,
std::vector< unsigned int > &  active_surfaces 
)
protectedinherited

Returns the external surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model.

Parameters
modelthe model number
activearray with entries being 'true' if the surface is active
[out]active_surfacesthe output

Definition at line 800 of file MultiPlasticityRawComponentAssembler.C.

803 {
804  active_surfaces.resize(0);
805  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
806  if (active[_surfaces_given_model[model][model_surface]])
807  active_surfaces.push_back(_surfaces_given_model[model][model_surface]);
808 }

Referenced by MultiPlasticityLinearSystem::calculateConstraints().

◆ anyActiveSurfaces()

bool MultiPlasticityRawComponentAssembler::anyActiveSurfaces ( int  model,
const std::vector< bool > &  active 
)
protectedinherited

returns true if any internal surfaces of the given model are active according to 'active'

Definition at line 791 of file MultiPlasticityRawComponentAssembler.C.

792 {
793  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
794  if (active[_surfaces_given_model[model][model_surface]])
795  return true;
796  return false;
797 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityLinearSystem::calculateRHS(), MultiPlasticityDebugger::checkSolution(), MultiPlasticityDebugger::dof_included(), MultiPlasticityLinearSystem::nrStep(), and residual2().

◆ applyKuhnTucker()

void ComputeMultiPlasticityStress::applyKuhnTucker ( const std::vector< Real > &  f,
const std::vector< Real > &  pm,
std::vector< bool > &  active 
)
protectedvirtual

Checks Kuhn-Tucker conditions, and alters "active" if appropriate.

Do not let the simplicity of this routine fool you! Explicitly: (1) checks that pm = 0 for all the f < 0. If not, then active is set to false for that constraint. This may be triggered if upon exit of the NR loops a constraint got deactivated due to linear dependence, and then f<0 and its pm>0. (2) checks that pm = 0 for all inactive constraints. This should always be true unless someone has screwed with the code. (3) if any pm < 0, then active is set to false for that constraint. This may be triggered if _deactivation_scheme!="optimized".

Parameters
fvalues of the active yield functions
pmvalues of all the plastic multipliers
activethe active constraints (true if active)
Returns
return false if any of the Kuhn-Tucker conditions were violated (and hence the set of active constraints was changed)

Definition at line 1313 of file ComputeMultiPlasticityStress.C.

1316 {
1317  bool turned_off = false;
1318  unsigned ind = 0;
1319 
1320  // turn off all active surfaces that have f<0 and pm!=0
1321  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1322  {
1323  if (active[surface])
1324  {
1325  if (f[ind++] < -_f[modelNumber(surface)]->_f_tol)
1326  if (pm[surface] != 0)
1327  {
1328  turned_off = true;
1329  active[surface] = false;
1330  }
1331  }
1332  else if (pm[surface] != 0)
1333  mooseError("Crash due to plastic multiplier not being zero. This occurred because of poor "
1334  "coding!!");
1335  }
1336 
1337  // if didn't turn off anything yet, turn off surface with minimum pm
1338  if (!turned_off)
1339  {
1340  int surface_to_turn_off = -1;
1341  Real min_pm = 0;
1342  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1343  if (pm[surface] < min_pm)
1344  {
1345  min_pm = pm[surface];
1346  surface_to_turn_off = surface;
1347  }
1348  if (surface_to_turn_off >= 0)
1349  active[surface_to_turn_off] = false;
1350  }
1351 }

Referenced by returnMap().

◆ buildActiveConstraints()

void MultiPlasticityRawComponentAssembler::buildActiveConstraints ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act 
)
protectedvirtualinherited

Constructs a set of active constraints, given the yield functions, f.

This uses TensorMechanicsPlasticModel::activeConstraints to identify the active constraints for each model.

Parameters
fyield functions (should be _num_surfaces of these)
stressstress tensor
intnlinternal parameters
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actthe set of active constraints (will be resized to _num_surfaces)

Definition at line 344 of file MultiPlasticityRawComponentAssembler.C.

349 {
350  mooseAssert(f.size() == _num_surfaces,
351  "buildActiveConstraints called with f.size = " << f.size() << " while there are "
352  << _num_surfaces << " surfaces");
353  mooseAssert(intnl.size() == _num_models,
354  "buildActiveConstraints called with intnl.size = "
355  << intnl.size() << " while there are " << _num_models << " models");
356 
357  if (_specialIC == "rock")
358  buildActiveConstraintsRock(f, stress, intnl, Eijkl, act);
359  else if (_specialIC == "joint")
360  buildActiveConstraintsJoint(f, stress, intnl, Eijkl, act);
361  else // no specialIC
362  {
363  act.resize(0);
364  unsigned ind = 0;
365  for (unsigned model = 0; model < _num_models; ++model)
366  {
367  std::vector<Real> model_f(0);
368  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
369  model_f.push_back(f[ind++]);
370  std::vector<bool> model_act;
371  RankTwoTensor returned_stress;
372  _f[model]->activeConstraints(
373  model_f, stress, intnl[model], Eijkl, model_act, returned_stress);
374  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
375  act.push_back(model_act[model_surface]);
376  }
377  }
378 }

Referenced by returnMap().

◆ buildActiveConstraintsJoint()

void MultiPlasticityRawComponentAssembler::buildActiveConstraintsJoint ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act 
)
privateinherited

"Joint" version Constructs a set of active constraints, given the yield functions, f.

This uses TensorMechanicsPlasticModel::activeConstraints to identify the active constraints for each model.

Parameters
fyield functions (should be _num_surfaces of these)
stressstress tensor
intnlinternal parameters
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actthe set of active constraints (will be resized to _num_surfaces)

Definition at line 381 of file MultiPlasticityRawComponentAssembler.C.

386 {
387  act.assign(2, false);
388 
389  RankTwoTensor returned_stress;
390  std::vector<bool> active_tensile;
391  std::vector<bool> active_shear;
392  std::vector<Real> f_single;
393 
394  // first try tensile alone
395  f_single.assign(1, 0);
396  f_single[0] = f[0];
397  _f[0]->activeConstraints(f_single, stress, intnl[0], Eijkl, active_tensile, returned_stress);
398  _f[1]->yieldFunctionV(returned_stress, intnl[1], f_single);
399  if (f_single[0] <= _f[1]->_f_tol)
400  {
401  act[0] = active_tensile[0];
402  return;
403  }
404 
405  // next try shear alone
406  f_single.assign(1, 0);
407  f_single[0] = f[1];
408  _f[1]->activeConstraints(f_single, stress, intnl[1], Eijkl, active_shear, returned_stress);
409  _f[0]->yieldFunctionV(returned_stress, intnl[0], f_single);
410  if (f_single[0] <= _f[0]->_f_tol)
411  {
412  act[1] = active_shear[0];
413  return;
414  }
415 
416  // must be mixed
417  act[0] = act[1] = true;
418  return;
419 }

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints().

◆ buildActiveConstraintsRock()

void MultiPlasticityRawComponentAssembler::buildActiveConstraintsRock ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act 
)
privateinherited

"Rock" version Constructs a set of active constraints, given the yield functions, f.

This uses TensorMechanicsPlasticModel::activeConstraints to identify the active constraints for each model.

Parameters
fyield functions (should be _num_surfaces of these)
stressstress tensor
intnlinternal parameters
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actthe set of active constraints (will be resized to _num_surfaces)

Definition at line 422 of file MultiPlasticityRawComponentAssembler.C.

427 {
428  act.assign(9, false);
429 
430  RankTwoTensor returned_stress;
431  std::vector<bool> active_tensile;
432  std::vector<bool> active_MC;
433  std::vector<Real> f_single;
434 
435  // first try tensile alone
436  f_single.assign(3, 0);
437  f_single[0] = f[0];
438  f_single[1] = f[1];
439  f_single[2] = f[2];
440  _f[0]->activeConstraints(f_single, stress, intnl[0], Eijkl, active_tensile, returned_stress);
441  _f[1]->yieldFunctionV(returned_stress, intnl[1], f_single);
442  if (f_single[0] <= _f[1]->_f_tol && f_single[1] <= _f[1]->_f_tol &&
443  f_single[2] <= _f[1]->_f_tol && f_single[3] <= _f[1]->_f_tol &&
444  f_single[4] <= _f[1]->_f_tol && f_single[5] <= _f[1]->_f_tol)
445  {
446  act[0] = active_tensile[0];
447  act[1] = active_tensile[1];
448  act[2] = active_tensile[2];
449  return;
450  }
451 
452  // next try MC alone
453  f_single.assign(6, 0);
454  f_single[0] = f[3];
455  f_single[1] = f[4];
456  f_single[2] = f[5];
457  f_single[3] = f[6];
458  f_single[4] = f[7];
459  f_single[5] = f[8];
460  _f[1]->activeConstraints(f_single, stress, intnl[1], Eijkl, active_MC, returned_stress);
461  _f[0]->yieldFunctionV(returned_stress, intnl[0], f_single);
462  if (f_single[0] <= _f[0]->_f_tol && f_single[1] <= _f[0]->_f_tol && f_single[2] <= _f[0]->_f_tol)
463  {
464  act[3] = active_MC[0];
465  act[4] = active_MC[1];
466  act[5] = active_MC[2];
467  act[6] = active_MC[3];
468  act[7] = active_MC[4];
469  act[8] = active_MC[5];
470  return;
471  }
472 
473  // must be a mix.
474  // The possibilities are enumerated below.
475 
476  // tensile=edge, MC=tip (two possibilities)
477  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
478  active_MC[0] == true && active_MC[1] == true && active_MC[2] == false &&
479  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
480  {
481  act[1] = act[2] = act[6] = true;
482  act[4] = true;
483  return;
484  }
485  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
486  active_MC[0] == false && active_MC[1] == true && active_MC[2] == false &&
487  active_MC[3] == true && active_MC[4] == false && active_MC[5] == true)
488  {
489  act[1] = act[2] = act[6] = true; // i don't think act[4] is necessary, is it?!
490  return;
491  }
492 
493  // tensile = edge, MC=edge (two possibilities)
494  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
495  active_MC[0] == false && active_MC[1] == true && active_MC[2] == false &&
496  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
497  {
498  act[1] = act[2] = act[4] = act[6] = true;
499  return;
500  }
501  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
502  active_MC[0] == false && active_MC[1] == false && active_MC[2] == false &&
503  active_MC[3] == true && active_MC[4] == false && active_MC[5] == true)
504  {
505  act[1] = act[2] = act[4] = act[6] = true;
506  return;
507  }
508 
509  // tensile = edge, MC=face
510  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
511  active_MC[0] == false && active_MC[1] == false && active_MC[2] == false &&
512  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
513  {
514  act[1] = act[2] = act[6] = true;
515  return;
516  }
517 
518  // tensile = face, MC=tip (two possibilities)
519  if (active_tensile[0] == false && active_tensile[1] == false && active_tensile[2] == true &&
520  active_MC[0] == true && active_MC[1] == true && active_MC[2] == false &&
521  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
522  {
523  act[2] = act[6] = true;
524  act[4] = true;
525  act[8] = true;
526  return;
527  }
528  if (active_tensile[0] == false && active_tensile[1] == false && active_tensile[2] == true &&
529  active_MC[0] == false && active_MC[1] == true && active_MC[2] == false &&
530  active_MC[3] == true && active_MC[4] == false && active_MC[5] == true)
531  {
532  act[2] = act[6] = true;
533  act[8] = true;
534  return;
535  }
536 
537  // tensile = face, MC=face
538  if (active_tensile[0] == false && active_tensile[1] == false && active_tensile[2] == true &&
539  active_MC[0] == false && active_MC[1] == false && active_MC[2] == false &&
540  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
541  {
542  act[1] = act[2] = act[6] = true;
543  return;
544  }
545 
546  // tensile = face, MC=edge (two possibilites).
547  act[2] = true; // tensile face
548  act[3] = active_MC[0];
549  act[4] = active_MC[1];
550  act[5] = active_MC[2];
551  act[6] = active_MC[3];
552  act[7] = active_MC[4];
553  act[8] = active_MC[5];
554  return;
555 }

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints().

◆ buildDumbOrder()

void ComputeMultiPlasticityStress::buildDumbOrder ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< unsigned int > &  dumb_order 
)
protected

Builds the order which "dumb" activation will take.

Parameters
stressstress to evaluate yield functions and derivatives at
intnlinternal parameters to evaluate yield functions and derivatives at
[out]dumb_orderdumb_order[0] will be the yield surface furthest away from (stress, intnl), dumb_order[1] will be the next yield surface, etc. The distance measure used is f/|df_dstress|. This array can then be fed into incrementDumb in order to first try the yield surfaces which are farthest away from the (stress, intnl).

Definition at line 1524 of file ComputeMultiPlasticityStress.C.

1527 {
1528  if (dumb_order.size() != 0)
1529  return;
1530 
1531  std::vector<bool> act;
1532  act.assign(_num_surfaces, true);
1533 
1534  std::vector<Real> f;
1535  yieldFunction(stress, intnl, act, f);
1536  std::vector<RankTwoTensor> df_dstress;
1537  dyieldFunction_dstress(stress, intnl, act, df_dstress);
1538 
1539  typedef std::pair<Real, unsigned> pair_for_sorting;
1540  std::vector<pair_for_sorting> dist(_num_surfaces);
1541  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1542  {
1543  dist[surface].first = f[surface] / df_dstress[surface].L2norm();
1544  dist[surface].second = surface;
1545  }
1546  std::sort(dist.begin(), dist.end()); // sorted in ascending order of f/df_dstress
1547 
1548  dumb_order.resize(_num_surfaces);
1549  for (unsigned i = 0; i < _num_surfaces; ++i)
1550  dumb_order[i] = dist[_num_surfaces - 1 - i].second;
1551  // now dumb_order[0] is the surface with the greatest f/df_dstress
1552 }

Referenced by changeScheme(), and returnMap().

◆ calculateConstraints()

void MultiPlasticityLinearSystem::calculateConstraints ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankTwoTensor &  delta_dp,
std::vector< Real > &  f,
std::vector< RankTwoTensor > &  r,
RankTwoTensor &  epp,
std::vector< Real > &  ic,
const std::vector< bool > &  active 
)
protectedvirtualinherited

The constraints.

These are set to zero (or <=0 in the case of the yield functions) by the Newton-Raphson process, except in the case of linear-dependence which complicates things.

Parameters
stressThe stress
intnl_oldold values of the internal parameters
intnlinternal parameters
pmCurrent value(s) of the plasticity multiplier(s) (consistency parameters)
delta_dpChange in plastic strain incurred so far during the return
[out]fActive yield function(s)
[out]rActive flow directions
[out]eppPlastic-strain increment constraint
[out]icActive internal-parameter constraint
activeThe active constraints.

Definition at line 228 of file MultiPlasticityLinearSystem.C.

238 {
239  // see comments at the start of .h file
240 
241  mooseAssert(intnl_old.size() == _num_models,
242  "Size of intnl_old is " << intnl_old.size()
243  << " which is incorrect in calculateConstraints");
244  mooseAssert(intnl.size() == _num_models,
245  "Size of intnl is " << intnl.size() << " which is incorrect in calculateConstraints");
246  mooseAssert(pm.size() == _num_surfaces,
247  "Size of pm is " << pm.size() << " which is incorrect in calculateConstraints");
248  mooseAssert(active.size() == _num_surfaces,
249  "Size of active is " << active.size()
250  << " which is incorrect in calculateConstraints");
251 
252  // yield functions
253  yieldFunction(stress, intnl, active, f);
254 
255  // flow directions and "epp"
256  flowPotential(stress, intnl, active, r);
257  epp = RankTwoTensor();
258  unsigned ind = 0;
259  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
260  if (active[surface])
261  epp += pm[surface] * r[ind++]; // note, even the deactivated_due_to_ld must get added in
262  epp -= delta_dp;
263 
264  // internal constraints
265  std::vector<Real> h;
266  hardPotential(stress, intnl, active, h);
267  ic.resize(0);
268  ind = 0;
269  std::vector<unsigned int> active_surfaces;
270  std::vector<unsigned int>::iterator active_surface;
271  for (unsigned model = 0; model < _num_models; ++model)
272  {
273  activeSurfaces(model, active, active_surfaces);
274  if (active_surfaces.size() > 0)
275  {
276  // some surfaces are active in this model, so must form an internal constraint
277  ic.push_back(intnl[model] - intnl_old[model]);
278  for (active_surface = active_surfaces.begin(); active_surface != active_surfaces.end();
279  ++active_surface)
280  ic[ic.size() - 1] += pm[*active_surface] * h[ind++]; // we know the correct one is h[ind]
281  // since it was constructed in the same
282  // manner
283  }
284  }
285 }

Referenced by MultiPlasticityLinearSystem::calculateRHS(), and lineSearch().

◆ calculateJacobian()

void MultiPlasticityLinearSystem::calculateJacobian ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankFourTensor &  E_inv,
const std::vector< bool > &  active,
const std::vector< bool > &  deactivated_due_to_ld,
std::vector< std::vector< Real >> &  jac 
)
protectedvirtualinherited

d(rhs)/d(dof)

Definition at line 367 of file MultiPlasticityLinearSystem.C.

374 {
375  // see comments at the start of .h file
376 
377  mooseAssert(intnl.size() == _num_models,
378  "Size of intnl is " << intnl.size() << " which is incorrect in calculateJacobian");
379  mooseAssert(pm.size() == _num_surfaces,
380  "Size of pm is " << pm.size() << " which is incorrect in calculateJacobian");
381  mooseAssert(active.size() == _num_surfaces,
382  "Size of active is " << active.size() << " which is incorrect in calculateJacobian");
383  mooseAssert(deactivated_due_to_ld.size() == _num_surfaces,
384  "Size of deactivated_due_to_ld is " << deactivated_due_to_ld.size()
385  << " which is incorrect in calculateJacobian");
386 
387  unsigned ind = 0;
388  unsigned active_surface_ind = 0;
389 
390  std::vector<bool> active_surface(_num_surfaces); // active and not deactivated_due_to_ld
391  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
392  active_surface[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
393  unsigned num_active_surface = 0;
394  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
395  if (active_surface[surface])
396  num_active_surface++;
397 
398  std::vector<bool> active_model(
399  _num_models); // whether a model has surfaces that are active and not deactivated_due_to_ld
400  for (unsigned model = 0; model < _num_models; ++model)
401  active_model[model] = anyActiveSurfaces(model, active_surface);
402 
403  unsigned num_active_model = 0;
404  for (unsigned model = 0; model < _num_models; ++model)
405  if (active_model[model])
406  num_active_model++;
407 
408  ind = 0;
409  std::vector<unsigned int> active_model_index(_num_models);
410  for (unsigned model = 0; model < _num_models; ++model)
411  if (active_model[model])
412  active_model_index[model] = ind++;
413  else
414  active_model_index[model] =
415  _num_models + 1; // just a dummy, that will probably cause a crash if something goes wrong
416 
417  std::vector<RankTwoTensor> df_dstress;
418  dyieldFunction_dstress(stress, intnl, active_surface, df_dstress);
419 
420  std::vector<Real> df_dintnl;
421  dyieldFunction_dintnl(stress, intnl, active_surface, df_dintnl);
422 
423  std::vector<RankTwoTensor> r;
424  flowPotential(stress, intnl, active, r);
425 
426  std::vector<RankFourTensor> dr_dstress;
427  dflowPotential_dstress(stress, intnl, active, dr_dstress);
428 
429  std::vector<RankTwoTensor> dr_dintnl;
430  dflowPotential_dintnl(stress, intnl, active, dr_dintnl);
431 
432  std::vector<Real> h;
433  hardPotential(stress, intnl, active, h);
434 
435  std::vector<RankTwoTensor> dh_dstress;
436  dhardPotential_dstress(stress, intnl, active, dh_dstress);
437 
438  std::vector<Real> dh_dintnl;
439  dhardPotential_dintnl(stress, intnl, active, dh_dintnl);
440 
441  // d(epp)/dstress = sum_{active alpha} pm[alpha]*dr_dstress
442  RankFourTensor depp_dstress;
443  ind = 0;
444  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
445  if (active[surface]) // includes deactivated_due_to_ld
446  depp_dstress += pm[surface] * dr_dstress[ind++];
447  depp_dstress += E_inv;
448 
449  // d(epp)/dpm_{active_surface_index} = r_{active_surface_index}
450  std::vector<RankTwoTensor> depp_dpm;
451  depp_dpm.resize(num_active_surface);
452  ind = 0;
453  active_surface_ind = 0;
454  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
455  {
456  if (active[surface])
457  {
458  if (active_surface[surface]) // do not include the deactived_due_to_ld, since their pm are not
459  // dofs in the NR
460  depp_dpm[active_surface_ind++] = r[ind];
461  ind++;
462  }
463  }
464 
465  // d(epp)/dintnl_{active_model_index} = sum(pm[asdf]*dr_dintnl[fdsa])
466  std::vector<RankTwoTensor> depp_dintnl;
467  depp_dintnl.assign(num_active_model, RankTwoTensor());
468  ind = 0;
469  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
470  {
471  if (active[surface])
472  {
473  unsigned int model_num = modelNumber(surface);
474  if (active_model[model_num]) // only include models with surfaces which are still active after
475  // deactivated_due_to_ld
476  depp_dintnl[active_model_index[model_num]] += pm[surface] * dr_dintnl[ind];
477  ind++;
478  }
479  }
480 
481  // df_dstress has been calculated above
482  // df_dpm is always zero
483  // df_dintnl has been calculated above, but only the active_surface+active_model stuff needs to be
484  // included in Jacobian: see below
485 
486  std::vector<RankTwoTensor> dic_dstress;
487  dic_dstress.assign(num_active_model, RankTwoTensor());
488  ind = 0;
489  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
490  {
491  if (active[surface])
492  {
493  unsigned int model_num = modelNumber(surface);
494  if (active_model[model_num]) // only include ic for models with active_surface (ie, if model
495  // only contains deactivated_due_to_ld don't include it)
496  dic_dstress[active_model_index[model_num]] += pm[surface] * dh_dstress[ind];
497  ind++;
498  }
499  }
500 
501  std::vector<std::vector<Real>> dic_dpm;
502  dic_dpm.resize(num_active_model);
503  ind = 0;
504  active_surface_ind = 0;
505  for (unsigned model = 0; model < num_active_model; ++model)
506  dic_dpm[model].assign(num_active_surface, 0);
507  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
508  {
509  if (active[surface])
510  {
511  if (active_surface[surface]) // only take derivs wrt active-but-not-deactivated_due_to_ld pm
512  {
513  unsigned int model_num = modelNumber(surface);
514  // if (active_model[model_num]) // do not need this check as if the surface has
515  // active_surface, the model must be deemed active!
516  dic_dpm[active_model_index[model_num]][active_surface_ind] = h[ind];
517  active_surface_ind++;
518  }
519  ind++;
520  }
521  }
522 
523  std::vector<std::vector<Real>> dic_dintnl;
524  dic_dintnl.resize(num_active_model);
525  for (unsigned model = 0; model < num_active_model; ++model)
526  {
527  dic_dintnl[model].assign(num_active_model, 0);
528  dic_dintnl[model][model] = 1; // deriv wrt internal parameter
529  }
530  ind = 0;
531  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
532  {
533  if (active[surface])
534  {
535  unsigned int model_num = modelNumber(surface);
536  if (active_model[model_num]) // only the models that contain surfaces that are still active
537  // after deactivation_due_to_ld
538  dic_dintnl[active_model_index[model_num]][active_model_index[model_num]] +=
539  pm[surface] * dh_dintnl[ind];
540  ind++;
541  }
542  }
543 
544  unsigned int dim = 3;
545  unsigned int system_size =
546  6 + num_active_surface + num_active_model; // "6" comes from symmeterizing epp
547  jac.resize(system_size);
548  for (unsigned i = 0; i < system_size; ++i)
549  jac[i].assign(system_size, 0);
550 
551  unsigned int row_num = 0;
552  unsigned int col_num = 0;
553  for (unsigned i = 0; i < dim; ++i)
554  for (unsigned j = 0; j <= i; ++j)
555  {
556  for (unsigned k = 0; k < dim; ++k)
557  for (unsigned l = 0; l <= k; ++l)
558  jac[col_num][row_num++] =
559  depp_dstress(i, j, k, l) +
560  (k != l ? depp_dstress(i, j, l, k)
561  : 0); // extra part is needed because i assume dstress(i, j) = dstress(j, i)
562  for (unsigned surface = 0; surface < num_active_surface; ++surface)
563  jac[col_num][row_num++] = depp_dpm[surface](i, j);
564  for (unsigned a = 0; a < num_active_model; ++a)
565  jac[col_num][row_num++] = depp_dintnl[a](i, j);
566  row_num = 0;
567  col_num++;
568  }
569 
570  ind = 0;
571  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
572  if (active_surface[surface])
573  {
574  for (unsigned k = 0; k < dim; ++k)
575  for (unsigned l = 0; l <= k; ++l)
576  jac[col_num][row_num++] =
577  df_dstress[ind](k, l) +
578  (k != l ? df_dstress[ind](l, k)
579  : 0); // extra part is needed because i assume dstress(i, j) = dstress(j, i)
580  for (unsigned beta = 0; beta < num_active_surface; ++beta)
581  jac[col_num][row_num++] = 0; // df_dpm
582  for (unsigned model = 0; model < _num_models; ++model)
583  if (active_model[model]) // only use df_dintnl for models in active_model
584  {
585  if (modelNumber(surface) == model)
586  jac[col_num][row_num++] = df_dintnl[ind];
587  else
588  jac[col_num][row_num++] = 0;
589  }
590  ind++;
591  row_num = 0;
592  col_num++;
593  }
594 
595  for (unsigned a = 0; a < num_active_model; ++a)
596  {
597  for (unsigned k = 0; k < dim; ++k)
598  for (unsigned l = 0; l <= k; ++l)
599  jac[col_num][row_num++] =
600  dic_dstress[a](k, l) +
601  (k != l ? dic_dstress[a](l, k)
602  : 0); // extra part is needed because i assume dstress(i, j) = dstress(j, i)
603  for (unsigned alpha = 0; alpha < num_active_surface; ++alpha)
604  jac[col_num][row_num++] = dic_dpm[a][alpha];
605  for (unsigned b = 0; b < num_active_model; ++b)
606  jac[col_num][row_num++] = dic_dintnl[a][b];
607  row_num = 0;
608  col_num++;
609  }
610 
611  mooseAssert(col_num == system_size, "Incorrect filling of cols in Jacobian");
612 }

Referenced by MultiPlasticityDebugger::checkJacobian(), MultiPlasticityDebugger::checkSolution(), and MultiPlasticityLinearSystem::nrStep().

◆ calculateRHS()

void MultiPlasticityLinearSystem::calculateRHS ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankTwoTensor &  delta_dp,
std::vector< Real > &  rhs,
const std::vector< bool > &  active,
bool  eliminate_ld,
std::vector< bool > &  deactivated_due_to_ld 
)
protectedvirtualinherited

Calculate the RHS which is rhs = -(epp(0,0), epp(1,0), epp(1,1), epp(2,0), epp(2,1), epp(2,2), f[0], f[1], ..., f[num_f], ic[0], ic[1], ..., ic[num_ic])

Note that the 'epp' components only contain the upper diagonal. These contain flow directions and plasticity-multipliers for all active surfaces, even the deactivated_due_to_ld surfaces. Note that the 'f' components only contain the active and not deactivated_due_to_ld surfaces Note that the 'ic' components only contain the internal constraints for models which contain active and not deactivated_due_to_ld surfaces. They contain hardening-potentials and plasticity-multipliers for the active surfaces, even the deactivated_due_to_ld surfaces

Parameters
stressThe stress
intnl_oldold values of the internal parameters
intnlinternal parameters
pmCurrent value(s) of the plasticity multiplier(s) (consistency parameters)
delta_dpChange in plastic strain incurred so far during the return
[out]rhsthe rhs
activeThe active constraints.
eliminate_ldCheck for linear dependence of constraints and put the results into deactivated_due_to_ld. Usually this should be true, but for certain debug operations it should be false
[out]deactivated_due_to_ldconstraints deactivated due to linear-dependence of flow directions

Definition at line 288 of file MultiPlasticityLinearSystem.C.

297 {
298  // see comments at the start of .h file
299 
300  mooseAssert(intnl_old.size() == _num_models,
301  "Size of intnl_old is " << intnl_old.size() << " which is incorrect in calculateRHS");
302  mooseAssert(intnl.size() == _num_models,
303  "Size of intnl is " << intnl.size() << " which is incorrect in calculateRHS");
304  mooseAssert(pm.size() == _num_surfaces,
305  "Size of pm is " << pm.size() << " which is incorrect in calculateRHS");
306  mooseAssert(active.size() == _num_surfaces,
307  "Size of active is " << active.size() << " which is incorrect in calculateRHS");
308 
309  std::vector<Real> f; // the yield functions
310  RankTwoTensor epp; // the plastic-strain constraint ("direction constraint")
311  std::vector<Real> ic; // the "internal constraints"
312 
313  std::vector<RankTwoTensor> r;
314  calculateConstraints(stress, intnl_old, intnl, pm, delta_dp, f, r, epp, ic, active);
315 
316  if (eliminate_ld)
317  eliminateLinearDependence(stress, intnl, f, r, active, deactivated_due_to_ld);
318  else
319  deactivated_due_to_ld.assign(_num_surfaces, false);
320 
321  std::vector<bool> active_not_deact(_num_surfaces);
322  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
323  active_not_deact[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
324 
325  unsigned num_active_f = 0;
326  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
327  if (active_not_deact[surface])
328  num_active_f++;
329 
330  unsigned num_active_ic = 0;
331  for (unsigned model = 0; model < _num_models; ++model)
332  if (anyActiveSurfaces(model, active_not_deact))
333  num_active_ic++;
334 
335  unsigned int dim = 3;
336  unsigned int system_size = 6 + num_active_f + num_active_ic; // "6" comes from symmeterizing epp,
337  // num_active_f comes from "f",
338  // num_active_f comes from "ic"
339 
340  rhs.resize(system_size);
341 
342  unsigned ind = 0;
343  for (unsigned i = 0; i < dim; ++i)
344  for (unsigned j = 0; j <= i; ++j)
345  rhs[ind++] = -epp(i, j);
346  unsigned active_surface = 0;
347  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
348  if (active[surface])
349  {
350  if (!deactivated_due_to_ld[surface])
351  rhs[ind++] = -f[active_surface];
352  active_surface++;
353  }
354  unsigned active_model = 0;
355  for (unsigned model = 0; model < _num_models; ++model)
356  if (anyActiveSurfaces(model, active))
357  {
358  if (anyActiveSurfaces(model, active_not_deact))
359  rhs[ind++] = -ic[active_model];
360  active_model++;
361  }
362 
363  mooseAssert(ind == system_size, "Incorrect filling of the rhs in calculateRHS");
364 }

Referenced by MultiPlasticityDebugger::checkSolution(), MultiPlasticityDebugger::fdJacobian(), and MultiPlasticityLinearSystem::nrStep().

◆ canAddConstraints()

bool ComputeMultiPlasticityStress::canAddConstraints ( const std::vector< bool > &  act,
const std::vector< Real > &  all_f 
)
protected

Definition at line 1015 of file ComputeMultiPlasticityStress.C.

1017 {
1018  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1019  if (!act[surface] && (all_f[surface] > _f[modelNumber(surface)]->_f_tol))
1020  return true;
1021  return false;
1022 }

Referenced by returnMap().

◆ canChangeScheme()

bool ComputeMultiPlasticityStress::canChangeScheme ( DeactivationSchemeEnum  current_deactivation_scheme,
bool  can_revert_to_dumb 
)
protected

Definition at line 1025 of file ComputeMultiPlasticityStress.C.

1027 {
1028  if (current_deactivation_scheme == optimized && _deactivation_scheme == optimized_to_safe)
1029  return true;
1030 
1031  if (current_deactivation_scheme == optimized && _deactivation_scheme == optimized_to_safe_to_dumb)
1032  return true;
1033 
1034  if (current_deactivation_scheme == safe && _deactivation_scheme == safe_to_dumb &&
1035  can_revert_to_dumb)
1036  return true;
1037 
1038  if (current_deactivation_scheme == safe && _deactivation_scheme == optimized_to_safe_to_dumb &&
1039  can_revert_to_dumb)
1040  return true;
1041 
1042  if (current_deactivation_scheme == optimized && _deactivation_scheme == optimized_to_dumb &&
1043  can_revert_to_dumb)
1044  return true;
1045 
1046  return false;
1047 }

Referenced by returnMap().

◆ canIncrementDumb()

bool ComputeMultiPlasticityStress::canIncrementDumb ( int  dumb_iteration)
protected

Definition at line 1567 of file ComputeMultiPlasticityStress.C.

1568 {
1569  // (1 << _num_surfaces) = 2^_num_surfaces
1570  return ((dumb_iteration + 1) < (1 << _num_surfaces));
1571 }

Referenced by returnMap().

◆ changeScheme()

void ComputeMultiPlasticityStress::changeScheme ( const std::vector< bool > &  initial_act,
bool  can_revert_to_dumb,
const RankTwoTensor &  initial_stress,
const std::vector< Real > &  intnl_old,
DeactivationSchemeEnum &  current_deactivation_scheme,
std::vector< bool > &  act,
int &  dumb_iteration,
std::vector< unsigned int > &  dumb_order 
)
protected

Definition at line 1050 of file ComputeMultiPlasticityStress.C.

1058 {
1059  if (current_deactivation_scheme == optimized &&
1060  (_deactivation_scheme == optimized_to_safe ||
1061  _deactivation_scheme == optimized_to_safe_to_dumb))
1062  {
1063  current_deactivation_scheme = safe;
1064  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1065  act[surface] = initial_act[surface];
1066  }
1067  else if ((current_deactivation_scheme == safe &&
1068  (_deactivation_scheme == optimized_to_safe_to_dumb ||
1069  _deactivation_scheme == safe_to_dumb) &&
1070  can_revert_to_dumb) ||
1071  (current_deactivation_scheme == optimized && _deactivation_scheme == optimized_to_dumb &&
1072  can_revert_to_dumb))
1073  {
1074  current_deactivation_scheme = dumb;
1075  dumb_iteration = 0;
1076  buildDumbOrder(initial_stress, intnl_old, dumb_order);
1077  incrementDumb(dumb_iteration, dumb_order, act);
1078  }
1079 }

Referenced by returnMap().

◆ checkAdmissible()

bool ComputeMultiPlasticityStress::checkAdmissible ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< Real > &  all_f 
)
protectedvirtual

Checks whether the yield functions are in the admissible region.

Parameters
stressstress
intnlinternal parameters
[out]all_fthe values of all the yield functions
Returns
return false if any yield functions exceed their tolerance

Definition at line 1271 of file ComputeMultiPlasticityStress.C.

1274 {
1275  std::vector<bool> act;
1276  act.assign(_num_surfaces, true);
1277 
1278  yieldFunction(stress, intnl, act, all_f);
1279 
1280  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1281  if (all_f[surface] > _f[modelNumber(surface)]->_f_tol)
1282  return false;
1283 
1284  return true;
1285 }

Referenced by returnMap().

◆ checkDerivatives()

void MultiPlasticityDebugger::checkDerivatives ( )
inherited

Checks the derivatives, eg dyieldFunction_dstress by using finite difference approximations.

Definition at line 85 of file MultiPlasticityDebugger.C.

86 {
87  Moose::err
88  << "\n\n++++++++++++++++++++++++\nChecking the derivatives\n++++++++++++++++++++++++\n";
89  outputAndCheckDebugParameters();
90 
91  std::vector<bool> act;
92  act.assign(_num_surfaces, true);
93 
94  Moose::err << "\ndyieldFunction_dstress. Relative L2 norms.\n";
95  std::vector<RankTwoTensor> df_dstress;
96  std::vector<RankTwoTensor> fddf_dstress;
97  dyieldFunction_dstress(_fspb_debug_stress, _fspb_debug_intnl, act, df_dstress);
98  fddyieldFunction_dstress(_fspb_debug_stress, _fspb_debug_intnl, fddf_dstress);
99  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
100  {
101  Moose::err << "surface = " << surface << " Relative L2norm = "
102  << 2 * (df_dstress[surface] - fddf_dstress[surface]).L2norm() /
103  (df_dstress[surface] + fddf_dstress[surface]).L2norm()
104  << "\n";
105  Moose::err << "Coded:\n";
106  df_dstress[surface].print();
107  Moose::err << "Finite difference:\n";
108  fddf_dstress[surface].print();
109  }
110 
111  Moose::err << "\ndyieldFunction_dintnl.\n";
112  std::vector<Real> df_dintnl;
113  dyieldFunction_dintnl(_fspb_debug_stress, _fspb_debug_intnl, act, df_dintnl);
114  Moose::err << "Coded:\n";
115  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
116  Moose::err << df_dintnl[surface] << " ";
117  Moose::err << "\n";
118  std::vector<Real> fddf_dintnl;
119  fddyieldFunction_dintnl(_fspb_debug_stress, _fspb_debug_intnl, fddf_dintnl);
120  Moose::err << "Finite difference:\n";
121  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
122  Moose::err << fddf_dintnl[surface] << " ";
123  Moose::err << "\n";
124 
125  Moose::err << "\ndflowPotential_dstress. Relative L2 norms.\n";
126  std::vector<RankFourTensor> dr_dstress;
127  std::vector<RankFourTensor> fddr_dstress;
128  dflowPotential_dstress(_fspb_debug_stress, _fspb_debug_intnl, act, dr_dstress);
129  fddflowPotential_dstress(_fspb_debug_stress, _fspb_debug_intnl, fddr_dstress);
130  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
131  {
132  Moose::err << "surface = " << surface << " Relative L2norm = "
133  << 2 * (dr_dstress[surface] - fddr_dstress[surface]).L2norm() /
134  (dr_dstress[surface] + fddr_dstress[surface]).L2norm()
135  << "\n";
136  Moose::err << "Coded:\n";
137  dr_dstress[surface].print();
138  Moose::err << "Finite difference:\n";
139  fddr_dstress[surface].print();
140  }
141 
142  Moose::err << "\ndflowPotential_dintnl. Relative L2 norms.\n";
143  std::vector<RankTwoTensor> dr_dintnl;
144  std::vector<RankTwoTensor> fddr_dintnl;
145  dflowPotential_dintnl(_fspb_debug_stress, _fspb_debug_intnl, act, dr_dintnl);
146  fddflowPotential_dintnl(_fspb_debug_stress, _fspb_debug_intnl, fddr_dintnl);
147  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
148  {
149  Moose::err << "surface = " << surface << " Relative L2norm = "
150  << 2 * (dr_dintnl[surface] - fddr_dintnl[surface]).L2norm() /
151  (dr_dintnl[surface] + fddr_dintnl[surface]).L2norm()
152  << "\n";
153  Moose::err << "Coded:\n";
154  dr_dintnl[surface].print();
155  Moose::err << "Finite difference:\n";
156  fddr_dintnl[surface].print();
157  }
158 }

Referenced by initQpStatefulProperties(), and plasticStep().

◆ checkJacobian()

void MultiPlasticityDebugger::checkJacobian ( const RankFourTensor &  E_inv,
const std::vector< Real > &  intnl_old 
)
inherited

Checks the full Jacobian, which is just certain linear combinations of the dyieldFunction_dstress, etc, by using finite difference approximations.

Definition at line 161 of file MultiPlasticityDebugger.C.

163 {
164  Moose::err << "\n\n+++++++++++++++++++++\nChecking the Jacobian\n+++++++++++++++++++++\n";
165  outputAndCheckDebugParameters();
166 
167  std::vector<bool> act;
168  act.assign(_num_surfaces, true);
169  std::vector<bool> deactivated_due_to_ld;
170  deactivated_due_to_ld.assign(_num_surfaces, false);
171 
172  RankTwoTensor delta_dp = -E_inv * _fspb_debug_stress;
173 
174  std::vector<std::vector<Real>> jac;
175  calculateJacobian(_fspb_debug_stress,
176  _fspb_debug_intnl,
177  _fspb_debug_pm,
178  E_inv,
179  act,
180  deactivated_due_to_ld,
181  jac);
182 
183  std::vector<std::vector<Real>> fdjac;
184  fdJacobian(_fspb_debug_stress,
185  intnl_old,
186  _fspb_debug_intnl,
187  _fspb_debug_pm,
188  delta_dp,
189  E_inv,
190  false,
191  fdjac);
192 
193  Real L2_numer = 0;
194  Real L2_denom = 0;
195  for (unsigned row = 0; row < jac.size(); ++row)
196  for (unsigned col = 0; col < jac.size(); ++col)
197  {
198  L2_numer += Utility::pow<2>(jac[row][col] - fdjac[row][col]);
199  L2_denom += Utility::pow<2>(jac[row][col] + fdjac[row][col]);
200  }
201  Moose::err << "\nRelative L2norm = " << std::sqrt(L2_numer / L2_denom) / 0.5 << "\n";
202 
203  Moose::err << "\nHand-coded Jacobian:\n";
204  for (unsigned row = 0; row < jac.size(); ++row)
205  {
206  for (unsigned col = 0; col < jac.size(); ++col)
207  Moose::err << jac[row][col] << " ";
208  Moose::err << "\n";
209  }
210 
211  Moose::err << "Finite difference Jacobian:\n";
212  for (unsigned row = 0; row < fdjac.size(); ++row)
213  {
214  for (unsigned col = 0; col < fdjac.size(); ++col)
215  Moose::err << fdjac[row][col] << " ";
216  Moose::err << "\n";
217  }
218 }

Referenced by computeQpStress(), and plasticStep().

◆ checkKuhnTucker()

bool ComputeMultiPlasticityStress::checkKuhnTucker ( const std::vector< Real > &  f,
const std::vector< Real > &  pm,
const std::vector< bool > &  active 
)
protectedvirtual

Checks Kuhn-Tucker conditions, and alters "active" if appropriate.

Do not let the simplicity of this routine fool you! Explicitly: (1) checks that pm = 0 for all the f < 0. If not, then active is set to false for that constraint. This may be triggered if upon exit of the NR loops a constraint got deactivated due to linear dependence, and then f<0 and its pm>0. (2) checks that pm = 0 for all inactive constraints. This should always be true unless someone has screwed with the code. (3) if any pm < 0, then active is set to false for that constraint. This may be triggered if _deactivation_scheme!="optimized".

Parameters
fvalues of the active yield functions
pmvalues of all the plastic multipliers
activethe active constraints (true if active)
Returns
return false if any of the Kuhn-Tucker conditions were violated (and hence the set of active constraints was changed)

Definition at line 1288 of file ComputeMultiPlasticityStress.C.

1291 {
1292  unsigned ind = 0;
1293  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1294  {
1295  if (active[surface])
1296  {
1297  if (f[ind++] < -_f[modelNumber(surface)]->_f_tol)
1298  if (pm[surface] != 0)
1299  return false;
1300  }
1301  else if (pm[surface] != 0)
1302  mooseError("Crash due to plastic multiplier not being zero. This occurred because of poor "
1303  "coding!!");
1304  }
1305  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1306  if (pm[surface] < 0)
1307  return false;
1308 
1309  return true;
1310 }

Referenced by returnMap().

◆ checkSolution()

void MultiPlasticityDebugger::checkSolution ( const RankFourTensor &  E_inv)
inherited

Checks that Ax does equal b in the NR procedure.

Definition at line 367 of file MultiPlasticityDebugger.C.

368 {
369  Moose::err << "\n\n+++++++++++++++++++++\nChecking the Solution\n";
370  Moose::err << "(Ie, checking Ax = b)\n+++++++++++++++++++++\n";
371  outputAndCheckDebugParameters();
372 
373  std::vector<bool> act;
374  act.assign(_num_surfaces, true);
375  std::vector<bool> deactivated_due_to_ld;
376  deactivated_due_to_ld.assign(_num_surfaces, false);
377 
378  RankTwoTensor delta_dp = -E_inv * _fspb_debug_stress;
379 
380  std::vector<Real> orig_rhs;
381  calculateRHS(_fspb_debug_stress,
382  _fspb_debug_intnl,
383  _fspb_debug_intnl,
384  _fspb_debug_pm,
385  delta_dp,
386  orig_rhs,
387  act,
388  true,
389  deactivated_due_to_ld);
390 
391  Moose::err << "\nb = ";
392  for (unsigned i = 0; i < orig_rhs.size(); ++i)
393  Moose::err << orig_rhs[i] << " ";
394  Moose::err << "\n\n";
395 
396  std::vector<std::vector<Real>> jac_coded;
397  calculateJacobian(_fspb_debug_stress,
398  _fspb_debug_intnl,
399  _fspb_debug_pm,
400  E_inv,
401  act,
402  deactivated_due_to_ld,
403  jac_coded);
404 
405  Moose::err
406  << "Before checking Ax=b is correct, check that the Jacobians given below are equal.\n";
407  Moose::err
408  << "The hand-coded Jacobian is used in calculating the solution 'x', given 'b' above.\n";
409  Moose::err << "Note that this only includes degrees of freedom that aren't deactivated due to "
410  "linear dependence.\n";
411  Moose::err << "Hand-coded Jacobian:\n";
412  for (unsigned row = 0; row < jac_coded.size(); ++row)
413  {
414  for (unsigned col = 0; col < jac_coded.size(); ++col)
415  Moose::err << jac_coded[row][col] << " ";
416  Moose::err << "\n";
417  }
418 
419  deactivated_due_to_ld.assign(_num_surfaces,
420  false); // this potentially gets changed by nrStep, below
421  RankTwoTensor dstress;
422  std::vector<Real> dpm;
423  std::vector<Real> dintnl;
424  nrStep(_fspb_debug_stress,
425  _fspb_debug_intnl,
426  _fspb_debug_intnl,
427  _fspb_debug_pm,
428  E_inv,
429  delta_dp,
430  dstress,
431  dpm,
432  dintnl,
433  act,
434  deactivated_due_to_ld);
435 
436  std::vector<bool> active_not_deact(_num_surfaces);
437  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
438  active_not_deact[surface] = !deactivated_due_to_ld[surface];
439 
440  std::vector<Real> x;
441  x.assign(orig_rhs.size(), 0);
442  unsigned ind = 0;
443  for (unsigned i = 0; i < 3; ++i)
444  for (unsigned j = 0; j <= i; ++j)
445  x[ind++] = dstress(i, j);
446  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
447  if (active_not_deact[surface])
448  x[ind++] = dpm[surface];
449  for (unsigned model = 0; model < _num_models; ++model)
450  if (anyActiveSurfaces(model, active_not_deact))
451  x[ind++] = dintnl[model];
452 
453  mooseAssert(ind == orig_rhs.size(),
454  "Incorrect extracting of changes from NR solution in the "
455  "finite-difference checking of nrStep");
456 
457  Moose::err << "\nThis yields x =";
458  for (unsigned i = 0; i < orig_rhs.size(); ++i)
459  Moose::err << x[i] << " ";
460  Moose::err << "\n";
461 
462  std::vector<std::vector<Real>> jac_fd;
463  fdJacobian(_fspb_debug_stress,
464  _fspb_debug_intnl,
465  _fspb_debug_intnl,
466  _fspb_debug_pm,
467  delta_dp,
468  E_inv,
469  true,
470  jac_fd);
471 
472  Moose::err << "\nThe finite-difference Jacobian is used to multiply by this 'x',\n";
473  Moose::err << "in order to check that the solution is correct\n";
474  Moose::err << "Finite-difference Jacobian:\n";
475  for (unsigned row = 0; row < jac_fd.size(); ++row)
476  {
477  for (unsigned col = 0; col < jac_fd.size(); ++col)
478  Moose::err << jac_fd[row][col] << " ";
479  Moose::err << "\n";
480  }
481 
482  Real L2_numer = 0;
483  Real L2_denom = 0;
484  for (unsigned row = 0; row < jac_coded.size(); ++row)
485  for (unsigned col = 0; col < jac_coded.size(); ++col)
486  {
487  L2_numer += Utility::pow<2>(jac_coded[row][col] - jac_fd[row][col]);
488  L2_denom += Utility::pow<2>(jac_coded[row][col] + jac_fd[row][col]);
489  }
490  Moose::err << "Relative L2norm of the hand-coded and finite-difference Jacobian is "
491  << std::sqrt(L2_numer / L2_denom) / 0.5 << "\n";
492 
493  std::vector<Real> fd_times_x;
494  fd_times_x.assign(orig_rhs.size(), 0);
495  for (unsigned row = 0; row < orig_rhs.size(); ++row)
496  for (unsigned col = 0; col < orig_rhs.size(); ++col)
497  fd_times_x[row] += jac_fd[row][col] * x[col];
498 
499  Moose::err << "\n(Finite-difference Jacobian)*x =\n";
500  for (unsigned i = 0; i < orig_rhs.size(); ++i)
501  Moose::err << fd_times_x[i] << " ";
502  Moose::err << "\n";
503  Moose::err << "Recall that b = \n";
504  for (unsigned i = 0; i < orig_rhs.size(); ++i)
505  Moose::err << orig_rhs[i] << " ";
506  Moose::err << "\n";
507 
508  L2_numer = 0;
509  L2_denom = 0;
510  for (unsigned i = 0; i < orig_rhs.size(); ++i)
511  {
512  L2_numer += Utility::pow<2>(orig_rhs[i] - fd_times_x[i]);
513  L2_denom += Utility::pow<2>(orig_rhs[i] + fd_times_x[i]);
514  }
515  Moose::err << "\nRelative L2norm of these is " << std::sqrt(L2_numer / L2_denom) / 0.5 << "\n";
516 }

Referenced by computeQpStress(), and plasticStep().

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 50 of file ComputeStressBase.C.

51 {
52  computeQpStress();
53 
54  // Add in extra stress
55  _stress[_qp] += _extra_stress[_qp];
56 }

◆ computeQpStress()

void ComputeMultiPlasticityStress::computeQpStress ( )
protectedvirtual

Compute the stress and store it in the _stress material property for the current quadrature point.

Implements ComputeStressBase.

Definition at line 221 of file ComputeMultiPlasticityStress.C.

222 {
223  // the following "_my" variables can get rotated by preReturnMap and postReturnMap
224  _my_elasticity_tensor = _elasticity_tensor[_qp];
225  _my_strain_increment = _strain_increment[_qp];
226  if (_cosserat)
227  {
228  _my_flexural_rigidity_tensor = (*_elastic_flexural_rigidity_tensor)[_qp];
229  _my_curvature = (*_curvature)[_qp];
230  }
231 
232  if (_fspb_debug == "jacobian_and_linear_system")
233  {
234  // cannot do this at initQpStatefulProperties level since E_ijkl is not defined
235  checkJacobian(_elasticity_tensor[_qp].invSymm(), _intnl_old[_qp]);
236  checkSolution(_elasticity_tensor[_qp].invSymm());
237  mooseError("Finite-differencing completed. Exiting with no error");
238  }
239 
240  preReturnMap(); // do rotations to new frame if necessary
241 
242  unsigned int number_iterations;
243  bool linesearch_needed = false;
244  bool ld_encountered = false;
245  bool constraints_added = false;
246 
247  _cumulative_pm.assign(_num_surfaces, 0);
248  // try a "quick" return first - this can be purely elastic, or a customised plastic return defined
249  // by a TensorMechanicsPlasticXXXX UserObject
250  const bool found_solution = quickStep(rot(_stress_old[_qp]),
251  _stress[_qp],
252  _intnl_old[_qp],
253  _intnl[_qp],
254  _dummy_pm,
255  _cumulative_pm,
256  rot(_plastic_strain_old[_qp]),
257  _plastic_strain[_qp],
258  _my_elasticity_tensor,
259  _my_strain_increment,
260  _yf[_qp],
261  number_iterations,
262  _Jacobian_mult[_qp],
263  computeQpStress_function,
264  true);
265 
266  // if not purely elastic or the customised stuff failed, do some plastic return
267  if (!found_solution)
268  plasticStep(rot(_stress_old[_qp]),
269  _stress[_qp],
270  _intnl_old[_qp],
271  _intnl[_qp],
272  rot(_plastic_strain_old[_qp]),
273  _plastic_strain[_qp],
274  _my_elasticity_tensor,
275  _my_strain_increment,
276  _yf[_qp],
277  number_iterations,
278  linesearch_needed,
279  ld_encountered,
280  constraints_added,
281  _Jacobian_mult[_qp]);
282 
283  if (_cosserat)
284  {
285  (*_couple_stress)[_qp] = (*_elastic_flexural_rigidity_tensor)[_qp] * _my_curvature;
286  (*_Jacobian_mult_couple)[_qp] = _my_flexural_rigidity_tensor;
287  }
288 
289  postReturnMap(); // rotate back from new frame if necessary
290 
291  _iter[_qp] = 1.0 * number_iterations;
292  _linesearch_needed[_qp] = linesearch_needed;
293  _ld_encountered[_qp] = ld_encountered;
294  _constraints_added[_qp] = constraints_added;
295 
296  // Update measures of strain
297  _elastic_strain[_qp] = _elastic_strain_old[_qp] + _my_strain_increment -
298  (_plastic_strain[_qp] - _plastic_strain_old[_qp]);
299 
300  // Rotate the tensors to the current configuration
301  if (_perform_finite_strain_rotations)
302  {
303  _stress[_qp] = _rotation_increment[_qp] * _stress[_qp] * _rotation_increment[_qp].transpose();
304  _elastic_strain[_qp] =
305  _rotation_increment[_qp] * _elastic_strain[_qp] * _rotation_increment[_qp].transpose();
306  _plastic_strain[_qp] =
307  _rotation_increment[_qp] * _plastic_strain[_qp] * _rotation_increment[_qp].transpose();
308  }
309 }

◆ consistentTangentOperator()

RankFourTensor ComputeMultiPlasticityStress::consistentTangentOperator ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  E_ijkl,
const std::vector< Real > &  pm_this_step,
const std::vector< Real > &  cumulative_pm 
)
protected

Computes the consistent tangent operator (another name for the jacobian = d(stress_rate)/d(strain_rate)

The computations performed depend upon _tangent_operator_type

Parameters
stressThe value of stress after the return map algorithm has converged
intnlThe internal parameters after the return map has converged
E_ijklThe elasticity tensor (in the case of no plasticity this is the jacobian)
pm_this_stepThe plastic multipliers coming from the final strain increment. In many cases these will be equal to cumulative_pm, but in the case where the returnMap algorithm had to be performed in multiple substeps of smaller applied strain increments, pm_this_step are just the plastic multipliers for the final application of the strain incrment
cumulative_pmThe plastic multipliers needed for this current Return (this is the sum of the plastic multipliers over all substeps if the strain increment was applied in small substeps)

Definition at line 1585 of file ComputeMultiPlasticityStress.C.

1590 {
1591 
1592  if (_tangent_operator_type == elastic)
1593  return E_ijkl;
1594 
1595  // Typically act_at_some_step = act, but it is possible
1596  // that when subdividing a strain increment, a surface
1597  // is only active for one sub-step
1598  std::vector<bool> act_at_some_step(_num_surfaces);
1599  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1600  act_at_some_step[surface] = (cumulative_pm[surface] > 0);
1601 
1602  // "act" might contain surfaces that are linearly dependent
1603  // with others. Only the plastic multipliers that are > 0
1604  // for this strain increment need to be varied to find
1605  // the consistent tangent operator
1606  std::vector<bool> act_vary(_num_surfaces);
1607  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1608  act_vary[surface] = (pm_this_step[surface] > 0);
1609 
1610  std::vector<RankTwoTensor> df_dstress;
1611  dyieldFunction_dstress(stress, intnl, act_vary, df_dstress);
1612  std::vector<Real> df_dintnl;
1613  dyieldFunction_dintnl(stress, intnl, act_vary, df_dintnl);
1614  std::vector<RankTwoTensor> r;
1615  flowPotential(stress, intnl, act_vary, r);
1616  std::vector<RankFourTensor> dr_dstress_at_some_step;
1617  dflowPotential_dstress(stress, intnl, act_at_some_step, dr_dstress_at_some_step);
1618  std::vector<RankTwoTensor> dr_dintnl_at_some_step;
1619  dflowPotential_dintnl(stress, intnl, act_at_some_step, dr_dintnl_at_some_step);
1620  std::vector<Real> h;
1621  hardPotential(stress, intnl, act_vary, h);
1622 
1623  unsigned ind1;
1624  unsigned ind2;
1625 
1626  // r_minus_stuff[alpha] = r[alpha] -
1627  // pm_cumulatve[gamma]*dr[gamma]_dintnl[a]_at_some_step*h[a][alpha], with alpha only being in
1628  // act_vary, but gamma being act_at_some_step
1629  std::vector<RankTwoTensor> r_minus_stuff;
1630  ind1 = 0;
1631  for (unsigned surface1 = 0; surface1 < _num_surfaces; ++surface1)
1632  if (act_vary[surface1])
1633  {
1634  r_minus_stuff.push_back(r[ind1]);
1635  ind2 = 0;
1636  for (unsigned surface2 = 0; surface2 < _num_surfaces; ++surface2)
1637  if (act_at_some_step[surface2])
1638  {
1639  if (modelNumber(surface1) == modelNumber(surface2))
1640  {
1641  r_minus_stuff.back() -=
1642  cumulative_pm[surface2] * dr_dintnl_at_some_step[ind2] * h[ind1];
1643  }
1644  ind2++;
1645  }
1646  ind1++;
1647  }
1648 
1649  unsigned int num_currently_active = 0;
1650  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1651  if (act_vary[surface])
1652  num_currently_active += 1;
1653 
1654  // zzz is a matrix in the form that can be easily
1655  // inverted by MatrixTools::inverse
1656  // Eg for num_currently_active = 3
1657  // (zzz[0] zzz[1] zzz[2])
1658  // (zzz[3] zzz[4] zzz[5])
1659  // (zzz[6] zzz[7] zzz[8])
1660  std::vector<PetscScalar> zzz;
1661  zzz.assign(num_currently_active * num_currently_active, 0.0);
1662 
1663  ind1 = 0;
1664  RankTwoTensor r2;
1665  for (unsigned surface1 = 0; surface1 < _num_surfaces; ++surface1)
1666  if (act_vary[surface1])
1667  {
1668  ind2 = 0;
1669  for (unsigned surface2 = 0; surface2 < _num_surfaces; ++surface2)
1670  if (act_vary[surface2])
1671  {
1672  r2 = df_dstress[ind1] * (E_ijkl * r_minus_stuff[ind2]);
1673  zzz[ind1 * num_currently_active + ind2] += r2(0, 0) + r2(1, 1) + r2(2, 2);
1674  if (modelNumber(surface1) == modelNumber(surface2))
1675  zzz[ind1 * num_currently_active + ind2] += df_dintnl[ind1] * h[ind2];
1676  ind2++;
1677  }
1678  ind1++;
1679  }
1680 
1681  if (num_currently_active > 0)
1682  {
1683  // invert zzz, in place. if num_currently_active = 0 then zzz is not needed.
1684  try
1685  {
1686  MatrixTools::inverse(zzz, num_currently_active);
1687  }
1688  catch (const MooseException & e)
1689  {
1690  // in the very rare case of zzz being singular, just return the "elastic" tangent operator
1691  return E_ijkl;
1692  }
1693  }
1694 
1695  RankFourTensor strain_coeff = E_ijkl;
1696  ind1 = 0;
1697  for (unsigned surface1 = 0; surface1 < _num_surfaces; ++surface1)
1698  if (act_vary[surface1])
1699  {
1700  RankTwoTensor part1 = E_ijkl * r_minus_stuff[ind1];
1701  ind2 = 0;
1702  for (unsigned surface2 = 0; surface2 < _num_surfaces; ++surface2)
1703  if (act_vary[surface2])
1704  {
1705  RankTwoTensor part2 = E_ijkl * df_dstress[ind2];
1706  for (unsigned i = 0; i < 3; i++)
1707  for (unsigned j = 0; j < 3; j++)
1708  for (unsigned k = 0; k < 3; k++)
1709  for (unsigned l = 0; l < 3; l++)
1710  strain_coeff(i, j, k, l) -=
1711  part1(i, j) * part2(k, l) * zzz[ind1 * num_currently_active + ind2];
1712  ind2++;
1713  }
1714  ind1++;
1715  }
1716 
1717  if (_tangent_operator_type == linear)
1718  return strain_coeff;
1719 
1720  RankFourTensor stress_coeff(RankFourTensor::initIdentitySymmetricFour);
1721 
1722  RankFourTensor part3;
1723  ind1 = 0;
1724  for (unsigned surface1 = 0; surface1 < _num_surfaces; ++surface1)
1725  if (act_at_some_step[surface1])
1726  {
1727  part3 += cumulative_pm[surface1] * E_ijkl * dr_dstress_at_some_step[ind1];
1728  ind1++;
1729  }
1730 
1731  stress_coeff += part3;
1732 
1733  part3 = part3.transposeMajor(); // this is because below i want df_dstress[ind2]*part3, and this
1734  // equals (part3.transposeMajor())*df_dstress[ind2]
1735 
1736  ind1 = 0;
1737  for (unsigned surface1 = 0; surface1 < _num_surfaces; ++surface1)
1738  if (act_vary[surface1])
1739  {
1740  RankTwoTensor part1 = E_ijkl * r_minus_stuff[ind1];
1741  ind2 = 0;
1742  for (unsigned surface2 = 0; surface2 < _num_surfaces; ++surface2)
1743  if (act_vary[surface2])
1744  {
1745  RankTwoTensor part2 = part3 * df_dstress[ind2];
1746  for (unsigned i = 0; i < 3; i++)
1747  for (unsigned j = 0; j < 3; j++)
1748  for (unsigned k = 0; k < 3; k++)
1749  for (unsigned l = 0; l < 3; l++)
1750  stress_coeff(i, j, k, l) -=
1751  part1(i, j) * part2(k, l) * zzz[ind1 * num_currently_active + ind2];
1752  ind2++;
1753  }
1754  ind1++;
1755  }
1756 
1757  // need to find the inverse of stress_coeff, but remember
1758  // stress_coeff does not have the symmetries commonly found
1759  // in tensor mechanics:
1760  // stress_coeff(i, j, k, l) = stress_coeff(j, i, k, l) = stress_coeff(i, j, l, k) !=
1761  // stress_coeff(k, l, i, j)
1762  // (note the final not-equals). We want s_inv, such that
1763  // s_inv(i, j, m, n)*stress_coeff(m, n, k, l) = (de_ik*de_jl + de_il*de_jk)/2
1764  // where de_ij = 1 if i=j, and 0 otherwise.
1765  RankFourTensor s_inv;
1766  try
1767  {
1768  s_inv = stress_coeff.invSymm();
1769  }
1770  catch (const MooseException & e)
1771  {
1772  return strain_coeff; // when stress_coeff is singular (perhaps for incompressible plasticity?)
1773  // return the "linear" tangent operator
1774  }
1775 
1776  return s_inv * strain_coeff;
1777 }

Referenced by quickStep(), and returnMap().

◆ dflowPotential_dintnl()

void MultiPlasticityRawComponentAssembler::dflowPotential_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  dr_dintnl 
)
protectedvirtualinherited

The derivative of the active flow potentials with respect to the active internal parameters The UserObjects explicitly assume that r[alpha] is only dependent on intnl[alpha].

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dr_dintnl"
[out]dr_dintnlthe derivatives. dr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl[alpha]

Definition at line 235 of file MultiPlasticityRawComponentAssembler.C.

239 {
240  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
241  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
242 
243  dr_dintnl.resize(0);
244  std::vector<unsigned int> active_surfaces_of_model;
245  std::vector<unsigned int>::iterator active_surface;
246  std::vector<RankTwoTensor> model_dr_dintnl;
247  for (unsigned model = 0; model < _num_models; ++model)
248  {
249  activeModelSurfaces(model, active, active_surfaces_of_model);
250  if (active_surfaces_of_model.size() > 0)
251  {
252  _f[model]->dflowPotential_dintnlV(stress, intnl[model], model_dr_dintnl);
253  for (active_surface = active_surfaces_of_model.begin();
254  active_surface != active_surfaces_of_model.end();
255  ++active_surface)
256  dr_dintnl.push_back(model_dr_dintnl[*active_surface]);
257  }
258  }
259 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityDebugger::checkDerivatives(), and consistentTangentOperator().

◆ dflowPotential_dstress()

void MultiPlasticityRawComponentAssembler::dflowPotential_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankFourTensor > &  dr_dstress 
)
protectedvirtualinherited

The derivative of the active flow potential(s) with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dr_dstress"
[out]dr_dstressthe derivative. dr_dstress[alpha](i, j, k, l) = dr[alpha](i, j)/dstress(k, l)

Definition at line 207 of file MultiPlasticityRawComponentAssembler.C.

212 {
213  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
214  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
215 
216  dr_dstress.resize(0);
217  std::vector<unsigned int> active_surfaces_of_model;
218  std::vector<unsigned int>::iterator active_surface;
219  std::vector<RankFourTensor> model_dr_dstress;
220  for (unsigned model = 0; model < _num_models; ++model)
221  {
222  activeModelSurfaces(model, active, active_surfaces_of_model);
223  if (active_surfaces_of_model.size() > 0)
224  {
225  _f[model]->dflowPotential_dstressV(stress, intnl[model], model_dr_dstress);
226  for (active_surface = active_surfaces_of_model.begin();
227  active_surface != active_surfaces_of_model.end();
228  ++active_surface)
229  dr_dstress.push_back(model_dr_dstress[*active_surface]);
230  }
231  }
232 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityDebugger::checkDerivatives(), and consistentTangentOperator().

◆ dhardPotential_dintnl()

void MultiPlasticityRawComponentAssembler::dhardPotential_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  dh_dintnl 
)
protectedvirtualinherited

The derivative of the active hardening potentials with respect to the active internal parameters.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dh_dintnl"
[out]dh_dintnlthe derivatives. dh_dintnl[a][alpha][b] = dh[a][alpha]/dintnl[b]. Note that the userobjects assume that there is exactly one internal parameter per yield function, so the derivative is only nonzero for a=alpha=b, so that is all we calculate

Definition at line 317 of file MultiPlasticityRawComponentAssembler.C.

321 {
322  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
323  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
324 
325  dh_dintnl.resize(0);
326  std::vector<unsigned int> active_surfaces_of_model;
327  std::vector<unsigned int>::iterator active_surface;
328  std::vector<Real> model_dh_dintnl;
329  for (unsigned model = 0; model < _num_models; ++model)
330  {
331  activeModelSurfaces(model, active, active_surfaces_of_model);
332  if (active_surfaces_of_model.size() > 0)
333  {
334  _f[model]->dhardPotential_dintnlV(stress, intnl[model], model_dh_dintnl);
335  for (active_surface = active_surfaces_of_model.begin();
336  active_surface != active_surfaces_of_model.end();
337  ++active_surface)
338  dh_dintnl.push_back(model_dh_dintnl[*active_surface]);
339  }
340  }
341 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian().

◆ dhardPotential_dstress()

void MultiPlasticityRawComponentAssembler::dhardPotential_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  dh_dstress 
)
protectedvirtualinherited

The derivative of the active hardening potentials with respect to stress By assumption in the Userobjects, the h[a][alpha] is nonzero only for a = alpha, so we only calculate those here.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dh_dstress"
[out]dh_dstressthe derivative. dh_dstress[a](i, j) = dh[a]/dstress(k, l)

Definition at line 289 of file MultiPlasticityRawComponentAssembler.C.

294 {
295  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
296  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
297 
298  dh_dstress.resize(0);
299  std::vector<unsigned int> active_surfaces_of_model;
300  std::vector<unsigned int>::iterator active_surface;
301  std::vector<RankTwoTensor> model_dh_dstress;
302  for (unsigned model = 0; model < _num_models; ++model)
303  {
304  activeModelSurfaces(model, active, active_surfaces_of_model);
305  if (active_surfaces_of_model.size() > 0)
306  {
307  _f[model]->dhardPotential_dstressV(stress, intnl[model], model_dh_dstress);
308  for (active_surface = active_surfaces_of_model.begin();
309  active_surface != active_surfaces_of_model.end();
310  ++active_surface)
311  dh_dstress.push_back(model_dh_dstress[*active_surface]);
312  }
313  }
314 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian().

◆ dof_included()

bool MultiPlasticityDebugger::dof_included ( unsigned int  dof,
const std::vector< bool > &  deactivated_due_to_ld 
)
privateinherited

Definition at line 349 of file MultiPlasticityDebugger.C.

351 {
352  if (dof < unsigned(6))
353  // these are the stress components
354  return true;
355  unsigned eff_dof = dof - 6;
356  if (eff_dof < _num_surfaces)
357  // these are the plastic multipliers, pm
358  return !deactivated_due_to_ld[eff_dof];
359  eff_dof -= _num_surfaces; // now we know the dof is an intnl parameter
360  std::vector<bool> active_surface(_num_surfaces);
361  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
362  active_surface[surface] = !deactivated_due_to_ld[surface];
363  return anyActiveSurfaces(eff_dof, active_surface);
364 }

Referenced by MultiPlasticityDebugger::fdJacobian().

◆ dyieldFunction_dintnl()

void MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  df_dintnl 
)
protectedvirtualinherited

The derivative of active yield function(s) with respect to their internal parameters (the user objects assume there is exactly one internal param per yield function)

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "df_dintnl"
[out]df_dintnlthe derivatives. df_dstress[alpha] = dyieldFunction[alpha]/dintnl[alpha]

Definition at line 153 of file MultiPlasticityRawComponentAssembler.C.

157 {
158  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
159  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
160 
161  df_dintnl.resize(0);
162  std::vector<unsigned int> active_surfaces_of_model;
163  std::vector<unsigned int>::iterator active_surface;
164  std::vector<Real> model_df_dintnl;
165  for (unsigned model = 0; model < _num_models; ++model)
166  {
167  activeModelSurfaces(model, active, active_surfaces_of_model);
168  if (active_surfaces_of_model.size() > 0)
169  {
170  _f[model]->dyieldFunction_dintnlV(stress, intnl[model], model_df_dintnl);
171  for (active_surface = active_surfaces_of_model.begin();
172  active_surface != active_surfaces_of_model.end();
173  ++active_surface)
174  df_dintnl.push_back(model_df_dintnl[*active_surface]);
175  }
176  }
177 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityDebugger::checkDerivatives(), and consistentTangentOperator().

◆ dyieldFunction_dstress()

void MultiPlasticityRawComponentAssembler::dyieldFunction_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  df_dstress 
)
protectedvirtualinherited

The derivative of the active yield function(s) with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "df_dstress"
[out]df_dstressthe derivative (or derivatives in the case of multisurface plasticity). df_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)

Definition at line 125 of file MultiPlasticityRawComponentAssembler.C.

130 {
131  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
132  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
133 
134  df_dstress.resize(0);
135  std::vector<unsigned int> active_surfaces_of_model;
136  std::vector<unsigned int>::iterator active_surface;
137  std::vector<RankTwoTensor> model_df_dstress;
138  for (unsigned model = 0; model < _num_models; ++model)
139  {
140  activeModelSurfaces(model, active, active_surfaces_of_model);
141  if (active_surfaces_of_model.size() > 0)
142  {
143  _f[model]->dyieldFunction_dstressV(stress, intnl[model], model_df_dstress);
144  for (active_surface = active_surfaces_of_model.begin();
145  active_surface != active_surfaces_of_model.end();
146  ++active_surface)
147  df_dstress.push_back(model_df_dstress[*active_surface]);
148  }
149  }
150 }

Referenced by buildDumbOrder(), MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityDebugger::checkDerivatives(), consistentTangentOperator(), and MultiPlasticityLinearSystem::eliminateLinearDependence().

◆ eliminateLinearDependence()

void MultiPlasticityLinearSystem::eliminateLinearDependence ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< Real > &  f,
const std::vector< RankTwoTensor > &  r,
const std::vector< bool > &  active,
std::vector< bool > &  deactivated_due_to_ld 
)
privatevirtualinherited

Performs a number of singular-value decompositions to check for linear-dependence of the active directions "r" If linear dependence is found, then deactivated_due_to_ld will contain 'true' entries where surfaces need to be deactivated_due_to_ld.

Parameters
stressthe current stress
intnlthe current values of internal parameters
fActive yield function values
rthe flow directions that for those yield functions that are active upon entry to this function
activetrue if active
[out]deactivated_due_to_ldYield functions deactivated due to linearly-dependent flow directions

Definition at line 107 of file MultiPlasticityLinearSystem.C.

113 {
114  deactivated_due_to_ld.resize(_num_surfaces, false);
115 
116  unsigned num_active = r.size();
117 
118  if (num_active <= 1)
119  return;
120 
121  std::vector<double> s;
122  int info = singularValuesOfR(r, s);
123  if (info != 0)
124  mooseError("In finding the SVD in the return-map algorithm, the PETSC LAPACK gesvd routine "
125  "returned with error code ",
126  info);
127 
128  // num_lin_dep are the number of linearly dependent
129  // "r vectors", if num_active <= 6
130  unsigned int num_lin_dep = 0;
131 
132  unsigned i = s.size();
133  while (i-- > 0)
134  if (s[i] < _svd_tol * s[0])
135  num_lin_dep++;
136  else
137  break;
138 
139  if (num_lin_dep == 0 && num_active <= 6)
140  return;
141 
142  // From here on, some flow directions are linearly dependent
143 
144  // Find the signed "distance" of the current (stress, internal) configuration
145  // from the yield surfaces. This distance will not be precise, but
146  // i want to preferentially deactivate yield surfaces that are close
147  // to the current stress point.
148  std::vector<RankTwoTensor> df_dstress;
149  dyieldFunction_dstress(stress, intnl, active, df_dstress);
150 
151  typedef std::pair<Real, unsigned> pair_for_sorting;
152  std::vector<pair_for_sorting> dist(num_active);
153  for (unsigned i = 0; i < num_active; ++i)
154  {
155  dist[i].first = f[i] / df_dstress[i].L2norm();
156  dist[i].second = i;
157  }
158  std::sort(dist.begin(), dist.end()); // sorted in ascending order
159 
160  // There is a potential problem when we have equal f[i], for it can give oscillations
161  bool equals_detected = false;
162  for (unsigned i = 0; i < num_active - 1; ++i)
163  if (std::abs(dist[i].first - dist[i + 1].first) < _min_f_tol)
164  {
165  equals_detected = true;
166  dist[i].first += _min_f_tol * (MooseRandom::rand() - 0.5);
167  }
168  if (equals_detected)
169  std::sort(dist.begin(), dist.end()); // sorted in ascending order
170 
171  std::vector<bool> scheduled_for_deactivation;
172  scheduled_for_deactivation.assign(num_active, false);
173 
174  // In the following loop we go through all the flow directions, from
175  // the one with the largest dist, to the one with the smallest dist,
176  // adding them one-by-one into r_tmp. Upon each addition we check
177  // for linear-dependence. if LD is found, we schedule the most
178  // recently added flow direction for deactivation, and pop it
179  // back off r_tmp
180  unsigned current_yf;
181  current_yf = dist[num_active - 1].second;
182  // the one with largest dist
183  std::vector<RankTwoTensor> r_tmp = {r[current_yf]};
184 
185  unsigned num_kept_active = 1;
186  for (unsigned yf_to_try = 2; yf_to_try <= num_active; ++yf_to_try)
187  {
188  current_yf = dist[num_active - yf_to_try].second;
189  if (num_active == 2) // shortcut to we don't have to singularValuesOfR
190  scheduled_for_deactivation[current_yf] = true;
191  else if (num_kept_active >= 6) // shortcut to we don't have to singularValuesOfR: there can
192  // never be > 6 linearly-independent r vectors
193  scheduled_for_deactivation[current_yf] = true;
194  else
195  {
196  r_tmp.push_back(r[current_yf]);
197  info = singularValuesOfR(r_tmp, s);
198  if (info != 0)
199  mooseError("In finding the SVD in the return-map algorithm, the PETSC LAPACK gesvd routine "
200  "returned with error code ",
201  info);
202  if (s[s.size() - 1] < _svd_tol * s[0])
203  {
204  scheduled_for_deactivation[current_yf] = true;
205  r_tmp.pop_back();
206  num_lin_dep--;
207  }
208  else
209  num_kept_active++;
210  if (num_lin_dep == 0 && num_active <= 6)
211  // have taken out all the vectors that were linearly dependent
212  // so no point continuing
213  break;
214  }
215  }
216 
217  unsigned int old_active_number = 0;
218  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
219  if (active[surface])
220  {
221  if (scheduled_for_deactivation[old_active_number])
222  deactivated_due_to_ld[surface] = true;
223  old_active_number++;
224  }
225 }

Referenced by MultiPlasticityLinearSystem::calculateRHS().

◆ fddflowPotential_dintnl()

void MultiPlasticityDebugger::fddflowPotential_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< RankTwoTensor > &  dr_dintnl 
)
privatevirtualinherited

The finite-difference derivative of the flow potentials with respect to internal parameters.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
[out]dr_dintnlthe derivatives. dr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl[alpha]

Definition at line 607 of file MultiPlasticityDebugger.C.

610 {
611  dr_dintnl.resize(_num_surfaces);
612 
613  std::vector<bool> act;
614  act.assign(_num_surfaces, true);
615 
616  std::vector<RankTwoTensor> origr;
617  flowPotential(stress, intnl, act, origr);
618 
619  std::vector<Real> intnlep;
620  intnlep.resize(_num_models);
621  for (unsigned model = 0; model < _num_models; ++model)
622  intnlep[model] = intnl[model];
623  Real ep;
624  std::vector<RankTwoTensor> rep;
625  unsigned int model;
626  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
627  {
628  model = modelNumber(surface);
629  ep = _fspb_debug_intnl_change[model];
630  intnlep[model] += ep;
631  flowPotential(stress, intnlep, act, rep);
632  dr_dintnl[surface] = (rep[surface] - origr[surface]) / ep;
633  intnlep[model] -= ep;
634  }
635 }

Referenced by MultiPlasticityDebugger::checkDerivatives().

◆ fddflowPotential_dstress()

void MultiPlasticityDebugger::fddflowPotential_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< RankFourTensor > &  dr_dstress 
)
privatevirtualinherited

The finite-difference derivative of the flow potential(s) with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
[out]dr_dstressthe derivative. dr_dstress[alpha](i, j, k, l) = dr[alpha](i, j)/dstress(k, l)

Definition at line 578 of file MultiPlasticityDebugger.C.

581 {
582  dr_dstress.assign(_num_surfaces, RankFourTensor());
583 
584  std::vector<bool> act;
585  act.assign(_num_surfaces, true);
586 
587  Real ep = _fspb_debug_stress_change;
588  RankTwoTensor stressep;
589  std::vector<RankTwoTensor> rep, rep_minus;
590  for (unsigned i = 0; i < 3; ++i)
591  for (unsigned j = 0; j < 3; ++j)
592  {
593  stressep = stress;
594  // do a central difference
595  stressep(i, j) += ep / 2.0;
596  flowPotential(stressep, intnl, act, rep);
597  stressep(i, j) -= ep;
598  flowPotential(stressep, intnl, act, rep_minus);
599  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
600  for (unsigned k = 0; k < 3; ++k)
601  for (unsigned l = 0; l < 3; ++l)
602  dr_dstress[surface](k, l, i, j) = (rep[surface](k, l) - rep_minus[surface](k, l)) / ep;
603  }
604 }

Referenced by MultiPlasticityDebugger::checkDerivatives().

◆ fddyieldFunction_dintnl()

void MultiPlasticityDebugger::fddyieldFunction_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< Real > &  df_dintnl 
)
privateinherited

The finite-difference derivative of yield function(s) with respect to internal parameter(s)

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
[out]df_dintnlthe derivative (or derivatives in the case of multisurface plasticity). df_dintnl[alpha] = dyieldFunction[alpha]/dintnl[alpha]

Definition at line 547 of file MultiPlasticityDebugger.C.

550 {
551  df_dintnl.resize(_num_surfaces);
552 
553  std::vector<bool> act;
554  act.assign(_num_surfaces, true);
555 
556  std::vector<Real> origf;
557  yieldFunction(stress, intnl, act, origf);
558 
559  std::vector<Real> intnlep;
560  intnlep.resize(_num_models);
561  for (unsigned model = 0; model < _num_models; ++model)
562  intnlep[model] = intnl[model];
563  Real ep;
564  std::vector<Real> fep;
565  unsigned int model;
566  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
567  {
568  model = modelNumber(surface);
569  ep = _fspb_debug_intnl_change[model];
570  intnlep[model] += ep;
571  yieldFunction(stress, intnlep, act, fep);
572  df_dintnl[surface] = (fep[surface] - origf[surface]) / ep;
573  intnlep[model] -= ep;
574  }
575 }

Referenced by MultiPlasticityDebugger::checkDerivatives().

◆ fddyieldFunction_dstress()

void MultiPlasticityDebugger::fddyieldFunction_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< RankTwoTensor > &  df_dstress 
)
privateinherited

The finite-difference derivative of yield function(s) with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
[out]df_dstressthe derivative (or derivatives in the case of multisurface plasticity). df_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)

Definition at line 519 of file MultiPlasticityDebugger.C.

522 {
523  df_dstress.assign(_num_surfaces, RankTwoTensor());
524 
525  std::vector<bool> act;
526  act.assign(_num_surfaces, true);
527 
528  Real ep = _fspb_debug_stress_change;
529  RankTwoTensor stressep;
530  std::vector<Real> fep, fep_minus;
531  for (unsigned i = 0; i < 3; ++i)
532  for (unsigned j = 0; j < 3; ++j)
533  {
534  stressep = stress;
535  // do a central difference to attempt to capture discontinuities
536  // such as those encountered in tensile and Mohr-Coulomb
537  stressep(i, j) += ep / 2.0;
538  yieldFunction(stressep, intnl, act, fep);
539  stressep(i, j) -= ep;
540  yieldFunction(stressep, intnl, act, fep_minus);
541  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
542  df_dstress[surface](i, j) = (fep[surface] - fep_minus[surface]) / ep;
543  }
544 }

Referenced by MultiPlasticityDebugger::checkDerivatives().

◆ fdJacobian()

void MultiPlasticityDebugger::fdJacobian ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankTwoTensor &  delta_dp,
const RankFourTensor &  E_inv,
bool  eliminate_ld,
std::vector< std::vector< Real >> &  jac 
)
privatevirtualinherited

The Jacobian calculated using finite differences.

The output should be equal to calculateJacobian(...) if everything is coded correctly.

Parameters
stressthe stress at which to calculate the Jacobian
intnl_oldthe old values of internal variables (jacobian is inependent of these, but they are needed to do the finite-differencing cleanly)
intnlthe vector of internal parameters at which to calculate the Jacobian
pmthe plasticity multipliers at which to calculate the Jacobian
delta_dpplastic_strain - plastic_strain_old (Jacobian is independent of this, but it is needed to do the finite-differencing cleanly)
E_invinverse of the elasticity tensor
eliminate_ldonly calculate the Jacobian for the linearly independent constraints
[out]jacthe finite-difference Jacobian

Definition at line 221 of file MultiPlasticityDebugger.C.

229 {
230  std::vector<bool> active;
231  active.assign(_num_surfaces, true);
232 
233  std::vector<bool> deactivated_due_to_ld;
234  std::vector<bool> deactivated_due_to_ld_ep;
235 
236  std::vector<Real> orig_rhs;
237  calculateRHS(stress,
238  intnl_old,
239  intnl,
240  pm,
241  delta_dp,
242  orig_rhs,
243  active,
244  eliminate_ld,
245  deactivated_due_to_ld); // this calculates RHS, and also set deactivated_due_to_ld.
246  // The latter stays fixed for the rest of this routine
247 
248  unsigned int whole_system_size = 6 + _num_surfaces + _num_models;
249  unsigned int system_size =
250  orig_rhs.size(); // will be = whole_system_size if eliminate_ld = false, since all active=true
251  jac.resize(system_size);
252  for (unsigned row = 0; row < system_size; ++row)
253  jac[row].assign(system_size, 0);
254 
255  std::vector<Real> rhs_ep;
256  unsigned col = 0;
257 
258  RankTwoTensor stressep;
259  RankTwoTensor delta_dpep;
260  Real ep = _fspb_debug_stress_change;
261  for (unsigned i = 0; i < 3; ++i)
262  for (unsigned j = 0; j <= i; ++j)
263  {
264  stressep = stress;
265  stressep(i, j) += ep;
266  if (i != j)
267  stressep(j, i) += ep;
268  delta_dpep = delta_dp;
269  for (unsigned k = 0; k < 3; ++k)
270  for (unsigned l = 0; l < 3; ++l)
271  {
272  delta_dpep(k, l) -= E_inv(k, l, i, j) * ep;
273  if (i != j)
274  delta_dpep(k, l) -= E_inv(k, l, j, i) * ep;
275  }
276  active.assign(_num_surfaces, true);
277  calculateRHS(stressep,
278  intnl_old,
279  intnl,
280  pm,
281  delta_dpep,
282  rhs_ep,
283  active,
284  false,
285  deactivated_due_to_ld_ep);
286  unsigned row = 0;
287  for (unsigned dof = 0; dof < whole_system_size; ++dof)
288  if (dof_included(dof, deactivated_due_to_ld))
289  {
290  jac[row][col] =
291  -(rhs_ep[dof] - orig_rhs[row]) / ep; // remember jacobian = -d(rhs)/d(something)
292  row++;
293  }
294  col++; // all of the first 6 columns are dof_included since they're stresses
295  }
296 
297  std::vector<Real> pmep;
298  pmep.resize(_num_surfaces);
299  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
300  pmep[surface] = pm[surface];
301  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
302  {
303  if (!dof_included(6 + surface, deactivated_due_to_ld))
304  continue;
305  ep = _fspb_debug_pm_change[surface];
306  pmep[surface] += ep;
307  active.assign(_num_surfaces, true);
308  calculateRHS(
309  stress, intnl_old, intnl, pmep, delta_dp, rhs_ep, active, false, deactivated_due_to_ld_ep);
310  unsigned row = 0;
311  for (unsigned dof = 0; dof < whole_system_size; ++dof)
312  if (dof_included(dof, deactivated_due_to_ld))
313  {
314  jac[row][col] =
315  -(rhs_ep[dof] - orig_rhs[row]) / ep; // remember jacobian = -d(rhs)/d(something)
316  row++;
317  }
318  pmep[surface] -= ep;
319  col++;
320  }
321 
322  std::vector<Real> intnlep;
323  intnlep.resize(_num_models);
324  for (unsigned model = 0; model < _num_models; ++model)
325  intnlep[model] = intnl[model];
326  for (unsigned model = 0; model < _num_models; ++model)
327  {
328  if (!dof_included(6 + _num_surfaces + model, deactivated_due_to_ld))
329  continue;
330  ep = _fspb_debug_intnl_change[model];
331  intnlep[model] += ep;
332  active.assign(_num_surfaces, true);
333  calculateRHS(
334  stress, intnl_old, intnlep, pm, delta_dp, rhs_ep, active, false, deactivated_due_to_ld_ep);
335  unsigned row = 0;
336  for (unsigned dof = 0; dof < whole_system_size; ++dof)
337  if (dof_included(dof, deactivated_due_to_ld))
338  {
339  jac[row][col] =
340  -(rhs_ep[dof] - orig_rhs[row]) / ep; // remember jacobian = -d(rhs)/d(something)
341  row++;
342  }
343  intnlep[model] -= ep;
344  col++;
345  }
346 }

Referenced by MultiPlasticityDebugger::checkJacobian(), and MultiPlasticityDebugger::checkSolution().

◆ flowPotential()

void MultiPlasticityRawComponentAssembler::flowPotential ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  r 
)
protectedvirtualinherited

The active flow potential(s) - one for each yield function.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
activeset of active constraints - only the active flow potentials are put into "r"
[out]rthe flow potential (flow potentials in the multi-surface case)

Definition at line 180 of file MultiPlasticityRawComponentAssembler.C.

184 {
185  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
186  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
187 
188  r.resize(0);
189  std::vector<unsigned int> active_surfaces_of_model;
190  std::vector<unsigned int>::iterator active_surface;
191  std::vector<RankTwoTensor> model_r;
192  for (unsigned model = 0; model < _num_models; ++model)
193  {
194  activeModelSurfaces(model, active, active_surfaces_of_model);
195  if (active_surfaces_of_model.size() > 0)
196  {
197  _f[model]->flowPotentialV(stress, intnl[model], model_r);
198  for (active_surface = active_surfaces_of_model.begin();
199  active_surface != active_surfaces_of_model.end();
200  ++active_surface)
201  r.push_back(model_r[*active_surface]);
202  }
203  }
204 }

Referenced by MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), consistentTangentOperator(), MultiPlasticityDebugger::fddflowPotential_dintnl(), and MultiPlasticityDebugger::fddflowPotential_dstress().

◆ hardPotential()

void MultiPlasticityRawComponentAssembler::hardPotential ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  h 
)
protectedvirtualinherited

The active hardening potentials (one for each internal parameter and for each yield function) by assumption in the Userobjects, the h[a][alpha] is nonzero only if the surface alpha is part of model a, so we only calculate those here.

Parameters
stressthe stress at which to calculate the hardening potential
intnlvector of internal parameters
activeset of active constraints - only the active hardening potentials are put into "h"
[out]hthe hardening potentials. h[alpha] = hardening potential for yield fcn alpha (and, by the above assumption we know which hardening parameter, a, this belongs to)

Definition at line 262 of file MultiPlasticityRawComponentAssembler.C.

266 {
267  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
268  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
269 
270  h.resize(0);
271  std::vector<unsigned int> active_surfaces_of_model;
272  std::vector<unsigned int>::iterator active_surface;
273  std::vector<Real> model_h;
274  for (unsigned model = 0; model < _num_models; ++model)
275  {
276  activeModelSurfaces(model, active, active_surfaces_of_model);
277  if (active_surfaces_of_model.size() > 0)
278  {
279  _f[model]->hardPotentialV(stress, intnl[model], model_h);
280  for (active_surface = active_surfaces_of_model.begin();
281  active_surface != active_surfaces_of_model.end();
282  ++active_surface)
283  h.push_back(model_h[*active_surface]);
284  }
285  }
286 }

Referenced by MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), and consistentTangentOperator().

◆ incrementDumb()

void ComputeMultiPlasticityStress::incrementDumb ( int &  dumb_iteration,
const std::vector< unsigned int > &  dumb_order,
std::vector< bool > &  act 
)
protectedvirtual

Increments "dumb_iteration" by 1, and sets "act" appropriately (act[alpha] = true iff alpha_th bit of dumb_iteration == 1)

Parameters
[in,out]dumb_iterationUsed to set act bitwise - the "dumb" scheme tries all possible combinations of act until a successful return
[in]dumb_orderdumb_order dumb_order[0] will be the yield surface furthest away from (stress, intnl), dumb_order[1] will be the next yield surface, etc. The distance measure used is f/|df_dstress|. This array can then be fed into incrementDumb in order to first try the yield surfaces which are farthest away from the (stress, intnl).
[out]actactive constraints

Definition at line 1555 of file ComputeMultiPlasticityStress.C.

1558 {
1559  dumb_iteration += 1;
1560  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1561  act[dumb_order[surface]] =
1562  (dumb_iteration &
1563  (1 << surface)); // returns true if the surface_th bit of dumb_iteration == 1
1564 }

Referenced by changeScheme(), and returnMap().

◆ initQpStatefulProperties()

void ComputeMultiPlasticityStress::initQpStatefulProperties ( )
protectedvirtual

Reimplemented from ComputeStressBase.

Definition at line 191 of file ComputeMultiPlasticityStress.C.

192 {
193  ComputeStressBase::initQpStatefulProperties();
194 
195  _plastic_strain[_qp].zero();
196 
197  _intnl[_qp].assign(_num_models, 0);
198 
199  _yf[_qp].assign(_num_surfaces, 0);
200 
201  _dummy_pm.assign(_num_surfaces, 0);
202 
203  _iter[_qp] = 0.0; // this is really an unsigned int, but for visualisation i convert it to Real
204  _linesearch_needed[_qp] = 0;
205  _ld_encountered[_qp] = 0;
206  _constraints_added[_qp] = 0;
207 
208  _n[_qp] = _n_input;
209 
210  if (_cosserat)
211  (*_couple_stress)[_qp].zero();
212 
213  if (_fspb_debug == "jacobian")
214  {
215  checkDerivatives();
216  mooseError("Finite-differencing completed. Exiting with no error");
217  }
218 }

◆ lineSearch()

bool ComputeMultiPlasticityStress::lineSearch ( Real &  nr_res2,
RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
std::vector< Real > &  intnl,
std::vector< Real > &  pm,
const RankFourTensor &  E_inv,
RankTwoTensor &  delta_dp,
const RankTwoTensor &  dstress,
const std::vector< Real > &  dpm,
const std::vector< Real > &  dintnl,
std::vector< Real > &  f,
RankTwoTensor &  epp,
std::vector< Real > &  ic,
const std::vector< bool > &  active,
const std::vector< bool > &  deactivated_due_to_ld,
bool &  linesearch_needed 
)
protectedvirtual

Performs a line search.

Algorithm is taken straight from "Numerical Recipes". Given the changes dstress, dpm and dintnl provided by the nrStep routine, a line-search looks for an appropriate under-relaxation that reduces the residual-squared (nr_res2).

Most variables are input/output variables: they enter the function with their values at the start of the Newton step, and they exit the function with values attained after applying the under-relaxation

Parameters
[in,out]nr_res2The residual-squared
[out]stressThe stress after returning to the yield surface
intnl_oldThe internal variables at the previous "time" step
[in,out]intnlThe internal variables
[in,out]pmThe plasticity multiplier(s) (consistency parameter(s))
E_invinverse of the elasticity tensor
[in,out]delta_dpChange in plastic strain from start of "time" step to current configuration (plastic_strain - plastic_strain_old)
dstressChange in stress for a full Newton step
dpmChange in plasticity multiplier for a full Newton step
dintnlchange in internal parameter(s) for a full Newton step
[in,out]fYield function(s). In this routine, only the active constraints that are not deactivated_due_to_ld are contained in f.
[in,out]eppPlastic strain increment constraint
[in,out]icInternal constraint. In this routine, only the active constraints that are not deactivated_due_to_ld are contained in ic.
activeThe active constraints.
deactivated_due_to_ldTrue if a constraint has temporarily been made deactive due to linear dependence.
[out]linesearch_neededTrue if the full Newton-Raphson step was cut by the linesearch
Returns
true if successfully found a step that reduces the residual-squared

Definition at line 1391 of file ComputeMultiPlasticityStress.C.

1407 {
1408  // Line search algorithm straight out of "Numerical Recipes"
1409 
1410  bool success =
1411  true; // return value: will be false if linesearch couldn't reduce the residual-squared
1412 
1413  // Aim is to decrease residual2
1414 
1415  Real lam = 1.0; // the line-search parameter: 1.0 is a full Newton step
1416  Real lam_min =
1417  1E-10; // minimum value of lam allowed - perhaps this should be dynamically calculated?
1418  Real f0 = nr_res2; // initial value of residual2
1419  Real slope = -2 * nr_res2; // "Numerical Recipes" uses -b*A*x, in order to check for roundoff, but
1420  // i hope the nrStep would warn if there were problems.
1421  Real tmp_lam; // cached value of lam used in quadratic & cubic line search
1422  Real f2 = nr_res2; // cached value of f = residual2 used in the cubic in the line search
1423  Real lam2 = lam; // cached value of lam used in the cubic in the line search
1424 
1425  // pm during the line-search
1426  std::vector<Real> ls_pm;
1427  ls_pm.resize(pm.size());
1428 
1429  // delta_dp during the line-search
1430  RankTwoTensor ls_delta_dp;
1431 
1432  // internal parameter during the line-search
1433  std::vector<Real> ls_intnl;
1434  ls_intnl.resize(intnl.size());
1435 
1436  // stress during the line-search
1437  RankTwoTensor ls_stress;
1438 
1439  // flow directions (not used in line search, but calculateConstraints returns this parameter)
1440  std::vector<RankTwoTensor> r;
1441 
1442  while (true)
1443  {
1444  // update the variables using this line-search parameter
1445  for (unsigned alpha = 0; alpha < pm.size(); ++alpha)
1446  ls_pm[alpha] = pm[alpha] + dpm[alpha] * lam;
1447  ls_delta_dp = delta_dp - E_inv * dstress * lam;
1448  for (unsigned a = 0; a < intnl.size(); ++a)
1449  ls_intnl[a] = intnl[a] + dintnl[a] * lam;
1450  ls_stress = stress + dstress * lam;
1451 
1452  // calculate the new active yield functions, epp and active internal constraints
1453  calculateConstraints(ls_stress, intnl_old, ls_intnl, ls_pm, ls_delta_dp, f, r, epp, ic, active);
1454 
1455  // calculate the new residual-squared
1456  nr_res2 = residual2(ls_pm, f, epp, ic, active, deactivated_due_to_ld);
1457 
1458  if (nr_res2 < f0 + 1E-4 * lam * slope)
1459  break;
1460  else if (lam < lam_min)
1461  {
1462  success = false;
1463  // restore plastic multipliers, yield functions, etc to original values
1464  for (unsigned alpha = 0; alpha < pm.size(); ++alpha)
1465  ls_pm[alpha] = pm[alpha];
1466  ls_delta_dp = delta_dp;
1467  for (unsigned a = 0; a < intnl.size(); ++a)
1468  ls_intnl[a] = intnl[a];
1469  ls_stress = stress;
1470  calculateConstraints(
1471  ls_stress, intnl_old, ls_intnl, ls_pm, ls_delta_dp, f, r, epp, ic, active);
1472  nr_res2 = residual2(ls_pm, f, epp, ic, active, deactivated_due_to_ld);
1473  break;
1474  }
1475  else if (lam == 1.0)
1476  {
1477  // model as a quadratic
1478  tmp_lam = -slope / 2.0 / (nr_res2 - f0 - slope);
1479  }
1480  else
1481  {
1482  // model as a cubic
1483  Real rhs1 = nr_res2 - f0 - lam * slope;
1484  Real rhs2 = f2 - f0 - lam2 * slope;
1485  Real a = (rhs1 / Utility::pow<2>(lam) - rhs2 / Utility::pow<2>(lam2)) / (lam - lam2);
1486  Real b =
1487  (-lam2 * rhs1 / Utility::pow<2>(lam) + lam * rhs2 / Utility::pow<2>(lam2)) / (lam - lam2);
1488  if (a == 0)
1489  tmp_lam = -slope / (2.0 * b);
1490  else
1491  {
1492  Real disc = Utility::pow<2>(b) - 3 * a * slope;
1493  if (disc < 0)
1494  tmp_lam = 0.5 * lam;
1495  else if (b <= 0)
1496  tmp_lam = (-b + std::sqrt(disc)) / (3.0 * a);
1497  else
1498  tmp_lam = -slope / (b + std::sqrt(disc));
1499  }
1500  if (tmp_lam > 0.5 * lam)
1501  tmp_lam = 0.5 * lam;
1502  }
1503  lam2 = lam;
1504  f2 = nr_res2;
1505  lam = std::max(tmp_lam, 0.1 * lam);
1506  }
1507 
1508  if (lam < 1.0)
1509  linesearch_needed = true;
1510 
1511  // assign the quantities found in the line-search
1512  // back to the originals
1513  for (unsigned alpha = 0; alpha < pm.size(); ++alpha)
1514  pm[alpha] = ls_pm[alpha];
1515  delta_dp = ls_delta_dp;
1516  for (unsigned a = 0; a < intnl.size(); ++a)
1517  intnl[a] = ls_intnl[a];
1518  stress = ls_stress;
1519 
1520  return success;
1521 }

Referenced by singleStep().

◆ modelNumber()

unsigned int MultiPlasticityRawComponentAssembler::modelNumber ( unsigned int  surface)
protectedinherited

returns the model number, given the surface number

Definition at line 785 of file MultiPlasticityRawComponentAssembler.C.

786 {
787  return _model_given_surface[surface];
788 }

Referenced by applyKuhnTucker(), MultiPlasticityLinearSystem::calculateJacobian(), canAddConstraints(), checkAdmissible(), checkKuhnTucker(), consistentTangentOperator(), MultiPlasticityDebugger::fddflowPotential_dintnl(), MultiPlasticityDebugger::fddyieldFunction_dintnl(), residual2(), returnMap(), and singleStep().

◆ nrStep()

void MultiPlasticityLinearSystem::nrStep ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankFourTensor &  E_inv,
const RankTwoTensor &  delta_dp,
RankTwoTensor &  dstress,
std::vector< Real > &  dpm,
std::vector< Real > &  dintnl,
const std::vector< bool > &  active,
std::vector< bool > &  deactivated_due_to_ld 
)
protectedvirtualinherited

Performs one Newton-Raphson step.

The purpose here is to find the changes, dstress, dpm and dintnl according to the Newton-Raphson procedure

Parameters
stressCurrent value of stress
intnl_oldThe internal variables at the previous "time" step
intnlCurrent value of the internal variables
pmCurrent value of the plasticity multipliers (consistency parameters)
E_invinverse of the elasticity tensor
delta_dpCurrent value of the plastic-strain increment (ie plastic_strain - plastic_strain_old)
[out]dstressThe change in stress for a full Newton step
[out]dpmThe change in all plasticity multipliers for a full Newton step
[out]dintnlThe change in all internal variables for a full Newton step
activeThe active constraints
[out]deactivated_due_to_ldThe constraints deactivated due to linear-dependence of the flow directions

Definition at line 615 of file MultiPlasticityLinearSystem.C.

626 {
627  // Calculate RHS and Jacobian
628  std::vector<Real> rhs;
629  calculateRHS(stress, intnl_old, intnl, pm, delta_dp, rhs, active, true, deactivated_due_to_ld);
630 
631  std::vector<std::vector<Real>> jac;
632  calculateJacobian(stress, intnl, pm, E_inv, active, deactivated_due_to_ld, jac);
633 
634  // prepare for LAPACKgesv_ routine provided by PETSc
635  int system_size = rhs.size();
636 
637  std::vector<double> a(system_size * system_size);
638  // Fill in the a "matrix" by going down columns
639  unsigned ind = 0;
640  for (int col = 0; col < system_size; ++col)
641  for (int row = 0; row < system_size; ++row)
642  a[ind++] = jac[row][col];
643 
644  int nrhs = 1;
645  std::vector<int> ipiv(system_size);
646  int info;
647  LAPACKgesv_(&system_size, &nrhs, &a[0], &system_size, &ipiv[0], &rhs[0], &system_size, &info);
648 
649  if (info != 0)
650  mooseError("In solving the linear system in a Newton-Raphson process, the PETSC LAPACK gsev "
651  "routine returned with error code ",
652  info);
653 
654  // Extract the results back to dstress, dpm and dintnl
655  std::vector<bool> active_not_deact(_num_surfaces);
656  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
657  active_not_deact[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
658 
659  unsigned int dim = 3;
660  ind = 0;
661 
662  for (unsigned i = 0; i < dim; ++i)
663  for (unsigned j = 0; j <= i; ++j)
664  dstress(i, j) = dstress(j, i) = rhs[ind++];
665  dpm.assign(_num_surfaces, 0);
666  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
667  if (active_not_deact[surface])
668  dpm[surface] = rhs[ind++];
669  dintnl.assign(_num_models, 0);
670  for (unsigned model = 0; model < _num_models; ++model)
671  if (anyActiveSurfaces(model, active_not_deact))
672  dintnl[model] = rhs[ind++];
673 
674  mooseAssert(static_cast<int>(ind) == system_size,
675  "Incorrect extracting of changes from NR solution in nrStep");
676 }

Referenced by MultiPlasticityDebugger::checkSolution(), and singleStep().

◆ numberActive()

unsigned ComputeMultiPlasticityStress::numberActive ( const std::vector< bool > &  active)
protectedvirtual

counts the number of active constraints

Definition at line 1261 of file ComputeMultiPlasticityStress.C.

1262 {
1263  unsigned num_active = 0;
1264  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1265  if (active[surface])
1266  num_active++;
1267  return num_active;
1268 }

Referenced by returnMap(), and singleStep().

◆ outputAndCheckDebugParameters()

void MultiPlasticityDebugger::outputAndCheckDebugParameters ( )
inherited

Outputs the debug parameters: _fspb_debug_stress, _fspd_debug_pm, etc and checks that they are sized correctly.

Definition at line 55 of file MultiPlasticityDebugger.C.

56 {
57  Moose::err << "Debug Parameters are as follows\n";
58  Moose::err << "stress = \n";
59  _fspb_debug_stress.print();
60 
61  if (_fspb_debug_pm.size() != _num_surfaces || _fspb_debug_intnl.size() != _num_models ||
62  _fspb_debug_pm_change.size() != _num_surfaces ||
63  _fspb_debug_intnl_change.size() != _num_models)
64  mooseError("The debug parameters have the wrong size\n");
65 
66  Moose::err << "plastic multipliers =\n";
67  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
68  Moose::err << _fspb_debug_pm[surface] << "\n";
69 
70  Moose::err << "internal parameters =\n";
71  for (unsigned model = 0; model < _num_models; ++model)
72  Moose::err << _fspb_debug_intnl[model] << "\n";
73 
74  Moose::err << "finite-differencing parameter for stress-changes:\n"
75  << _fspb_debug_stress_change << "\n";
76  Moose::err << "finite-differencing parameter(s) for plastic-multiplier(s):\n";
77  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
78  Moose::err << _fspb_debug_pm_change[surface] << "\n";
79  Moose::err << "finite-differencing parameter(s) for internal-parameter(s):\n";
80  for (unsigned model = 0; model < _num_models; ++model)
81  Moose::err << _fspb_debug_intnl_change[model] << "\n";
82 }

Referenced by MultiPlasticityDebugger::checkDerivatives(), MultiPlasticityDebugger::checkJacobian(), and MultiPlasticityDebugger::checkSolution().

◆ plasticStep()

bool ComputeMultiPlasticityStress::plasticStep ( const RankTwoTensor &  stress_old,
RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
std::vector< Real > &  intnl,
const RankTwoTensor &  plastic_strain_old,
RankTwoTensor &  plastic_strain,
const RankFourTensor &  E_ijkl,
const RankTwoTensor &  strain_increment,
std::vector< Real > &  yf,
unsigned int &  iterations,
bool &  linesearch_needed,
bool &  ld_encountered,
bool &  constraints_added,
RankFourTensor &  consistent_tangent_operator 
)
protectedvirtual

performs a plastic step

Parameters
stress_oldThe value of stress at the previous "time" step
[out]stressstress after returning to the yield surface
intnl_oldThe internal variables at the previous "time" step
[out]intnlinternal variables after returning to the yield surface
plastic_strain_oldThe value of plastic strain at the previous "time" step
[out]plastic_strainplastic_strain after returning to the yield surface
E_ijklThe elasticity tensor.
strain_incrementThe applied strain increment
[out]yfAll the yield functions at (stress, intnl)
[out]iterationsThe total number of Newton-Raphson iterations used
[out]linesearch_neededTrue if a linesearch was needed at any stage during the Newton-Raphson proceedure
[out]ld_encounteredTrue if a linear-dependence of the flow directions was encountered at any stage during the Newton-Raphson proceedure
[out]constraints_addedTrue if constraints were added into the active set at any stage during the Newton-Raphson proceedure
[out]consistent_tangent_operatorThe consistent tangent operator d(stress_rate)/d(strain_rate)
Returns
true if the (stress, intnl) are admissible. Otherwise, if _ignore_failures==true, the output variables will be the best admissible ones found during the return-map. Otherwise, if _ignore_failures==false, this routine will perform some finite-diference checks and call mooseError

the idea in the following is to potentially subdivide the strain increment into smaller fractions, of size "step_size". First step_size = 1 is tried, and if that fails then 0.5 is tried, then 0.25, etc. As soon as a step is successful, its results are put into the "good" variables, which are used if a subsequent step fails. If >= 2 consecutive steps are successful, the step_size is increased by 1.2

The point of all this is that I hope that the time-step for the entire mesh need not be cut if there are only a few "bad" quadpoints where the return-map is difficult

Definition at line 468 of file ComputeMultiPlasticityStress.C.

482 {
498  bool return_successful = false;
499 
500  // total number of Newton-Raphson iterations used
501  unsigned int iter = 0;
502 
503  Real step_size = 1.0;
504  Real time_simulated = 0.0;
505 
506  // the "good" variables hold the latest admissible stress
507  // and internal parameters.
508  RankTwoTensor stress_good = stress_old;
509  RankTwoTensor plastic_strain_good = plastic_strain_old;
510  std::vector<Real> intnl_good(_num_models);
511  for (unsigned model = 0; model < _num_models; ++model)
512  intnl_good[model] = intnl_old[model];
513  std::vector<Real> yf_good(_num_surfaces);
514 
515  // Following is necessary because I want strain_increment to be "const"
516  // but I also want to be able to subdivide an initial_stress
517  RankTwoTensor this_strain_increment = strain_increment;
518 
519  RankTwoTensor dep = step_size * this_strain_increment;
520 
521  _cumulative_pm.assign(_num_surfaces, 0);
522 
523  unsigned int num_consecutive_successes = 0;
524  while (time_simulated < 1.0 && step_size >= _min_stepsize)
525  {
526  iter = 0;
527  return_successful = returnMap(stress_good,
528  stress,
529  intnl_good,
530  intnl,
531  plastic_strain_good,
532  plastic_strain,
533  E_ijkl,
534  dep,
535  yf,
536  iter,
537  step_size <= _max_stepsize_for_dumb,
538  linesearch_needed,
539  ld_encountered,
540  constraints_added,
541  time_simulated + step_size >= 1,
542  consistent_tangent_operator,
543  _cumulative_pm);
544  iterations += iter;
545 
546  if (return_successful)
547  {
548  num_consecutive_successes += 1;
549  time_simulated += step_size;
550 
551  if (time_simulated < 1.0) // this condition is just for optimization: if time_simulated=1 then
552  // the "good" quantities are no longer needed
553  {
554  stress_good = stress;
555  plastic_strain_good = plastic_strain;
556  for (unsigned model = 0; model < _num_models; ++model)
557  intnl_good[model] = intnl[model];
558  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
559  yf_good[surface] = yf[surface];
560  if (num_consecutive_successes >= 2)
561  step_size *= 1.2;
562  }
563  step_size = std::min(step_size, 1.0 - time_simulated); // avoid overshoots
564  }
565  else
566  {
567  step_size *= 0.5;
568  num_consecutive_successes = 0;
569  stress = stress_good;
570  plastic_strain = plastic_strain_good;
571  for (unsigned model = 0; model < _num_models; ++model)
572  intnl[model] = intnl_good[model];
573  yf.resize(_num_surfaces); // might have excited with junk
574  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
575  yf[surface] = yf_good[surface];
576  dep = step_size * this_strain_increment;
577  }
578  }
579 
580  if (!return_successful)
581  {
582  if (_ignore_failures)
583  {
584  stress = stress_good;
585  plastic_strain = plastic_strain_good;
586  for (unsigned model = 0; model < _num_models; ++model)
587  intnl[model] = intnl_good[model];
588  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
589  yf[surface] = yf_good[surface];
590  }
591  else
592  {
593  Moose::out << "After reducing the stepsize to " << step_size
594  << " with original strain increment with L2norm " << this_strain_increment.L2norm()
595  << " the returnMap algorithm failed\n";
596 
597  _fspb_debug_stress = stress_good + E_ijkl * dep;
598  _fspb_debug_pm.assign(
599  _num_surfaces,
600  1); // this is chosen arbitrarily - please change if a more suitable value occurs to you!
601  _fspb_debug_intnl.resize(_num_models);
602  for (unsigned model = 0; model < _num_models; ++model)
603  _fspb_debug_intnl[model] = intnl_good[model];
604  checkDerivatives();
605  checkJacobian(_my_elasticity_tensor.invSymm(), _intnl_old[_qp]);
606  checkSolution(_my_elasticity_tensor.invSymm());
607  mooseError("Exiting\n");
608  }
609  }
610 
611  return return_successful;
612 }

Referenced by computeQpStress().

◆ postReturnMap()

void ComputeMultiPlasticityStress::postReturnMap ( )
protectedvirtual

Definition at line 339 of file ComputeMultiPlasticityStress.C.

340 {
341  if (_n_supplied)
342  {
343  // this is a rotation matrix that will rotate "z" axis back to _n
344  _rot = _rot.transpose();
345 
346  // rotate the tensors back to original frame where _n is correctly oriented
347  _my_elasticity_tensor.rotate(_rot);
348  _Jacobian_mult[_qp].rotate(_rot);
349  _my_strain_increment.rotate(_rot);
350  _stress[_qp].rotate(_rot);
351  _plastic_strain[_qp].rotate(_rot);
352  if (_cosserat)
353  {
354  _my_flexural_rigidity_tensor.rotate(_rot);
355  (*_Jacobian_mult_couple)[_qp].rotate(_rot);
356  _my_curvature.rotate(_rot);
357  (*_couple_stress)[_qp].rotate(_rot);
358  }
359 
360  // Rotate n by _rotation_increment
361  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
362  {
363  _n[_qp](i) = 0;
364  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
365  _n[_qp](i) += _rotation_increment[_qp](i, j) * _n_old[_qp](j);
366  }
367  }
368 }

Referenced by computeQpStress().

◆ preReturnMap()

void ComputeMultiPlasticityStress::preReturnMap ( )
protectedvirtual

Definition at line 320 of file ComputeMultiPlasticityStress.C.

321 {
322  if (_n_supplied)
323  {
324  // this is a rotation matrix that will rotate _n to the "z" axis
325  _rot = RotationMatrix::rotVecToZ(_n[_qp]);
326 
327  // rotate the tensors to this frame
328  _my_elasticity_tensor.rotate(_rot);
329  _my_strain_increment.rotate(_rot);
330  if (_cosserat)
331  {
332  _my_flexural_rigidity_tensor.rotate(_rot);
333  _my_curvature.rotate(_rot);
334  }
335  }
336 }

Referenced by computeQpStress().

◆ quickStep()

bool ComputeMultiPlasticityStress::quickStep ( const RankTwoTensor &  stress_old,
RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
std::vector< Real > &  intnl,
std::vector< Real > &  pm,
std::vector< Real > &  cumulative_pm,
const RankTwoTensor &  plastic_strain_old,
RankTwoTensor &  plastic_strain,
const RankFourTensor &  E_ijkl,
const RankTwoTensor &  strain_increment,
std::vector< Real > &  yf,
unsigned int &  iterations,
RankFourTensor &  consistent_tangent_operator,
const quickStep_called_from_t  called_from,
bool  final_step 
)
protectedvirtual

Attempts to find an admissible (stress, intnl) by using the customized return-map algorithms defined through the TensorMechanicsPlasticXXXX.returnMap functions.

Parameters
stress_oldThe value of stress at the previous "time" step
[out]stressIf returnvalue=true then this is the returned value of stress. Otherwise, this is undefined
intnl_oldThe internal variables at the previous "time" step
[out]intnlIf returnvalue=true then this is the value of the internal parameters after returning. Otherwise, this is undefined
[out]pmIf returnvalue=true, this is the plastic multipliers needed to bring aout the return. Otherwise, this is undefined
[in/out]cumulative_pm If returnvalue=true, this is cumulative plastic multipliers, updated with pm. Otherwise, this is untouched by the algorithm
plastic_strain_oldThe value of plastic strain at the previous "time" step
[out]plastic_strainIf returnvalue=true, this is the new plastic strain. Otherwise it is set to plastic_strain_old
E_ijklThe elasticity tensor.
strain_incrementThe applied strain increment
[out]yfIf returnvalue=true, then all the yield functions at (stress, intnl). Otherwise, all the yield functions at (stress_old, intnl_old)
[out]iterationsNumber of NR iterations used, which is always zero in the current implementation.
called_fromThis can be called from computeQpStress, in which case it can actually provde an answer to the returnmap algorithm, or from returnMap in which case it is probably only providing an answer to a particular subdivision of the returnmap algorithm. The consistent tangent operator is calculated idfferently in each case
final_stepThe consistent tangent operator is calculated if this is true
[out]consistent_tangent_operatorIf final_step==true and returnvalue=true, then this is the consistent tangent operator d(stress_rate)/d(strain_rate). Otherwise it is undefined.
Returns
true if the (stress, intnl) are admissible, in which case the (stress_old, intnl_old) could have been admissible, or exactly one of the plastic models successfully used its custom returnMap function to provide the returned (stress, intnl) values and all other plastic models are admissible at that configuration. Or, false, then (stress_old, intnl_old) is not admissible according to >=1 plastic model and the custom returnMap functions failed in some way.

Definition at line 371 of file ComputeMultiPlasticityStress.C.

386 {
387  iterations = 0;
388 
389  unsigned num_plastic_returns;
390  RankTwoTensor delta_dp;
391 
392  // the following does the customized returnMap algorithm
393  // for all the plastic models.
394  unsigned custom_model = 0;
395  bool successful_return = returnMapAll(stress_old + E_ijkl * strain_increment,
396  intnl_old,
397  E_ijkl,
398  _epp_tol,
399  stress,
400  intnl,
401  pm,
402  cumulative_pm,
403  delta_dp,
404  yf,
405  num_plastic_returns,
406  custom_model);
407 
408  // the following updates the plastic_strain, when necessary
409  // and calculates the consistent_tangent_operator, when necessary
410  if (num_plastic_returns == 0)
411  {
412  // if successful_return = true, then a purely elastic step
413  // if successful_return = false, then >=1 plastic model is in
414  // inadmissible zone and failed to return using its customized
415  // returnMap function.
416  // In either case:
417  plastic_strain = plastic_strain_old;
418  if (successful_return && final_step)
419  {
420  if (called_from == computeQpStress_function)
421  consistent_tangent_operator = E_ijkl;
422  else // cannot necessarily use E_ijkl since different plastic models may have been active
423  // during other substeps
424  consistent_tangent_operator =
425  consistentTangentOperator(stress, intnl, E_ijkl, pm, cumulative_pm);
426  }
427  return successful_return;
428  }
429  else if (num_plastic_returns == 1 && successful_return)
430  {
431  // one model has successfully completed its custom returnMap algorithm
432  // and the other models have signalled they are elastic at
433  // the trial stress
434  plastic_strain = plastic_strain_old + delta_dp;
435  if (final_step)
436  {
437  if (called_from == computeQpStress_function && _f[custom_model]->useCustomCTO())
438  {
439  if (_tangent_operator_type == elastic)
440  consistent_tangent_operator = E_ijkl;
441  else
442  {
443  std::vector<Real> custom_model_pm;
444  for (unsigned surface = 0; surface < _f[custom_model]->numberSurfaces(); ++surface)
445  custom_model_pm.push_back(cumulative_pm[_surfaces_given_model[custom_model][surface]]);
446  consistent_tangent_operator =
447  _f[custom_model]->consistentTangentOperator(stress_old + E_ijkl * strain_increment,
448  intnl_old[custom_model],
449  stress,
450  intnl[custom_model],
451  E_ijkl,
452  custom_model_pm);
453  }
454  }
455  else // cannot necessarily use the custom consistentTangentOperator since different plastic
456  // models may have been active during other substeps or the custom model says not to use
457  // its custom CTO algorithm
458  consistent_tangent_operator =
459  consistentTangentOperator(stress, intnl, E_ijkl, pm, cumulative_pm);
460  }
461  return true;
462  }
463  else // presumably returnMapAll is incorrectly coded!
464  mooseError("ComputeMultiPlasticityStress::quickStep should not get here!");
465 }

Referenced by computeQpStress(), and returnMap().

◆ reinstateLinearDependentConstraints()

bool ComputeMultiPlasticityStress::reinstateLinearDependentConstraints ( std::vector< bool > &  deactivated_due_to_ld)
protectedvirtual

makes all deactivated_due_to_ld false, and if >0 of them were initially true, returns true

Definition at line 1248 of file ComputeMultiPlasticityStress.C.

1250 {
1251  bool reinstated_actives = false;
1252  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1253  if (deactivated_due_to_ld[surface])
1254  reinstated_actives = true;
1255 
1256  deactivated_due_to_ld.assign(_num_surfaces, false);
1257  return reinstated_actives;
1258 }

Referenced by singleStep().

◆ residual2()

Real ComputeMultiPlasticityStress::residual2 ( const std::vector< Real > &  pm,
const std::vector< Real > &  f,
const RankTwoTensor &  epp,
const std::vector< Real > &  ic,
const std::vector< bool > &  active,
const std::vector< bool > &  deactivated_due_to_ld 
)
protectedvirtual

The residual-squared.

Parameters
pmthe plastic multipliers for all constraints
fthe active yield function(s) (not including the ones that are deactivated_due_to_ld)
eppthe plastic strain increment constraint
icthe active internal constraint(s) (not including the ones that are deactivated_due_to_ld)
activetrue if constraint is active
deactivated_due_to_ldtrue if constraint has been temporarily deactivated due to linear dependence of flow directions

Definition at line 1354 of file ComputeMultiPlasticityStress.C.

1360 {
1361  Real nr_res2 = 0;
1362  unsigned ind = 0;
1363 
1364  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1365  if (active[surface])
1366  {
1367  if (!deactivated_due_to_ld[surface])
1368  {
1369  if (!(pm[surface] == 0 && f[ind] <= 0))
1370  nr_res2 += 0.5 * Utility::pow<2>(f[ind] / _f[modelNumber(surface)]->_f_tol);
1371  }
1372  else if (deactivated_due_to_ld[surface] && f[ind] > 0)
1373  nr_res2 += 0.5 * Utility::pow<2>(f[ind] / _f[modelNumber(surface)]->_f_tol);
1374  ind++;
1375  }
1376 
1377  nr_res2 += 0.5 * Utility::pow<2>(epp.L2norm() / _epp_tol);
1378 
1379  std::vector<bool> active_not_deact(_num_surfaces);
1380  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1381  active_not_deact[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
1382  ind = 0;
1383  for (unsigned model = 0; model < _num_models; ++model)
1384  if (anyActiveSurfaces(model, active_not_deact))
1385  nr_res2 += 0.5 * Utility::pow<2>(ic[ind++] / _f[model]->_ic_tol);
1386 
1387  return nr_res2;
1388 }

Referenced by lineSearch(), and singleStep().

◆ returnMap()

bool ComputeMultiPlasticityStress::returnMap ( const RankTwoTensor &  stress_old,
RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
std::vector< Real > &  intnl,
const RankTwoTensor &  plastic_strain_old,
RankTwoTensor &  plastic_strain,
const RankFourTensor &  E_ijkl,
const RankTwoTensor &  strain_increment,
std::vector< Real > &  f,
unsigned int &  iter,
bool  can_revert_to_dumb,
bool &  linesearch_needed,
bool &  ld_encountered,
bool &  constraints_added,
bool  final_step,
RankFourTensor &  consistent_tangent_operator,
std::vector< Real > &  cumulative_pm 
)
protectedvirtual

Implements the return map.

Note that this algorithm doesn't do any rotations. In order to find the final stress and plastic_strain must be rotated using _rotation_increment. This is usually done in computeQpStress

Parameters
stress_oldThe value of stress at the previous "time" step
[out]stressThe stress after returning to the yield surface
intnl_oldThe internal variables at the previous "time" step
[out]intnlAll the internal variables after returning to the yield surface
plastic_strain_oldThe value of plastic strain at the previous "time" step
[out]plastic_strainThe value of plastic strain after returning to the yield surface
E_ijklThe elasticity tensor. If no plasiticity then stress = stress_old + E_ijkl*strain_increment
strain_incrementThe applied strain increment
[out]fAll the yield functions after returning to the yield surface
[out]iterThe number of Newton-Raphson iterations used
can_revert_to_dumbIf the _deactivation_scheme is set to revert to dumb, it will only be allowed to do so if this parameter is true
[out]linesearch_neededTrue if a linesearch was needed at any stage during the Newton-Raphson proceedure
[out]ld_encounteredTrue if a linear-dependence of the flow directions was encountered at any stage during the Newton-Raphson proceedure
[out]constraints_addedTrue if constraints were added into the active set at any stage during the Newton-Raphson proceedure
final_stepEach strain increment may be decomposed into a sum of smaller increments if the return-map algorithm fails. This flag indicates whether this is the last application of incremental strain
[out]consistent_tangent_operatorThe consistent tangent operator d(stress_rate)/d(strain_rate). This is only output if final_step=true, and the return value of returnMap is also true.
[in,out]cumulative_pmUpon input: the plastic multipliers before the return map. Upon output: the plastic multipliers after this return map, if the return map was successful
Returns
true if the stress was successfully returned to the yield surface

Definition at line 615 of file ComputeMultiPlasticityStress.C.

632 {
633 
634  // The "consistency parameters" (plastic multipliers)
635  // Change in plastic strain in this timestep = pm*flowPotential
636  // Each pm must be non-negative
637  std::vector<Real> pm;
638  pm.assign(_num_surfaces, 0.0);
639 
640  bool successful_return = quickStep(stress_old,
641  stress,
642  intnl_old,
643  intnl,
644  pm,
645  cumulative_pm,
646  plastic_strain_old,
647  plastic_strain,
648  E_ijkl,
649  strain_increment,
650  f,
651  iter,
652  consistent_tangent_operator,
653  returnMap_function,
654  final_step);
655 
656  if (successful_return)
657  return successful_return;
658 
659  // Here we know that the trial stress and intnl_old
660  // is inadmissible, and we have to return from those values
661  // value to the yield surface. There are three
662  // types of constraints we have to satisfy, listed
663  // below, and calculated in calculateConstraints(...)
664  // f<=0, epp=0, ic=0 (up to tolerances), and these are
665  // guaranteed to hold if nr_res2<=0.5
666  //
667  // Kuhn-Tucker conditions must also be satisfied
668  // These are:
669  // pm>=0, which may not hold upon exit of the NR loops
670  // due to _deactivation_scheme!=optimized;
671  // pm*f=0 (up to tolerances), which may not hold upon exit
672  // of the NR loops if a constraint got deactivated
673  // due to linear dependence, and then f<0, and its pm>0
674 
675  // Plastic strain constraint, L2 norm must be zero (up to a tolerance)
676  RankTwoTensor epp;
677 
678  // Yield function constraint passed to this function as
679  // std::vector<Real> & f
680  // Each yield function must be <= 0 (up to tolerance)
681  // Note that only the constraints that are active will be
682  // contained in f until the final few lines of returnMap
683 
684  // Internal constraint(s), must be zero (up to a tolerance)
685  // Note that only the constraints that are active will be
686  // contained in ic.
687  std::vector<Real> ic;
688 
689  // Record the stress before Newton-Raphson in case of failure-and-restart
690  RankTwoTensor initial_stress = stress;
691 
692  iter = 0;
693 
694  // Initialize the set of active constraints
695  // At this stage, the active constraints are
696  // those that exceed their _f_tol
697  // active constraints.
698  std::vector<bool> act;
699  buildActiveConstraints(f, stress, intnl, E_ijkl, act);
700 
701  // Inverse of E_ijkl (assuming symmetric)
702  RankFourTensor E_inv = E_ijkl.invSymm();
703 
704  // convenience variable that holds the change in plastic strain incurred during the return
705  // delta_dp = plastic_strain - plastic_strain_old
706  // delta_dp = E^{-1}*(initial_stress - stress), where initial_stress = E*(strain -
707  // plastic_strain_old)
708  RankTwoTensor delta_dp = RankTwoTensor();
709 
710  // whether single step was successful (whether line search was successful, and whether turning off
711  // constraints was successful)
712  bool single_step_success = true;
713 
714  // deactivation scheme
715  DeactivationSchemeEnum deact_scheme = _deactivation_scheme;
716 
717  // For complicated deactivation schemes we have to record the initial active set
718  std::vector<bool> initial_act;
719  initial_act.resize(_num_surfaces);
720  if (_deactivation_scheme == optimized_to_safe ||
721  _deactivation_scheme == optimized_to_safe_to_dumb ||
722  _deactivation_scheme == optimized_to_dumb)
723  {
724  // if "optimized" fails we can change the deactivation scheme to "safe", etc
725  deact_scheme = optimized;
726  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
727  initial_act[surface] = act[surface];
728  }
729 
730  if (_deactivation_scheme == safe_to_dumb)
731  deact_scheme = safe;
732 
733  // For "dumb" deactivation, the active set takes all combinations until a solution is found
734  int dumb_iteration = 0;
735  std::vector<unsigned int> dumb_order;
736 
737  if (_deactivation_scheme == dumb ||
738  (_deactivation_scheme == optimized_to_safe_to_dumb && can_revert_to_dumb) ||
739  (_deactivation_scheme == safe_to_dumb && can_revert_to_dumb) ||
740  (_deactivation_scheme == optimized_to_dumb && can_revert_to_dumb))
741  buildDumbOrder(stress, intnl, dumb_order);
742 
743  if (_deactivation_scheme == dumb)
744  {
745  incrementDumb(dumb_iteration, dumb_order, act);
746  yieldFunction(stress, intnl, act, f);
747  }
748 
749  // To avoid any re-trials of "act" combinations that
750  // we've already tried and rejected, i record the
751  // combinations in actives_tried
752  std::set<unsigned int> actives_tried;
753  actives_tried.insert(activeCombinationNumber(act));
754 
755  // The residual-squared that the line-search will reduce
756  // Later it will get contributions from epp and ic, but
757  // at present these are zero
758  Real nr_res2 = 0;
759  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
760  if (act[surface])
761  nr_res2 += 0.5 * Utility::pow<2>(f[surface] / _f[modelNumber(surface)]->_f_tol);
762 
763  successful_return = false;
764 
765  bool still_finding_solution = true;
766  while (still_finding_solution)
767  {
768  single_step_success = true;
769  unsigned int local_iter = 0;
770 
771  // The Newton-Raphson loops
772  while (nr_res2 > 0.5 && local_iter++ < _max_iter && single_step_success)
773  single_step_success = singleStep(nr_res2,
774  stress,
775  intnl_old,
776  intnl,
777  pm,
778  delta_dp,
779  E_inv,
780  f,
781  epp,
782  ic,
783  act,
784  deact_scheme,
785  linesearch_needed,
786  ld_encountered);
787 
788  bool nr_good = (nr_res2 <= 0.5 && local_iter <= _max_iter && single_step_success);
789 
790  iter += local_iter;
791 
792  // 'act' might have changed due to using deact_scheme = optimized, so
793  actives_tried.insert(activeCombinationNumber(act));
794 
795  if (!nr_good)
796  {
797  // failure of NR routine.
798  // We might be able to change the deactivation_scheme and
799  // then re-try the NR starting from the initial values
800  // Or, if deact_scheme == "dumb", we just increarse the
801  // dumb_iteration number and re-try
802  bool change_scheme = false;
803  bool increment_dumb = false;
804  change_scheme = canChangeScheme(deact_scheme, can_revert_to_dumb);
805  if (!change_scheme && deact_scheme == dumb)
806  increment_dumb = canIncrementDumb(dumb_iteration);
807 
808  still_finding_solution = (change_scheme || increment_dumb);
809 
810  if (change_scheme)
811  changeScheme(initial_act,
812  can_revert_to_dumb,
813  initial_stress,
814  intnl_old,
815  deact_scheme,
816  act,
817  dumb_iteration,
818  dumb_order);
819 
820  if (increment_dumb)
821  incrementDumb(dumb_iteration, dumb_order, act);
822 
823  if (!still_finding_solution)
824  {
825  // we cannot change the scheme, or have run out of "dumb" options
826  successful_return = false;
827  break;
828  }
829  }
830 
831  bool kt_good = false;
832  if (nr_good)
833  {
834  // check Kuhn-Tucker
835  kt_good = checkKuhnTucker(f, pm, act);
836  if (!kt_good)
837  {
838  if (deact_scheme != dumb)
839  {
840  applyKuhnTucker(f, pm, act);
841 
842  // true if we haven't tried this active set before
843  still_finding_solution =
844  (actives_tried.find(activeCombinationNumber(act)) == actives_tried.end());
845  if (!still_finding_solution)
846  {
847  // must have tried turning off the constraints already.
848  // so try changing the scheme
849  if (canChangeScheme(deact_scheme, can_revert_to_dumb))
850  {
851  still_finding_solution = true;
852  changeScheme(initial_act,
853  can_revert_to_dumb,
854  initial_stress,
855  intnl_old,
856  deact_scheme,
857  act,
858  dumb_iteration,
859  dumb_order);
860  }
861  }
862  }
863  else
864  {
865  bool increment_dumb = false;
866  increment_dumb = canIncrementDumb(dumb_iteration);
867  still_finding_solution = increment_dumb;
868  if (increment_dumb)
869  incrementDumb(dumb_iteration, dumb_order, act);
870  }
871 
872  if (!still_finding_solution)
873  {
874  // have tried turning off the constraints already,
875  // or have run out of "dumb" options
876  successful_return = false;
877  break;
878  }
879  }
880  }
881 
882  bool admissible = false;
883  if (nr_good && kt_good)
884  {
885  // check admissible
886  std::vector<Real> all_f;
887  if (_num_surfaces == 1)
888  admissible = true; // for a single surface if NR has exited successfully then (stress,
889  // intnl) must be admissible
890  else
891  admissible = checkAdmissible(stress, intnl, all_f);
892 
893  if (!admissible)
894  {
895  // Not admissible.
896  // We can try adding constraints back in
897  // We can try changing the deactivation scheme
898  // Or, if deact_scheme == dumb, just increase dumb_iteration
899  bool add_constraints = canAddConstraints(act, all_f);
900  if (add_constraints)
901  {
902  constraints_added = true;
903  std::vector<bool> act_plus(_num_surfaces,
904  false); // "act" with the positive constraints added in
905  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
906  if (act[surface] ||
907  (!act[surface] && (all_f[surface] > _f[modelNumber(surface)]->_f_tol)))
908  act_plus[surface] = true;
909  if (actives_tried.find(activeCombinationNumber(act_plus)) == actives_tried.end())
910  {
911  // haven't tried this combination of actives yet
912  constraints_added = true;
913  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
914  act[surface] = act_plus[surface];
915  }
916  else
917  add_constraints = false; // haven't managed to add a new combination
918  }
919 
920  bool change_scheme = false;
921  bool increment_dumb = false;
922 
923  if (!add_constraints)
924  change_scheme = canChangeScheme(deact_scheme, can_revert_to_dumb);
925 
926  if (!add_constraints && !change_scheme && deact_scheme == dumb)
927  increment_dumb = canIncrementDumb(dumb_iteration);
928 
929  still_finding_solution = (add_constraints || change_scheme || increment_dumb);
930 
931  if (change_scheme)
932  changeScheme(initial_act,
933  can_revert_to_dumb,
934  initial_stress,
935  intnl_old,
936  deact_scheme,
937  act,
938  dumb_iteration,
939  dumb_order);
940 
941  if (increment_dumb)
942  incrementDumb(dumb_iteration, dumb_order, act);
943 
944  if (!still_finding_solution)
945  {
946  // we cannot change the scheme, or have run out of "dumb" options
947  successful_return = false;
948  break;
949  }
950  }
951  }
952 
953  successful_return = (nr_good && admissible && kt_good);
954  if (successful_return)
955  break;
956 
957  if (still_finding_solution)
958  {
959  stress = initial_stress;
960  delta_dp = RankTwoTensor(); // back to zero change in plastic strain
961  for (unsigned model = 0; model < _num_models; ++model)
962  intnl[model] = intnl_old[model]; // back to old internal params
963  pm.assign(_num_surfaces, 0.0); // back to zero plastic multipliers
964 
965  unsigned num_active = numberActive(act);
966  if (num_active == 0)
967  {
968  successful_return = false;
969  break; // failure
970  }
971 
972  actives_tried.insert(activeCombinationNumber(act));
973 
974  // Since "act" set has changed, either by changing deact_scheme, or by KT failing, so need to
975  // re-calculate nr_res2
976  yieldFunction(stress, intnl, act, f);
977 
978  nr_res2 = 0;
979  unsigned ind = 0;
980  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
981  if (act[surface])
982  {
983  if (f[ind] > _f[modelNumber(surface)]->_f_tol)
984  nr_res2 += 0.5 * Utility::pow<2>(f[ind] / _f[modelNumber(surface)]->_f_tol);
985  ind++;
986  }
987  }
988  }
989 
990  // returned, with either success or failure
991  if (successful_return)
992  {
993  plastic_strain += delta_dp;
994 
995  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
996  cumulative_pm[surface] += pm[surface];
997 
998  if (final_step)
999  consistent_tangent_operator =
1000  consistentTangentOperator(stress, intnl, E_ijkl, pm, cumulative_pm);
1001 
1002  if (f.size() != _num_surfaces)
1003  {
1004  // at this stage f.size() = num_active, but we need to return with all the yield functions
1005  // evaluated, so:
1006  act.assign(_num_surfaces, true);
1007  yieldFunction(stress, intnl, act, f);
1008  }
1009  }
1010 
1011  return successful_return;
1012 }

Referenced by plasticStep().

◆ returnMapAll()

bool MultiPlasticityRawComponentAssembler::returnMapAll ( const RankTwoTensor &  trial_stress,
const std::vector< Real > &  intnl_old,
const RankFourTensor &  E_ijkl,
Real  ep_plastic_tolerance,
RankTwoTensor &  stress,
std::vector< Real > &  intnl,
std::vector< Real > &  pm,
std::vector< Real > &  cumulative_pm,
RankTwoTensor &  delta_dp,
std::vector< Real > &  yf,
unsigned &  num_successful_plastic_returns,
unsigned &  custom_model 
)
protectedinherited

Performs a returnMap for each plastic model using their inbuilt returnMap functions.

Performs a returnMap for each plastic model.

This may be used to quickly ascertain whether a (trial_stress, intnl_old) configuration is admissible, or whether a single model's customized returnMap function can provide a solution to the return-map problem, or whether a full Newton-Raphson approach such as implemented in ComputeMultiPlasticityStress is needed.

There are three cases mentioned below: (A) The (trial_stress, intnl_old) configuration is admissible according to all plastic models (B) The (trial_stress, intnl_old) configuration is inadmissible to exactly one plastic model, and that model can successfully use its customized returnMap function to provide a returned (stress, intnl) configuration, and that configuration is admissible according to all plastic models (C) All other cases. This includes customized returnMap functions failing, or more than one plastic_model being inadmissible, etc

Parameters
trial_stressthe trial stress
intnl_oldthe old values of the internal parameters
E_ijklthe elasticity tensor
ep_plastic_tolerancethe tolerance on the plastic strain
[out]stressis set to trial_stress in case (A) or (C), and the returned value of stress in case (B).
[out]intnlis set to intnl_old in case (A) or (C), and the returned value of intnl in case (B)
[out]pmZero in case (A) or (C), otherwise the plastic multipliers needed to bring about the returnMap in case (B)
[in/out]cumulative_pm cumulative plastic multipliers, updated in case (B), otherwise left untouched
[out]delta_dpis unchanged in case (A) or (C), and is set to the change in plastic strain in case(B)
[out]yfwill contain the yield function values at (stress, intnl)
[out]num_successful_plastic_returnswill be 0 for (A) and (C), and 1 for (B)
Returns
true in case (A) and (B), and false in case (C)

If all models actually signal "elastic" by returning true from their returnMap, and by returning model_plastically_active=0, then yf will contain the yield function values num_successful_plastic_returns will be zero intnl = intnl_old delta_dp will be unchanged from its input value stress will be set to trial_stress pm will be zero cumulative_pm will be unchanged return value will be true num_successful_plastic_returns = 0

If only one model signals "plastically active" by returning true from its returnMap, and by returning model_plastically_active=1, then yf will contain the yield function values num_successful_plastic_returns will be one intnl will be set by the returnMap algorithm delta_dp will be set by the returnMap algorithm stress will be set by the returnMap algorithm pm will be nonzero for the single model, and zero for other models cumulative_pm will be updated return value will be true num_successful_plastic_returns = 1

If >1 model signals "plastically active" or if >=1 model's returnMap fails, then yf will contain the yield function values num_successful_plastic_returns will be set appropriately intnl = intnl_old delta_dp will be unchanged from its input value stress will be set to trial_stress pm will be zero cumulative_pm will be unchanged return value will be true if all returnMap functions returned true, otherwise it will be false num_successful_plastic_returns is set appropriately

Definition at line 599 of file MultiPlasticityRawComponentAssembler.C.

611 {
612  mooseAssert(intnl_old.size() == _num_models,
613  "returnMapAll: Incorrect size of internal parameters");
614  mooseAssert(intnl.size() == _num_models, "returnMapAll: Incorrect size of internal parameters");
615  mooseAssert(pm.size() == _num_surfaces, "returnMapAll: Incorrect size of pm");
616 
617  num_successful_plastic_returns = 0;
618  yf.resize(0);
619  pm.assign(_num_surfaces, 0.0);
620 
621  RankTwoTensor returned_stress; // each model will give a returned_stress. if only one model is
622  // plastically active, i set stress=returned_stress, so as to
623  // record this returned value
624  std::vector<Real> model_f;
625  RankTwoTensor model_delta_dp;
626  std::vector<Real> model_pm;
627  bool trial_stress_inadmissible;
628  bool successful_return = true;
629  unsigned the_single_plastic_model = 0;
630  bool using_custom_return_map = true;
631 
632  // run through all the plastic models, performing their
633  // returnMap algorithms.
634  // If one finds (trial_stress, intnl) inadmissible and
635  // successfully returns, break from the loop to evaluate
636  // all the others at that returned stress
637  for (unsigned model = 0; model < _num_models; ++model)
638  {
639  if (using_custom_return_map)
640  {
641  model_pm.assign(_f[model]->numberSurfaces(), 0.0);
642  bool model_returned = _f[model]->returnMap(trial_stress,
643  intnl_old[model],
644  E_ijkl,
645  ep_plastic_tolerance,
646  returned_stress,
647  intnl[model],
648  model_pm,
649  model_delta_dp,
650  model_f,
651  trial_stress_inadmissible);
652  if (!trial_stress_inadmissible)
653  {
654  // in the elastic zone: record the yield-function values (returned_stress, intnl, model_pm
655  // and model_delta_dp are undefined)
656  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces();
657  ++model_surface)
658  yf.push_back(model_f[model_surface]);
659  }
660  else if (trial_stress_inadmissible && !model_returned)
661  {
662  // in the plastic zone, and the customized returnMap failed
663  // for some reason (or wasn't implemented). The coder
664  // should have correctly returned model_f(trial_stress, intnl_old)
665  // so record them
666  // (returned_stress, intnl, model_pm and model_delta_dp are undefined)
667  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces();
668  ++model_surface)
669  yf.push_back(model_f[model_surface]);
670  // now there's almost zero point in using the custom
671  // returnMap functions
672  using_custom_return_map = false;
673  successful_return = false;
674  }
675  else
676  {
677  // in the plastic zone, and the customized returnMap
678  // succeeded.
679  // record the first returned_stress and delta_dp if everything is going OK
680  // as they could be the actual answer
681  if (trial_stress_inadmissible)
682  num_successful_plastic_returns++;
683  the_single_plastic_model = model;
684  stress = returned_stress;
685  // note that i break here, and don't push_back
686  // model_f to yf. So now yf contains only the values of
687  // model_f from previous models to the_single_plastic_model
688  // also i don't set delta_dp = model_delta_dp yet, because
689  // i might find problems later on
690  // also, don't increment cumulative_pm for the same reason
691 
692  break;
693  }
694  }
695  else
696  {
697  // not using custom returnMap functions because one
698  // has already failed and that one said trial_stress
699  // was inadmissible. So now calculate the yield functions
700  // at the trial stress
701  _f[model]->yieldFunctionV(trial_stress, intnl_old[model], model_f);
702  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
703  yf.push_back(model_f[model_surface]);
704  }
705  }
706 
707  if (num_successful_plastic_returns == 0)
708  {
709  // here either all the models were elastic (successful_return=true),
710  // or some were plastic and either the customized returnMap failed
711  // or wasn't implemented (successful_return=false).
712  // In either case, have to set the following:
713  stress = trial_stress;
714  for (unsigned model = 0; model < _num_models; ++model)
715  intnl[model] = intnl_old[model];
716  return successful_return;
717  }
718 
719  // Now we know that num_successful_plastic_returns == 1 and all the other
720  // models (with model number < the_single_plastic_model) must have been
721  // admissible at (trial_stress, intnl). However, all models might
722  // not be admissible at (trial_stress, intnl), so must check that
723  std::vector<Real> yf_at_returned_stress(0);
724  bool all_admissible = true;
725  for (unsigned model = 0; model < _num_models; ++model)
726  {
727  if (model == the_single_plastic_model)
728  {
729  // no need to spend time calculating the yield function: we know its admissible
730  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
731  yf_at_returned_stress.push_back(model_f[model_surface]);
732  continue;
733  }
734  _f[model]->yieldFunctionV(stress, intnl_old[model], model_f);
735  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
736  {
737  if (model_f[model_surface] > _f[model]->_f_tol)
738  // bummer, this model is not admissible at the returned_stress
739  all_admissible = false;
740  yf_at_returned_stress.push_back(model_f[model_surface]);
741  }
742  if (!all_admissible)
743  // no point in continuing computing yield functions
744  break;
745  }
746 
747  if (!all_admissible)
748  {
749  // we tried using the returned value of stress predicted by
750  // the_single_plastic_model, but it wasn't admissible according
751  // to other plastic models. We need to set:
752  stress = trial_stress;
753  for (unsigned model = 0; model < _num_models; ++model)
754  intnl[model] = intnl_old[model];
755  // and calculate the remainder of the yield functions at trial_stress
756  for (unsigned model = the_single_plastic_model; model < _num_models; ++model)
757  {
758  _f[model]->yieldFunctionV(trial_stress, intnl[model], model_f);
759  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
760  yf.push_back(model_f[model_surface]);
761  }
762  num_successful_plastic_returns = 0;
763  return false;
764  }
765 
766  // So the customized returnMap algorithm can provide a returned
767  // (stress, intnl) configuration, and that is admissible according
768  // to all plastic models
769  yf.resize(0);
770  for (unsigned surface = 0; surface < yf_at_returned_stress.size(); ++surface)
771  yf.push_back(yf_at_returned_stress[surface]);
772  delta_dp = model_delta_dp;
773  for (unsigned model_surface = 0; model_surface < _f[the_single_plastic_model]->numberSurfaces();
774  ++model_surface)
775  {
776  cumulative_pm[_surfaces_given_model[the_single_plastic_model][model_surface]] +=
777  model_pm[model_surface];
778  pm[_surfaces_given_model[the_single_plastic_model][model_surface]] = model_pm[model_surface];
779  }
780  custom_model = the_single_plastic_model;
781  return true;
782 }

Referenced by quickStep().

◆ rot()

RankTwoTensor ComputeMultiPlasticityStress::rot ( const RankTwoTensor &  tens)
private

Definition at line 312 of file ComputeMultiPlasticityStress.C.

313 {
314  if (!_n_supplied)
315  return tens;
316  return tens.rotated(_rot);
317 }

Referenced by computeQpStress().

◆ singleStep()

bool ComputeMultiPlasticityStress::singleStep ( Real &  nr_res2,
RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
std::vector< Real > &  intnl,
std::vector< Real > &  pm,
RankTwoTensor &  delta_dp,
const RankFourTensor &  E_inv,
std::vector< Real > &  f,
RankTwoTensor &  epp,
std::vector< Real > &  ic,
std::vector< bool > &  active,
DeactivationSchemeEnum  deactivation_scheme,
bool &  linesearch_needed,
bool &  ld_encountered 
)
protectedvirtual

Performs a single Newton-Raphson + linesearch step Constraints are deactivated and the step is re-done if deactivation_scheme is set appropriately.

Parameters
[in,out]nr_res2Residual-squared that the line-search will reduce
[in,out]stressstress
[in]intnl_oldold values of the internal parameters
[in,out]intnlinternal parameters
[in,out]pmplastic multipliers
[in,out]delta_dpChange in plastic strain from start of "time" step to current configuration (plastic_strain - plastic_strain_old)
[in]E_invInverse of the elasticity tensor
[in,out]fYield function(s). Upon successful exit only the active constraints are contained in f
[in,out]eppPlastic strain increment constraint
[in,out]icInternal constraint. Upon successful exit only the active constraints are contained in ic
activeThe active constraints. This is may be modified, depending upon deactivation_scheme
deactivation_schemeThe scheme used for deactivating constraints
[out]linesearch_neededTrue if a linesearch was employed during this Newton-Raphson step
[out]ld_encounteredTrue if a linear-dependence of the flow directions was encountered at any stage during the Newton-Raphson proceedure
Returns
true if the step was successful, ie, if the linesearch was successful and the number of constraints wasn't reduced to zero via deactivation

Definition at line 1082 of file ComputeMultiPlasticityStress.C.

1096 {
1097  bool successful_step; // return value
1098 
1099  Real nr_res2_before_step = nr_res2;
1100  RankTwoTensor stress_before_step;
1101  std::vector<Real> intnl_before_step;
1102  std::vector<Real> pm_before_step;
1103  RankTwoTensor delta_dp_before_step;
1104 
1105  if (deactivation_scheme == optimized)
1106  {
1107  // we potentially use the "before_step" quantities, so record them here
1108  stress_before_step = stress;
1109  intnl_before_step.resize(_num_models);
1110  for (unsigned model = 0; model < _num_models; ++model)
1111  intnl_before_step[model] = intnl[model];
1112  pm_before_step.resize(_num_surfaces);
1113  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1114  pm_before_step[surface] = pm[surface];
1115  delta_dp_before_step = delta_dp;
1116  }
1117 
1118  // During the Newton-Raphson procedure, we'll be
1119  // changing the following parameters in order to
1120  // (attempt to) satisfy the constraints.
1121  RankTwoTensor dstress; // change in stress
1122  std::vector<Real> dpm; // change in plasticity multipliers ("consistency parameters"). For ALL
1123  // contraints (active and deactive)
1124  std::vector<Real>
1125  dintnl; // change in internal parameters. For ALL internal params (active and deactive)
1126 
1127  // The constraints that have been deactivated for this NR step
1128  // due to the flow directions being linearly dependent
1129  std::vector<bool> deact_ld;
1130  deact_ld.assign(_num_surfaces, false);
1131 
1132  /* After NR and linesearch, if _deactivation_scheme == "optimized", the
1133  * active plasticity multipliers are checked for non-negativity. If some
1134  * are negative then they are deactivated forever, and the NR step is
1135  * re-done starting from the _before_step quantities recorded above
1136  */
1137  bool constraints_changing = true;
1138  bool reinstated_actives;
1139  while (constraints_changing)
1140  {
1141  // calculate dstress, dpm and dintnl for one full Newton-Raphson step
1142  nrStep(stress, intnl_old, intnl, pm, E_inv, delta_dp, dstress, dpm, dintnl, active, deact_ld);
1143 
1144  for (unsigned surface = 0; surface < deact_ld.size(); ++surface)
1145  if (deact_ld[surface])
1146  {
1147  ld_encountered = true;
1148  break;
1149  }
1150 
1151  // perform a line search
1152  // The line-search will exit with updated values
1153  successful_step = lineSearch(nr_res2,
1154  stress,
1155  intnl_old,
1156  intnl,
1157  pm,
1158  E_inv,
1159  delta_dp,
1160  dstress,
1161  dpm,
1162  dintnl,
1163  f,
1164  epp,
1165  ic,
1166  active,
1167  deact_ld,
1168  linesearch_needed);
1169 
1170  if (!successful_step)
1171  // completely bomb out
1172  return successful_step;
1173 
1174  // See if any active constraints need to be removed, and the step re-done
1175  constraints_changing = false;
1176  if (deactivation_scheme == optimized)
1177  {
1178  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1179  if (active[surface] && pm[surface] < 0.0)
1180  constraints_changing = true;
1181  }
1182 
1183  if (constraints_changing)
1184  {
1185  stress = stress_before_step;
1186  delta_dp = delta_dp_before_step;
1187  nr_res2 = nr_res2_before_step;
1188  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1189  {
1190  if (active[surface] && pm[surface] < 0.0)
1191  {
1192  // turn off the constraint forever
1193  active[surface] = false;
1194  pm_before_step[surface] = 0.0;
1195  intnl_before_step[modelNumber(surface)] =
1196  intnl_old[modelNumber(surface)]; // don't want to muck-up hardening!
1197  }
1198  intnl[modelNumber(surface)] = intnl_before_step[modelNumber(surface)];
1199  pm[surface] = pm_before_step[surface];
1200  }
1201  if (numberActive(active) == 0)
1202  {
1203  // completely bomb out
1204  successful_step = false;
1205  return successful_step;
1206  }
1207  }
1208 
1209  // reinstate any active values that have been turned off due to linear-dependence
1210  reinstated_actives = reinstateLinearDependentConstraints(deact_ld);
1211  } // ends "constraints_changing" loop
1212 
1213  // if active constraints were reinstated then nr_res2 needs to be re-calculated so it is correct
1214  // upson returning
1215  if (reinstated_actives)
1216  {
1217  bool completely_converged = true;
1218  if (successful_step && nr_res2 < 0.5)
1219  {
1220  // Here we have converged to the correct solution if
1221  // all the yield functions are < 0. Excellent!
1222  //
1223  // This is quite tricky - perhaps i can refactor to make it more obvious.
1224  //
1225  // Because actives are now reinstated, the residual2
1226  // calculation below will give nr_res2 > 0.5, because it won't
1227  // realise that we only had to set the active-but-not-deactivated f = 0,
1228  // and not the entire active set. If we pass that nr_res2 back from
1229  // this function then the calling function will not realise we've converged!
1230  // Therefore, check for this case
1231  unsigned ind = 0;
1232  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
1233  if (active[surface])
1234  if (f[ind++] > _f[modelNumber(surface)]->_f_tol)
1235  completely_converged = false;
1236  }
1237  else
1238  completely_converged = false;
1239 
1240  if (!completely_converged)
1241  nr_res2 = residual2(pm, f, epp, ic, active, deact_ld);
1242  }
1243 
1244  return successful_step;
1245 }

Referenced by returnMap().

◆ singularValuesOfR()

int MultiPlasticityLinearSystem::singularValuesOfR ( const std::vector< RankTwoTensor > &  r,
std::vector< Real > &  s 
)
privatevirtualinherited

Performs a singular-value decomposition of r and returns the singular values.

Example: If r has size 5 then the singular values of the following matrix are returned: ( r[0](0,0) r[0](0,1) r[0](0,2) r[0](1,1) r[0](1,2) r[0](2,2) ) ( r[1](0,0) r[1](0,1) r[1](0,2) r[1](1,1) r[1](1,2) r[1](2,2) ) a = ( r[2](0,0) r[2](0,1) r[2](0,2) r[2](1,1) r[2](1,2) r[2](2,2) ) ( r[3](0,0) r[3](0,1) r[3](0,2) r[3](1,1) r[3](1,2) r[3](2,2) ) ( r[4](0,0) r[4](0,1) r[4](0,2) r[4](1,1) r[4](1,2) r[4](2,2) )

Parameters
rThe flow directions
[out]sThe singular values
Returns
The return value from the PETSc LAPACK gesvd reoutine

Definition at line 47 of file MultiPlasticityLinearSystem.C.

49 {
50  int bm = r.size();
51  int bn = 6;
52 
53  s.resize(std::min(bm, bn));
54 
55  // prepare for gesvd or gesdd routine provided by PETSc
56  // Want to find the singular values of matrix
57  // ( r[0](0,0) r[0](0,1) r[0](0,2) r[0](1,1) r[0](1,2) r[0](2,2) )
58  // ( r[1](0,0) r[1](0,1) r[1](0,2) r[1](1,1) r[1](1,2) r[1](2,2) )
59  // a = ( r[2](0,0) r[2](0,1) r[2](0,2) r[2](1,1) r[2](1,2) r[2](2,2) )
60  // ( r[3](0,0) r[3](0,1) r[3](0,2) r[3](1,1) r[3](1,2) r[3](2,2) )
61  // ( r[4](0,0) r[4](0,1) r[4](0,2) r[4](1,1) r[4](1,2) r[4](2,2) )
62  // bm = 5
63 
64  std::vector<double> a(bm * 6);
65  // Fill in the a "matrix" by going down columns
66  unsigned ind = 0;
67  for (int col = 0; col < 3; ++col)
68  for (int row = 0; row < bm; ++row)
69  a[ind++] = r[row](0, col);
70  for (int col = 3; col < 5; ++col)
71  for (int row = 0; row < bm; ++row)
72  a[ind++] = r[row](1, col - 2);
73  for (int row = 0; row < bm; ++row)
74  a[ind++] = r[row](2, 2);
75 
76  // u and vt are dummy variables because they won't
77  // get referenced due to the "N" and "N" choices
78  int sizeu = 1;
79  std::vector<double> u(sizeu);
80  int sizevt = 1;
81  std::vector<double> vt(sizevt);
82 
83  int sizework = 16 * (bm + 6); // this is above the lowerbound specified in the LAPACK doco
84  std::vector<double> work(sizework);
85 
86  int info;
87 
88  LAPACKgesvd_("N",
89  "N",
90  &bm,
91  &bn,
92  &a[0],
93  &bm,
94  &s[0],
95  &u[0],
96  &sizeu,
97  &vt[0],
98  &sizevt,
99  &work[0],
100  &sizework,
101  &info);
102 
103  return info;
104 }

Referenced by MultiPlasticityLinearSystem::eliminateLinearDependence().

◆ validParams()

InputParameters ComputeMultiPlasticityStress::validParams ( )
static

Definition at line 24 of file ComputeMultiPlasticityStress.C.

25 {
26  InputParameters params = ComputeStressBase::validParams();
27  params += MultiPlasticityDebugger::validParams();
28  params.addClassDescription("Base class for multi-surface finite-strain plasticity");
29  params.addRangeCheckedParam<unsigned int>("max_NR_iterations",
30  20,
31  "max_NR_iterations>0",
32  "Maximum number of Newton-Raphson iterations allowed");
33  params.addRequiredParam<Real>("ep_plastic_tolerance",
34  "The Newton-Raphson process is only deemed "
35  "converged if the plastic strain increment "
36  "constraints have L2 norm less than this.");
37  params.addRangeCheckedParam<Real>(
38  "min_stepsize",
39  0.01,
40  "min_stepsize>0 & min_stepsize<=1",
41  "If ordinary Newton-Raphson + line-search fails, then the applied strain increment is "
42  "subdivided, and the return-map is tried again. This parameter is the minimum fraction of "
43  "applied strain increment that may be applied before the algorithm gives up entirely");
44  params.addRangeCheckedParam<Real>("max_stepsize_for_dumb",
45  0.01,
46  "max_stepsize_for_dumb>0 & max_stepsize_for_dumb<=1",
47  "If your deactivation_scheme is 'something_to_dumb', then "
48  "'dumb' will only be used if the stepsize falls below this "
49  "value. This parameter is useful because the 'dumb' scheme is "
50  "computationally expensive");
51  MooseEnum deactivation_scheme("optimized safe dumb optimized_to_safe safe_to_dumb "
52  "optimized_to_safe_to_dumb optimized_to_dumb",
53  "optimized");
54  params.addParam<MooseEnum>(
55  "deactivation_scheme",
56  deactivation_scheme,
57  "Scheme by which constraints are deactivated. (NOTE: This is irrelevant if there is only "
58  "one yield surface.) safe: return to the yield surface and then deactivate constraints with "
59  "negative plasticity multipliers. optimized: deactivate a constraint as soon as its "
60  "plasticity multiplier becomes negative. dumb: iteratively try all combinations of active "
61  "constraints until the solution is found. You may specify fall-back options. Eg "
62  "optimized_to_safe: first use 'optimized', and if that fails, try the return with 'safe'.");
63  params.addParam<RealVectorValue>(
64  "transverse_direction",
65  "If this parameter is provided, before the return-map algorithm is "
66  "called a rotation is performed so that the 'z' axis in the new "
67  "frame lies along the transverse_direction in the original frame. "
68  "After returning, the inverse rotation is performed. The "
69  "transverse_direction will itself rotate with large strains. This "
70  "is so that transversely-isotropic plasticity models may be easily "
71  "defined in the frame where the isotropy holds in the x-y plane.");
72  params.addParam<bool>("ignore_failures",
73  false,
74  "The return-map algorithm will return with the best admissible "
75  "stresses and internal parameters that it can, even if they don't "
76  "fully correspond to the applied strain increment. To speed "
77  "computations, this flag can be set to true, the max_NR_iterations "
78  "set small, and the min_stepsize large.");
79  MooseEnum tangent_operator("elastic linear nonlinear", "nonlinear");
80  params.addParam<MooseEnum>("tangent_operator",
81  tangent_operator,
82  "Type of tangent operator to return. 'elastic': return the "
83  "elasticity tensor. 'linear': return the consistent tangent operator "
84  "that is correct for plasticity with yield functions linear in "
85  "stress. 'nonlinear': return the full, general consistent tangent "
86  "operator. The calculations assume the hardening potentials are "
87  "independent of stress and hardening parameters.");
88  params.addParam<bool>("perform_finite_strain_rotations",
89  true,
90  "Tensors are correctly rotated in "
91  "finite-strain simulations. For "
92  "optimal performance you can set "
93  "this to 'false' if you are only "
94  "ever using small strains");
95  params.addClassDescription("Material for multi-surface finite-strain plasticity");
96  return params;
97 }

◆ yieldFunction()

void MultiPlasticityRawComponentAssembler::yieldFunction ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  f 
)
protectedvirtualinherited

The active yield function(s)

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
activeset of active constraints - only the active yield functions are put into "f"
[out]fthe yield function (or functions in the case of multisurface plasticity)

Definition at line 98 of file MultiPlasticityRawComponentAssembler.C.

102 {
103  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
104  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
105 
106  f.resize(0);
107  std::vector<unsigned int> active_surfaces_of_model;
108  std::vector<unsigned int>::iterator active_surface;
109  std::vector<Real> model_f;
110  for (unsigned model = 0; model < _num_models; ++model)
111  {
112  activeModelSurfaces(model, active, active_surfaces_of_model);
113  if (active_surfaces_of_model.size() > 0)
114  {
115  _f[model]->yieldFunctionV(stress, intnl[model], model_f);
116  for (active_surface = active_surfaces_of_model.begin();
117  active_surface != active_surfaces_of_model.end();
118  ++active_surface)
119  f.push_back(model_f[*active_surface]);
120  }
121  }
122 }

Referenced by buildDumbOrder(), MultiPlasticityLinearSystem::calculateConstraints(), checkAdmissible(), MultiPlasticityDebugger::fddyieldFunction_dintnl(), MultiPlasticityDebugger::fddyieldFunction_dstress(), and returnMap().

Member Data Documentation

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

Base name prepended to all material property names to allow for multi-material systems.

Definition at line 45 of file ComputeStressBase.h.

Referenced by ComputeLinearElasticStress::initialSetup(), and ComputeCosseratLinearElasticStress::initialSetup().

◆ _constraints_added

MaterialProperty<Real>& ComputeMultiPlasticityStress::_constraints_added
protected

Whether constraints were added in during the latest Newton-Raphson process (1 if true, 0 otherwise)

Definition at line 132 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _cosserat

bool ComputeMultiPlasticityStress::_cosserat
protected

whether Cosserat mechanics should be used

Definition at line 156 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), initQpStatefulProperties(), postReturnMap(), and preReturnMap().

◆ _couple_stress

MaterialProperty<RankTwoTensor>* ComputeMultiPlasticityStress::_couple_stress
protected

the Cosserat couple-stress

Definition at line 165 of file ComputeMultiPlasticityStress.h.

◆ _couple_stress_old

const MaterialProperty<RankTwoTensor>* ComputeMultiPlasticityStress::_couple_stress_old
protected

the old value of Cosserat couple-stress

Definition at line 168 of file ComputeMultiPlasticityStress.h.

◆ _cumulative_pm

std::vector<Real> ComputeMultiPlasticityStress::_cumulative_pm
protected

the sum of the plastic multipliers over all the sub-steps.

This is used for calculating the consistent tangent operator

Definition at line 71 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and plasticStep().

◆ _curvature

const MaterialProperty<RankTwoTensor>* ComputeMultiPlasticityStress::_curvature
protected

The Cosserat curvature strain.

Definition at line 159 of file ComputeMultiPlasticityStress.h.

◆ _deactivation_scheme

enum ComputeMultiPlasticityStress::DeactivationSchemeEnum ComputeMultiPlasticityStress::_deactivation_scheme
protected

Referenced by canChangeScheme(), changeScheme(), ComputeMultiPlasticityStress(), and returnMap().

◆ _dummy_pm

std::vector<Real> ComputeMultiPlasticityStress::_dummy_pm
protected

dummy "consistency parameters" (plastic multipliers) used in quickStep when called from computeQpStress

Definition at line 65 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _elastic_flexural_rigidity_tensor

const MaterialProperty<RankFourTensor>* ComputeMultiPlasticityStress::_elastic_flexural_rigidity_tensor
protected

The Cosserat elastic flexural rigidity tensor.

Definition at line 162 of file ComputeMultiPlasticityStress.h.

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

Elastic strain material property.

Definition at line 52 of file ComputeStressBase.h.

Referenced by ComputeSmearedCrackingStress::computeCrackStrainAndOrientation(), ComputeLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeMultipleInelasticStress::finiteStrainRotation(), and ComputeStressBase::initQpStatefulProperties().

◆ _elastic_strain_old

const MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_elastic_strain_old
protected

Old value of elastic strain.

Definition at line 153 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress().

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& ComputeMultiPlasticityStress::_elasticity_tensor
protected

Elasticity tensor material property.

Definition at line 105 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress().

◆ _elasticity_tensor_name

const std::string ComputeMultiPlasticityStress::_elasticity_tensor_name
protected

Name of the elasticity tensor material property.

Definition at line 103 of file ComputeMultiPlasticityStress.h.

◆ _epp_tol

Real ComputeMultiPlasticityStress::_epp_tol
protected

Tolerance on the plastic strain increment ("direction") constraint.

Definition at line 62 of file ComputeMultiPlasticityStress.h.

Referenced by ComputeMultiPlasticityStress(), quickStep(), and residual2().

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 55 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _f

std::vector<const TensorMechanicsPlasticModel *> MultiPlasticityRawComponentAssembler::_f
protectedinherited

User objects that define the yield functions, flow potentials, etc.

Definition at line 75 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::activeModelSurfaces(), MultiPlasticityRawComponentAssembler::activeSurfaces(), MultiPlasticityRawComponentAssembler::anyActiveSurfaces(), applyKuhnTucker(), MultiPlasticityRawComponentAssembler::buildActiveConstraints(), MultiPlasticityRawComponentAssembler::buildActiveConstraintsJoint(), MultiPlasticityRawComponentAssembler::buildActiveConstraintsRock(), canAddConstraints(), checkAdmissible(), checkKuhnTucker(), MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), MultiPlasticityLinearSystem::MultiPlasticityLinearSystem(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), quickStep(), residual2(), returnMap(), MultiPlasticityRawComponentAssembler::returnMapAll(), singleStep(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ _fspb_debug

MooseEnum MultiPlasticityDebugger::_fspb_debug
protectedinherited

none - don't do any debugging crash - currently inactive jacobian - check the jacobian entries jacobian_and_linear_system - check entire jacobian and check that Ax=b

Definition at line 62 of file MultiPlasticityDebugger.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _fspb_debug_intnl

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_intnl
protectedinherited

Debug the Jacobian entires at these internal parameters.

Definition at line 71 of file MultiPlasticityDebugger.h.

Referenced by MultiPlasticityDebugger::checkDerivatives(), MultiPlasticityDebugger::checkJacobian(), MultiPlasticityDebugger::checkSolution(), MultiPlasticityDebugger::outputAndCheckDebugParameters(), and plasticStep().

◆ _fspb_debug_intnl_change

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_intnl_change
protectedinherited

Debug finite-differencing parameters for the internal parameters.

Definition at line 80 of file MultiPlasticityDebugger.h.

Referenced by MultiPlasticityDebugger::fddflowPotential_dintnl(), MultiPlasticityDebugger::fddyieldFunction_dintnl(), MultiPlasticityDebugger::fdJacobian(), and MultiPlasticityDebugger::outputAndCheckDebugParameters().

◆ _fspb_debug_pm

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_pm
protectedinherited

Debug the Jacobian entires at these plastic multipliers.

Definition at line 68 of file MultiPlasticityDebugger.h.

Referenced by MultiPlasticityDebugger::checkJacobian(), MultiPlasticityDebugger::checkSolution(), MultiPlasticityDebugger::outputAndCheckDebugParameters(), and plasticStep().

◆ _fspb_debug_pm_change

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_pm_change
protectedinherited

Debug finite-differencing parameters for the plastic multipliers.

Definition at line 77 of file MultiPlasticityDebugger.h.

Referenced by MultiPlasticityDebugger::fdJacobian(), and MultiPlasticityDebugger::outputAndCheckDebugParameters().

◆ _fspb_debug_stress

RankTwoTensor MultiPlasticityDebugger::_fspb_debug_stress
protectedinherited

Debug the Jacobian entries at this stress.

Definition at line 65 of file MultiPlasticityDebugger.h.

Referenced by MultiPlasticityDebugger::checkDerivatives(), MultiPlasticityDebugger::checkJacobian(), MultiPlasticityDebugger::checkSolution(), MultiPlasticityDebugger::outputAndCheckDebugParameters(), and plasticStep().

◆ _fspb_debug_stress_change

Real MultiPlasticityDebugger::_fspb_debug_stress_change
protectedinherited

Debug finite-differencing parameter for the stress.

Definition at line 74 of file MultiPlasticityDebugger.h.

Referenced by MultiPlasticityDebugger::fddflowPotential_dstress(), MultiPlasticityDebugger::fddyieldFunction_dstress(), MultiPlasticityDebugger::fdJacobian(), and MultiPlasticityDebugger::outputAndCheckDebugParameters().

◆ _ignore_failures

bool ComputeMultiPlasticityStress::_ignore_failures
protected

Even if the returnMap fails, return the best values found for stress and internal parameters.

Definition at line 51 of file ComputeMultiPlasticityStress.h.

Referenced by plasticStep().

◆ _initial_stress_fcn

std::vector<const Function *> ComputeStressBase::_initial_stress_fcn
protectedinherited

initial stress components

Definition at line 58 of file ComputeStressBase.h.

◆ _intnl

MaterialProperty<std::vector<Real> >& ComputeMultiPlasticityStress::_intnl
protected

internal parameters

Definition at line 114 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _intnl_old

const MaterialProperty<std::vector<Real> >& ComputeMultiPlasticityStress::_intnl_old
protected

old values of internal parameters

Definition at line 117 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and plasticStep().

◆ _iter

MaterialProperty<Real>& ComputeMultiPlasticityStress::_iter
protected

Number of Newton-Raphson iterations used in the return-map.

Definition at line 123 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

derivative of stress w.r.t. strain (_dstress_dstrain)

Definition at line 61 of file ComputeStressBase.h.

Referenced by ComputeStrainIncrementBasedStress::computeQpJacobian(), FiniteStrainHyperElasticViscoPlastic::computeQpJacobian(), ComputeMultipleInelasticCosseratStress::computeQpJacobianMult(), ComputeMultipleInelasticStress::computeQpJacobianMult(), ComputeLinearElasticStress::computeQpStress(), ComputeDamageStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), ComputeLinearElasticPFFractureStress::computeStressSpectral(), FiniteStrainUObasedCP::elasticTangentModuli(), FiniteStrainUObasedCP::elastoPlasticTangentModuli(), ComputeMultipleInelasticStress::finiteStrainRotation(), postReturnMap(), FiniteStrainCrystalPlasticity::postSolveQp(), FiniteStrainCrystalPlasticity::preSolveQp(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

◆ _Jacobian_mult_couple

MaterialProperty<RankFourTensor>* ComputeMultiPlasticityStress::_Jacobian_mult_couple
protected

derivative of couple-stress w.r.t. curvature

Definition at line 171 of file ComputeMultiPlasticityStress.h.

◆ _ld_encountered

MaterialProperty<Real>& ComputeMultiPlasticityStress::_ld_encountered
protected

Whether linear-dependence was encountered in the latest Newton-Raphson process (1 if true, 0 otherwise)

Definition at line 129 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _linesearch_needed

MaterialProperty<Real>& ComputeMultiPlasticityStress::_linesearch_needed
protected

Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise)

Definition at line 126 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _max_iter

unsigned int ComputeMultiPlasticityStress::_max_iter
protected

Maximum number of Newton-Raphson iterations allowed.

Definition at line 42 of file ComputeMultiPlasticityStress.h.

Referenced by returnMap().

◆ _max_stepsize_for_dumb

Real ComputeMultiPlasticityStress::_max_stepsize_for_dumb
protected

"dumb" deactivation will only be used if the stepsize falls below this quantity

Definition at line 48 of file ComputeMultiPlasticityStress.h.

Referenced by plasticStep().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

Mechanical strain material property.

Definition at line 48 of file ComputeStressBase.h.

Referenced by ComputeLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), and ComputeLinearElasticPFFractureStress::computeStressSpectral().

◆ _min_f_tol

Real MultiPlasticityLinearSystem::_min_f_tol
protectedinherited

Minimum value of the _f_tol parameters for the Yield Function User Objects.

Definition at line 136 of file MultiPlasticityLinearSystem.h.

Referenced by MultiPlasticityLinearSystem::eliminateLinearDependence(), and MultiPlasticityLinearSystem::MultiPlasticityLinearSystem().

◆ _min_stepsize

Real ComputeMultiPlasticityStress::_min_stepsize
protected

Minimum fraction of applied strain that may be applied during adaptive stepsizing.

Definition at line 45 of file ComputeMultiPlasticityStress.h.

Referenced by plasticStep().

◆ _model_given_surface

std::vector<unsigned int> MultiPlasticityRawComponentAssembler::_model_given_surface
privateinherited

given a surface number, this returns the model number

Definition at line 295 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::modelNumber(), and MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _model_surface_given_surface

std::vector<unsigned int> MultiPlasticityRawComponentAssembler::_model_surface_given_surface
privateinherited

given a surface number, this returns the corresponding-model's internal surface number

Definition at line 298 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _my_curvature

RankTwoTensor ComputeMultiPlasticityStress::_my_curvature
protected

Curvature that can be rotated by this class, and split into multiple increments (ie, its not const)

Definition at line 183 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), postReturnMap(), and preReturnMap().

◆ _my_elasticity_tensor

RankFourTensor ComputeMultiPlasticityStress::_my_elasticity_tensor
protected

Elasticity tensor that can be rotated by this class (ie, its not const)

Definition at line 174 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), plasticStep(), postReturnMap(), and preReturnMap().

◆ _my_flexural_rigidity_tensor

RankFourTensor ComputeMultiPlasticityStress::_my_flexural_rigidity_tensor
protected

Flexual rigidity tensor that can be rotated by this class (ie, its not const)

Definition at line 180 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), postReturnMap(), and preReturnMap().

◆ _my_strain_increment

RankTwoTensor ComputeMultiPlasticityStress::_my_strain_increment
protected

Strain increment that can be rotated by this class, and split into multiple increments (ie, its not const)

Definition at line 177 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), postReturnMap(), and preReturnMap().

◆ _n

MaterialProperty<RealVectorValue>& ComputeMultiPlasticityStress::_n
protected

current value of transverse direction

Definition at line 135 of file ComputeMultiPlasticityStress.h.

Referenced by initQpStatefulProperties(), postReturnMap(), and preReturnMap().

◆ _n_input

RealVectorValue ComputeMultiPlasticityStress::_n_input
protected

the supplied transverse direction vector

Definition at line 94 of file ComputeMultiPlasticityStress.h.

Referenced by ComputeMultiPlasticityStress(), and initQpStatefulProperties().

◆ _n_old

const MaterialProperty<RealVectorValue>& ComputeMultiPlasticityStress::_n_old
protected

old value of transverse direction

Definition at line 138 of file ComputeMultiPlasticityStress.h.

Referenced by postReturnMap().

◆ _n_supplied

bool ComputeMultiPlasticityStress::_n_supplied
protected

User supplied the transverse direction vector.

Definition at line 91 of file ComputeMultiPlasticityStress.h.

Referenced by ComputeMultiPlasticityStress(), postReturnMap(), preReturnMap(), and rot().

◆ _num_models

unsigned int MultiPlasticityRawComponentAssembler::_num_models
protectedinherited

Number of plastic models for this material.

Definition at line 57 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints(), MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityLinearSystem::calculateRHS(), MultiPlasticityDebugger::checkSolution(), MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityDebugger::fddflowPotential_dintnl(), MultiPlasticityDebugger::fddyieldFunction_dintnl(), MultiPlasticityDebugger::fdJacobian(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), initQpStatefulProperties(), MultiPlasticityLinearSystem::MultiPlasticityLinearSystem(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), MultiPlasticityLinearSystem::nrStep(), MultiPlasticityDebugger::outputAndCheckDebugParameters(), plasticStep(), residual2(), returnMap(), MultiPlasticityRawComponentAssembler::returnMapAll(), singleStep(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ _num_surfaces

unsigned int MultiPlasticityRawComponentAssembler::_num_surfaces
protectedinherited

Number of surfaces within the plastic models.

For many situations this will be = _num_models since each model will contain just one surface. More generally it is >= _num_models. For instance, Mohr-Coulomb is a single model with 6 surfaces

Definition at line 66 of file MultiPlasticityRawComponentAssembler.h.

Referenced by activeCombinationNumber(), applyKuhnTucker(), MultiPlasticityRawComponentAssembler::buildActiveConstraints(), buildDumbOrder(), MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityLinearSystem::calculateRHS(), canAddConstraints(), canIncrementDumb(), changeScheme(), checkAdmissible(), MultiPlasticityDebugger::checkDerivatives(), MultiPlasticityDebugger::checkJacobian(), checkKuhnTucker(), MultiPlasticityDebugger::checkSolution(), ComputeMultiPlasticityStress(), computeQpStress(), consistentTangentOperator(), MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityDebugger::dof_included(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityLinearSystem::eliminateLinearDependence(), MultiPlasticityDebugger::fddflowPotential_dintnl(), MultiPlasticityDebugger::fddflowPotential_dstress(), MultiPlasticityDebugger::fddyieldFunction_dintnl(), MultiPlasticityDebugger::fddyieldFunction_dstress(), MultiPlasticityDebugger::fdJacobian(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), incrementDumb(), initQpStatefulProperties(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), MultiPlasticityLinearSystem::nrStep(), numberActive(), MultiPlasticityDebugger::outputAndCheckDebugParameters(), plasticStep(), reinstateLinearDependentConstraints(), residual2(), returnMap(), MultiPlasticityRawComponentAssembler::returnMapAll(), singleStep(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ _params

const InputParameters& MultiPlasticityRawComponentAssembler::_params
protectedinherited

Definition at line 54 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _perform_finite_strain_rotations

bool ComputeMultiPlasticityStress::_perform_finite_strain_rotations
protected

whether to perform the rotations necessary in finite-strain simulations

Definition at line 100 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress().

◆ _plastic_strain

MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_plastic_strain
protected

plastic strain

Definition at line 108 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), initQpStatefulProperties(), and postReturnMap().

◆ _plastic_strain_old

const MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_plastic_strain_old
protected

Old value of plastic strain.

Definition at line 111 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress().

◆ _rot

RealTensorValue ComputeMultiPlasticityStress::_rot
protected

rotation matrix that takes _n to (0, 0, 1)

Definition at line 97 of file ComputeMultiPlasticityStress.h.

Referenced by postReturnMap(), preReturnMap(), and rot().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_rotation_increment
protected

Rotation increment (coming from ComputeIncrementalSmallStrain, for example)

Definition at line 147 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and postReturnMap().

◆ _specialIC

MooseEnum MultiPlasticityRawComponentAssembler::_specialIC
protectedinherited

Allows initial set of active constraints to be chosen optimally.

Definition at line 72 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints(), and MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _strain_increment

const MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_strain_increment
protected

strain increment (coming from ComputeIncrementalSmallStrain, for example)

Definition at line 141 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress().

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Stress material property.

Definition at line 50 of file ComputeStressBase.h.

Referenced by ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), ComputeMultipleInelasticStress::computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticStress::computeQpStress(), ComputeDamageStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), ComputeLinearElasticPFFractureStress::computeStressSpectral(), ComputeMultipleInelasticStress::finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeSmearedCrackingStress::updateCrackingStateAndStress(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

◆ _stress_old

const MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_stress_old
protected

Old value of stress.

Definition at line 150 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress().

◆ _surfaces_given_model

std::vector<std::vector<unsigned int> > MultiPlasticityRawComponentAssembler::_surfaces_given_model
protectedinherited

_surfaces_given_model[model_number] = vector of surface numbers for this model

Definition at line 69 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::activeModelSurfaces(), MultiPlasticityRawComponentAssembler::activeSurfaces(), MultiPlasticityRawComponentAssembler::anyActiveSurfaces(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), quickStep(), and MultiPlasticityRawComponentAssembler::returnMapAll().

◆ _svd_tol

Real MultiPlasticityLinearSystem::_svd_tol
protectedinherited

Tolerance on the minimum ratio of singular values before flow-directions are deemed linearly dependent.

Definition at line 133 of file MultiPlasticityLinearSystem.h.

Referenced by MultiPlasticityLinearSystem::eliminateLinearDependence().

◆ _tangent_operator_type

enum ComputeMultiPlasticityStress::TangentOperatorEnum ComputeMultiPlasticityStress::_tangent_operator_type
protected

Referenced by consistentTangentOperator(), and quickStep().

◆ _total_strain_old

const MaterialProperty<RankTwoTensor>& ComputeMultiPlasticityStress::_total_strain_old
protected

Old value of total strain (coming from ComputeIncrementalSmallStrain, for example)

Definition at line 144 of file ComputeMultiPlasticityStress.h.

◆ _yf

MaterialProperty<std::vector<Real> >& ComputeMultiPlasticityStress::_yf
protected

yield functions

Definition at line 120 of file ComputeMultiPlasticityStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

The documentation for this class was generated from the following files:
ComputeMultiPlasticityStress::_n
MaterialProperty< RealVectorValue > & _n
current value of transverse direction
Definition: ComputeMultiPlasticityStress.h:135
ComputeMultiPlasticityStress::safe
Definition: ComputeMultiPlasticityStress.h:82
MultiPlasticityRawComponentAssembler::dflowPotential_dintnl
virtual void dflowPotential_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dr_dintnl)
The derivative of the active flow potentials with respect to the active internal parameters The UserO...
Definition: MultiPlasticityRawComponentAssembler.C:235
MultiPlasticityDebugger::validParams
static InputParameters validParams()
Definition: MultiPlasticityDebugger.C:17
ComputeMultiPlasticityStress::_linesearch_needed
MaterialProperty< Real > & _linesearch_needed
Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise)
Definition: ComputeMultiPlasticityStress.h:126
ComputeMultiPlasticityStress::_dummy_pm
std::vector< Real > _dummy_pm
dummy "consistency parameters" (plastic multipliers) used in quickStep when called from computeQpStre...
Definition: ComputeMultiPlasticityStress.h:65
ComputeMultiPlasticityStress::elastic
Definition: ComputeMultiPlasticityStress.h:56
MultiPlasticityRawComponentAssembler::modelNumber
unsigned int modelNumber(unsigned int surface)
returns the model number, given the surface number
Definition: MultiPlasticityRawComponentAssembler.C:785
ComputeMultiPlasticityStress::_ld_encountered
MaterialProperty< Real > & _ld_encountered
Whether linear-dependence was encountered in the latest Newton-Raphson process (1 if true,...
Definition: ComputeMultiPlasticityStress.h:129
ComputeMultiPlasticityStress::changeScheme
void changeScheme(const std::vector< bool > &initial_act, bool can_revert_to_dumb, const RankTwoTensor &initial_stress, const std::vector< Real > &intnl_old, DeactivationSchemeEnum &current_deactivation_scheme, std::vector< bool > &act, int &dumb_iteration, std::vector< unsigned int > &dumb_order)
Definition: ComputeMultiPlasticityStress.C:1050
ComputeMultiPlasticityStress::checkKuhnTucker
virtual bool checkKuhnTucker(const std::vector< Real > &f, const std::vector< Real > &pm, const std::vector< bool > &active)
Checks Kuhn-Tucker conditions, and alters "active" if appropriate.
Definition: ComputeMultiPlasticityStress.C:1288
ComputeMultiPlasticityStress::lineSearch
virtual bool lineSearch(Real &nr_res2, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, const RankFourTensor &E_inv, RankTwoTensor &delta_dp, const RankTwoTensor &dstress, const std::vector< Real > &dpm, const std::vector< Real > &dintnl, std::vector< Real > &f, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, bool &linesearch_needed)
Performs a line search.
Definition: ComputeMultiPlasticityStress.C:1391
MultiPlasticityLinearSystem::nrStep
virtual void nrStep(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const RankTwoTensor &delta_dp, RankTwoTensor &dstress, std::vector< Real > &dpm, std::vector< Real > &dintnl, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
Performs one Newton-Raphson step.
Definition: MultiPlasticityLinearSystem.C:615
ComputeMultiPlasticityStress::rot
RankTwoTensor rot(const RankTwoTensor &tens)
Definition: ComputeMultiPlasticityStress.C:312
ComputeMultiPlasticityStress::canIncrementDumb
bool canIncrementDumb(int dumb_iteration)
Definition: ComputeMultiPlasticityStress.C:1567
ComputeStressBase::_stress
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
Definition: ComputeStressBase.h:50
ComputeMultiPlasticityStress::_tangent_operator_type
enum ComputeMultiPlasticityStress::TangentOperatorEnum _tangent_operator_type
MultiPlasticityRawComponentAssembler::dhardPotential_dintnl
virtual void dhardPotential_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &dh_dintnl)
The derivative of the active hardening potentials with respect to the active internal parameters.
Definition: MultiPlasticityRawComponentAssembler.C:317
MultiPlasticityRawComponentAssembler::_model_given_surface
std::vector< unsigned int > _model_given_surface
given a surface number, this returns the model number
Definition: MultiPlasticityRawComponentAssembler.h:295
ComputeMultiPlasticityStress::returnMap
virtual bool returnMap(const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &f, unsigned int &iter, bool can_revert_to_dumb, bool &linesearch_needed, bool &ld_encountered, bool &constraints_added, bool final_step, RankFourTensor &consistent_tangent_operator, std::vector< Real > &cumulative_pm)
Implements the return map.
Definition: ComputeMultiPlasticityStress.C:615
MultiPlasticityRawComponentAssembler::flowPotential
virtual void flowPotential(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &r)
The active flow potential(s) - one for each yield function.
Definition: MultiPlasticityRawComponentAssembler.C:180
ComputeStressBase::_extra_stress
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.
Definition: ComputeStressBase.h:55
ComputeMultiPlasticityStress::incrementDumb
virtual void incrementDumb(int &dumb_iteration, const std::vector< unsigned int > &dumb_order, std::vector< bool > &act)
Increments "dumb_iteration" by 1, and sets "act" appropriately (act[alpha] = true iff alpha_th bit of...
Definition: ComputeMultiPlasticityStress.C:1555
ComputeStressBase::_Jacobian_mult
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
Definition: ComputeStressBase.h:61
ComputeMultiPlasticityStress::DeactivationSchemeEnum
DeactivationSchemeEnum
Definition: ComputeMultiPlasticityStress.h:79
ComputeMultiPlasticityStress::_n_input
RealVectorValue _n_input
the supplied transverse direction vector
Definition: ComputeMultiPlasticityStress.h:94
ComputeMultiPlasticityStress::_max_stepsize_for_dumb
Real _max_stepsize_for_dumb
"dumb" deactivation will only be used if the stepsize falls below this quantity
Definition: ComputeMultiPlasticityStress.h:48
ComputeMultiPlasticityStress::_couple_stress_old
const MaterialProperty< RankTwoTensor > * _couple_stress_old
the old value of Cosserat couple-stress
Definition: ComputeMultiPlasticityStress.h:168
MultiPlasticityDebugger::_fspb_debug_intnl
std::vector< Real > _fspb_debug_intnl
Debug the Jacobian entires at these internal parameters.
Definition: MultiPlasticityDebugger.h:71
MultiPlasticityRawComponentAssembler::_specialIC
MooseEnum _specialIC
Allows initial set of active constraints to be chosen optimally.
Definition: MultiPlasticityRawComponentAssembler.h:72
ComputeMultiPlasticityStress::_epp_tol
Real _epp_tol
Tolerance on the plastic strain increment ("direction") constraint.
Definition: ComputeMultiPlasticityStress.h:62
ComputeMultiPlasticityStress::_intnl_old
const MaterialProperty< std::vector< Real > > & _intnl_old
old values of internal parameters
Definition: ComputeMultiPlasticityStress.h:117
ComputeMultiPlasticityStress::checkAdmissible
virtual bool checkAdmissible(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &all_f)
Checks whether the yield functions are in the admissible region.
Definition: ComputeMultiPlasticityStress.C:1271
MultiPlasticityDebugger::fdJacobian
virtual void fdJacobian(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, const RankFourTensor &E_inv, bool eliminate_ld, std::vector< std::vector< Real >> &jac)
The Jacobian calculated using finite differences.
Definition: MultiPlasticityDebugger.C:221
ComputeMultiPlasticityStress::_min_stepsize
Real _min_stepsize
Minimum fraction of applied strain that may be applied during adaptive stepsizing.
Definition: ComputeMultiPlasticityStress.h:45
MultiPlasticityLinearSystem::calculateConstraints
virtual void calculateConstraints(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &f, std::vector< RankTwoTensor > &r, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active)
The constraints.
Definition: MultiPlasticityLinearSystem.C:228
ComputeMultiPlasticityStress::residual2
virtual Real residual2(const std::vector< Real > &pm, const std::vector< Real > &f, const RankTwoTensor &epp, const std::vector< Real > &ic, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld)
The residual-squared.
Definition: ComputeMultiPlasticityStress.C:1354
ComputeMultiPlasticityStress::dumb
Definition: ComputeMultiPlasticityStress.h:83
MultiPlasticityRawComponentAssembler::_num_surfaces
unsigned int _num_surfaces
Number of surfaces within the plastic models.
Definition: MultiPlasticityRawComponentAssembler.h:66
ComputeMultiPlasticityStress::_cumulative_pm
std::vector< Real > _cumulative_pm
the sum of the plastic multipliers over all the sub-steps.
Definition: ComputeMultiPlasticityStress.h:71
ComputeMultiPlasticityStress::_plastic_strain
MaterialProperty< RankTwoTensor > & _plastic_strain
plastic strain
Definition: ComputeMultiPlasticityStress.h:108
MultiPlasticityLinearSystem::eliminateLinearDependence
virtual void eliminateLinearDependence(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &f, const std::vector< RankTwoTensor > &r, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
Performs a number of singular-value decompositions to check for linear-dependence of the active direc...
Definition: MultiPlasticityLinearSystem.C:107
MultiPlasticityRawComponentAssembler::hardPotential
virtual void hardPotential(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &h)
The active hardening potentials (one for each internal parameter and for each yield function) by assu...
Definition: MultiPlasticityRawComponentAssembler.C:262
MultiPlasticityDebugger::_fspb_debug
MooseEnum _fspb_debug
none - don't do any debugging crash - currently inactive jacobian - check the jacobian entries jacobi...
Definition: MultiPlasticityDebugger.h:62
ComputeMultiPlasticityStress::_my_curvature
RankTwoTensor _my_curvature
Curvature that can be rotated by this class, and split into multiple increments (ie,...
Definition: ComputeMultiPlasticityStress.h:183
MultiPlasticityRawComponentAssembler::returnMapAll
bool returnMapAll(const RankTwoTensor &trial_stress, const std::vector< Real > &intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &stress, std::vector< Real > &intnl, std::vector< Real > &pm, std::vector< Real > &cumulative_pm, RankTwoTensor &delta_dp, std::vector< Real > &yf, unsigned &num_successful_plastic_returns, unsigned &custom_model)
Performs a returnMap for each plastic model using their inbuilt returnMap functions.
Definition: MultiPlasticityRawComponentAssembler.C:599
ComputeMultiPlasticityStress::returnMap_function
Definition: ComputeMultiPlasticityStress.h:448
ComputeMultiPlasticityStress::_curvature
const MaterialProperty< RankTwoTensor > * _curvature
The Cosserat curvature strain.
Definition: ComputeMultiPlasticityStress.h:159
ComputeMultiPlasticityStress::consistentTangentOperator
RankFourTensor consistentTangentOperator(const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &pm_this_step, const std::vector< Real > &cumulative_pm)
Computes the consistent tangent operator (another name for the jacobian = d(stress_rate)/d(strain_rat...
Definition: ComputeMultiPlasticityStress.C:1585
ComputeMultiPlasticityStress::quickStep
virtual bool quickStep(const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, std::vector< Real > &cumulative_pm, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &yf, unsigned int &iterations, RankFourTensor &consistent_tangent_operator, const quickStep_called_from_t called_from, bool final_step)
Attempts to find an admissible (stress, intnl) by using the customized return-map algorithms defined ...
Definition: ComputeMultiPlasticityStress.C:371
MultiPlasticityDebugger::fddyieldFunction_dstress
void fddyieldFunction_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &df_dstress)
The finite-difference derivative of yield function(s) with respect to stress.
Definition: MultiPlasticityDebugger.C:519
ComputeStressBase::computeQpStress
virtual void computeQpStress()=0
Compute the stress and store it in the _stress material property for the current quadrature point.
ComputeMultiPlasticityStress::_deactivation_scheme
enum ComputeMultiPlasticityStress::DeactivationSchemeEnum _deactivation_scheme
MultiPlasticityLinearSystem::calculateJacobian
virtual void calculateJacobian(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, std::vector< std::vector< Real >> &jac)
d(rhs)/d(dof)
Definition: MultiPlasticityLinearSystem.C:367
ComputeMultiPlasticityStress::_elasticity_tensor
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
Definition: ComputeMultiPlasticityStress.h:105
ComputeMultiPlasticityStress::_plastic_strain_old
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
Old value of plastic strain.
Definition: ComputeMultiPlasticityStress.h:111
ComputeMultiPlasticityStress::plasticStep
virtual bool plasticStep(const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &yf, unsigned int &iterations, bool &linesearch_needed, bool &ld_encountered, bool &constraints_added, RankFourTensor &consistent_tangent_operator)
performs a plastic step
Definition: ComputeMultiPlasticityStress.C:468
ComputeMultiPlasticityStress::_couple_stress
MaterialProperty< RankTwoTensor > * _couple_stress
the Cosserat couple-stress
Definition: ComputeMultiPlasticityStress.h:165
ComputeMultiPlasticityStress::_rot
RealTensorValue _rot
rotation matrix that takes _n to (0, 0, 1)
Definition: ComputeMultiPlasticityStress.h:97
MultiPlasticityLinearSystem::calculateRHS
virtual void calculateRHS(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &rhs, const std::vector< bool > &active, bool eliminate_ld, std::vector< bool > &deactivated_due_to_ld)
Calculate the RHS which is rhs = -(epp(0,0), epp(1,0), epp(1,1), epp(2,0), epp(2,1),...
Definition: MultiPlasticityLinearSystem.C:288
MultiPlasticityDebugger::fddflowPotential_dintnl
virtual void fddflowPotential_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &dr_dintnl)
The finite-difference derivative of the flow potentials with respect to internal parameters.
Definition: MultiPlasticityDebugger.C:607
ComputeStressBase::validParams
static InputParameters validParams()
Definition: ComputeStressBase.C:17
MultiPlasticityDebugger::_fspb_debug_stress
RankTwoTensor _fspb_debug_stress
Debug the Jacobian entries at this stress.
Definition: MultiPlasticityDebugger.h:65
ComputeMultiPlasticityStress::buildDumbOrder
void buildDumbOrder(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< unsigned int > &dumb_order)
Builds the order which "dumb" activation will take.
Definition: ComputeMultiPlasticityStress.C:1524
ComputeMultiPlasticityStress::_Jacobian_mult_couple
MaterialProperty< RankFourTensor > * _Jacobian_mult_couple
derivative of couple-stress w.r.t. curvature
Definition: ComputeMultiPlasticityStress.h:171
ComputeMultiPlasticityStress::numberActive
virtual unsigned int numberActive(const std::vector< bool > &active)
counts the number of active constraints
Definition: ComputeMultiPlasticityStress.C:1261
MultiPlasticityDebugger::dof_included
bool dof_included(unsigned int dof, const std::vector< bool > &deactivated_due_to_ld)
Definition: MultiPlasticityDebugger.C:349
MultiPlasticityRawComponentAssembler::_f
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
Definition: MultiPlasticityRawComponentAssembler.h:75
ComputeMultiPlasticityStress::_cosserat
bool _cosserat
whether Cosserat mechanics should be used
Definition: ComputeMultiPlasticityStress.h:156
MultiPlasticityDebugger::_fspb_debug_stress_change
Real _fspb_debug_stress_change
Debug finite-differencing parameter for the stress.
Definition: MultiPlasticityDebugger.h:74
ComputeMultiPlasticityStress::_total_strain_old
const MaterialProperty< RankTwoTensor > & _total_strain_old
Old value of total strain (coming from ComputeIncrementalSmallStrain, for example)
Definition: ComputeMultiPlasticityStress.h:144
ComputeMultiPlasticityStress::_elastic_flexural_rigidity_tensor
const MaterialProperty< RankFourTensor > * _elastic_flexural_rigidity_tensor
The Cosserat elastic flexural rigidity tensor.
Definition: ComputeMultiPlasticityStress.h:162
ComputeMultiPlasticityStress::_perform_finite_strain_rotations
bool _perform_finite_strain_rotations
whether to perform the rotations necessary in finite-strain simulations
Definition: ComputeMultiPlasticityStress.h:100
MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl
virtual void dyieldFunction_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &df_dintnl)
The derivative of active yield function(s) with respect to their internal parameters (the user object...
Definition: MultiPlasticityRawComponentAssembler.C:153
ComputeStressBase::initQpStatefulProperties
virtual void initQpStatefulProperties() override
Definition: ComputeStressBase.C:43
ComputeStressBase::_base_name
const std::string _base_name
Base name prepended to all material property names to allow for multi-material systems.
Definition: ComputeStressBase.h:45
MultiPlasticityRawComponentAssembler::buildActiveConstraints
virtual void buildActiveConstraints(const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
Constructs a set of active constraints, given the yield functions, f.
Definition: MultiPlasticityRawComponentAssembler.C:344
MultiPlasticityRawComponentAssembler::dhardPotential_dstress
virtual void dhardPotential_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dh_dstress)
The derivative of the active hardening potentials with respect to stress By assumption in the Userobj...
Definition: MultiPlasticityRawComponentAssembler.C:289
ComputeMultiPlasticityStress::optimized_to_dumb
Definition: ComputeMultiPlasticityStress.h:87
ComputeMultiPlasticityStress::optimized_to_safe
Definition: ComputeMultiPlasticityStress.h:84
MultiPlasticityDebugger::_fspb_debug_intnl_change
std::vector< Real > _fspb_debug_intnl_change
Debug finite-differencing parameters for the internal parameters.
Definition: MultiPlasticityDebugger.h:80
ComputeMultiPlasticityStress::reinstateLinearDependentConstraints
virtual bool reinstateLinearDependentConstraints(std::vector< bool > &deactivated_due_to_ld)
makes all deactivated_due_to_ld false, and if >0 of them were initially true, returns true
Definition: ComputeMultiPlasticityStress.C:1248
ComputeMultiPlasticityStress::_rotation_increment
const MaterialProperty< RankTwoTensor > & _rotation_increment
Rotation increment (coming from ComputeIncrementalSmallStrain, for example)
Definition: ComputeMultiPlasticityStress.h:147
MultiPlasticityDebugger::fddflowPotential_dstress
virtual void fddflowPotential_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankFourTensor > &dr_dstress)
The finite-difference derivative of the flow potential(s) with respect to stress.
Definition: MultiPlasticityDebugger.C:578
ComputeMultiPlasticityStress::singleStep
virtual bool singleStep(Real &nr_res2, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, RankTwoTensor &delta_dp, const RankFourTensor &E_inv, std::vector< Real > &f, RankTwoTensor &epp, std::vector< Real > &ic, std::vector< bool > &active, DeactivationSchemeEnum deactivation_scheme, bool &linesearch_needed, bool &ld_encountered)
Performs a single Newton-Raphson + linesearch step Constraints are deactivated and the step is re-don...
Definition: ComputeMultiPlasticityStress.C:1082
ComputeMultiPlasticityStress::_constraints_added
MaterialProperty< Real > & _constraints_added
Whether constraints were added in during the latest Newton-Raphson process (1 if true,...
Definition: ComputeMultiPlasticityStress.h:132
MultiPlasticityRawComponentAssembler::dyieldFunction_dstress
virtual void dyieldFunction_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &df_dstress)
The derivative of the active yield function(s) with respect to stress.
Definition: MultiPlasticityRawComponentAssembler.C:125
ComputeMultiPlasticityStress::preReturnMap
virtual void preReturnMap()
Definition: ComputeMultiPlasticityStress.C:320
ComputeMultiPlasticityStress::optimized
Definition: ComputeMultiPlasticityStress.h:81
ComputeMultiPlasticityStress::computeQpStress_function
Definition: ComputeMultiPlasticityStress.h:447
ComputeMultiPlasticityStress::safe_to_dumb
Definition: ComputeMultiPlasticityStress.h:85
MultiPlasticityDebugger::MultiPlasticityDebugger
MultiPlasticityDebugger(const MooseObject *moose_object)
Definition: MultiPlasticityDebugger.C:42
MultiPlasticityDebugger::outputAndCheckDebugParameters
void outputAndCheckDebugParameters()
Outputs the debug parameters: _fspb_debug_stress, _fspd_debug_pm, etc and checks that they are sized ...
Definition: MultiPlasticityDebugger.C:55
MultiPlasticityRawComponentAssembler::activeModelSurfaces
void activeModelSurfaces(int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces_of_model)
Returns the internal surface number(s) of the active surfaces of the given model This may be of size=...
Definition: MultiPlasticityRawComponentAssembler.C:811
ComputeMultiPlasticityStress::canAddConstraints
bool canAddConstraints(const std::vector< bool > &act, const std::vector< Real > &all_f)
Definition: ComputeMultiPlasticityStress.C:1015
MultiPlasticityDebugger::_fspb_debug_pm
std::vector< Real > _fspb_debug_pm
Debug the Jacobian entires at these plastic multipliers.
Definition: MultiPlasticityDebugger.h:68
RankTwoTensor
RankTwoTensorTempl< Real > RankTwoTensor
Definition: ACGrGrElasticDrivingForce.h:17
ComputeStressBase::_elastic_strain
MaterialProperty< RankTwoTensor > & _elastic_strain
Elastic strain material property.
Definition: ComputeStressBase.h:52
MultiPlasticityDebugger::checkDerivatives
void checkDerivatives()
Checks the derivatives, eg dyieldFunction_dstress by using finite difference approximations.
Definition: MultiPlasticityDebugger.C:85
ComputeMultiPlasticityStress::_strain_increment
const MaterialProperty< RankTwoTensor > & _strain_increment
strain increment (coming from ComputeIncrementalSmallStrain, for example)
Definition: ComputeMultiPlasticityStress.h:141
MultiPlasticityRawComponentAssembler::_num_models
unsigned int _num_models
Number of plastic models for this material.
Definition: MultiPlasticityRawComponentAssembler.h:57
MultiPlasticityRawComponentAssembler::anyActiveSurfaces
bool anyActiveSurfaces(int model, const std::vector< bool > &active)
returns true if any internal surfaces of the given model are active according to 'active'
Definition: MultiPlasticityRawComponentAssembler.C:791
ComputeMultiPlasticityStress::_n_old
const MaterialProperty< RealVectorValue > & _n_old
old value of transverse direction
Definition: ComputeMultiPlasticityStress.h:138
ComputeMultiPlasticityStress::_stress_old
const MaterialProperty< RankTwoTensor > & _stress_old
Old value of stress.
Definition: ComputeMultiPlasticityStress.h:150
RankFourTensorTempl
Definition: ACGrGrElasticDrivingForce.h:20
ComputeMultiPlasticityStress::_intnl
MaterialProperty< std::vector< Real > > & _intnl
internal parameters
Definition: ComputeMultiPlasticityStress.h:114
ComputeMultiPlasticityStress::_max_iter
unsigned int _max_iter
Maximum number of Newton-Raphson iterations allowed.
Definition: ComputeMultiPlasticityStress.h:42
MultiPlasticityDebugger::fddyieldFunction_dintnl
void fddyieldFunction_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &df_dintnl)
The finite-difference derivative of yield function(s) with respect to internal parameter(s)
Definition: MultiPlasticityDebugger.C:547
ComputeMultiPlasticityStress::canChangeScheme
bool canChangeScheme(DeactivationSchemeEnum current_deactivation_scheme, bool can_revert_to_dumb)
Definition: ComputeMultiPlasticityStress.C:1025
ComputeMultiPlasticityStress::activeCombinationNumber
unsigned int activeCombinationNumber(const std::vector< bool > &act)
Definition: ComputeMultiPlasticityStress.C:1574
MultiPlasticityDebugger::checkJacobian
void checkJacobian(const RankFourTensor &E_inv, const std::vector< Real > &intnl_old)
Checks the full Jacobian, which is just certain linear combinations of the dyieldFunction_dstress,...
Definition: MultiPlasticityDebugger.C:161
ComputeMultiPlasticityStress::_ignore_failures
bool _ignore_failures
Even if the returnMap fails, return the best values found for stress and internal parameters.
Definition: ComputeMultiPlasticityStress.h:51
MultiPlasticityRawComponentAssembler::buildActiveConstraintsJoint
void buildActiveConstraintsJoint(const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
"Joint" version Constructs a set of active constraints, given the yield functions,...
Definition: MultiPlasticityRawComponentAssembler.C:381
MultiPlasticityDebugger::checkSolution
void checkSolution(const RankFourTensor &E_inv)
Checks that Ax does equal b in the NR procedure.
Definition: MultiPlasticityDebugger.C:367
ComputeMultiPlasticityStress::_my_flexural_rigidity_tensor
RankFourTensor _my_flexural_rigidity_tensor
Flexual rigidity tensor that can be rotated by this class (ie, its not const)
Definition: ComputeMultiPlasticityStress.h:180
ComputeMultiPlasticityStress::TangentOperatorEnum
TangentOperatorEnum
The type of tangent operator to return. tangent operator = d(stress_rate)/d(strain_rate).
Definition: ComputeMultiPlasticityStress.h:54
MultiPlasticityLinearSystem::_min_f_tol
Real _min_f_tol
Minimum value of the _f_tol parameters for the Yield Function User Objects.
Definition: MultiPlasticityLinearSystem.h:136
ComputeMultiPlasticityStress::_my_elasticity_tensor
RankFourTensor _my_elasticity_tensor
Elasticity tensor that can be rotated by this class (ie, its not const)
Definition: ComputeMultiPlasticityStress.h:174
MultiPlasticityRawComponentAssembler::yieldFunction
virtual void yieldFunction(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &f)
The active yield function(s)
Definition: MultiPlasticityRawComponentAssembler.C:98
MultiPlasticityLinearSystem::singularValuesOfR
virtual int singularValuesOfR(const std::vector< RankTwoTensor > &r, std::vector< Real > &s)
Performs a singular-value decomposition of r and returns the singular values.
Definition: MultiPlasticityLinearSystem.C:47
MultiPlasticityDebugger::_fspb_debug_pm_change
std::vector< Real > _fspb_debug_pm_change
Debug finite-differencing parameters for the plastic multipliers.
Definition: MultiPlasticityDebugger.h:77
ComputeMultiPlasticityStress::_iter
MaterialProperty< Real > & _iter
Number of Newton-Raphson iterations used in the return-map.
Definition: ComputeMultiPlasticityStress.h:123
RankTwoTensorTempl< Real >
MultiPlasticityLinearSystem::_svd_tol
Real _svd_tol
Tolerance on the minimum ratio of singular values before flow-directions are deemed linearly dependen...
Definition: MultiPlasticityLinearSystem.h:133
ComputeMultiPlasticityStress::linear
Definition: ComputeMultiPlasticityStress.h:57
ComputeMultiPlasticityStress::_my_strain_increment
RankTwoTensor _my_strain_increment
Strain increment that can be rotated by this class, and split into multiple increments (ie,...
Definition: ComputeMultiPlasticityStress.h:177
ComputeMultiPlasticityStress::_elasticity_tensor_name
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
Definition: ComputeMultiPlasticityStress.h:103
ComputeMultiPlasticityStress::applyKuhnTucker
virtual void applyKuhnTucker(const std::vector< Real > &f, const std::vector< Real > &pm, std::vector< bool > &active)
Checks Kuhn-Tucker conditions, and alters "active" if appropriate.
Definition: ComputeMultiPlasticityStress.C:1313
ComputeMultiPlasticityStress::_elastic_strain_old
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Old value of elastic strain.
Definition: ComputeMultiPlasticityStress.h:153
MultiPlasticityRawComponentAssembler::activeSurfaces
void activeSurfaces(int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces)
Returns the external surface number(s) of the active surfaces of the given model This may be of size=...
Definition: MultiPlasticityRawComponentAssembler.C:800
ComputeMultiPlasticityStress::nonlinear
Definition: ComputeMultiPlasticityStress.h:58
MultiPlasticityRawComponentAssembler::dflowPotential_dstress
virtual void dflowPotential_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankFourTensor > &dr_dstress)
The derivative of the active flow potential(s) with respect to stress.
Definition: MultiPlasticityRawComponentAssembler.C:207
RankFourTensor
RankFourTensorTempl< Real > RankFourTensor
Definition: ACGrGrElasticDrivingForce.h:20
RankTwoScalarTools::L2norm
T L2norm(const RankTwoTensorTempl< T > &r2tensor)
Definition: RankTwoScalarTools.h:98
ComputeMultiPlasticityStress::postReturnMap
virtual void postReturnMap()
Definition: ComputeMultiPlasticityStress.C:339
MultiPlasticityRawComponentAssembler::_surfaces_given_model
std::vector< std::vector< unsigned int > > _surfaces_given_model
_surfaces_given_model[model_number] = vector of surface numbers for this model
Definition: MultiPlasticityRawComponentAssembler.h:69
MultiPlasticityRawComponentAssembler::buildActiveConstraintsRock
void buildActiveConstraintsRock(const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
"Rock" version Constructs a set of active constraints, given the yield functions, f.
Definition: MultiPlasticityRawComponentAssembler.C:422
ComputeMultiPlasticityStress::optimized_to_safe_to_dumb
Definition: ComputeMultiPlasticityStress.h:86
ComputeStressBase::ComputeStressBase
ComputeStressBase(const InputParameters &parameters)
Definition: ComputeStressBase.C:28
ComputeMultiPlasticityStress::_n_supplied
bool _n_supplied
User supplied the transverse direction vector.
Definition: ComputeMultiPlasticityStress.h:91
ComputeMultiPlasticityStress::_yf
MaterialProperty< std::vector< Real > > & _yf
yield functions
Definition: ComputeMultiPlasticityStress.h:120
Generated on Sat Jan 25 2020 12:04:58 for www.mooseframework.org by   doxygen 1.8.16