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 std::string _elasticity_tensor_name
 Name of the elasticity tensor material property. More...
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 Elasticity tensor material property. More...
 
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
 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< 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  _elasticity_tensor_name(_base_name + "elasticity_tensor"),
70  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
71  _deformation_gradient(getMaterialProperty<RankTwoTensor>(_base_name + "deformation_gradient")),
73  getMaterialPropertyOld<RankTwoTensor>(_base_name + "deformation_gradient")),
74  _rotation_increment(getMaterialProperty<RankTwoTensor>(_base_name + "rotation_increment"))
75 {
77 
79 
82  _resid.resize(_num_flow_rate_uos);
83 
85 }
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.
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
std::vector< UserObjectName > _int_var_rate_uo_names
Names of internal variable rate user objects.
const std::string _base_name
Base name prepended to all material property names to allow for multi-material systems.
ComputeStressBase(const InputParameters &parameters)
const MaterialProperty< RankTwoTensor > & _rotation_increment
const MaterialProperty< RankTwoTensor > & _deformation_gradient
Real _resid_rel_tol
Relative tolerance for residual convergence check.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.

Member Function Documentation

◆ computeDeeDce()

void FiniteStrainHyperElasticViscoPlastic::computeDeeDce ( )
protectedvirtual

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

Definition at line 488 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initJacobianVariables().

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

◆ computeDpk2Dfpinv()

void FiniteStrainHyperElasticViscoPlastic::computeDpk2Dfpinv ( )
protectedvirtual

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

Definition at line 518 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

519 {
520  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
521  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
522  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
523  _dfe_dfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
524 
526 }

◆ computeElasticPlasticDeformGrad()

void FiniteStrainHyperElasticViscoPlastic::computeElasticPlasticDeformGrad ( )
protectedvirtual

Computes elastic and plastic deformation gradients.

Definition at line 504 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

505 {
506  RankTwoTensor iden(RankTwoTensor::initIdentity);
507 
508  RankTwoTensor val;
509  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
510  val += _flow_rate(i) * _flow_dirn[i] * _dt_substep;
511 
512  _fp_tmp_inv = _fp_tmp_old_inv * (iden - val);
513  _fp_tmp_inv = std::pow(_fp_tmp_inv.det(), -1.0 / 3.0) * _fp_tmp_inv;
515 }
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 498 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

499 {
500  _ce[_qp] = _fe.transpose() * _fe;
501 }

◆ computeElasticStrain()

void FiniteStrainHyperElasticViscoPlastic::computeElasticStrain ( )
protectedvirtual

Computes elastic Lagrangian strain.

Definition at line 481 of file FiniteStrainHyperElasticViscoPlastic.C.

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

482 {
483  RankTwoTensor iden(RankTwoTensor::initIdentity);
484  _ee = 0.5 * (_ce[_qp] - iden);
485 }

◆ computeFlowDirection()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowDirection ( )
protectedvirtual

Calls user objects to compute flow directions.

Definition at line 438 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

439 {
440  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
441  {
442  if (!_flow_rate_uo[i]->computeDirection(_qp, _flow_dirn[i]))
443  return false;
444  }
445  return true;
446 }
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 449 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

450 {
451  Real val = 0;
452  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
453  {
454  if (_flow_rate_uo[i]->computeValue(_qp, val))
455  _resid(i) = -val;
456  else
457  return false;
458  }
459  return true;
460 }
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 400 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

401 {
405 
406  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
407  for (unsigned int j = 0; j < _num_strength_uos; ++j)
408  _flow_rate_uo[i]->computeDerivative(_qp, _strength_uo_names[j], _dflowrate_dstrength(i, j));
409 
410  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
411  _flow_rate_uo[i]->computeTensorDerivative(_qp, _pk2_prop_name, _dflowrate_dpk2[i]);
412 
414 
415  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
416  {
419  }
420 
421  DenseMatrix<Real> dflowrate_dflowrate;
422  dflowrate_dflowrate = _dflowrate_dstrength;
423  dflowrate_dflowrate.right_multiply(_dstrength_dintvar);
424  dflowrate_dflowrate.right_multiply(_dintvar_dflowrate);
425 
426  _jac.zero();
427  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
428  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
429  {
430  if (i == j)
431  _jac(i, j) = 1;
432  _jac(i, j) -= dflowrate_dflowrate(i, j);
433  _jac(i, j) -= _dflowrate_dpk2[i].doubleContraction(_dpk2_dflowrate[j]);
434  }
435 }
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 374 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

375 {
376  if (!computeIntVarRates())
377  return false;
378 
379  if (!computeIntVar())
380  return false;
381 
382  if (!computeStrength())
383  return false;
384 
387 
389  return false;
390 
391  if (!computeFlowDirection())
392  return false;
393 
394  _resid += _flow_rate;
395 
396  return true;
397 }
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 594 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

595 {
596  Real val = 0;
597  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
598  {
599  if (_int_var_uo[i]->computeValue(_qp, _dt_substep, val))
600  (*_int_var_stateful_prop[i])[_qp] = val;
601  else
602  return false;
603  }
604  return true;
605 }
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 635 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

636 {
637  Real val = 0;
638 
639  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
640  for (unsigned int j = 0; j < _num_int_var_rate_uos; ++j)
641  {
642  _int_var_uo[i]->computeDerivative(_qp, _dt_substep, _int_var_rate_uo_names[j], val);
643  _dintvar_dintvarrate(i, j) = val;
644  }
645 
646  _dintvar_dintvar.zero();
647 
648  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
649  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
650  {
651  if (i == j)
652  _dintvar_dintvar(i, j) = 1;
653  for (unsigned int k = 0; k < _num_int_var_rate_uos; ++k)
655  }
656 
657  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
659 
660  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
661  {
662  _dintvar_dintvar_x.zero();
664  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
666  }
667 }
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 622 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

623 {
624  Real val = 0;
625 
626  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
627  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
628  {
629  _int_var_rate_uo[i]->computeDerivative(_qp, _flow_rate_uo_names[j], val);
630  _dintvarrate_dflowrate[j](i) = val;
631  }
632 }
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 580 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

581 {
582  Real val = 0;
583  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
584  {
585  if (_int_var_rate_uo[i]->computeValue(_qp, val))
586  (*_int_var_rate_prop[i])[_qp] = val;
587  else
588  return false;
589  }
590  return true;
591 }
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 529 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

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

◆ computePK2StressAndDerivative()

void FiniteStrainHyperElasticViscoPlastic::computePK2StressAndDerivative ( )
protectedvirtual

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

Reimplemented in HyperElasticPhaseFieldIsoDamage.

Definition at line 463 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

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

◆ computeQpJacobian()

void FiniteStrainHyperElasticViscoPlastic::computeQpJacobian ( )
protectedvirtual

This function computes the Jacobian.

Reimplemented in HyperElasticPhaseFieldIsoDamage.

Definition at line 548 of file FiniteStrainHyperElasticViscoPlastic.C.

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

549 {
550  _tan_mod = _fe.mixedProductIkJl(_fe) * _dpk2_dfe;
551  _pk2_fet = _pk2[_qp] * _fe.transpose();
552  _fe_pk2 = _fe * _pk2[_qp];
553 
554  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
555  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
556  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
557  {
558  _tan_mod(i, j, i, l) += _pk2_fet(l, j);
559  _tan_mod(i, j, j, l) += _fe_pk2(i, l);
560  }
561 
562  _tan_mod /= _fe.det();
563 
564  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
565  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
566  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
567  _dfe_df(i, j, i, l) = _fp_tmp_inv(l, j);
568 
569  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
570  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
571  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
572  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
573  _df_dstretch_inc(i, j, k, l) =
574  _rotation_increment[_qp](i, k) * _deformation_gradient_old[_qp](l, j);
575 
577 }
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 49 of file ComputeStressBase.C.

50 {
52 
53  // Add in extra stress
54  _stress[_qp] += _extra_stress[_qp];
55 }
virtual void computeQpStress()=0
Compute the stress and store it in the _stress material property for the current quadrature point...
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.

◆ computeQpStress()

void FiniteStrainHyperElasticViscoPlastic::computeQpStress ( )
protectedvirtual

This function computes the Cauchy stress.

Implements ComputeStressBase.

Definition at line 209 of file FiniteStrainHyperElasticViscoPlastic.C.

210 {
211  bool converge;
213  unsigned int num_substep = 1;
214  unsigned int substep_iter = 1;
215 
216  saveOldState();
217 
218  do
219  {
220  preSolveQp();
221 
222  converge = true;
223  _dt_substep = _dt / num_substep;
224 
225  for (unsigned int istep = 0; istep < num_substep; ++istep)
226  {
227  _dfgrd_tmp = (istep + 1.0) * delta_dfgrd / num_substep + _deformation_gradient_old[_qp];
228  if (!solveQp())
229  {
230  converge = false;
231  substep_iter++;
232  num_substep *= 2;
233  break;
234  }
235  }
236 
237  if (substep_iter > _max_substep_iter)
238  mooseError("Constitutive failure with substepping at quadrature point ",
239  _q_point[_qp](0),
240  " ",
241  _q_point[_qp](1),
242  " ",
243  _q_point[_qp](2));
244  } while (!converge);
245 
246  postSolveQp();
247 }
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 608 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

609 {
610  Real val = 0;
611  for (unsigned int i = 0; i < _num_strength_uos; ++i)
612  {
613  if (_strength_uo[i]->computeValue(_qp, val))
614  (*_strength_prop[i])[_qp] = val;
615  else
616  return false;
617  }
618  return true;
619 }
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 670 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

671 {
672  Real val = 0;
673 
674  for (unsigned int i = 0; i < _num_strength_uos; ++i)
675  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
676  {
677  _strength_uo[i]->computeDerivative(_qp, _int_var_uo_names[j], val);
678  _dstrength_dintvar(i, j) = val;
679  }
680 }
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 156 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic().

157 {
159  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
161 
163  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
165 
173 
177 
179 
180  computeDeeDce();
181 }
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 111 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

113 {
114  uo_num = uo_names.size();
115 }

◆ 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 119 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

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

◆ 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 130 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

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

◆ initQpStatefulProperties()

void FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties ( )
protectedvirtual

Initializes state.

Reimplemented from ComputeStressBase.

Definition at line 184 of file FiniteStrainHyperElasticViscoPlastic.C.

185 {
186  _stress[_qp].zero();
187  _ce[_qp].zero();
188  _pk2[_qp].zero();
189  _fp[_qp].setToIdentity();
190 
191  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
192  (*_flow_rate_prop[i])[_qp] = 0.0;
193 
194  for (unsigned int i = 0; i < _num_strength_uos; ++i)
195  (*_strength_prop[i])[_qp] = 0.0;
196 
197  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
198  {
199  (*_int_var_stateful_prop[i])[_qp] = 0.0;
200  // TODO: remove this nasty const_cast if you can figure out how
201  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = 0.0;
202  }
203 
204  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
205  (*_int_var_rate_prop[i])[_qp] = 0.0;
206 }
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
Stress material property.
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 88 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic().

89 {
94 
99 
101 
106 
107  _int_var_old.resize(_num_int_var_uos, 0.0);
108 }
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 142 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

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

◆ postSolveFlowrate()

void FiniteStrainHyperElasticViscoPlastic::postSolveFlowrate ( )
protectedvirtual

Update state for output (Inside substepping)

Definition at line 363 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

364 {
366 
367  // TODO: remove this nasty const_cast if you can figure out how
368  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
369  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] =
370  (*_int_var_stateful_prop[i])[_qp];
371 }
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 281 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

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

◆ preSolveFlowrate()

void FiniteStrainHyperElasticViscoPlastic::preSolveFlowrate ( )
protectedvirtual

Sets state for solve (Inside substepping)

Definition at line 300 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

301 {
302  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
303  {
304  _flow_rate(i) = 0.0;
305  (*_flow_rate_prop[i])[_qp] = 0.0;
306  }
307 
308  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
309  (*_int_var_stateful_prop[i])[_qp] = (*_int_var_stateful_prop_old[i])[_qp];
310 
313 }
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 257 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

258 {
259  _fp_tmp_old_inv = _fp_old[_qp].inverse();
260 
261  // TODO: remove this nasty const_cast if you can figure out how
262  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
263  (*_int_var_stateful_prop[i])[_qp] =
264  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
265 
267 }
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
Elasticity tensor material property.

◆ recoverOldState()

void FiniteStrainHyperElasticViscoPlastic::recoverOldState ( )
protectedvirtual

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

Definition at line 292 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by postSolveQp().

293 {
294  // TODO: remove this nasty const_cast if you can figure out how
295  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
296  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
297 }
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 250 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

251 {
252  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
254 }
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 316 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

317 {
318  Real resid0, rnorm;
319  unsigned int iter = 0;
320 
321 #ifdef DEBUG
322  std::vector<Real> rnormst(_maxiters + 1), flowratest(_maxiters + 1);
323 #endif
324 
326  return false;
327 
328  rnorm = computeNorm(_resid.get_values());
329  resid0 = rnorm;
330 
331 #ifdef DEBUG
332  rnormst[iter] = rnorm;
333  flowratest[iter] = computeNorm(_flow_rate.get_values());
334 #endif
335 
336  while (rnorm > _resid_abs_tol && rnorm > _resid_rel_tol * resid0 && iter < _maxiters)
337  {
339 
340  updateFlowRate();
341 
343 
345  return false;
346 
347  rnorm = computeNorm(_resid.get_values());
348  iter++;
349 
350 #ifdef DEBUG
351  rnormst[iter] = rnorm;
352  flowratest[iter] = computeNorm(_flow_rate.get_values());
353 #endif
354  }
355 
356  if (iter == _maxiters && rnorm > _resid_abs_tol && rnorm > _resid_rel_tol * resid0)
357  return false;
358 
359  return true;
360 }
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 270 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

271 {
273  if (!solveFlowrate())
274  return false;
276 
277  return true;
278 }
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 538 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

539 {
540  _jac.lu_solve(_resid, _dflow_rate);
542 
543  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
544  (*_flow_rate_prop[i])[_qp] = _flow_rate(i);
545 }
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

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

Definition at line 44 of file ComputeStressBase.h.

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

◆ _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 204 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpStress().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient_old
protected

Definition at line 205 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 219 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>& FiniteStrainHyperElasticViscoPlastic::_elasticity_tensor
protected

Elasticity tensor material property.

Definition at line 203 of file FiniteStrainHyperElasticViscoPlastic.h.

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

◆ _elasticity_tensor_name

const std::string FiniteStrainHyperElasticViscoPlastic::_elasticity_tensor_name
protected

Name of the elasticity tensor material property.

Definition at line 201 of file FiniteStrainHyperElasticViscoPlastic.h.

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 54 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _fe

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe
protected

◆ _fe_pk2

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe_pk2
protected

Definition at line 218 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 57 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 218 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 206 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

Stress material property.

Definition at line 49 of file ComputeStressBase.h.

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

◆ _tan_mod

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_tan_mod
protected

Definition at line 221 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().


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