www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
FiniteStrainHyperElasticViscoPlastic Class Reference

This class solves the viscoplastic flow rate equations in the total form Involves 4 different types of user objects that calculates: Internal variable rates - functions of internal variables and flow rates Internal variables - functions of internal variables Strengths - functions of internal variables Flow rates - functions of strengths and PK2 stress Flow directions - functions of strengths and PK2 stress The associated derivatives from user objects are assembled and the system is solved using NR. More...

#include <FiniteStrainHyperElasticViscoPlastic.h>

Inheritance diagram for FiniteStrainHyperElasticViscoPlastic:
[legend]

Public Member Functions

 FiniteStrainHyperElasticViscoPlastic (const InputParameters &parameters)
 

Protected Member Functions

virtual void initUOVariables ()
 This function initializes the properties, stateful properties and user objects The properties and stateful properties associated with user objects are only initialized here The properties have the same name as the user object name. More...
 
void initNumUserObjects (const std::vector< UserObjectName > &, unsigned int &)
 This function calculates the number of each user object type. More...
 
template<typename T >
void initProp (const std::vector< UserObjectName > &, unsigned int, std::vector< MaterialProperty< T > *> &)
 This function initializes properties for each user object. More...
 
template<typename T >
void initPropOld (const std::vector< UserObjectName > &, unsigned int, std::vector< const MaterialProperty< T > *> &)
 This function initializes old for stateful properties associated with user object Only user objects that update internal variables have an associated old property. More...
 
template<typename T >
void initUserObjects (const std::vector< UserObjectName > &, unsigned int, std::vector< const T *> &)
 This function initializes user objects. More...
 
virtual void initJacobianVariables ()
 This function initialize variables required for Jacobian calculation. More...
 
virtual void initQpStatefulProperties ()
 Initializes state. More...
 
virtual void computeQpStress ()
 This function computes the Cauchy stress. More...
 
virtual void computeQpJacobian ()
 This function computes the Jacobian. More...
 
virtual void saveOldState ()
 This function saves the old stateful properties that is modified during sub stepping. More...
 
virtual void preSolveQp ()
 Sets state for solve. More...
 
virtual bool solveQp ()
 Solve state. More...
 
virtual void postSolveQp ()
 Update state for output (Outside substepping) More...
 
virtual void recoverOldState ()
 This function restores the the old stateful properties after a successful solve. More...
 
virtual void preSolveFlowrate ()
 Sets state for solve (Inside substepping) More...
 
virtual bool solveFlowrate ()
 Solve for flow rate and state. More...
 
virtual void postSolveFlowrate ()
 Update state for output (Inside substepping) More...
 
virtual bool computeFlowRateFunction ()
 Calls user objects to compute flow rates. More...
 
virtual bool computeFlowDirection ()
 Calls user objects to compute flow directions. More...
 
virtual void computeElasticRightCauchyGreenTensor ()
 Computes elastic Right Cauchy Green Tensor. More...
 
virtual void computePK2StressAndDerivative ()
 Computes PK2 stress and derivative w.r.t elastic Right Cauchy Green Tensor. More...
 
virtual void computeElasticStrain ()
 Computes elastic Lagrangian strain. More...
 
virtual void computeDeeDce ()
 Computes derivative of elastic strain w.r.t elastic Right Cauchy Green Tensor. More...
 
virtual bool computeFlowRateResidual ()
 Computes flow rate residual vector. More...
 
virtual void computeFlowRateJacobian ()
 Computes flow rate Jacobian matrix. More...
 
virtual void computeElasticPlasticDeformGrad ()
 Computes elastic and plastic deformation gradients. More...
 
virtual Real computeNorm (const std::vector< Real > &)
 Computes norm of residual vector. More...
 
virtual void updateFlowRate ()
 Update flow rate. More...
 
virtual void computeDpk2Dfpinv ()
 Computes derivative of PK2 stress wrt inverse of plastic deformation gradient. More...
 
virtual bool computeIntVarRates ()
 This function call user objects to calculate rate of internal variables. More...
 
virtual bool computeIntVar ()
 This function call user objects to integrate internal variables. More...
 
virtual bool computeStrength ()
 This function call user objects to compute strength. More...
 
virtual void computeIntVarRateDerivatives ()
 This function call user objects to compute dintvar_rate/dintvar and dintvarrate/dflowrate. More...
 
virtual void computeIntVarDerivatives ()
 This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate. More...
 
void computeStrengthDerivatives ()
 This function call user objects to compute dstrength/dintvar. More...
 
virtual void computeQpProperties () override
 

Protected Attributes

Real _resid_abs_tol
 Absolute tolerance for residual convergence check. More...
 
Real _resid_rel_tol
 Relative tolerance for residual convergence check. More...
 
unsigned int _maxiters
 Maximum number of iterations. More...
 
unsigned int _max_substep_iter
 Maximum number of substep iterations. More...
 
std::vector< UserObjectName > _flow_rate_uo_names
 Names of flow rate user objects. More...
 
std::vector< UserObjectName > _strength_uo_names
 Names of strength user objects. More...
 
std::vector< UserObjectName > _int_var_uo_names
 Names of internal variable user objects. More...
 
std::vector< UserObjectName > _int_var_rate_uo_names
 Names of internal variable rate user objects. More...
 
unsigned int _num_flow_rate_uos
 Number of flow rate user objects. More...
 
unsigned int _num_strength_uos
 Number of strength user objects. More...
 
unsigned int _num_int_var_uos
 Number of internal variable user objects. More...
 
unsigned int _num_int_var_rate_uos
 Number of internal variable rate user objects. More...
 
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
 Flow rate user objects. More...
 
std::vector< const HEVPStrengthUOBase * > _strength_uo
 Strength user objects. More...
 
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
 Internal variable user objects. More...
 
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
 Internal variable rate user objects. More...
 
std::string _pk2_prop_name
 
MaterialProperty< RankTwoTensor > & _pk2
 
MaterialProperty< RankTwoTensor > & _fp
 
const MaterialProperty< RankTwoTensor > & _fp_old
 
MaterialProperty< RankTwoTensor > & _ce
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
 
const MaterialProperty< RankTwoTensor > & _rotation_increment
 
std::vector< MaterialProperty< Real > * > _flow_rate_prop
 
std::vector< MaterialProperty< Real > * > _strength_prop
 
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
 
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
 
std::vector< MaterialProperty< Real > * > _int_var_rate_prop
 
std::vector< Real > _int_var_old
 
RankTwoTensor _dfgrd_tmp
 
RankTwoTensor _fp_tmp_inv
 
RankTwoTensor _fp_tmp_old_inv
 
RankTwoTensor _fe
 
RankTwoTensor _ee
 
RankTwoTensor _pk2_fet
 
RankTwoTensor _fe_pk2
 
RankFourTensor _dpk2_dce
 
RankFourTensor _dpk2_dfe
 
RankFourTensor _dfe_dfpinv
 
RankFourTensor _dpk2_dfpinv
 
RankFourTensor _dee_dce
 
RankFourTensor _dce_dfe
 
RankFourTensor _dfe_df
 
RankFourTensor _tan_mod
 
RankFourTensor _df_dstretch_inc
 
std::vector< RankTwoTensor > _flow_dirn
 
std::vector< RankTwoTensor > _dflowrate_dpk2
 
std::vector< RankTwoTensor > _dpk2_dflowrate
 
std::vector< RankTwoTensor > _dfpinv_dflowrate
 
DenseVector< Real > _dflow_rate
 
DenseVector< Real > _flow_rate
 
DenseVector< Real > _resid
 
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
 Jacobian variables. More...
 
std::vector< DenseVector< Real > > _dintvar_dflowrate_tmp
 
DenseMatrix< Real > _dintvarrate_dintvar
 
DenseMatrix< Real > _dintvar_dintvarrate
 
DenseMatrix< Real > _dintvar_dintvar
 
DenseMatrix< Real > _dintvar_dflowrate
 
DenseMatrix< Real > _dstrength_dintvar
 
DenseMatrix< Real > _dflowrate_dstrength
 
DenseVector< Real > _dintvar_dintvar_x
 
DenseMatrix< Real > _jac
 
Real _dt_substep
 
const std::string _base_name
 
const std::string _elasticity_tensor_name
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< RankTwoTensor > & _stress
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 
const MaterialProperty< RankTwoTensor > & _extra_stress
 Extra stress tensor. More...
 
std::vector< Function * > _initial_stress_fcn
 initial stress components More...
 
MaterialProperty< RankFourTensor > & _Jacobian_mult
 derivative of stress w.r.t. strain (_dstress_dstrain) More...
 

Detailed Description

This class solves the viscoplastic flow rate equations in the total form Involves 4 different types of user objects that calculates: Internal variable rates - functions of internal variables and flow rates Internal variables - functions of internal variables Strengths - functions of internal variables Flow rates - functions of strengths and PK2 stress Flow directions - functions of strengths and PK2 stress The associated derivatives from user objects are assembled and the system is solved using NR.

Definition at line 34 of file FiniteStrainHyperElasticViscoPlastic.h.

Constructor & Destructor Documentation

◆ FiniteStrainHyperElasticViscoPlastic()

FiniteStrainHyperElasticViscoPlastic::FiniteStrainHyperElasticViscoPlastic ( const InputParameters &  parameters)

Definition at line 44 of file FiniteStrainHyperElasticViscoPlastic.C.

46  : ComputeStressBase(parameters),
47  _resid_abs_tol(getParam<Real>("resid_abs_tol")),
48  _resid_rel_tol(getParam<Real>("resid_rel_tol")),
49  _maxiters(getParam<unsigned int>("maxiters")),
50  _max_substep_iter(getParam<unsigned int>("max_substep_iteration")),
51  _flow_rate_uo_names(isParamValid("flow_rate_user_objects")
52  ? getParam<std::vector<UserObjectName>>("flow_rate_user_objects")
53  : std::vector<UserObjectName>(0)),
54  _strength_uo_names(isParamValid("strength_user_objects")
55  ? getParam<std::vector<UserObjectName>>("strength_user_objects")
56  : std::vector<UserObjectName>(0)),
57  _int_var_uo_names(isParamValid("internal_var_user_objects")
58  ? getParam<std::vector<UserObjectName>>("internal_var_user_objects")
59  : std::vector<UserObjectName>(0)),
61  isParamValid("internal_var_rate_user_objects")
62  ? getParam<std::vector<UserObjectName>>("internal_var_rate_user_objects")
63  : std::vector<UserObjectName>(0)),
64  _pk2_prop_name(_base_name + "pk2"),
65  _pk2(declareProperty<RankTwoTensor>(_pk2_prop_name)),
66  _fp(declareProperty<RankTwoTensor>(_base_name + "fp")),
67  _fp_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "fp")),
68  _ce(declareProperty<RankTwoTensor>(_base_name + "ce")),
69  _deformation_gradient(getMaterialProperty<RankTwoTensor>(_base_name + "deformation_gradient")),
71  getMaterialPropertyOld<RankTwoTensor>(_base_name + "deformation_gradient")),
72  _rotation_increment(getMaterialProperty<RankTwoTensor>(_base_name + "rotation_increment"))
73 {
75 
77 
80  _resid.resize(_num_flow_rate_uos);
81 
83 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _maxiters
Maximum number of iterations.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
std::vector< UserObjectName > _strength_uo_names
Names of strength user objects.
virtual void initUOVariables()
This function initializes the properties, stateful properties and user objects The properties and sta...
const MaterialProperty< RankTwoTensor > & _fp_old
Real _resid_abs_tol
Absolute tolerance for residual convergence check.
std::vector< UserObjectName > _flow_rate_uo_names
Names of flow rate user objects.
virtual void initJacobianVariables()
This function initialize variables required for Jacobian calculation.
std::vector< UserObjectName > _int_var_uo_names
Names of internal variable user objects.
std::vector< UserObjectName > _int_var_rate_uo_names
Names of internal variable rate user objects.
const std::string _base_name
ComputeStressBase(const InputParameters &parameters)
const MaterialProperty< RankTwoTensor > & _rotation_increment
const MaterialProperty< RankTwoTensor > & _deformation_gradient
Real _resid_rel_tol
Relative tolerance for residual convergence check.

Member Function Documentation

◆ computeDeeDce()

void FiniteStrainHyperElasticViscoPlastic::computeDeeDce ( )
protectedvirtual

Computes derivative of elastic strain w.r.t elastic Right Cauchy Green Tensor.

Definition at line 490 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initJacobianVariables().

491 {
492  _dee_dce.zero();
493 
494  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
495  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
496  _dee_dce(i, j, i, j) = 0.5;
497 }

◆ computeDpk2Dfpinv()

void FiniteStrainHyperElasticViscoPlastic::computeDpk2Dfpinv ( )
protectedvirtual

Computes derivative of PK2 stress wrt inverse of plastic deformation gradient.

Definition at line 521 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

522 {
523  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
524  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
525  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
526  _dfe_dfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
527 
529 }

◆ computeElasticPlasticDeformGrad()

void FiniteStrainHyperElasticViscoPlastic::computeElasticPlasticDeformGrad ( )
protectedvirtual

Computes elastic and plastic deformation gradients.

Definition at line 506 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

507 {
508  RankTwoTensor iden;
509  iden.addIa(1.0);
510 
511  RankTwoTensor val;
512  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
513  val += _flow_rate(i) * _flow_dirn[i] * _dt_substep;
514 
515  _fp_tmp_inv = _fp_tmp_old_inv * (iden - val);
516  _fp_tmp_inv = std::pow(_fp_tmp_inv.det(), -1.0 / 3.0) * _fp_tmp_inv;
518 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ computeElasticRightCauchyGreenTensor()

void FiniteStrainHyperElasticViscoPlastic::computeElasticRightCauchyGreenTensor ( )
protectedvirtual

Computes elastic Right Cauchy Green Tensor.

Definition at line 500 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

501 {
502  _ce[_qp] = _fe.transpose() * _fe;
503 }

◆ computeElasticStrain()

void FiniteStrainHyperElasticViscoPlastic::computeElasticStrain ( )
protectedvirtual

Computes elastic Lagrangian strain.

Definition at line 482 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by HyperElasticPhaseFieldIsoDamage::computePK2StressAndDerivative(), and computePK2StressAndDerivative().

483 {
484  RankTwoTensor iden;
485  iden.addIa(1.0);
486  _ee = 0.5 * (_ce[_qp] - iden);
487 }

◆ computeFlowDirection()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowDirection ( )
protectedvirtual

Calls user objects to compute flow directions.

Definition at line 439 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

440 {
441  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
442  {
443  if (!_flow_rate_uo[i]->computeDirection(_qp, _flow_dirn[i]))
444  return false;
445  }
446  return true;
447 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.

◆ computeFlowRateFunction()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowRateFunction ( )
protectedvirtual

Calls user objects to compute flow rates.

Definition at line 450 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

451 {
452  Real val = 0;
453  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
454  {
455  if (_flow_rate_uo[i]->computeValue(_qp, val))
456  _resid(i) = -val;
457  else
458  return false;
459  }
460  return true;
461 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.

◆ computeFlowRateJacobian()

void FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian ( )
protectedvirtual

Computes flow rate Jacobian matrix.

Definition at line 401 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

402 {
406 
407  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
408  for (unsigned int j = 0; j < _num_strength_uos; ++j)
409  _flow_rate_uo[i]->computeDerivative(_qp, _strength_uo_names[j], _dflowrate_dstrength(i, j));
410 
411  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
412  _flow_rate_uo[i]->computeTensorDerivative(_qp, _pk2_prop_name, _dflowrate_dpk2[i]);
413 
415 
416  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
417  {
420  }
421 
422  DenseMatrix<Real> dflowrate_dflowrate;
423  dflowrate_dflowrate = _dflowrate_dstrength;
424  dflowrate_dflowrate.right_multiply(_dstrength_dintvar);
425  dflowrate_dflowrate.right_multiply(_dintvar_dflowrate);
426 
427  _jac.zero();
428  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
429  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
430  {
431  if (i == j)
432  _jac(i, j) = 1;
433  _jac(i, j) -= dflowrate_dflowrate(i, j);
434  _jac(i, j) -= _dflowrate_dpk2[i].doubleContraction(_dpk2_dflowrate[j]);
435  }
436 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< UserObjectName > _strength_uo_names
Names of strength user objects.
virtual void computeDpk2Dfpinv()
Computes derivative of PK2 stress wrt inverse of plastic deformation gradient.
void computeStrengthDerivatives()
This function call user objects to compute dstrength/dintvar.
virtual void computeIntVarRateDerivatives()
This function call user objects to compute dintvar_rate/dintvar and dintvarrate/dflowrate.
virtual void computeIntVarDerivatives()
This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate.
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.

◆ computeFlowRateResidual()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual ( )
protectedvirtual

Computes flow rate residual vector.

Definition at line 375 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

376 {
377  if (!computeIntVarRates())
378  return false;
379 
380  if (!computeIntVar())
381  return false;
382 
383  if (!computeStrength())
384  return false;
385 
388 
390  return false;
391 
392  if (!computeFlowDirection())
393  return false;
394 
395  _resid += _flow_rate;
396 
397  return true;
398 }
virtual bool computeFlowDirection()
Calls user objects to compute flow directions.
virtual void computeElasticRightCauchyGreenTensor()
Computes elastic Right Cauchy Green Tensor.
virtual void computePK2StressAndDerivative()
Computes PK2 stress and derivative w.r.t elastic Right Cauchy Green Tensor.
virtual bool computeFlowRateFunction()
Calls user objects to compute flow rates.
virtual bool computeIntVar()
This function call user objects to integrate internal variables.
virtual bool computeStrength()
This function call user objects to compute strength.
virtual bool computeIntVarRates()
This function call user objects to calculate rate of internal variables.

◆ computeIntVar()

bool FiniteStrainHyperElasticViscoPlastic::computeIntVar ( )
protectedvirtual

This function call user objects to integrate internal variables.

Definition at line 597 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

598 {
599  Real val = 0;
600  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
601  {
602  if (_int_var_uo[i]->computeValue(_qp, _dt_substep, val))
603  (*_int_var_stateful_prop[i])[_qp] = val;
604  else
605  return false;
606  }
607  return true;
608 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop

◆ computeIntVarDerivatives()

void FiniteStrainHyperElasticViscoPlastic::computeIntVarDerivatives ( )
protectedvirtual

This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate.

Definition at line 638 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

639 {
640  Real val = 0;
641 
642  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
643  for (unsigned int j = 0; j < _num_int_var_rate_uos; ++j)
644  {
645  _int_var_uo[i]->computeDerivative(_qp, _dt_substep, _int_var_rate_uo_names[j], val);
646  _dintvar_dintvarrate(i, j) = val;
647  }
648 
649  _dintvar_dintvar.zero();
650 
651  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
652  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
653  {
654  if (i == j)
655  _dintvar_dintvar(i, j) = 1;
656  for (unsigned int k = 0; k < _num_int_var_rate_uos; ++k)
658  }
659 
660  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
662 
663  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
664  {
665  _dintvar_dintvar_x.zero();
667  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
669  }
670 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< DenseVector< Real > > _dintvar_dflowrate_tmp
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
Jacobian variables.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
std::vector< UserObjectName > _int_var_rate_uo_names
Names of internal variable rate user objects.

◆ computeIntVarRateDerivatives()

void FiniteStrainHyperElasticViscoPlastic::computeIntVarRateDerivatives ( )
protectedvirtual

This function call user objects to compute dintvar_rate/dintvar and dintvarrate/dflowrate.

Definition at line 625 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

626 {
627  Real val = 0;
628 
629  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
630  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
631  {
632  _int_var_rate_uo[i]->computeDerivative(_qp, _flow_rate_uo_names[j], val);
633  _dintvarrate_dflowrate[j](i) = val;
634  }
635 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
Jacobian variables.
std::vector< UserObjectName > _flow_rate_uo_names
Names of flow rate user objects.
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
Internal variable rate user objects.

◆ computeIntVarRates()

bool FiniteStrainHyperElasticViscoPlastic::computeIntVarRates ( )
protectedvirtual

This function call user objects to calculate rate of internal variables.

Definition at line 583 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

584 {
585  Real val = 0;
586  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
587  {
588  if (_int_var_rate_uo[i]->computeValue(_qp, val))
589  (*_int_var_rate_prop[i])[_qp] = val;
590  else
591  return false;
592  }
593  return true;
594 }
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
Internal variable rate user objects.
std::vector< MaterialProperty< Real > * > _int_var_rate_prop

◆ computeNorm()

Real FiniteStrainHyperElasticViscoPlastic::computeNorm ( const std::vector< Real > &  var)
protectedvirtual

Computes norm of residual vector.

Definition at line 532 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

533 {
534  Real val = 0.0;
535  for (unsigned int i = 0; i < var.size(); ++i)
536  val += Utility::pow<2>(var[i]);
537  return std::sqrt(val);
538 }

◆ computePK2StressAndDerivative()

void FiniteStrainHyperElasticViscoPlastic::computePK2StressAndDerivative ( )
protectedvirtual

Computes PK2 stress and derivative w.r.t elastic Right Cauchy Green Tensor.

Reimplemented in HyperElasticPhaseFieldIsoDamage.

Definition at line 464 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

465 {
467  _pk2[_qp] = _elasticity_tensor[_qp] * _ee;
468 
469  _dce_dfe.zero();
470  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
471  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
472  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
473  {
474  _dce_dfe(i, j, k, i) = _dce_dfe(i, j, k, i) + _fe(k, j);
475  _dce_dfe(i, j, k, j) = _dce_dfe(i, j, k, j) + _fe(k, i);
476  }
477 
479 }
virtual void computeElasticStrain()
Computes elastic Lagrangian strain.
const MaterialProperty< RankFourTensor > & _elasticity_tensor

◆ computeQpJacobian()

void FiniteStrainHyperElasticViscoPlastic::computeQpJacobian ( )
protectedvirtual

This function computes the Jacobian.

Reimplemented in HyperElasticPhaseFieldIsoDamage.

Definition at line 551 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by HyperElasticPhaseFieldIsoDamage::computeQpJacobian(), and postSolveQp().

552 {
553  _tan_mod = _fe.mixedProductIkJl(_fe) * _dpk2_dfe;
554  _pk2_fet = _pk2[_qp] * _fe.transpose();
555  _fe_pk2 = _fe * _pk2[_qp];
556 
557  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
558  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
559  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
560  {
561  _tan_mod(i, j, i, l) += _pk2_fet(l, j);
562  _tan_mod(i, j, j, l) += _fe_pk2(i, l);
563  }
564 
565  _tan_mod /= _fe.det();
566 
567  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
568  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
569  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
570  _dfe_df(i, j, i, l) = _fp_tmp_inv(l, j);
571 
572  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
573  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
574  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
575  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
576  _df_dstretch_inc(i, j, k, l) =
577  _rotation_increment[_qp](i, k) * _deformation_gradient_old[_qp](l, j);
578 
580 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
const MaterialProperty< RankTwoTensor > & _rotation_increment

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 51 of file ComputeStressBase.C.

52 {
54 
55  // Add in extra stress
56  _stress[_qp] += _extra_stress[_qp];
57 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.

◆ computeQpStress()

void FiniteStrainHyperElasticViscoPlastic::computeQpStress ( )
protectedvirtual

This function computes the Cauchy stress.

Implements ComputeStressBase.

Definition at line 210 of file FiniteStrainHyperElasticViscoPlastic.C.

211 {
212  bool converge;
213  RankTwoTensor delta_dfgrd = _deformation_gradient[_qp] - _deformation_gradient_old[_qp];
214  unsigned int num_substep = 1;
215  unsigned int substep_iter = 1;
216 
217  saveOldState();
218 
219  do
220  {
221  preSolveQp();
222 
223  converge = true;
224  _dt_substep = _dt / num_substep;
225 
226  for (unsigned int istep = 0; istep < num_substep; ++istep)
227  {
228  _dfgrd_tmp = (istep + 1.0) * delta_dfgrd / num_substep + _deformation_gradient_old[_qp];
229  if (!solveQp())
230  {
231  converge = false;
232  substep_iter++;
233  num_substep *= 2;
234  break;
235  }
236  }
237 
238  if (substep_iter > _max_substep_iter)
239  mooseError("Constitutive failure with substepping at quadrature point ",
240  _q_point[_qp](0),
241  " ",
242  _q_point[_qp](1),
243  " ",
244  _q_point[_qp](2));
245  } while (!converge);
246 
247  postSolveQp();
248 }
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
virtual void postSolveQp()
Update state for output (Outside substepping)
virtual void saveOldState()
This function saves the old stateful properties that is modified during sub stepping.
const MaterialProperty< RankTwoTensor > & _deformation_gradient

◆ computeStrength()

bool FiniteStrainHyperElasticViscoPlastic::computeStrength ( )
protectedvirtual

This function call user objects to compute strength.

Definition at line 611 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

612 {
613  Real val = 0;
614  for (unsigned int i = 0; i < _num_strength_uos; ++i)
615  {
616  if (_strength_uo[i]->computeValue(_qp, val))
617  (*_strength_prop[i])[_qp] = val;
618  else
619  return false;
620  }
621  return true;
622 }
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< MaterialProperty< Real > * > _strength_prop

◆ computeStrengthDerivatives()

void FiniteStrainHyperElasticViscoPlastic::computeStrengthDerivatives ( )
protected

This function call user objects to compute dstrength/dintvar.

Definition at line 673 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

674 {
675  Real val = 0;
676 
677  for (unsigned int i = 0; i < _num_strength_uos; ++i)
678  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
679  {
680  _strength_uo[i]->computeDerivative(_qp, _int_var_uo_names[j], val);
681  _dstrength_dintvar(i, j) = val;
682  }
683 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< UserObjectName > _int_var_uo_names
Names of internal variable user objects.
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
unsigned int _num_strength_uos
Number of strength user objects.

◆ initJacobianVariables()

void FiniteStrainHyperElasticViscoPlastic::initJacobianVariables ( )
protectedvirtual

This function initialize variables required for Jacobian calculation.

Definition at line 154 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic().

155 {
157  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
159 
161  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
163 
171 
175 
177 
178  computeDeeDce();
179 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< DenseVector< Real > > _dintvar_dflowrate_tmp
virtual void computeDeeDce()
Computes derivative of elastic strain w.r.t elastic Right Cauchy Green Tensor.
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
Jacobian variables.
unsigned int _num_strength_uos
Number of strength user objects.

◆ initNumUserObjects()

void FiniteStrainHyperElasticViscoPlastic::initNumUserObjects ( const std::vector< UserObjectName > &  uo_names,
unsigned int &  uo_num 
)
protected

This function calculates the number of each user object type.

Definition at line 109 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

111 {
112  uo_num = uo_names.size();
113 }

◆ initProp()

template<typename T >
void FiniteStrainHyperElasticViscoPlastic::initProp ( const std::vector< UserObjectName > &  uo_names,
unsigned int  uo_num,
std::vector< MaterialProperty< T > *> &  uo_prop 
)
protected

This function initializes properties for each user object.

Definition at line 117 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

120 {
121  uo_prop.resize(uo_num);
122  for (unsigned int i = 0; i < uo_num; ++i)
123  uo_prop[i] = &declareProperty<T>(uo_names[i]);
124 }

◆ initPropOld()

template<typename T >
void FiniteStrainHyperElasticViscoPlastic::initPropOld ( const std::vector< UserObjectName > &  uo_names,
unsigned int  uo_num,
std::vector< const MaterialProperty< T > *> &  uo_prop_old 
)
protected

This function initializes old for stateful properties associated with user object Only user objects that update internal variables have an associated old property.

Definition at line 128 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

132 {
133  uo_prop_old.resize(uo_num);
134  for (unsigned int i = 0; i < uo_num; ++i)
135  uo_prop_old[i] = &getMaterialPropertyOld<T>(uo_names[i]);
136 }

◆ initQpStatefulProperties()

void FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties ( )
protectedvirtual

Initializes state.

Reimplemented from ComputeStressBase.

Definition at line 182 of file FiniteStrainHyperElasticViscoPlastic.C.

183 {
184  _stress[_qp].zero();
185 
186  _fp[_qp].zero();
187  _fp[_qp].addIa(1.0);
188  _ce[_qp].zero();
189 
190  _pk2[_qp].zero();
191 
192  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
193  (*_flow_rate_prop[i])[_qp] = 0.0;
194 
195  for (unsigned int i = 0; i < _num_strength_uos; ++i)
196  (*_strength_prop[i])[_qp] = 0.0;
197 
198  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
199  {
200  (*_int_var_stateful_prop[i])[_qp] = 0.0;
201  // TODO: remove this nasty const_cast if you can figure out how
202  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = 0.0;
203  }
204 
205  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
206  (*_int_var_rate_prop[i])[_qp] = 0.0;
207 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
MaterialProperty< RankTwoTensor > & _stress
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
std::vector< MaterialProperty< Real > * > _int_var_rate_prop
std::vector< MaterialProperty< Real > * > _strength_prop

◆ initUOVariables()

void FiniteStrainHyperElasticViscoPlastic::initUOVariables ( )
protectedvirtual

This function initializes the properties, stateful properties and user objects The properties and stateful properties associated with user objects are only initialized here The properties have the same name as the user object name.

Definition at line 86 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic().

87 {
92 
97 
99 
104 
105  _int_var_old.resize(_num_int_var_uos, 0.0);
106 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
void initPropOld(const std::vector< UserObjectName > &, unsigned int, std::vector< const MaterialProperty< T > *> &)
This function initializes old for stateful properties associated with user object Only user objects t...
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
void initNumUserObjects(const std::vector< UserObjectName > &, unsigned int &)
This function calculates the number of each user object type.
std::vector< UserObjectName > _strength_uo_names
Names of strength user objects.
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< UserObjectName > _flow_rate_uo_names
Names of flow rate user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< UserObjectName > _int_var_uo_names
Names of internal variable user objects.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< UserObjectName > _int_var_rate_uo_names
Names of internal variable rate user objects.
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
Internal variable rate user objects.
void initProp(const std::vector< UserObjectName > &, unsigned int, std::vector< MaterialProperty< T > *> &)
This function initializes properties for each user object.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.
void initUserObjects(const std::vector< UserObjectName > &, unsigned int, std::vector< const T *> &)
This function initializes user objects.
std::vector< MaterialProperty< Real > * > _int_var_rate_prop
std::vector< MaterialProperty< Real > * > _strength_prop

◆ initUserObjects()

template<typename T >
void FiniteStrainHyperElasticViscoPlastic::initUserObjects ( const std::vector< UserObjectName > &  uo_names,
unsigned int  uo_num,
std::vector< const T *> &  uo 
)
protected

This function initializes user objects.

Definition at line 140 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

143 {
144  uo.resize(uo_num);
145 
146  if (uo_num == 0)
147  mooseError("Specify atleast one user object of type", typeid(T).name());
148 
149  for (unsigned int i = 0; i < uo_num; ++i)
150  uo[i] = &getUserObjectByName<T>(uo_names[i]);
151 }
const std::string name
Definition: Setup.h:22

◆ postSolveFlowrate()

void FiniteStrainHyperElasticViscoPlastic::postSolveFlowrate ( )
protectedvirtual

Update state for output (Inside substepping)

Definition at line 364 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

365 {
367 
368  // TODO: remove this nasty const_cast if you can figure out how
369  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
370  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] =
371  (*_int_var_stateful_prop[i])[_qp];
372 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop

◆ postSolveQp()

void FiniteStrainHyperElasticViscoPlastic::postSolveQp ( )
protectedvirtual

Update state for output (Outside substepping)

Definition at line 282 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

283 {
284  recoverOldState();
285 
286  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
287  _fp[_qp] = _fp_tmp_inv.inverse();
288 
290 }
virtual void recoverOldState()
This function restores the the old stateful properties after a successful solve.
MaterialProperty< RankTwoTensor > & _stress
virtual void computeQpJacobian()
This function computes the Jacobian.

◆ preSolveFlowrate()

void FiniteStrainHyperElasticViscoPlastic::preSolveFlowrate ( )
protectedvirtual

Sets state for solve (Inside substepping)

Definition at line 301 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

302 {
303  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
304  {
305  _flow_rate(i) = 0.0;
306  (*_flow_rate_prop[i])[_qp] = 0.0;
307  }
308 
309  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
310  (*_int_var_stateful_prop[i])[_qp] = (*_int_var_stateful_prop_old[i])[_qp];
311 
314 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop

◆ preSolveQp()

void FiniteStrainHyperElasticViscoPlastic::preSolveQp ( )
protectedvirtual

Sets state for solve.

Definition at line 258 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

259 {
260  _fp_tmp_old_inv = _fp_old[_qp].inverse();
261 
262  // TODO: remove this nasty const_cast if you can figure out how
263  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
264  (*_int_var_stateful_prop[i])[_qp] =
265  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
266 
268 }
const MaterialProperty< RankTwoTensor > & _fp_old
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
const MaterialProperty< RankFourTensor > & _elasticity_tensor

◆ recoverOldState()

void FiniteStrainHyperElasticViscoPlastic::recoverOldState ( )
protectedvirtual

This function restores the the old stateful properties after a successful solve.

Definition at line 293 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by postSolveQp().

294 {
295  // TODO: remove this nasty const_cast if you can figure out how
296  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
297  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
298 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old

◆ saveOldState()

void FiniteStrainHyperElasticViscoPlastic::saveOldState ( )
protectedvirtual

This function saves the old stateful properties that is modified during sub stepping.

Definition at line 251 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

252 {
253  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
255 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old

◆ solveFlowrate()

bool FiniteStrainHyperElasticViscoPlastic::solveFlowrate ( )
protectedvirtual

Solve for flow rate and state.

Definition at line 317 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

318 {
319  Real resid0, rnorm;
320  unsigned int iter = 0;
321 
322 #ifdef DEBUG
323  std::vector<Real> rnormst(_maxiters + 1), flowratest(_maxiters + 1);
324 #endif
325 
327  return false;
328 
329  rnorm = computeNorm(_resid.get_values());
330  resid0 = rnorm;
331 
332 #ifdef DEBUG
333  rnormst[iter] = rnorm;
334  flowratest[iter] = computeNorm(_flow_rate.get_values());
335 #endif
336 
337  while (rnorm > _resid_abs_tol && rnorm > _resid_rel_tol * resid0 && iter < _maxiters)
338  {
340 
341  updateFlowRate();
342 
344 
346  return false;
347 
348  rnorm = computeNorm(_resid.get_values());
349  iter++;
350 
351 #ifdef DEBUG
352  rnormst[iter] = rnorm;
353  flowratest[iter] = computeNorm(_flow_rate.get_values());
354 #endif
355  }
356 
357  if (iter == _maxiters && rnorm > _resid_abs_tol && rnorm > _resid_rel_tol * resid0)
358  return false;
359 
360  return true;
361 }
unsigned int _maxiters
Maximum number of iterations.
Real _resid_abs_tol
Absolute tolerance for residual convergence check.
virtual bool computeFlowRateResidual()
Computes flow rate residual vector.
virtual void computeFlowRateJacobian()
Computes flow rate Jacobian matrix.
virtual Real computeNorm(const std::vector< Real > &)
Computes norm of residual vector.
virtual void computeElasticPlasticDeformGrad()
Computes elastic and plastic deformation gradients.
Real _resid_rel_tol
Relative tolerance for residual convergence check.

◆ solveQp()

bool FiniteStrainHyperElasticViscoPlastic::solveQp ( )
protectedvirtual

Solve state.

Definition at line 271 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

272 {
274  if (!solveFlowrate())
275  return false;
277 
278  return true;
279 }
virtual bool solveFlowrate()
Solve for flow rate and state.
virtual void preSolveFlowrate()
Sets state for solve (Inside substepping)
virtual void postSolveFlowrate()
Update state for output (Inside substepping)

◆ updateFlowRate()

void FiniteStrainHyperElasticViscoPlastic::updateFlowRate ( )
protectedvirtual

Update flow rate.

Definition at line 541 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

542 {
543  _jac.lu_solve(_resid, _dflow_rate);
545 
546  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
547  (*_flow_rate_prop[i])[_qp] = _flow_rate(i);
548 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop

Member Data Documentation

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

◆ _ce

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_ce
protected

◆ _dce_dfe

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dce_dfe
protected

◆ _dee_dce

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dee_dce
protected

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient
protected

Definition at line 200 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpStress().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient_old
protected

Definition at line 201 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian(), and computeQpStress().

◆ _df_dstretch_inc

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_df_dstretch_inc
protected

◆ _dfe_df

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dfe_df
protected

◆ _dfe_dfpinv

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dfe_dfpinv
protected

Definition at line 215 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeDpk2Dfpinv().

◆ _dfgrd_tmp

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_dfgrd_tmp
protected

◆ _dflow_rate

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_dflow_rate
protected

◆ _dflowrate_dpk2

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dflowrate_dpk2
protected

◆ _dflowrate_dstrength

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dflowrate_dstrength
protected

◆ _dfpinv_dflowrate

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dfpinv_dflowrate
protected

◆ _dintvar_dflowrate

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dflowrate
protected

◆ _dintvar_dflowrate_tmp

std::vector<DenseVector<Real> > FiniteStrainHyperElasticViscoPlastic::_dintvar_dflowrate_tmp
protected

◆ _dintvar_dintvar

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvar
protected

◆ _dintvar_dintvar_x

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvar_x
protected

◆ _dintvar_dintvarrate

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvarrate
protected

◆ _dintvarrate_dflowrate

std::vector<DenseVector<Real> > FiniteStrainHyperElasticViscoPlastic::_dintvarrate_dflowrate
protected

◆ _dintvarrate_dintvar

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvarrate_dintvar
protected

◆ _dpk2_dce

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dce
protected

◆ _dpk2_dfe

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dfe
protected

◆ _dpk2_dflowrate

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dpk2_dflowrate
protected

◆ _dpk2_dfpinv

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dfpinv
protected

◆ _dstrength_dintvar

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dstrength_dintvar
protected

◆ _dt_substep

Real FiniteStrainHyperElasticViscoPlastic::_dt_substep
protected

◆ _ee

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_ee
protected

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& ComputeStressBase::_elasticity_tensor
protectedinherited

◆ _elasticity_tensor_name

const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 47 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _fe

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe
protected

◆ _fe_pk2

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe_pk2
protected

Definition at line 214 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().

◆ _flow_dirn

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_flow_dirn
protected

◆ _flow_rate

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_flow_rate
protected

◆ _flow_rate_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_flow_rate_prop
protected

◆ _flow_rate_uo

std::vector<const HEVPFlowRateUOBase *> FiniteStrainHyperElasticViscoPlastic::_flow_rate_uo
protected

◆ _flow_rate_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_flow_rate_uo_names
protected

Names of flow rate user objects.

Definition at line 168 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeIntVarRateDerivatives(), and initUOVariables().

◆ _fp

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_fp
protected

◆ _fp_old

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_fp_old
protected

Definition at line 197 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by preSolveQp().

◆ _fp_tmp_inv

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fp_tmp_inv
protected

◆ _fp_tmp_old_inv

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fp_tmp_old_inv
protected

◆ _initial_stress_fcn

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

initial stress components

Definition at line 50 of file ComputeStressBase.h.

◆ _int_var_old

std::vector<Real> FiniteStrainHyperElasticViscoPlastic::_int_var_old
protected

◆ _int_var_rate_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_prop
protected

◆ _int_var_rate_uo

std::vector<const HEVPInternalVarRateUOBase *> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_uo
protected

Internal variable rate user objects.

Definition at line 192 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeIntVarRateDerivatives(), computeIntVarRates(), and initUOVariables().

◆ _int_var_rate_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_uo_names
protected

Names of internal variable rate user objects.

Definition at line 174 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeIntVarDerivatives(), and initUOVariables().

◆ _int_var_stateful_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_stateful_prop
protected

◆ _int_var_stateful_prop_old

std::vector<const MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_stateful_prop_old
protected

◆ _int_var_uo

std::vector<const HEVPInternalVarUOBase *> FiniteStrainHyperElasticViscoPlastic::_int_var_uo
protected

Internal variable user objects.

Definition at line 190 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeIntVar(), computeIntVarDerivatives(), and initUOVariables().

◆ _int_var_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_int_var_uo_names
protected

Names of internal variable user objects.

Definition at line 172 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeStrengthDerivatives(), and initUOVariables().

◆ _jac

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_jac
protected

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _max_substep_iter

unsigned int FiniteStrainHyperElasticViscoPlastic::_max_substep_iter
protected

Maximum number of substep iterations.

Definition at line 165 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpStress().

◆ _maxiters

unsigned int FiniteStrainHyperElasticViscoPlastic::_maxiters
protected

Maximum number of iterations.

Definition at line 163 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by solveFlowrate().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _num_flow_rate_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_flow_rate_uos
protected

◆ _num_int_var_rate_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_int_var_rate_uos
protected

◆ _num_int_var_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_int_var_uos
protected

◆ _num_strength_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_strength_uos
protected

◆ _pk2

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_pk2
protected

◆ _pk2_fet

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_pk2_fet
protected

Definition at line 214 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().

◆ _pk2_prop_name

std::string FiniteStrainHyperElasticViscoPlastic::_pk2_prop_name
protected

Definition at line 194 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeFlowRateJacobian().

◆ _resid

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_resid
protected

◆ _resid_abs_tol

Real FiniteStrainHyperElasticViscoPlastic::_resid_abs_tol
protected

Absolute tolerance for residual convergence check.

Definition at line 159 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by solveFlowrate().

◆ _resid_rel_tol

Real FiniteStrainHyperElasticViscoPlastic::_resid_rel_tol
protected

Relative tolerance for residual convergence check.

Definition at line 161 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by solveFlowrate().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_rotation_increment
protected

Definition at line 202 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().

◆ _strength_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_strength_prop
protected

◆ _strength_uo

std::vector<const HEVPStrengthUOBase *> FiniteStrainHyperElasticViscoPlastic::_strength_uo
protected

Strength user objects.

Definition at line 188 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeStrength(), computeStrengthDerivatives(), and initUOVariables().

◆ _strength_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_strength_uo_names
protected

Names of strength user objects.

Definition at line 170 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeFlowRateJacobian(), and initUOVariables().

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

◆ _tan_mod

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_tan_mod
protected

Definition at line 217 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().


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