15 #include "libmesh/utility.h" 18 #include <petscblaslapack.h> 24 params.
addClassDescription(
"Return-map and Jacobian algorithms for plastic models where the " 25 "yield function and flow directions depend on multiple functions of " 30 "max_NR_iterations>0",
31 "Maximum number of Newton-Raphson iterations allowed during the return-map algorithm");
32 params.
addParam<
bool>(
"perform_finite_strain_rotations",
34 "Tensors are correctly rotated " 35 "in finite-strain simulations. " 36 "For optimal performance you can " 37 "set this to 'false' if you are " 38 "only ever using small strains");
41 "Intersections of the yield surfaces will be smoothed by this amount (this " 42 "is measured in units of stress). Often this is related to other physical " 43 "parameters (eg, 0.1*cohesion) but it is important to set this small enough " 44 "so that the individual yield surfaces do not mix together in the smoothing " 45 "process to produce a result where no stress is admissible (for example, " 46 "mixing together tensile and compressive failure envelopes).");
48 "The return-map process will be deemed to have converged if all " 49 "yield functions are within yield_function_tol of zero. If this " 50 "is set very low then precision-loss might be encountered: if the " 51 "code detects precision loss then it also deems the return-map " 52 "process has converged.");
53 MooseEnum tangent_operator(
"elastic nonlinear",
"nonlinear");
56 "In order to help the Newton-Raphson procedure, the applied strain " 57 "increment may be applied in sub-increments of size greater than this " 58 "value. Usually it is better for Moose's nonlinear convergence to " 59 "increase max_NR_iterations rather than decrease this parameter.");
60 params.
addParam<
bool>(
"warn_about_precision_loss",
62 "Output a message to the console every " 63 "time precision-loss is encountered " 64 "during the Newton-Raphson process");
65 params.
addParam<std::vector<Real>>(
"admissible_stress",
66 "A single admissible value of the value of the stress " 67 "parameters for internal parameters = 0. This is used " 68 "to initialize the return-mapping algorithm during the first " 69 "nonlinear iteration. If not given then it is assumed that " 70 "stress parameters = 0 is admissible.");
71 MooseEnum smoother_fcn_enum(
"cos poly1 poly2 poly3",
"cos");
74 "Type of smoother function to use. 'cos' means (-a/pi)cos(pi x/2/a), " 75 "'polyN' means a polynomial of degree 2N+2");
81 const InputParameters & parameters,
unsigned num_sp,
unsigned num_yf,
unsigned num_intnl)
84 _definitely_ok_sp(isParamValid(
"admissible_stress")
85 ? getParam<
std::vector<
Real>>(
"admissible_stress")
86 :
std::vector<
Real>(_num_sp, 0.0)),
87 _Eij(num_sp,
std::vector<
Real>(num_sp)),
89 _Cij(num_sp,
std::vector<
Real>(num_sp)),
91 _num_intnl(num_intnl),
92 _max_nr_its(getParam<unsigned>(
"max_NR_iterations")),
93 _perform_finite_strain_rotations(getParam<bool>(
"perform_finite_strain_rotations")),
94 _smoothing_tol(getParam<
Real>(
"smoothing_tol")),
95 _smoothing_tol2(Utility::
pow<2>(getParam<
Real>(
"smoothing_tol"))),
96 _f_tol(getParam<
Real>(
"yield_function_tol")),
97 _f_tol2(Utility::
pow<2>(getParam<
Real>(
"yield_function_tol"))),
98 _min_step_size(getParam<
Real>(
"min_step_size")),
99 _step_one(declareRestartableData<bool>(
"step_one", true)),
100 _warn_about_precision_loss(getParam<bool>(
"warn_about_precision_loss")),
102 _plastic_strain(declareProperty<
RankTwoTensor>(_base_name +
"plastic_strain")),
103 _plastic_strain_old(getMaterialPropertyOld<
RankTwoTensor>(_base_name +
"plastic_strain")),
104 _intnl(declareProperty<
std::vector<
Real>>(_base_name +
"plastic_internal_parameter")),
106 getMaterialPropertyOld<
std::vector<
Real>>(_base_name +
"plastic_internal_parameter")),
107 _yf(declareProperty<
std::vector<
Real>>(_base_name +
"plastic_yield_function")),
108 _iter(declareProperty<
Real>(_base_name +
109 "plastic_NR_iterations")),
111 _max_iter_used(declareProperty<
Real>(
112 _base_name +
"max_plastic_NR_iterations")),
114 _max_iter_used_old(getMaterialPropertyOld<
Real>(_base_name +
"max_plastic_NR_iterations")),
116 declareProperty<
Real>(_base_name +
"plastic_linesearch_needed")),
123 _dvar_dtrial(num_sp + 1,
std::vector<
Real>(num_sp, 0.0)),
125 _ok_intnl(num_intnl),
126 _del_stress_params(num_sp),
128 _current_intnl(num_intnl),
129 _smoother_function_type(
133 mooseError(
"MultiParameterPlasticityStressUpdate: admissible_stress parameter must consist of ",
165 bool compute_full_tangent_operator,
194 inelastic_strain_increment.
zero();
235 for (
unsigned i = 0; i <
_num_sp; ++i)
238 Real step_taken = 0.0;
240 Real step_size = 1.0;
241 Real gaE_total = 0.0;
249 bool smoothed_q_calculated =
false;
253 if (1.0 - step_taken < step_size)
255 step_size = 1.0 - step_taken;
258 for (
unsigned i = 0; i <
_num_sp; ++i)
270 unsigned step_iter = 0;
279 smoothed_q_calculated =
false;
288 smoothed_q_calculated =
true;
293 while (res2 >
_f_tol2 && step_iter <
_max_nr_its && nr_failure == 0 && ls_failure == 0)
305 Moose::err <<
"MultiParameterPlasticityStressUpdate: precision-loss in element " 307 for (
unsigned i = 0; i <
_num_sp; ++i)
309 Moose::err <<
" gaE = " << gaE <<
"\n";
328 if (res2 <=
_f_tol2 && step_iter <
_max_nr_its && nr_failure == 0 && ls_failure == 0 &&
334 step_taken += step_size;
345 compute_full_tangent_operator,
347 if (static_cast<Real>(step_iter) >
_iter[
_qp])
361 errorHandler(
"MultiParameterPlasticityStressUpdate: Minimum step-size violated");
367 if (!smoothed_q_calculated)
378 inelastic_strain_increment);
380 strain_increment = strain_increment - inelastic_strain_increment;
392 compute_full_tangent_operator,
399 const std::vector<Real> & intnl)
const 428 std::size_t res_f = 0;
430 for (std::size_t
a = 1;
a < all_q.size(); ++
a)
444 const Real f_diff = all_q[res_f].f - all_q[
a].f;
450 for (
unsigned i = 0; i <
_num_sp; ++i)
453 all_q[res_f].d2g[i][
j] =
454 0.5 * (all_q[res_f].d2g[i][
j] + all_q[
a].d2g[i][
j]) +
455 dsm * (all_q[res_f].df[
j] - all_q[
a].df[
j]) * (all_q[res_f].dg[i] - all_q[
a].dg[i]) +
456 sm * (all_q[res_f].d2g[i][
j] - all_q[
a].d2g[i][
j]);
458 all_q[res_f].d2g_di[i][
j] = 0.5 * (all_q[res_f].d2g_di[i][
j] + all_q[
a].d2g_di[i][
j]) +
459 dsm * (all_q[res_f].df_di[
j] - all_q[
a].df_di[
j]) *
460 (all_q[res_f].dg[i] - all_q[
a].dg[i]) +
461 sm * (all_q[res_f].d2g_di[i][
j] - all_q[
a].d2g_di[i][
j]);
463 for (
unsigned i = 0; i <
_num_sp; ++i)
465 all_q[res_f].df[i] = 0.5 * (all_q[res_f].df[i] + all_q[
a].df[i]) +
466 sm * (all_q[res_f].df[i] - all_q[
a].df[i]);
469 all_q[res_f].dg[i] = 0.5 * (all_q[res_f].dg[i] + all_q[
a].dg[i]) +
470 sm * (all_q[res_f].dg[i] - all_q[
a].dg[i]);
473 all_q[res_f].df_di[i] = 0.5 * (all_q[res_f].df_di[i] + all_q[
a].df_di[i]) +
474 sm * (all_q[res_f].df_di[i] - all_q[
a].df_di[i]);
475 all_q[res_f].f = 0.5 * (all_q[res_f].f + all_q[
a].f +
_smoothing_tol) + ism;
552 std::vector<Real> & stress_params,
554 const std::vector<Real> & trial_stress_params,
556 const std::vector<Real> & intnl_ok,
557 std::vector<Real> & intnl,
558 std::vector<Real> & rhs,
559 Real & linesearch_needed)
const 561 const Real res2_old = res2;
562 const std::vector<Real> sp_params_old = stress_params;
563 const Real gaE_old = gaE;
564 const std::vector<Real> delta_nr_params = rhs;
567 const Real lam_min = 1E-10;
568 const Real slope = -2.0 * res2_old;
578 for (
unsigned i = 0; i <
_num_sp; ++i)
579 stress_params[i] = sp_params_old[i] + lam * delta_nr_params[i];
580 gaE = gaE_old + lam * delta_nr_params[
_num_sp];
588 calculateRHS(trial_stress_params, stress_params, gaE, smoothed_q, rhs);
592 if (res2 < res2_old + 1E-4 * lam * slope)
594 else if (lam < lam_min)
599 tmp_lam = -0.5 * slope / (res2 - res2_old - slope);
604 const Real rhs1 = res2 - res2_old - lam * slope;
605 const Real rhs2 = f2 - res2_old - lam2 * slope;
606 const Real a = (rhs1 / Utility::pow<2>(lam) - rhs2 / Utility::pow<2>(lam2)) / (lam - lam2);
608 (-lam2 * rhs1 / Utility::pow<2>(lam) + lam * rhs2 / Utility::pow<2>(lam2)) / (lam - lam2);
610 tmp_lam = -slope / (2.0 *
b);
613 const Real disc = Utility::pow<2>(
b) - 3.0 *
a * slope;
621 if (tmp_lam > 0.5 * lam)
626 lam = std::max(tmp_lam, 0.1 * lam);
630 linesearch_needed = 1.0;
636 const std::vector<Real> & trial_stress_params,
637 const std::vector<Real> & stress_params,
638 const std::vector<Real> & intnl,
640 std::vector<Real> & rhs)
const 646 dnRHSdVar(smoothed_q, dintnl, stress_params, gaE, jac);
649 const PetscBLASInt nrhs = 1;
650 std::vector<PetscBLASInt> ipiv(
_num_sp + 1);
652 const PetscBLASInt gesv_num_rhs =
_num_sp + 1;
654 &gesv_num_rhs, &nrhs, &jac[0], &gesv_num_rhs, &ipiv[0], &rhs[0], &gesv_num_rhs, &
info);
677 const std::vector<Real> & ,
679 const std::vector<Real> & ,
680 const std::vector<Real> & ,
687 const std::vector<Real> & intnl)
const 689 std::vector<Real> yfs(
_num_yf);
698 for (std::size_t i = 1; i < yfs.size(); ++i)
712 const std::vector<Real> & intnl_old,
713 std::vector<Real> & stress_params,
715 std::vector<Real> & intnl)
const 718 std::copy(trial_stress_params.begin(), trial_stress_params.end(), stress_params.begin());
719 std::copy(intnl_old.begin(), intnl_old.end(), intnl.begin());
725 const std::vector<Real> & ,
727 const std::vector<Real> & ,
731 bool compute_full_tangent_operator,
732 const std::vector<std::vector<Real>> & dvar_dtrial,
736 if (!compute_full_tangent_operator)
751 s -= dsp[
b] *
_Cij[
b][
c] * dvar_dtrial[
c][
a];
753 cto -= s.outerProduct(t);
759 for (
unsigned i = 0; i <
_num_sp; ++i)
760 Tijab += smoothed_q.
dg[i] * d2sp[i];
772 mooseWarning(
"MultiParameterPlasticityStressUpdate: Cannot invert 1+T in consistent tangent " 773 "operator computation at quadpoint ",
794 for (
unsigned i = 0; i <
_num_sp; ++i)
795 inelastic_strain_increment += smoothed_q.
dg[i] * dsp_dstress[i];
796 inelastic_strain_increment *= gaE /
_En;
803 for (
const auto & r : rhs)
810 const std::vector<Real> & stress_params,
813 std::vector<Real> & rhs)
const 816 for (
unsigned i = 0; i <
_num_sp; ++i)
818 rhs[i] = stress_params[i] - trial_stress_params[i];
820 rhs[i] += ga *
_Eij[i][
j] * smoothed_q.
dg[
j];
827 const std::vector<std::vector<Real>> & dintnl,
828 const std::vector<Real> & ,
830 std::vector<double> & jac)
const 832 for (
auto & jac_entry : jac)
838 for (
unsigned var = 0; var <
_num_sp; ++var)
840 for (
unsigned rhs = 0; rhs <
_num_sp; ++rhs)
846 jac[ind] -= ga *
_Eij[rhs][
j] * smoothed_q.
d2g[
j][var];
848 jac[ind] -= ga *
_Eij[rhs][
j] * smoothed_q.
d2g_di[
j][
k] * dintnl[
k][var];
853 jac[ind] -= smoothed_q.
df[var];
855 jac[ind] -= smoothed_q.
df_di[
k] * dintnl[
k][var];
860 for (
unsigned rhs = 0; rhs <
_num_sp; ++rhs)
863 jac[ind] -= (1.0 /
_En) *
_Eij[rhs][
j] * smoothed_q.
dg[
j];
872 const std::vector<Real> & trial_stress_params,
873 const std::vector<Real> & stress_params,
875 const std::vector<Real> & intnl,
878 bool compute_full_tangent_operator,
879 std::vector<std::vector<Real>> & dvar_dtrial)
const 884 if (!compute_full_tangent_operator)
891 dvar_dtrial[
v][
v] += step_size;
920 rhs_cto[ind] -= smoothed_q.
df_di[
k] * dintnl[
k][
a];
926 dnRHSdVar(smoothed_q, dintnl, stress_params, gaE, jac);
928 std::vector<PetscBLASInt> ipiv(
_num_sp + 1);
930 const PetscBLASInt gesv_num_rhs =
_num_sp + 1;
931 const PetscBLASInt gesv_num_pq =
_num_sp;
932 LAPACKgesv_(&gesv_num_rhs,
941 errorHandler(
"MultiParameterPlasticityStressUpdate: PETSC LAPACK gsev routine returned with " 946 std::vector<std::vector<Real>> dvarn_dtrialn(
_num_sp + 1, std::vector<Real>(
_num_sp, 0.0));
947 for (
unsigned spt = 0; spt <
_num_sp; ++spt)
951 dvarn_dtrialn[
v][spt] = rhs_cto[ind];
955 dvarn_dtrialn[
_num_sp][spt] = rhs_cto[ind];
959 const std::vector<std::vector<Real>> dvar_dtrial_old = dvar_dtrial;
963 for (
unsigned spt = 0; spt <
_num_sp; ++spt)
965 dvar_dtrial[
v][spt] = step_size * dvarn_dtrialn[
v][spt];
967 dvar_dtrial[
v][spt] += dvarn_dtrialn[
v][
a] * dvar_dtrial_old[
a][spt];
972 for (
unsigned spt = 0; spt <
_num_sp; ++spt)
974 dvar_dtrial[
v][spt] += step_size * dvarn_dtrialn[
v][spt];
976 dvar_dtrial[
v][spt] += dvarn_dtrialn[
v][
a] * dvar_dtrial_old[
a][spt];
982 const std::vector<Real> & stress_params,
987 for (
unsigned i = 0; i <
_num_sp; ++i)
const Real _smoothing_tol2
Square of the smoothing tolerance.
std::vector< Real > df_di
FEProblemBase & _fe_problem
MaterialProperty< std::vector< Real > > & _yf
yield functions
std::vector< Real > _trial_sp
"Trial" value of stress_params that initializes the return-map process This is derived from stress = ...
std::vector< std::vector< Real > > d2g
SmootherFunctionType
The type of smoother function.
virtual void initQpStatefulProperties() override
void dnRHSdVar(const yieldAndFlow &smoothed_q, const std::vector< std::vector< Real >> &dintnl, const std::vector< Real > &stress_params, Real gaE, std::vector< double > &jac) const
Derivative of -RHS with respect to the stress_params and gaE, placed into an array ready for solving ...
MultiParameterPlasticityStressUpdate(const InputParameters ¶meters, unsigned num_sp, unsigned num_yf, unsigned num_intnl)
virtual void setEffectiveElasticity(const RankFourTensor &Eijkl)=0
Sets _Eij and _En and _Cij.
const bool _warn_about_precision_loss
Output a warning message if precision loss is encountered during the return-map process.
virtual void errorHandler(const std::string &message) const
Performs any necessary cleaning-up, then throw MooseException(message)
virtual void initializeVarsV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &intnl_old, std::vector< Real > &stress_params, Real &gaE, std::vector< Real > &intnl) const
Sets (stress_params, intnl) at "good guesses" of the solution to the Return-Map algorithm.
int lineSearch(Real &res2, std::vector< Real > &stress_params, Real &gaE, const std::vector< Real > &trial_stress_params, yieldAndFlow &smoothed_q, const std::vector< Real > &intnl_ok, std::vector< Real > &intnl, std::vector< Real > &rhs, Real &linesearch_needed) const
Performs a line-search to find stress_params and gaE Upon entry:
const unsigned _num_intnl
Number of internal parameters.
Real _En
normalising factor
Struct designed to hold info about a single yield function and its derivatives, as well as the flow d...
std::vector< std::vector< Real > > _Cij
_Cij[i, j] * _Eij[j, k] = 1 iff j == k
virtual void consistentTangentOperatorV(const RankTwoTensor &stress_trial, const std::vector< Real > &trial_stress_params, const RankTwoTensor &stress, const std::vector< Real > &stress_params, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, bool compute_full_tangent_operator, const std::vector< std::vector< Real >> &dvar_dtrial, RankFourTensor &cto)
Calculates the consistent tangent operator.
static InputParameters validParams()
virtual void setIntnlValuesV(const std::vector< Real > &trial_stress_params, const std::vector< Real > ¤t_stress_params, const std::vector< Real > &intnl_old, std::vector< Real > &intnl) const =0
Sets the internal parameters based on the trial values of stress_params, their current values...
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
std::vector< Real > _ok_intnl
The state (ok_sp, ok_intnl) is known to be admissible.
const Real _f_tol2
Square of the yield-function tolerance.
const Real _min_step_size
In order to help the Newton-Raphson procedure, the applied strain increment may be applied in sub-inc...
Real smoother(Real f_diff) const
Derivative of ismoother.
int nrStep(const yieldAndFlow &smoothed_q, const std::vector< Real > &trial_stress_params, const std::vector< Real > &stress_params, const std::vector< Real > &intnl, Real gaE, std::vector< Real > &rhs) const
Performs a Newton-Raphson step to attempt to zero rhs Upon return, rhs will contain the solution...
void mooseWarning(Args &&... args) const
void calculateRHS(const std::vector< Real > &trial_stress_params, const std::vector< Real > &stress_params, Real gaE, const yieldAndFlow &smoothed_q, std::vector< Real > &rhs) const
Calculates the RHS in the following 0 = rhs[0] = S[0] - S[0]^trial + ga * E[0, j] * dg/dS[j] 0 = rhs[...
virtual void updateState(RankTwoTensor &strain_increment, RankTwoTensor &inelastic_strain_increment, const RankTwoTensor &rotation_increment, RankTwoTensor &stress_new, const RankTwoTensor &stress_old, const RankFourTensor &elasticity_tensor, const RankTwoTensor &elastic_strain_old, bool compute_full_tangent_operator, RankFourTensor &tangent_operator) override
MaterialProperty< std::vector< Real > > & _intnl
internal parameters
yieldAndFlow smoothAllQuantities(const std::vector< Real > &stress_params, const std::vector< Real > &intnl) const
Calculates all yield functions and derivatives, and then performs the smoothing scheme.
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Real elasticity_tensor(unsigned int i, unsigned int j, unsigned int k, unsigned int l)
std::vector< std::vector< Real > > _Eij
E[i, j] in the system of equations to be solved.
const Real _f_tol
The yield-function tolerance.
const unsigned _num_yf
Number of yield functions.
virtual void computeAllQV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< yieldAndFlow > &all_q) const =0
Completely fills all_q with correct values.
virtual void setInelasticStrainIncrementAfterReturn(const RankTwoTensor &stress_trial, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &elasticity_tensor, const RankTwoTensor &returned_stress, RankTwoTensor &inelastic_strain_increment) const
Sets inelastic strain increment from the returned configuration This is called after the return-map p...
const std::vector< Real > _definitely_ok_sp
An admissible value of stress_params at the initial time.
const MaterialProperty< std::vector< Real > > & _intnl_old
old values of internal parameters
RankTwoTensor _stress_trial
"Trial" value of stress that is set at the beginning of the return-map process.
Real deriv(unsigned n, unsigned alpha, unsigned beta, Real x)
virtual void propagateQpStatefulProperties() override
If updateState is not called during a timestep, this will be.
Real f(Real x)
Test function for Brents method.
virtual void computeStressParams(const RankTwoTensor &stress, std::vector< Real > &stress_params) const =0
Computes stress_params, given stress.
const MaterialProperty< Real > & _max_iter_used_old
Old value of maximum number of Newton-Raphson iterations used in the return-map during the course of ...
Real calculateRes2(const std::vector< Real > &rhs) const
Calculates the residual-squared for the Newton-Raphson + line-search.
virtual void setIntnlDerivativesV(const std::vector< Real > &trial_stress_params, const std::vector< Real > ¤t_stress_params, const std::vector< Real > &intnl, std::vector< std::vector< Real >> &dintnl) const =0
Sets the derivatives of internal parameters, based on the trial values of stress_params, their current values, and the current values of the internal parameters.
RankFourTensorTempl< T > transposeMajor() const
virtual void yieldFunctionValuesV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< Real > &yf) const =0
Computes the values of the yield functions, given stress_params and intnl parameters.
StressUpdateBase is a material that is not called by MOOSE because of the compute=false flag set in t...
std::string stringify(const T &t)
Real dsmoother(Real f_diff) const
Derivative of smoother.
enum MultiParameterPlasticityStressUpdate::SmootherFunctionType _smoother_function_type
MaterialProperty< Real > & _iter
Number of Newton-Raphson iterations used in the return-map.
MaterialProperty< RankTwoTensor > & _plastic_strain
plastic strain
virtual std::vector< RankFourTensor > d2stress_param_dstress(const RankTwoTensor &stress) const =0
d2(stress_param[i])/d(stress)/d(stress) at given stress
const Real _smoothing_tol
Smoothing tolerance: edges of the yield surface get smoothed by this amount.
std::vector< Real > _current_intnl
The current values of the internal params during the Newton-Raphson.
virtual void finalizeReturnProcess(const RankTwoTensor &rotation_increment)
Derived classes may use this to perform calculations after the return-map process has completed succe...
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string v
static const std::string message
MaterialProperty< Real > & _linesearch_needed
Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise) ...
std::vector< std::vector< Real > > _dvar_dtrial
d({stress_param[i], gaE})/d(trial_stress_param[j])
std::vector< Real > _current_sp
The current values of the stress params during the Newton-Raphson.
bool _step_one
handles case of initial_stress that is inadmissible being supplied
void mooseError(Args &&... args) const
void dVardTrial(bool elastic_only, const std::vector< Real > &trial_stress_params, const std::vector< Real > &stress_params, Real gaE, const std::vector< Real > &intnl, const yieldAndFlow &smoothed_q, Real step_size, bool compute_full_tangent_operator, std::vector< std::vector< Real >> &dvar_dtrial) const
Calculates derivatives of the stress_params and gaE with repect to the trial values of the stress_par...
static InputParameters validParams()
virtual std::vector< RankTwoTensor > dstress_param_dstress(const RankTwoTensor &stress) const =0
d(stress_param[i])/d(stress) at given stress
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
virtual void setStressAfterReturnV(const RankTwoTensor &stress_trial, const std::vector< Real > &stress_params, Real gaE, const std::vector< Real > &intnl, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, RankTwoTensor &stress) const =0
Sets stress from the admissible parameters.
std::vector< std::vector< Real > > d2g_di
virtual void preReturnMapV(const std::vector< Real > &trial_stress_params, const RankTwoTensor &stress_trial, const std::vector< Real > &intnl_old, const std::vector< Real > &yf, const RankFourTensor &Eijkl)
Derived classes may employ this function to record stuff or do other computations prior to the return...
MaterialProperty< Real > & _max_iter_used
Maximum number of Newton-Raphson iterations used in the return-map during the course of the entire si...
bool precisionLoss(const std::vector< Real > &solution, const std::vector< Real > &stress_params, Real gaE) const
Check whether precision loss has occurred.
const bool & currentlyComputingJacobian() const
std::vector< Real > _del_stress_params
_del_stress_params = trial_stress_params - ok_sp This is fixed at the beginning of the return-map pro...
const unsigned _max_nr_its
Maximum number of Newton-Raphson iterations allowed in the return-map process.
virtual void initializeReturnProcess()
Derived classes may use this to perform calculations before any return-map process is performed...
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
Old value of plastic strain.
RankFourTensorTempl< T > invSymm() const
static const std::string k
Real yieldF(const std::vector< Real > &stress_params, const std::vector< Real > &intnl) const
Computes the smoothed yield function.
std::vector< Real > _rhs
0 = rhs[0] = S[0] - S[0]^trial + ga * E[0, i] * dg/dS[i] 0 = rhs[1] = S[1] - S[1]^trial + ga * E[1...
const Elem & get(const ElemType type_in)
const Elem *const & _current_elem
std::vector< Real > _ok_sp
The state (ok_sp, ok_intnl) is known to be admissible, so ok_sp are stress_params that are "OK"...
const unsigned _num_sp
Number of stress parameters.
Real ismoother(Real f_diff) const
Smooths yield functions.