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

Referenced by initJacobianVariables().

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

◆ computeDpk2Dfpinv()

void FiniteStrainHyperElasticViscoPlastic::computeDpk2Dfpinv ( )
protectedvirtual

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

Definition at line 516 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateJacobian().

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

◆ computeElasticPlasticDeformGrad()

void FiniteStrainHyperElasticViscoPlastic::computeElasticPlasticDeformGrad ( )
protectedvirtual

Computes elastic and plastic deformation gradients.

Definition at line 502 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

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

Referenced by computeFlowRateResidual().

497 {
498  _ce[_qp] = _fe.transpose() * _fe;
499 }

◆ computeElasticStrain()

void FiniteStrainHyperElasticViscoPlastic::computeElasticStrain ( )
protectedvirtual

Computes elastic Lagrangian strain.

Definition at line 479 of file FiniteStrainHyperElasticViscoPlastic.C.

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

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

◆ computeFlowDirection()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowDirection ( )
protectedvirtual

Calls user objects to compute flow directions.

Definition at line 436 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

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

Referenced by computeFlowRateResidual().

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

Referenced by solveFlowrate().

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

Referenced by solveFlowrate().

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

Referenced by computeFlowRateResidual().

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

Referenced by computeFlowRateJacobian().

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

Referenced by computeFlowRateJacobian().

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

Referenced by computeFlowRateResidual().

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

Referenced by solveFlowrate().

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

◆ computePK2StressAndDerivative()

void FiniteStrainHyperElasticViscoPlastic::computePK2StressAndDerivative ( )
protectedvirtual

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

Reimplemented in HyperElasticPhaseFieldIsoDamage.

Definition at line 461 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

462 {
464  _pk2[_qp] = _elasticity_tensor[_qp] * _ee;
465 
466  _dce_dfe.zero();
467  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
468  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
469  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
470  {
471  _dce_dfe(i, j, k, i) = _dce_dfe(i, j, k, i) + _fe(k, j);
472  _dce_dfe(i, j, k, j) = _dce_dfe(i, j, k, j) + _fe(k, i);
473  }
474 
476 }
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 546 of file FiniteStrainHyperElasticViscoPlastic.C.

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

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

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

Referenced by computeFlowRateResidual().

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

Referenced by computeFlowRateJacobian().

669 {
670  Real val = 0;
671 
672  for (unsigned int i = 0; i < _num_strength_uos; ++i)
673  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
674  {
675  _strength_uo[i]->computeDerivative(_qp, _int_var_uo_names[j], val);
676  _dstrength_dintvar(i, j) = val;
677  }
678 }
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  _ce[_qp].zero();
186  _pk2[_qp].zero();
187  _fp[_qp].setToIdentity();
188 
189  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
190  (*_flow_rate_prop[i])[_qp] = 0.0;
191 
192  for (unsigned int i = 0; i < _num_strength_uos; ++i)
193  (*_strength_prop[i])[_qp] = 0.0;
194 
195  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
196  {
197  (*_int_var_stateful_prop[i])[_qp] = 0.0;
198  // TODO: remove this nasty const_cast if you can figure out how
199  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = 0.0;
200  }
201 
202  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
203  (*_int_var_rate_prop[i])[_qp] = 0.0;
204 }
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 361 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

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

Referenced by computeQpStress().

280 {
281  recoverOldState();
282 
283  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
284  _fp[_qp] = _fp_tmp_inv.inverse();
285 
287 }
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 298 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

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

Referenced by computeQpStress().

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

Referenced by postSolveQp().

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

Referenced by computeQpStress().

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

Referenced by solveQp().

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

Referenced by computeQpStress().

269 {
271  if (!solveFlowrate())
272  return false;
274 
275  return true;
276 }
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 536 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

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

Definition at line 41 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 217 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().


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