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
 
void addQpInitialStress ()
 InitialStress Deprecation: remove this method. More...
 

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...
 
const bool _store_stress_old
 Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping. More...
 
const bool _initial_stress_provided
 Whether initial stress was provided. InitialStress Deprecation: remove this. More...
 
MaterialProperty< RankTwoTensor > * _initial_stress
 Initial stress, if provided. InitialStress Deprecation: remove this. More...
 
const MaterialProperty< RankTwoTensor > * _initial_stress_old
 Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this. 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 31 of file FiniteStrainHyperElasticViscoPlastic.h.

Constructor & Destructor Documentation

FiniteStrainHyperElasticViscoPlastic::FiniteStrainHyperElasticViscoPlastic ( const InputParameters &  parameters)

Definition at line 39 of file FiniteStrainHyperElasticViscoPlastic.C.

41  : ComputeStressBase(parameters),
42  _resid_abs_tol(getParam<Real>("resid_abs_tol")),
43  _resid_rel_tol(getParam<Real>("resid_rel_tol")),
44  _maxiters(getParam<unsigned int>("maxiters")),
45  _max_substep_iter(getParam<unsigned int>("max_substep_iteration")),
46  _flow_rate_uo_names(isParamValid("flow_rate_user_objects")
47  ? getParam<std::vector<UserObjectName>>("flow_rate_user_objects")
48  : std::vector<UserObjectName>(0)),
49  _strength_uo_names(isParamValid("strength_user_objects")
50  ? getParam<std::vector<UserObjectName>>("strength_user_objects")
51  : std::vector<UserObjectName>(0)),
52  _int_var_uo_names(isParamValid("internal_var_user_objects")
53  ? getParam<std::vector<UserObjectName>>("internal_var_user_objects")
54  : std::vector<UserObjectName>(0)),
56  isParamValid("internal_var_rate_user_objects")
57  ? getParam<std::vector<UserObjectName>>("internal_var_rate_user_objects")
58  : std::vector<UserObjectName>(0)),
59  _pk2_prop_name(_base_name + "pk2"),
60  _pk2(declareProperty<RankTwoTensor>(_pk2_prop_name)),
61  _fp(declareProperty<RankTwoTensor>(_base_name + "fp")),
62  _fp_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "fp")),
63  _ce(declareProperty<RankTwoTensor>(_base_name + "ce")),
64  _deformation_gradient(getMaterialProperty<RankTwoTensor>(_base_name + "deformation_gradient")),
66  getMaterialPropertyOld<RankTwoTensor>(_base_name + "deformation_gradient")),
67  _rotation_increment(getMaterialProperty<RankTwoTensor>(_base_name + "rotation_increment"))
68 {
70 
72 
75  _resid.resize(_num_flow_rate_uos);
76 
78 }
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

void ComputeStressBase::addQpInitialStress ( )
protectedinherited

InitialStress Deprecation: remove this method.

Adds initial stress, if it is provided, to _stress[_qp]. This function should NOT be used if you calculate stress using

stress = stress_old + elasticity * strain_increment

because stress_old will already include initial stress. However this function SHOULD be used if your code uses

stress = elasticity * (elastic_strain_old + strain_increment) or stress = elasticity * elastic_strain

since in these cases the elastic_strain and elastic_strain_old will not include any contribution from initial stress.

Definition at line 112 of file ComputeStressBase.C.

Referenced by ComputeLinearElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeStressBase::initQpStatefulProperties(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

113 {
115  _stress[_qp] += (*_initial_stress)[_qp];
116 }
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
void FiniteStrainHyperElasticViscoPlastic::computeDeeDce ( )
protectedvirtual

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

Definition at line 485 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initJacobianVariables().

486 {
487  _dee_dce.zero();
488 
489  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
490  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
491  _dee_dce(i, j, i, j) = 0.5;
492 }
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 }
void FiniteStrainHyperElasticViscoPlastic::computeElasticPlasticDeformGrad ( )
protectedvirtual

Computes elastic and plastic deformation gradients.

Definition at line 501 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

502 {
503  RankTwoTensor iden;
504  iden.addIa(1.0);
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)
void FiniteStrainHyperElasticViscoPlastic::computeElasticRightCauchyGreenTensor ( )
protectedvirtual

Computes elastic Right Cauchy Green Tensor.

Definition at line 495 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

496 {
497  _ce[_qp] = _fe.transpose() * _fe;
498 }
void FiniteStrainHyperElasticViscoPlastic::computeElasticStrain ( )
protectedvirtual

Computes elastic Lagrangian strain.

Definition at line 477 of file FiniteStrainHyperElasticViscoPlastic.C.

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

478 {
479  RankTwoTensor iden;
480  iden.addIa(1.0);
481  _ee = 0.5 * (_ce[_qp] - iden);
482 }
bool FiniteStrainHyperElasticViscoPlastic::computeFlowDirection ( )
protectedvirtual

Calls user objects to compute flow directions.

Definition at line 434 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

435 {
436  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
437  {
438  if (!_flow_rate_uo[i]->computeDirection(_qp, _flow_dirn[i]))
439  return false;
440  }
441  return true;
442 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.
bool FiniteStrainHyperElasticViscoPlastic::computeFlowRateFunction ( )
protectedvirtual

Calls user objects to compute flow rates.

Definition at line 445 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

446 {
447  Real val = 0;
448  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
449  {
450  if (_flow_rate_uo[i]->computeValue(_qp, val))
451  _resid(i) = -val;
452  else
453  return false;
454  }
455  return true;
456 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.
void FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian ( )
protectedvirtual

Computes flow rate Jacobian matrix.

Definition at line 396 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

397 {
401 
402  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
403  for (unsigned int j = 0; j < _num_strength_uos; ++j)
404  _flow_rate_uo[i]->computeDerivative(_qp, _strength_uo_names[j], _dflowrate_dstrength(i, j));
405 
406  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
407  _flow_rate_uo[i]->computeTensorDerivative(_qp, _pk2_prop_name, _dflowrate_dpk2[i]);
408 
410 
411  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
412  {
415  }
416 
417  DenseMatrix<Real> dflowrate_dflowrate;
418  dflowrate_dflowrate = _dflowrate_dstrength;
419  dflowrate_dflowrate.right_multiply(_dstrength_dintvar);
420  dflowrate_dflowrate.right_multiply(_dintvar_dflowrate);
421 
422  _jac.zero();
423  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
424  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
425  {
426  if (i == j)
427  _jac(i, j) = 1;
428  _jac(i, j) -= dflowrate_dflowrate(i, j);
429  _jac(i, j) -= _dflowrate_dpk2[i].doubleContraction(_dpk2_dflowrate[j]);
430  }
431 }
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.
bool FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual ( )
protectedvirtual

Computes flow rate residual vector.

Definition at line 370 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveFlowrate().

371 {
372  if (!computeIntVarRates())
373  return false;
374 
375  if (!computeIntVar())
376  return false;
377 
378  if (!computeStrength())
379  return false;
380 
383 
385  return false;
386 
387  if (!computeFlowDirection())
388  return false;
389 
390  _resid += _flow_rate;
391 
392  return true;
393 }
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.
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
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.
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.
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
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 }
void FiniteStrainHyperElasticViscoPlastic::computePK2StressAndDerivative ( )
protectedvirtual

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

Reimplemented in HyperElasticPhaseFieldIsoDamage.

Definition at line 459 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeFlowRateResidual().

460 {
462  _pk2[_qp] = _elasticity_tensor[_qp] * _ee;
463 
464  _dce_dfe.zero();
465  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
466  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
467  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
468  {
469  _dce_dfe(i, j, k, i) = _dce_dfe(i, j, k, i) + _fe(k, j);
470  _dce_dfe(i, j, k, j) = _dce_dfe(i, j, k, j) + _fe(k, i);
471  }
472 
474 }
virtual void computeElasticStrain()
Computes elastic Lagrangian strain.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
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
void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Reimplemented in ComputeBirchMurnaghanEquationOfStress, and ComputeStressEosBase.

Definition at line 99 of file ComputeStressBase.C.

100 {
101  // InitialStress Deprecation: remove the following 2 lines
103  (*_initial_stress)[_qp] = (*_initial_stress_old)[_qp];
104 
105  computeQpStress();
106 
107  // Add in extra stress
108  _stress[_qp] += _extra_stress[_qp];
109 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.
void FiniteStrainHyperElasticViscoPlastic::computeQpStress ( )
protectedvirtual

This function computes the Cauchy stress.

Implements ComputeStressBase.

Definition at line 205 of file FiniteStrainHyperElasticViscoPlastic.C.

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

This function initialize variables required for Jacobian calculation.

Definition at line 149 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic().

150 {
152  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
154 
155  _dintvar_dflowrate_tmp.resize(_num_flow_rate_uos);
156  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
158 
161  _dintvar_dflowrate.resize(_num_int_var_uos, _num_flow_rate_uos);
164  _dflowrate_dstrength.resize(_num_flow_rate_uos, _num_strength_uos);
166 
167  _dpk2_dflowrate.resize(_num_flow_rate_uos);
168  _dflowrate_dpk2.resize(_num_flow_rate_uos);
169  _dfpinv_dflowrate.resize(_num_flow_rate_uos);
170 
171  _jac.resize(_num_flow_rate_uos, _num_flow_rate_uos);
172 
173  computeDeeDce();
174 }
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.
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 104 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

106 {
107  uo_num = uo_names.size();
108 }
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 112 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

115 {
116  uo_prop.resize(uo_num);
117  for (unsigned int i = 0; i < uo_num; ++i)
118  uo_prop[i] = &declareProperty<T>(uo_names[i]);
119 }
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 123 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

127 {
128  uo_prop_old.resize(uo_num);
129  for (unsigned int i = 0; i < uo_num; ++i)
130  uo_prop_old[i] = &getMaterialPropertyOld<T>(uo_names[i]);
131 }
void FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties ( )
protectedvirtual

Initializes state.

Reimplemented from ComputeStressBase.

Definition at line 177 of file FiniteStrainHyperElasticViscoPlastic.C.

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

Referenced by FiniteStrainHyperElasticViscoPlastic().

82 {
87 
92 
94 
99 
100  _int_var_old.resize(_num_int_var_uos, 0.0);
101 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
void initUserObjects(const std::vector< UserObjectName > &, unsigned int, std::vector< const T * > &)
This function initializes user objects.
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.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
void initProp(const std::vector< UserObjectName > &, unsigned int, std::vector< MaterialProperty< T > * > &)
This function initializes properties for each user object.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.
std::vector< MaterialProperty< Real > * > _int_var_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...
std::vector< MaterialProperty< Real > * > _strength_prop
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 135 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by initUOVariables().

138 {
139  uo.resize(uo_num);
140 
141  if (uo_num == 0)
142  mooseError("Specify atleast one user object of type", typeid(T).name());
143 
144  for (unsigned int i = 0; i < uo_num; ++i)
145  uo[i] = &getUserObjectByName<T>(uo_names[i]);
146 }
void FiniteStrainHyperElasticViscoPlastic::postSolveFlowrate ( )
protectedvirtual

Update state for output (Inside substepping)

Definition at line 359 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

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

Update state for output (Outside substepping)

Definition at line 277 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

278 {
279  recoverOldState();
280 
281  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
282  _fp[_qp] = _fp_tmp_inv.inverse();
283 
285 }
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.
void FiniteStrainHyperElasticViscoPlastic::preSolveFlowrate ( )
protectedvirtual

Sets state for solve (Inside substepping)

Definition at line 296 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

297 {
298  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
299  {
300  _flow_rate(i) = 0.0;
301  (*_flow_rate_prop[i])[_qp] = 0.0;
302  }
303 
304  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
305  (*_int_var_stateful_prop[i])[_qp] = (*_int_var_stateful_prop_old[i])[_qp];
306 
309 }
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
void FiniteStrainHyperElasticViscoPlastic::preSolveQp ( )
protectedvirtual

Sets state for solve.

Definition at line 253 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

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

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

Definition at line 288 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by postSolveQp().

289 {
290  // TODO: remove this nasty const_cast if you can figure out how
291  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
292  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
293 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
void FiniteStrainHyperElasticViscoPlastic::saveOldState ( )
protectedvirtual

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

Definition at line 246 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

247 {
248  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
250 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
bool FiniteStrainHyperElasticViscoPlastic::solveFlowrate ( )
protectedvirtual

Solve for flow rate and state.

Definition at line 312 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by solveQp().

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

Solve state.

Definition at line 266 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computeQpStress().

267 {
269  if (!solveFlowrate())
270  return false;
272 
273  return true;
274 }
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)
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

const std::string ComputeStressBase::_base_name
protectedinherited
MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_ce
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dce_dfe
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dee_dce
protected
const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient
protected

Definition at line 197 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpStress().

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient_old
protected

Definition at line 198 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian(), and computeQpStress().

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_df_dstretch_inc
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dfe_df
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dfe_dfpinv
protected

Definition at line 212 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeDpk2Dfpinv().

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_dfgrd_tmp
protected
DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_dflow_rate
protected
std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dflowrate_dpk2
protected
DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dflowrate_dstrength
protected
std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dfpinv_dflowrate
protected
DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dflowrate
protected
std::vector<DenseVector<Real> > FiniteStrainHyperElasticViscoPlastic::_dintvar_dflowrate_tmp
protected
DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvar
protected
DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvar_x
protected
DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvarrate
protected
std::vector<DenseVector<Real> > FiniteStrainHyperElasticViscoPlastic::_dintvarrate_dflowrate
protected
DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvarrate_dintvar
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dce
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dfe
protected
std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dpk2_dflowrate
protected
RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dfpinv
protected
DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dstrength_dintvar
protected
Real FiniteStrainHyperElasticViscoPlastic::_dt_substep
protected
RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_ee
protected
MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited
const MaterialProperty<RankFourTensor>& ComputeStressBase::_elasticity_tensor
protectedinherited
const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited
RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe
protected
RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe_pk2
protected

Definition at line 211 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_flow_dirn
protected
DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_flow_rate
protected
std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_flow_rate_prop
protected
std::vector<const HEVPFlowRateUOBase *> FiniteStrainHyperElasticViscoPlastic::_flow_rate_uo
protected
std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_flow_rate_uo_names
protected

Names of flow rate user objects.

Definition at line 165 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeIntVarRateDerivatives(), and initUOVariables().

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_fp
protected
const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_fp_old
protected

Definition at line 194 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by preSolveQp().

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fp_tmp_inv
protected
RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fp_tmp_old_inv
protected
MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress
protectedinherited

Initial stress, if provided. InitialStress Deprecation: remove this.

Definition at line 74 of file ComputeStressBase.h.

Referenced by ComputeStressBase::initQpStatefulProperties().

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

initial stress components

Definition at line 62 of file ComputeStressBase.h.

Referenced by ComputeStressBase::ComputeStressBase(), and ComputeStressBase::initQpStatefulProperties().

const MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress_old
protectedinherited

Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this.

Definition at line 77 of file ComputeStressBase.h.

const bool ComputeStressBase::_initial_stress_provided
protectedinherited
std::vector<Real> FiniteStrainHyperElasticViscoPlastic::_int_var_old
protected
std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_prop
protected
std::vector<const HEVPInternalVarRateUOBase *> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_uo
protected

Internal variable rate user objects.

Definition at line 189 of file FiniteStrainHyperElasticViscoPlastic.h.

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

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

Names of internal variable rate user objects.

Definition at line 171 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeIntVarDerivatives(), and initUOVariables().

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

Internal variable user objects.

Definition at line 187 of file FiniteStrainHyperElasticViscoPlastic.h.

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

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

Names of internal variable user objects.

Definition at line 169 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeStrengthDerivatives(), and initUOVariables().

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_jac
protected
MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited
unsigned int FiniteStrainHyperElasticViscoPlastic::_max_substep_iter
protected

Maximum number of substep iterations.

Definition at line 162 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpStress().

unsigned int FiniteStrainHyperElasticViscoPlastic::_maxiters
protected

Maximum number of iterations.

Definition at line 160 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by solveFlowrate().

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited
unsigned int FiniteStrainHyperElasticViscoPlastic::_num_flow_rate_uos
protected
unsigned int FiniteStrainHyperElasticViscoPlastic::_num_int_var_rate_uos
protected
unsigned int FiniteStrainHyperElasticViscoPlastic::_num_int_var_uos
protected
unsigned int FiniteStrainHyperElasticViscoPlastic::_num_strength_uos
protected
MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_pk2
protected
RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_pk2_fet
protected

Definition at line 211 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().

std::string FiniteStrainHyperElasticViscoPlastic::_pk2_prop_name
protected

Definition at line 191 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeFlowRateJacobian().

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_resid
protected
Real FiniteStrainHyperElasticViscoPlastic::_resid_abs_tol
protected

Absolute tolerance for residual convergence check.

Definition at line 156 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by solveFlowrate().

Real FiniteStrainHyperElasticViscoPlastic::_resid_rel_tol
protected

Relative tolerance for residual convergence check.

Definition at line 158 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by solveFlowrate().

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_rotation_increment
protected

Definition at line 199 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().

const bool ComputeStressBase::_store_stress_old
protectedinherited

Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping.

Definition at line 68 of file ComputeStressBase.h.

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_strength_prop
protected
std::vector<const HEVPStrengthUOBase *> FiniteStrainHyperElasticViscoPlastic::_strength_uo
protected

Strength user objects.

Definition at line 185 of file FiniteStrainHyperElasticViscoPlastic.h.

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

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

Names of strength user objects.

Definition at line 167 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeFlowRateJacobian(), and initUOVariables().

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Definition at line 53 of file ComputeStressBase.h.

Referenced by ComputeStressBase::addQpInitialStress(), ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), ComputeMultipleInelasticStress::computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStressEosBase::computeQpProperties(), ComputeBirchMurnaghanEquationOfStress::computeQpProperties(), ComputeLinearElasticStress::computeQpStress(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeIsotropicLinearElasticPFFractureStress::computeQpStress(), ComputeFiniteStrainElasticStressBirchMurnaghan::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeSmearedCrackingStress::crackingStressRotation(), ComputeMultipleInelasticStress::finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeMultipleInelasticStress::updateQpState(), ComputeMultipleInelasticStress::updateQpStateSingleModel(), and LinearIsoElasticPFDamage::updateVar().

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_tan_mod
protected

Definition at line 214 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by computeQpJacobian().


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