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

MultiPlasticityDebugger computes various finite-difference things to help developers remove bugs in their derivatives, etc. More...

#include <MultiPlasticityDebugger.h>

Inheritance diagram for MultiPlasticityDebugger:
[legend]

Public Member Functions

 MultiPlasticityDebugger (const MooseObject *moose_object)
 
void outputAndCheckDebugParameters ()
 Outputs the debug parameters: _fspb_debug_stress, _fspd_debug_pm, etc and checks that they are sized correctly. More...
 
void checkDerivatives ()
 Checks the derivatives, eg dyieldFunction_dstress by using finite difference approximations. More...
 
void checkJacobian (const RankFourTensor &E_inv, const std::vector< Real > &intnl_old)
 Checks the full Jacobian, which is just certain linear combinations of the dyieldFunction_dstress, etc, by using finite difference approximations. More...
 
void checkSolution (const RankFourTensor &E_inv)
 Checks that Ax does equal b in the NR procedure. More...
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Member Functions

virtual void calculateConstraints (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &f, std::vector< RankTwoTensor > &r, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active)
 The constraints. More...
 
virtual void calculateRHS (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &rhs, const std::vector< bool > &active, bool eliminate_ld, std::vector< bool > &deactivated_due_to_ld)
 Calculate the RHS which is rhs = -(epp(0,0), epp(1,0), epp(1,1), epp(2,0), epp(2,1), epp(2,2), f[0], f[1], ..., f[num_f], ic[0], ic[1], ..., ic[num_ic]) More...
 
virtual void calculateJacobian (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, std::vector< std::vector< Real >> &jac)
 d(rhs)/d(dof) More...
 
virtual void nrStep (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const RankTwoTensor &delta_dp, RankTwoTensor &dstress, std::vector< Real > &dpm, std::vector< Real > &dintnl, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
 Performs one Newton-Raphson step. More...
 
virtual void yieldFunction (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &f)
 The active yield function(s) More...
 
virtual void dyieldFunction_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &df_dstress)
 The derivative of the active yield function(s) with respect to stress. More...
 
virtual void dyieldFunction_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &df_dintnl)
 The derivative of active yield function(s) with respect to their internal parameters (the user objects assume there is exactly one internal param per yield function) More...
 
virtual void flowPotential (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &r)
 The active flow potential(s) - one for each yield function. More...
 
virtual void dflowPotential_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankFourTensor > &dr_dstress)
 The derivative of the active flow potential(s) with respect to stress. More...
 
virtual void dflowPotential_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dr_dintnl)
 The derivative of the active flow potentials with respect to the active internal parameters The UserObjects explicitly assume that r[alpha] is only dependent on intnl[alpha]. More...
 
virtual void hardPotential (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &h)
 The active hardening potentials (one for each internal parameter and for each yield function) by assumption in the Userobjects, the h[a][alpha] is nonzero only if the surface alpha is part of model a, so we only calculate those here. More...
 
virtual void dhardPotential_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dh_dstress)
 The derivative of the active hardening potentials with respect to stress By assumption in the Userobjects, the h[a][alpha] is nonzero only for a = alpha, so we only calculate those here. More...
 
virtual void dhardPotential_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &dh_dintnl)
 The derivative of the active hardening potentials with respect to the active internal parameters. More...
 
virtual void buildActiveConstraints (const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
 Constructs a set of active constraints, given the yield functions, f. More...
 
unsigned int modelNumber (unsigned int surface)
 returns the model number, given the surface number More...
 
bool anyActiveSurfaces (int model, const std::vector< bool > &active)
 returns true if any internal surfaces of the given model are active according to 'active' More...
 
void activeModelSurfaces (int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces_of_model)
 Returns the internal surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model. More...
 
void activeSurfaces (int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces)
 Returns the external surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model. More...
 
bool returnMapAll (const RankTwoTensor &trial_stress, const std::vector< Real > &intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &stress, std::vector< Real > &intnl, std::vector< Real > &pm, std::vector< Real > &cumulative_pm, RankTwoTensor &delta_dp, std::vector< Real > &yf, unsigned &num_successful_plastic_returns, unsigned &custom_model)
 Performs a returnMap for each plastic model using their inbuilt returnMap functions. More...
 

Protected Attributes

MooseEnum _fspb_debug
 none - don't do any debugging crash - currently inactive jacobian - check the jacobian entries jacobian_and_linear_system - check entire jacobian and check that Ax=b More...
 
RankTwoTensor _fspb_debug_stress
 Debug the Jacobian entries at this stress. More...
 
std::vector< Real > _fspb_debug_pm
 Debug the Jacobian entires at these plastic multipliers. More...
 
std::vector< Real > _fspb_debug_intnl
 Debug the Jacobian entires at these internal parameters. More...
 
Real _fspb_debug_stress_change
 Debug finite-differencing parameter for the stress. More...
 
std::vector< Real > _fspb_debug_pm_change
 Debug finite-differencing parameters for the plastic multipliers. More...
 
std::vector< Real > _fspb_debug_intnl_change
 Debug finite-differencing parameters for the internal parameters. More...
 
Real _svd_tol
 Tolerance on the minimum ratio of singular values before flow-directions are deemed linearly dependent. More...
 
Real _min_f_tol
 Minimum value of the _f_tol parameters for the Yield Function User Objects. More...
 
const InputParameters & _params
 
unsigned int _num_models
 Number of plastic models for this material. More...
 
unsigned int _num_surfaces
 Number of surfaces within the plastic models. More...
 
std::vector< std::vector< unsigned int > > _surfaces_given_model
 _surfaces_given_model[model_number] = vector of surface numbers for this model More...
 
MooseEnum _specialIC
 Allows initial set of active constraints to be chosen optimally. More...
 
std::vector< const TensorMechanicsPlasticModel * > _f
 User objects that define the yield functions, flow potentials, etc. More...
 

Private Member Functions

void fddyieldFunction_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &df_dstress)
 The finite-difference derivative of yield function(s) with respect to stress. More...
 
void fddyieldFunction_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &df_dintnl)
 The finite-difference derivative of yield function(s) with respect to internal parameter(s) More...
 
virtual void fddflowPotential_dstress (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankFourTensor > &dr_dstress)
 The finite-difference derivative of the flow potential(s) with respect to stress. More...
 
virtual void fddflowPotential_dintnl (const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &dr_dintnl)
 The finite-difference derivative of the flow potentials with respect to internal parameters. More...
 
virtual void fdJacobian (const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, const RankFourTensor &E_inv, bool eliminate_ld, std::vector< std::vector< Real >> &jac)
 The Jacobian calculated using finite differences. More...
 
bool dof_included (unsigned int dof, const std::vector< bool > &deactivated_due_to_ld)
 
virtual int singularValuesOfR (const std::vector< RankTwoTensor > &r, std::vector< Real > &s)
 Performs a singular-value decomposition of r and returns the singular values. More...
 
virtual void eliminateLinearDependence (const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &f, const std::vector< RankTwoTensor > &r, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
 Performs a number of singular-value decompositions to check for linear-dependence of the active directions "r" If linear dependence is found, then deactivated_due_to_ld will contain 'true' entries where surfaces need to be deactivated_due_to_ld. More...
 
void buildActiveConstraintsRock (const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
 "Rock" version Constructs a set of active constraints, given the yield functions, f. More...
 
void buildActiveConstraintsJoint (const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
 "Joint" version Constructs a set of active constraints, given the yield functions, f. More...
 

Private Attributes

std::vector< unsigned int > _model_given_surface
 given a surface number, this returns the model number More...
 
std::vector< unsigned int > _model_surface_given_surface
 given a surface number, this returns the corresponding-model's internal surface number More...
 

Detailed Description

MultiPlasticityDebugger computes various finite-difference things to help developers remove bugs in their derivatives, etc.

Definition at line 24 of file MultiPlasticityDebugger.h.

Constructor & Destructor Documentation

◆ MultiPlasticityDebugger()

MultiPlasticityDebugger::MultiPlasticityDebugger ( const MooseObject *  moose_object)

Definition at line 42 of file MultiPlasticityDebugger.C.

43  : MultiPlasticityLinearSystem(moose_object),
44  _fspb_debug(_params.get<MooseEnum>("debug_fspb")),
45  _fspb_debug_stress(_params.get<RealTensorValue>("debug_jac_at_stress")),
46  _fspb_debug_pm(_params.get<std::vector<Real>>("debug_jac_at_pm")),
47  _fspb_debug_intnl(_params.get<std::vector<Real>>("debug_jac_at_intnl")),
48  _fspb_debug_stress_change(_params.get<Real>("debug_stress_change")),
49  _fspb_debug_pm_change(_params.get<std::vector<Real>>("debug_pm_change")),
50  _fspb_debug_intnl_change(_params.get<std::vector<Real>>("debug_intnl_change"))
51 {
52 }

Member Function Documentation

◆ activeModelSurfaces()

void MultiPlasticityRawComponentAssembler::activeModelSurfaces ( int  model,
const std::vector< bool > &  active,
std::vector< unsigned int > &  active_surfaces_of_model 
)
protectedinherited

Returns the internal surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model.

Parameters
modelthe model number
activearray with entries being 'true' if the surface is active
[out]active_surfaces_of_modelthe output

Definition at line 811 of file MultiPlasticityRawComponentAssembler.C.

815 {
816  active_surfaces_of_model.resize(0);
817  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
818  if (active[_surfaces_given_model[model][model_surface]])
819  active_surfaces_of_model.push_back(model_surface);
820 }

Referenced by MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ activeSurfaces()

void MultiPlasticityRawComponentAssembler::activeSurfaces ( int  model,
const std::vector< bool > &  active,
std::vector< unsigned int > &  active_surfaces 
)
protectedinherited

Returns the external surface number(s) of the active surfaces of the given model This may be of size=0 if there are no active surfaces of the given model.

Parameters
modelthe model number
activearray with entries being 'true' if the surface is active
[out]active_surfacesthe output

Definition at line 800 of file MultiPlasticityRawComponentAssembler.C.

803 {
804  active_surfaces.resize(0);
805  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
806  if (active[_surfaces_given_model[model][model_surface]])
807  active_surfaces.push_back(_surfaces_given_model[model][model_surface]);
808 }

Referenced by MultiPlasticityLinearSystem::calculateConstraints().

◆ anyActiveSurfaces()

bool MultiPlasticityRawComponentAssembler::anyActiveSurfaces ( int  model,
const std::vector< bool > &  active 
)
protectedinherited

returns true if any internal surfaces of the given model are active according to 'active'

Definition at line 791 of file MultiPlasticityRawComponentAssembler.C.

792 {
793  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
794  if (active[_surfaces_given_model[model][model_surface]])
795  return true;
796  return false;
797 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityLinearSystem::calculateRHS(), checkSolution(), dof_included(), MultiPlasticityLinearSystem::nrStep(), and ComputeMultiPlasticityStress::residual2().

◆ buildActiveConstraints()

void MultiPlasticityRawComponentAssembler::buildActiveConstraints ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act 
)
protectedvirtualinherited

Constructs a set of active constraints, given the yield functions, f.

This uses TensorMechanicsPlasticModel::activeConstraints to identify the active constraints for each model.

Parameters
fyield functions (should be _num_surfaces of these)
stressstress tensor
intnlinternal parameters
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actthe set of active constraints (will be resized to _num_surfaces)

Definition at line 344 of file MultiPlasticityRawComponentAssembler.C.

349 {
350  mooseAssert(f.size() == _num_surfaces,
351  "buildActiveConstraints called with f.size = " << f.size() << " while there are "
352  << _num_surfaces << " surfaces");
353  mooseAssert(intnl.size() == _num_models,
354  "buildActiveConstraints called with intnl.size = "
355  << intnl.size() << " while there are " << _num_models << " models");
356 
357  if (_specialIC == "rock")
358  buildActiveConstraintsRock(f, stress, intnl, Eijkl, act);
359  else if (_specialIC == "joint")
360  buildActiveConstraintsJoint(f, stress, intnl, Eijkl, act);
361  else // no specialIC
362  {
363  act.resize(0);
364  unsigned ind = 0;
365  for (unsigned model = 0; model < _num_models; ++model)
366  {
367  std::vector<Real> model_f(0);
368  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
369  model_f.push_back(f[ind++]);
370  std::vector<bool> model_act;
371  RankTwoTensor returned_stress;
372  _f[model]->activeConstraints(
373  model_f, stress, intnl[model], Eijkl, model_act, returned_stress);
374  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
375  act.push_back(model_act[model_surface]);
376  }
377  }
378 }

Referenced by ComputeMultiPlasticityStress::returnMap().

◆ buildActiveConstraintsJoint()

void MultiPlasticityRawComponentAssembler::buildActiveConstraintsJoint ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act 
)
privateinherited

"Joint" version Constructs a set of active constraints, given the yield functions, f.

This uses TensorMechanicsPlasticModel::activeConstraints to identify the active constraints for each model.

Parameters
fyield functions (should be _num_surfaces of these)
stressstress tensor
intnlinternal parameters
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actthe set of active constraints (will be resized to _num_surfaces)

Definition at line 381 of file MultiPlasticityRawComponentAssembler.C.

386 {
387  act.assign(2, false);
388 
389  RankTwoTensor returned_stress;
390  std::vector<bool> active_tensile;
391  std::vector<bool> active_shear;
392  std::vector<Real> f_single;
393 
394  // first try tensile alone
395  f_single.assign(1, 0);
396  f_single[0] = f[0];
397  _f[0]->activeConstraints(f_single, stress, intnl[0], Eijkl, active_tensile, returned_stress);
398  _f[1]->yieldFunctionV(returned_stress, intnl[1], f_single);
399  if (f_single[0] <= _f[1]->_f_tol)
400  {
401  act[0] = active_tensile[0];
402  return;
403  }
404 
405  // next try shear alone
406  f_single.assign(1, 0);
407  f_single[0] = f[1];
408  _f[1]->activeConstraints(f_single, stress, intnl[1], Eijkl, active_shear, returned_stress);
409  _f[0]->yieldFunctionV(returned_stress, intnl[0], f_single);
410  if (f_single[0] <= _f[0]->_f_tol)
411  {
412  act[1] = active_shear[0];
413  return;
414  }
415 
416  // must be mixed
417  act[0] = act[1] = true;
418  return;
419 }

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints().

◆ buildActiveConstraintsRock()

void MultiPlasticityRawComponentAssembler::buildActiveConstraintsRock ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act 
)
privateinherited

"Rock" version Constructs a set of active constraints, given the yield functions, f.

This uses TensorMechanicsPlasticModel::activeConstraints to identify the active constraints for each model.

Parameters
fyield functions (should be _num_surfaces of these)
stressstress tensor
intnlinternal parameters
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actthe set of active constraints (will be resized to _num_surfaces)

Definition at line 422 of file MultiPlasticityRawComponentAssembler.C.

427 {
428  act.assign(9, false);
429 
430  RankTwoTensor returned_stress;
431  std::vector<bool> active_tensile;
432  std::vector<bool> active_MC;
433  std::vector<Real> f_single;
434 
435  // first try tensile alone
436  f_single.assign(3, 0);
437  f_single[0] = f[0];
438  f_single[1] = f[1];
439  f_single[2] = f[2];
440  _f[0]->activeConstraints(f_single, stress, intnl[0], Eijkl, active_tensile, returned_stress);
441  _f[1]->yieldFunctionV(returned_stress, intnl[1], f_single);
442  if (f_single[0] <= _f[1]->_f_tol && f_single[1] <= _f[1]->_f_tol &&
443  f_single[2] <= _f[1]->_f_tol && f_single[3] <= _f[1]->_f_tol &&
444  f_single[4] <= _f[1]->_f_tol && f_single[5] <= _f[1]->_f_tol)
445  {
446  act[0] = active_tensile[0];
447  act[1] = active_tensile[1];
448  act[2] = active_tensile[2];
449  return;
450  }
451 
452  // next try MC alone
453  f_single.assign(6, 0);
454  f_single[0] = f[3];
455  f_single[1] = f[4];
456  f_single[2] = f[5];
457  f_single[3] = f[6];
458  f_single[4] = f[7];
459  f_single[5] = f[8];
460  _f[1]->activeConstraints(f_single, stress, intnl[1], Eijkl, active_MC, returned_stress);
461  _f[0]->yieldFunctionV(returned_stress, intnl[0], f_single);
462  if (f_single[0] <= _f[0]->_f_tol && f_single[1] <= _f[0]->_f_tol && f_single[2] <= _f[0]->_f_tol)
463  {
464  act[3] = active_MC[0];
465  act[4] = active_MC[1];
466  act[5] = active_MC[2];
467  act[6] = active_MC[3];
468  act[7] = active_MC[4];
469  act[8] = active_MC[5];
470  return;
471  }
472 
473  // must be a mix.
474  // The possibilities are enumerated below.
475 
476  // tensile=edge, MC=tip (two possibilities)
477  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
478  active_MC[0] == true && active_MC[1] == true && active_MC[2] == false &&
479  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
480  {
481  act[1] = act[2] = act[6] = true;
482  act[4] = true;
483  return;
484  }
485  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
486  active_MC[0] == false && active_MC[1] == true && active_MC[2] == false &&
487  active_MC[3] == true && active_MC[4] == false && active_MC[5] == true)
488  {
489  act[1] = act[2] = act[6] = true; // i don't think act[4] is necessary, is it?!
490  return;
491  }
492 
493  // tensile = edge, MC=edge (two possibilities)
494  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
495  active_MC[0] == false && active_MC[1] == true && active_MC[2] == false &&
496  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
497  {
498  act[1] = act[2] = act[4] = act[6] = true;
499  return;
500  }
501  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
502  active_MC[0] == false && active_MC[1] == false && active_MC[2] == false &&
503  active_MC[3] == true && active_MC[4] == false && active_MC[5] == true)
504  {
505  act[1] = act[2] = act[4] = act[6] = true;
506  return;
507  }
508 
509  // tensile = edge, MC=face
510  if (active_tensile[0] == false && active_tensile[1] == true && active_tensile[2] == true &&
511  active_MC[0] == false && active_MC[1] == false && active_MC[2] == false &&
512  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
513  {
514  act[1] = act[2] = act[6] = true;
515  return;
516  }
517 
518  // tensile = face, MC=tip (two possibilities)
519  if (active_tensile[0] == false && active_tensile[1] == false && active_tensile[2] == true &&
520  active_MC[0] == true && active_MC[1] == true && active_MC[2] == false &&
521  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
522  {
523  act[2] = act[6] = true;
524  act[4] = true;
525  act[8] = true;
526  return;
527  }
528  if (active_tensile[0] == false && active_tensile[1] == false && active_tensile[2] == true &&
529  active_MC[0] == false && active_MC[1] == true && active_MC[2] == false &&
530  active_MC[3] == true && active_MC[4] == false && active_MC[5] == true)
531  {
532  act[2] = act[6] = true;
533  act[8] = true;
534  return;
535  }
536 
537  // tensile = face, MC=face
538  if (active_tensile[0] == false && active_tensile[1] == false && active_tensile[2] == true &&
539  active_MC[0] == false && active_MC[1] == false && active_MC[2] == false &&
540  active_MC[3] == true && active_MC[4] == false && active_MC[5] == false)
541  {
542  act[1] = act[2] = act[6] = true;
543  return;
544  }
545 
546  // tensile = face, MC=edge (two possibilites).
547  act[2] = true; // tensile face
548  act[3] = active_MC[0];
549  act[4] = active_MC[1];
550  act[5] = active_MC[2];
551  act[6] = active_MC[3];
552  act[7] = active_MC[4];
553  act[8] = active_MC[5];
554  return;
555 }

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints().

◆ calculateConstraints()

void MultiPlasticityLinearSystem::calculateConstraints ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankTwoTensor &  delta_dp,
std::vector< Real > &  f,
std::vector< RankTwoTensor > &  r,
RankTwoTensor &  epp,
std::vector< Real > &  ic,
const std::vector< bool > &  active 
)
protectedvirtualinherited

The constraints.

These are set to zero (or <=0 in the case of the yield functions) by the Newton-Raphson process, except in the case of linear-dependence which complicates things.

Parameters
stressThe stress
intnl_oldold values of the internal parameters
intnlinternal parameters
pmCurrent value(s) of the plasticity multiplier(s) (consistency parameters)
delta_dpChange in plastic strain incurred so far during the return
[out]fActive yield function(s)
[out]rActive flow directions
[out]eppPlastic-strain increment constraint
[out]icActive internal-parameter constraint
activeThe active constraints.

Definition at line 228 of file MultiPlasticityLinearSystem.C.

238 {
239  // see comments at the start of .h file
240 
241  mooseAssert(intnl_old.size() == _num_models,
242  "Size of intnl_old is " << intnl_old.size()
243  << " which is incorrect in calculateConstraints");
244  mooseAssert(intnl.size() == _num_models,
245  "Size of intnl is " << intnl.size() << " which is incorrect in calculateConstraints");
246  mooseAssert(pm.size() == _num_surfaces,
247  "Size of pm is " << pm.size() << " which is incorrect in calculateConstraints");
248  mooseAssert(active.size() == _num_surfaces,
249  "Size of active is " << active.size()
250  << " which is incorrect in calculateConstraints");
251 
252  // yield functions
253  yieldFunction(stress, intnl, active, f);
254 
255  // flow directions and "epp"
256  flowPotential(stress, intnl, active, r);
257  epp = RankTwoTensor();
258  unsigned ind = 0;
259  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
260  if (active[surface])
261  epp += pm[surface] * r[ind++]; // note, even the deactivated_due_to_ld must get added in
262  epp -= delta_dp;
263 
264  // internal constraints
265  std::vector<Real> h;
266  hardPotential(stress, intnl, active, h);
267  ic.resize(0);
268  ind = 0;
269  std::vector<unsigned int> active_surfaces;
270  std::vector<unsigned int>::iterator active_surface;
271  for (unsigned model = 0; model < _num_models; ++model)
272  {
273  activeSurfaces(model, active, active_surfaces);
274  if (active_surfaces.size() > 0)
275  {
276  // some surfaces are active in this model, so must form an internal constraint
277  ic.push_back(intnl[model] - intnl_old[model]);
278  for (active_surface = active_surfaces.begin(); active_surface != active_surfaces.end();
279  ++active_surface)
280  ic[ic.size() - 1] += pm[*active_surface] * h[ind++]; // we know the correct one is h[ind]
281  // since it was constructed in the same
282  // manner
283  }
284  }
285 }

Referenced by MultiPlasticityLinearSystem::calculateRHS(), and ComputeMultiPlasticityStress::lineSearch().

◆ calculateJacobian()

void MultiPlasticityLinearSystem::calculateJacobian ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankFourTensor &  E_inv,
const std::vector< bool > &  active,
const std::vector< bool > &  deactivated_due_to_ld,
std::vector< std::vector< Real >> &  jac 
)
protectedvirtualinherited

d(rhs)/d(dof)

Definition at line 367 of file MultiPlasticityLinearSystem.C.

374 {
375  // see comments at the start of .h file
376 
377  mooseAssert(intnl.size() == _num_models,
378  "Size of intnl is " << intnl.size() << " which is incorrect in calculateJacobian");
379  mooseAssert(pm.size() == _num_surfaces,
380  "Size of pm is " << pm.size() << " which is incorrect in calculateJacobian");
381  mooseAssert(active.size() == _num_surfaces,
382  "Size of active is " << active.size() << " which is incorrect in calculateJacobian");
383  mooseAssert(deactivated_due_to_ld.size() == _num_surfaces,
384  "Size of deactivated_due_to_ld is " << deactivated_due_to_ld.size()
385  << " which is incorrect in calculateJacobian");
386 
387  unsigned ind = 0;
388  unsigned active_surface_ind = 0;
389 
390  std::vector<bool> active_surface(_num_surfaces); // active and not deactivated_due_to_ld
391  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
392  active_surface[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
393  unsigned num_active_surface = 0;
394  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
395  if (active_surface[surface])
396  num_active_surface++;
397 
398  std::vector<bool> active_model(
399  _num_models); // whether a model has surfaces that are active and not deactivated_due_to_ld
400  for (unsigned model = 0; model < _num_models; ++model)
401  active_model[model] = anyActiveSurfaces(model, active_surface);
402 
403  unsigned num_active_model = 0;
404  for (unsigned model = 0; model < _num_models; ++model)
405  if (active_model[model])
406  num_active_model++;
407 
408  ind = 0;
409  std::vector<unsigned int> active_model_index(_num_models);
410  for (unsigned model = 0; model < _num_models; ++model)
411  if (active_model[model])
412  active_model_index[model] = ind++;
413  else
414  active_model_index[model] =
415  _num_models + 1; // just a dummy, that will probably cause a crash if something goes wrong
416 
417  std::vector<RankTwoTensor> df_dstress;
418  dyieldFunction_dstress(stress, intnl, active_surface, df_dstress);
419 
420  std::vector<Real> df_dintnl;
421  dyieldFunction_dintnl(stress, intnl, active_surface, df_dintnl);
422 
423  std::vector<RankTwoTensor> r;
424  flowPotential(stress, intnl, active, r);
425 
426  std::vector<RankFourTensor> dr_dstress;
427  dflowPotential_dstress(stress, intnl, active, dr_dstress);
428 
429  std::vector<RankTwoTensor> dr_dintnl;
430  dflowPotential_dintnl(stress, intnl, active, dr_dintnl);
431 
432  std::vector<Real> h;
433  hardPotential(stress, intnl, active, h);
434 
435  std::vector<RankTwoTensor> dh_dstress;
436  dhardPotential_dstress(stress, intnl, active, dh_dstress);
437 
438  std::vector<Real> dh_dintnl;
439  dhardPotential_dintnl(stress, intnl, active, dh_dintnl);
440 
441  // d(epp)/dstress = sum_{active alpha} pm[alpha]*dr_dstress
442  RankFourTensor depp_dstress;
443  ind = 0;
444  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
445  if (active[surface]) // includes deactivated_due_to_ld
446  depp_dstress += pm[surface] * dr_dstress[ind++];
447  depp_dstress += E_inv;
448 
449  // d(epp)/dpm_{active_surface_index} = r_{active_surface_index}
450  std::vector<RankTwoTensor> depp_dpm;
451  depp_dpm.resize(num_active_surface);
452  ind = 0;
453  active_surface_ind = 0;
454  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
455  {
456  if (active[surface])
457  {
458  if (active_surface[surface]) // do not include the deactived_due_to_ld, since their pm are not
459  // dofs in the NR
460  depp_dpm[active_surface_ind++] = r[ind];
461  ind++;
462  }
463  }
464 
465  // d(epp)/dintnl_{active_model_index} = sum(pm[asdf]*dr_dintnl[fdsa])
466  std::vector<RankTwoTensor> depp_dintnl;
467  depp_dintnl.assign(num_active_model, RankTwoTensor());
468  ind = 0;
469  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
470  {
471  if (active[surface])
472  {
473  unsigned int model_num = modelNumber(surface);
474  if (active_model[model_num]) // only include models with surfaces which are still active after
475  // deactivated_due_to_ld
476  depp_dintnl[active_model_index[model_num]] += pm[surface] * dr_dintnl[ind];
477  ind++;
478  }
479  }
480 
481  // df_dstress has been calculated above
482  // df_dpm is always zero
483  // df_dintnl has been calculated above, but only the active_surface+active_model stuff needs to be
484  // included in Jacobian: see below
485 
486  std::vector<RankTwoTensor> dic_dstress;
487  dic_dstress.assign(num_active_model, RankTwoTensor());
488  ind = 0;
489  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
490  {
491  if (active[surface])
492  {
493  unsigned int model_num = modelNumber(surface);
494  if (active_model[model_num]) // only include ic for models with active_surface (ie, if model
495  // only contains deactivated_due_to_ld don't include it)
496  dic_dstress[active_model_index[model_num]] += pm[surface] * dh_dstress[ind];
497  ind++;
498  }
499  }
500 
501  std::vector<std::vector<Real>> dic_dpm;
502  dic_dpm.resize(num_active_model);
503  ind = 0;
504  active_surface_ind = 0;
505  for (unsigned model = 0; model < num_active_model; ++model)
506  dic_dpm[model].assign(num_active_surface, 0);
507  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
508  {
509  if (active[surface])
510  {
511  if (active_surface[surface]) // only take derivs wrt active-but-not-deactivated_due_to_ld pm
512  {
513  unsigned int model_num = modelNumber(surface);
514  // if (active_model[model_num]) // do not need this check as if the surface has
515  // active_surface, the model must be deemed active!
516  dic_dpm[active_model_index[model_num]][active_surface_ind] = h[ind];
517  active_surface_ind++;
518  }
519  ind++;
520  }
521  }
522 
523  std::vector<std::vector<Real>> dic_dintnl;
524  dic_dintnl.resize(num_active_model);
525  for (unsigned model = 0; model < num_active_model; ++model)
526  {
527  dic_dintnl[model].assign(num_active_model, 0);
528  dic_dintnl[model][model] = 1; // deriv wrt internal parameter
529  }
530  ind = 0;
531  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
532  {
533  if (active[surface])
534  {
535  unsigned int model_num = modelNumber(surface);
536  if (active_model[model_num]) // only the models that contain surfaces that are still active
537  // after deactivation_due_to_ld
538  dic_dintnl[active_model_index[model_num]][active_model_index[model_num]] +=
539  pm[surface] * dh_dintnl[ind];
540  ind++;
541  }
542  }
543 
544  unsigned int dim = 3;
545  unsigned int system_size =
546  6 + num_active_surface + num_active_model; // "6" comes from symmeterizing epp
547  jac.resize(system_size);
548  for (unsigned i = 0; i < system_size; ++i)
549  jac[i].assign(system_size, 0);
550 
551  unsigned int row_num = 0;
552  unsigned int col_num = 0;
553  for (unsigned i = 0; i < dim; ++i)
554  for (unsigned j = 0; j <= i; ++j)
555  {
556  for (unsigned k = 0; k < dim; ++k)
557  for (unsigned l = 0; l <= k; ++l)
558  jac[col_num][row_num++] =
559  depp_dstress(i, j, k, l) +
560  (k != l ? depp_dstress(i, j, l, k)
561  : 0); // extra part is needed because i assume dstress(i, j) = dstress(j, i)
562  for (unsigned surface = 0; surface < num_active_surface; ++surface)
563  jac[col_num][row_num++] = depp_dpm[surface](i, j);
564  for (unsigned a = 0; a < num_active_model; ++a)
565  jac[col_num][row_num++] = depp_dintnl[a](i, j);
566  row_num = 0;
567  col_num++;
568  }
569 
570  ind = 0;
571  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
572  if (active_surface[surface])
573  {
574  for (unsigned k = 0; k < dim; ++k)
575  for (unsigned l = 0; l <= k; ++l)
576  jac[col_num][row_num++] =
577  df_dstress[ind](k, l) +
578  (k != l ? df_dstress[ind](l, k)
579  : 0); // extra part is needed because i assume dstress(i, j) = dstress(j, i)
580  for (unsigned beta = 0; beta < num_active_surface; ++beta)
581  jac[col_num][row_num++] = 0; // df_dpm
582  for (unsigned model = 0; model < _num_models; ++model)
583  if (active_model[model]) // only use df_dintnl for models in active_model
584  {
585  if (modelNumber(surface) == model)
586  jac[col_num][row_num++] = df_dintnl[ind];
587  else
588  jac[col_num][row_num++] = 0;
589  }
590  ind++;
591  row_num = 0;
592  col_num++;
593  }
594 
595  for (unsigned a = 0; a < num_active_model; ++a)
596  {
597  for (unsigned k = 0; k < dim; ++k)
598  for (unsigned l = 0; l <= k; ++l)
599  jac[col_num][row_num++] =
600  dic_dstress[a](k, l) +
601  (k != l ? dic_dstress[a](l, k)
602  : 0); // extra part is needed because i assume dstress(i, j) = dstress(j, i)
603  for (unsigned alpha = 0; alpha < num_active_surface; ++alpha)
604  jac[col_num][row_num++] = dic_dpm[a][alpha];
605  for (unsigned b = 0; b < num_active_model; ++b)
606  jac[col_num][row_num++] = dic_dintnl[a][b];
607  row_num = 0;
608  col_num++;
609  }
610 
611  mooseAssert(col_num == system_size, "Incorrect filling of cols in Jacobian");
612 }

Referenced by checkJacobian(), checkSolution(), and MultiPlasticityLinearSystem::nrStep().

◆ calculateRHS()

void MultiPlasticityLinearSystem::calculateRHS ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankTwoTensor &  delta_dp,
std::vector< Real > &  rhs,
const std::vector< bool > &  active,
bool  eliminate_ld,
std::vector< bool > &  deactivated_due_to_ld 
)
protectedvirtualinherited

Calculate the RHS which is rhs = -(epp(0,0), epp(1,0), epp(1,1), epp(2,0), epp(2,1), epp(2,2), f[0], f[1], ..., f[num_f], ic[0], ic[1], ..., ic[num_ic])

Note that the 'epp' components only contain the upper diagonal. These contain flow directions and plasticity-multipliers for all active surfaces, even the deactivated_due_to_ld surfaces. Note that the 'f' components only contain the active and not deactivated_due_to_ld surfaces Note that the 'ic' components only contain the internal constraints for models which contain active and not deactivated_due_to_ld surfaces. They contain hardening-potentials and plasticity-multipliers for the active surfaces, even the deactivated_due_to_ld surfaces

Parameters
stressThe stress
intnl_oldold values of the internal parameters
intnlinternal parameters
pmCurrent value(s) of the plasticity multiplier(s) (consistency parameters)
delta_dpChange in plastic strain incurred so far during the return
[out]rhsthe rhs
activeThe active constraints.
eliminate_ldCheck for linear dependence of constraints and put the results into deactivated_due_to_ld. Usually this should be true, but for certain debug operations it should be false
[out]deactivated_due_to_ldconstraints deactivated due to linear-dependence of flow directions

Definition at line 288 of file MultiPlasticityLinearSystem.C.

297 {
298  // see comments at the start of .h file
299 
300  mooseAssert(intnl_old.size() == _num_models,
301  "Size of intnl_old is " << intnl_old.size() << " which is incorrect in calculateRHS");
302  mooseAssert(intnl.size() == _num_models,
303  "Size of intnl is " << intnl.size() << " which is incorrect in calculateRHS");
304  mooseAssert(pm.size() == _num_surfaces,
305  "Size of pm is " << pm.size() << " which is incorrect in calculateRHS");
306  mooseAssert(active.size() == _num_surfaces,
307  "Size of active is " << active.size() << " which is incorrect in calculateRHS");
308 
309  std::vector<Real> f; // the yield functions
310  RankTwoTensor epp; // the plastic-strain constraint ("direction constraint")
311  std::vector<Real> ic; // the "internal constraints"
312 
313  std::vector<RankTwoTensor> r;
314  calculateConstraints(stress, intnl_old, intnl, pm, delta_dp, f, r, epp, ic, active);
315 
316  if (eliminate_ld)
317  eliminateLinearDependence(stress, intnl, f, r, active, deactivated_due_to_ld);
318  else
319  deactivated_due_to_ld.assign(_num_surfaces, false);
320 
321  std::vector<bool> active_not_deact(_num_surfaces);
322  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
323  active_not_deact[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
324 
325  unsigned num_active_f = 0;
326  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
327  if (active_not_deact[surface])
328  num_active_f++;
329 
330  unsigned num_active_ic = 0;
331  for (unsigned model = 0; model < _num_models; ++model)
332  if (anyActiveSurfaces(model, active_not_deact))
333  num_active_ic++;
334 
335  unsigned int dim = 3;
336  unsigned int system_size = 6 + num_active_f + num_active_ic; // "6" comes from symmeterizing epp,
337  // num_active_f comes from "f",
338  // num_active_f comes from "ic"
339 
340  rhs.resize(system_size);
341 
342  unsigned ind = 0;
343  for (unsigned i = 0; i < dim; ++i)
344  for (unsigned j = 0; j <= i; ++j)
345  rhs[ind++] = -epp(i, j);
346  unsigned active_surface = 0;
347  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
348  if (active[surface])
349  {
350  if (!deactivated_due_to_ld[surface])
351  rhs[ind++] = -f[active_surface];
352  active_surface++;
353  }
354  unsigned active_model = 0;
355  for (unsigned model = 0; model < _num_models; ++model)
356  if (anyActiveSurfaces(model, active))
357  {
358  if (anyActiveSurfaces(model, active_not_deact))
359  rhs[ind++] = -ic[active_model];
360  active_model++;
361  }
362 
363  mooseAssert(ind == system_size, "Incorrect filling of the rhs in calculateRHS");
364 }

Referenced by checkSolution(), fdJacobian(), and MultiPlasticityLinearSystem::nrStep().

◆ checkDerivatives()

void MultiPlasticityDebugger::checkDerivatives ( )

Checks the derivatives, eg dyieldFunction_dstress by using finite difference approximations.

Definition at line 85 of file MultiPlasticityDebugger.C.

86 {
87  Moose::err
88  << "\n\n++++++++++++++++++++++++\nChecking the derivatives\n++++++++++++++++++++++++\n";
89  outputAndCheckDebugParameters();
90 
91  std::vector<bool> act;
92  act.assign(_num_surfaces, true);
93 
94  Moose::err << "\ndyieldFunction_dstress. Relative L2 norms.\n";
95  std::vector<RankTwoTensor> df_dstress;
96  std::vector<RankTwoTensor> fddf_dstress;
97  dyieldFunction_dstress(_fspb_debug_stress, _fspb_debug_intnl, act, df_dstress);
98  fddyieldFunction_dstress(_fspb_debug_stress, _fspb_debug_intnl, fddf_dstress);
99  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
100  {
101  Moose::err << "surface = " << surface << " Relative L2norm = "
102  << 2 * (df_dstress[surface] - fddf_dstress[surface]).L2norm() /
103  (df_dstress[surface] + fddf_dstress[surface]).L2norm()
104  << "\n";
105  Moose::err << "Coded:\n";
106  df_dstress[surface].print();
107  Moose::err << "Finite difference:\n";
108  fddf_dstress[surface].print();
109  }
110 
111  Moose::err << "\ndyieldFunction_dintnl.\n";
112  std::vector<Real> df_dintnl;
113  dyieldFunction_dintnl(_fspb_debug_stress, _fspb_debug_intnl, act, df_dintnl);
114  Moose::err << "Coded:\n";
115  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
116  Moose::err << df_dintnl[surface] << " ";
117  Moose::err << "\n";
118  std::vector<Real> fddf_dintnl;
119  fddyieldFunction_dintnl(_fspb_debug_stress, _fspb_debug_intnl, fddf_dintnl);
120  Moose::err << "Finite difference:\n";
121  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
122  Moose::err << fddf_dintnl[surface] << " ";
123  Moose::err << "\n";
124 
125  Moose::err << "\ndflowPotential_dstress. Relative L2 norms.\n";
126  std::vector<RankFourTensor> dr_dstress;
127  std::vector<RankFourTensor> fddr_dstress;
128  dflowPotential_dstress(_fspb_debug_stress, _fspb_debug_intnl, act, dr_dstress);
129  fddflowPotential_dstress(_fspb_debug_stress, _fspb_debug_intnl, fddr_dstress);
130  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
131  {
132  Moose::err << "surface = " << surface << " Relative L2norm = "
133  << 2 * (dr_dstress[surface] - fddr_dstress[surface]).L2norm() /
134  (dr_dstress[surface] + fddr_dstress[surface]).L2norm()
135  << "\n";
136  Moose::err << "Coded:\n";
137  dr_dstress[surface].print();
138  Moose::err << "Finite difference:\n";
139  fddr_dstress[surface].print();
140  }
141 
142  Moose::err << "\ndflowPotential_dintnl. Relative L2 norms.\n";
143  std::vector<RankTwoTensor> dr_dintnl;
144  std::vector<RankTwoTensor> fddr_dintnl;
145  dflowPotential_dintnl(_fspb_debug_stress, _fspb_debug_intnl, act, dr_dintnl);
146  fddflowPotential_dintnl(_fspb_debug_stress, _fspb_debug_intnl, fddr_dintnl);
147  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
148  {
149  Moose::err << "surface = " << surface << " Relative L2norm = "
150  << 2 * (dr_dintnl[surface] - fddr_dintnl[surface]).L2norm() /
151  (dr_dintnl[surface] + fddr_dintnl[surface]).L2norm()
152  << "\n";
153  Moose::err << "Coded:\n";
154  dr_dintnl[surface].print();
155  Moose::err << "Finite difference:\n";
156  fddr_dintnl[surface].print();
157  }
158 }

Referenced by ComputeMultiPlasticityStress::initQpStatefulProperties(), and ComputeMultiPlasticityStress::plasticStep().

◆ checkJacobian()

void MultiPlasticityDebugger::checkJacobian ( const RankFourTensor &  E_inv,
const std::vector< Real > &  intnl_old 
)

Checks the full Jacobian, which is just certain linear combinations of the dyieldFunction_dstress, etc, by using finite difference approximations.

Definition at line 161 of file MultiPlasticityDebugger.C.

163 {
164  Moose::err << "\n\n+++++++++++++++++++++\nChecking the Jacobian\n+++++++++++++++++++++\n";
165  outputAndCheckDebugParameters();
166 
167  std::vector<bool> act;
168  act.assign(_num_surfaces, true);
169  std::vector<bool> deactivated_due_to_ld;
170  deactivated_due_to_ld.assign(_num_surfaces, false);
171 
172  RankTwoTensor delta_dp = -E_inv * _fspb_debug_stress;
173 
174  std::vector<std::vector<Real>> jac;
175  calculateJacobian(_fspb_debug_stress,
176  _fspb_debug_intnl,
177  _fspb_debug_pm,
178  E_inv,
179  act,
180  deactivated_due_to_ld,
181  jac);
182 
183  std::vector<std::vector<Real>> fdjac;
184  fdJacobian(_fspb_debug_stress,
185  intnl_old,
186  _fspb_debug_intnl,
187  _fspb_debug_pm,
188  delta_dp,
189  E_inv,
190  false,
191  fdjac);
192 
193  Real L2_numer = 0;
194  Real L2_denom = 0;
195  for (unsigned row = 0; row < jac.size(); ++row)
196  for (unsigned col = 0; col < jac.size(); ++col)
197  {
198  L2_numer += Utility::pow<2>(jac[row][col] - fdjac[row][col]);
199  L2_denom += Utility::pow<2>(jac[row][col] + fdjac[row][col]);
200  }
201  Moose::err << "\nRelative L2norm = " << std::sqrt(L2_numer / L2_denom) / 0.5 << "\n";
202 
203  Moose::err << "\nHand-coded Jacobian:\n";
204  for (unsigned row = 0; row < jac.size(); ++row)
205  {
206  for (unsigned col = 0; col < jac.size(); ++col)
207  Moose::err << jac[row][col] << " ";
208  Moose::err << "\n";
209  }
210 
211  Moose::err << "Finite difference Jacobian:\n";
212  for (unsigned row = 0; row < fdjac.size(); ++row)
213  {
214  for (unsigned col = 0; col < fdjac.size(); ++col)
215  Moose::err << fdjac[row][col] << " ";
216  Moose::err << "\n";
217  }
218 }

Referenced by ComputeMultiPlasticityStress::computeQpStress(), and ComputeMultiPlasticityStress::plasticStep().

◆ checkSolution()

void MultiPlasticityDebugger::checkSolution ( const RankFourTensor &  E_inv)

Checks that Ax does equal b in the NR procedure.

Definition at line 367 of file MultiPlasticityDebugger.C.

368 {
369  Moose::err << "\n\n+++++++++++++++++++++\nChecking the Solution\n";
370  Moose::err << "(Ie, checking Ax = b)\n+++++++++++++++++++++\n";
371  outputAndCheckDebugParameters();
372 
373  std::vector<bool> act;
374  act.assign(_num_surfaces, true);
375  std::vector<bool> deactivated_due_to_ld;
376  deactivated_due_to_ld.assign(_num_surfaces, false);
377 
378  RankTwoTensor delta_dp = -E_inv * _fspb_debug_stress;
379 
380  std::vector<Real> orig_rhs;
381  calculateRHS(_fspb_debug_stress,
382  _fspb_debug_intnl,
383  _fspb_debug_intnl,
384  _fspb_debug_pm,
385  delta_dp,
386  orig_rhs,
387  act,
388  true,
389  deactivated_due_to_ld);
390 
391  Moose::err << "\nb = ";
392  for (unsigned i = 0; i < orig_rhs.size(); ++i)
393  Moose::err << orig_rhs[i] << " ";
394  Moose::err << "\n\n";
395 
396  std::vector<std::vector<Real>> jac_coded;
397  calculateJacobian(_fspb_debug_stress,
398  _fspb_debug_intnl,
399  _fspb_debug_pm,
400  E_inv,
401  act,
402  deactivated_due_to_ld,
403  jac_coded);
404 
405  Moose::err
406  << "Before checking Ax=b is correct, check that the Jacobians given below are equal.\n";
407  Moose::err
408  << "The hand-coded Jacobian is used in calculating the solution 'x', given 'b' above.\n";
409  Moose::err << "Note that this only includes degrees of freedom that aren't deactivated due to "
410  "linear dependence.\n";
411  Moose::err << "Hand-coded Jacobian:\n";
412  for (unsigned row = 0; row < jac_coded.size(); ++row)
413  {
414  for (unsigned col = 0; col < jac_coded.size(); ++col)
415  Moose::err << jac_coded[row][col] << " ";
416  Moose::err << "\n";
417  }
418 
419  deactivated_due_to_ld.assign(_num_surfaces,
420  false); // this potentially gets changed by nrStep, below
421  RankTwoTensor dstress;
422  std::vector<Real> dpm;
423  std::vector<Real> dintnl;
424  nrStep(_fspb_debug_stress,
425  _fspb_debug_intnl,
426  _fspb_debug_intnl,
427  _fspb_debug_pm,
428  E_inv,
429  delta_dp,
430  dstress,
431  dpm,
432  dintnl,
433  act,
434  deactivated_due_to_ld);
435 
436  std::vector<bool> active_not_deact(_num_surfaces);
437  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
438  active_not_deact[surface] = !deactivated_due_to_ld[surface];
439 
440  std::vector<Real> x;
441  x.assign(orig_rhs.size(), 0);
442  unsigned ind = 0;
443  for (unsigned i = 0; i < 3; ++i)
444  for (unsigned j = 0; j <= i; ++j)
445  x[ind++] = dstress(i, j);
446  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
447  if (active_not_deact[surface])
448  x[ind++] = dpm[surface];
449  for (unsigned model = 0; model < _num_models; ++model)
450  if (anyActiveSurfaces(model, active_not_deact))
451  x[ind++] = dintnl[model];
452 
453  mooseAssert(ind == orig_rhs.size(),
454  "Incorrect extracting of changes from NR solution in the "
455  "finite-difference checking of nrStep");
456 
457  Moose::err << "\nThis yields x =";
458  for (unsigned i = 0; i < orig_rhs.size(); ++i)
459  Moose::err << x[i] << " ";
460  Moose::err << "\n";
461 
462  std::vector<std::vector<Real>> jac_fd;
463  fdJacobian(_fspb_debug_stress,
464  _fspb_debug_intnl,
465  _fspb_debug_intnl,
466  _fspb_debug_pm,
467  delta_dp,
468  E_inv,
469  true,
470  jac_fd);
471 
472  Moose::err << "\nThe finite-difference Jacobian is used to multiply by this 'x',\n";
473  Moose::err << "in order to check that the solution is correct\n";
474  Moose::err << "Finite-difference Jacobian:\n";
475  for (unsigned row = 0; row < jac_fd.size(); ++row)
476  {
477  for (unsigned col = 0; col < jac_fd.size(); ++col)
478  Moose::err << jac_fd[row][col] << " ";
479  Moose::err << "\n";
480  }
481 
482  Real L2_numer = 0;
483  Real L2_denom = 0;
484  for (unsigned row = 0; row < jac_coded.size(); ++row)
485  for (unsigned col = 0; col < jac_coded.size(); ++col)
486  {
487  L2_numer += Utility::pow<2>(jac_coded[row][col] - jac_fd[row][col]);
488  L2_denom += Utility::pow<2>(jac_coded[row][col] + jac_fd[row][col]);
489  }
490  Moose::err << "Relative L2norm of the hand-coded and finite-difference Jacobian is "
491  << std::sqrt(L2_numer / L2_denom) / 0.5 << "\n";
492 
493  std::vector<Real> fd_times_x;
494  fd_times_x.assign(orig_rhs.size(), 0);
495  for (unsigned row = 0; row < orig_rhs.size(); ++row)
496  for (unsigned col = 0; col < orig_rhs.size(); ++col)
497  fd_times_x[row] += jac_fd[row][col] * x[col];
498 
499  Moose::err << "\n(Finite-difference Jacobian)*x =\n";
500  for (unsigned i = 0; i < orig_rhs.size(); ++i)
501  Moose::err << fd_times_x[i] << " ";
502  Moose::err << "\n";
503  Moose::err << "Recall that b = \n";
504  for (unsigned i = 0; i < orig_rhs.size(); ++i)
505  Moose::err << orig_rhs[i] << " ";
506  Moose::err << "\n";
507 
508  L2_numer = 0;
509  L2_denom = 0;
510  for (unsigned i = 0; i < orig_rhs.size(); ++i)
511  {
512  L2_numer += Utility::pow<2>(orig_rhs[i] - fd_times_x[i]);
513  L2_denom += Utility::pow<2>(orig_rhs[i] + fd_times_x[i]);
514  }
515  Moose::err << "\nRelative L2norm of these is " << std::sqrt(L2_numer / L2_denom) / 0.5 << "\n";
516 }

Referenced by ComputeMultiPlasticityStress::computeQpStress(), and ComputeMultiPlasticityStress::plasticStep().

◆ dflowPotential_dintnl()

void MultiPlasticityRawComponentAssembler::dflowPotential_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  dr_dintnl 
)
protectedvirtualinherited

The derivative of the active flow potentials with respect to the active internal parameters The UserObjects explicitly assume that r[alpha] is only dependent on intnl[alpha].

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dr_dintnl"
[out]dr_dintnlthe derivatives. dr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl[alpha]

Definition at line 235 of file MultiPlasticityRawComponentAssembler.C.

239 {
240  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
241  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
242 
243  dr_dintnl.resize(0);
244  std::vector<unsigned int> active_surfaces_of_model;
245  std::vector<unsigned int>::iterator active_surface;
246  std::vector<RankTwoTensor> model_dr_dintnl;
247  for (unsigned model = 0; model < _num_models; ++model)
248  {
249  activeModelSurfaces(model, active, active_surfaces_of_model);
250  if (active_surfaces_of_model.size() > 0)
251  {
252  _f[model]->dflowPotential_dintnlV(stress, intnl[model], model_dr_dintnl);
253  for (active_surface = active_surfaces_of_model.begin();
254  active_surface != active_surfaces_of_model.end();
255  ++active_surface)
256  dr_dintnl.push_back(model_dr_dintnl[*active_surface]);
257  }
258  }
259 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), checkDerivatives(), and ComputeMultiPlasticityStress::consistentTangentOperator().

◆ dflowPotential_dstress()

void MultiPlasticityRawComponentAssembler::dflowPotential_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankFourTensor > &  dr_dstress 
)
protectedvirtualinherited

The derivative of the active flow potential(s) with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dr_dstress"
[out]dr_dstressthe derivative. dr_dstress[alpha](i, j, k, l) = dr[alpha](i, j)/dstress(k, l)

Definition at line 207 of file MultiPlasticityRawComponentAssembler.C.

212 {
213  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
214  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
215 
216  dr_dstress.resize(0);
217  std::vector<unsigned int> active_surfaces_of_model;
218  std::vector<unsigned int>::iterator active_surface;
219  std::vector<RankFourTensor> model_dr_dstress;
220  for (unsigned model = 0; model < _num_models; ++model)
221  {
222  activeModelSurfaces(model, active, active_surfaces_of_model);
223  if (active_surfaces_of_model.size() > 0)
224  {
225  _f[model]->dflowPotential_dstressV(stress, intnl[model], model_dr_dstress);
226  for (active_surface = active_surfaces_of_model.begin();
227  active_surface != active_surfaces_of_model.end();
228  ++active_surface)
229  dr_dstress.push_back(model_dr_dstress[*active_surface]);
230  }
231  }
232 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), checkDerivatives(), and ComputeMultiPlasticityStress::consistentTangentOperator().

◆ dhardPotential_dintnl()

void MultiPlasticityRawComponentAssembler::dhardPotential_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  dh_dintnl 
)
protectedvirtualinherited

The derivative of the active hardening potentials with respect to the active internal parameters.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dh_dintnl"
[out]dh_dintnlthe derivatives. dh_dintnl[a][alpha][b] = dh[a][alpha]/dintnl[b]. Note that the userobjects assume that there is exactly one internal parameter per yield function, so the derivative is only nonzero for a=alpha=b, so that is all we calculate

Definition at line 317 of file MultiPlasticityRawComponentAssembler.C.

321 {
322  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
323  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
324 
325  dh_dintnl.resize(0);
326  std::vector<unsigned int> active_surfaces_of_model;
327  std::vector<unsigned int>::iterator active_surface;
328  std::vector<Real> model_dh_dintnl;
329  for (unsigned model = 0; model < _num_models; ++model)
330  {
331  activeModelSurfaces(model, active, active_surfaces_of_model);
332  if (active_surfaces_of_model.size() > 0)
333  {
334  _f[model]->dhardPotential_dintnlV(stress, intnl[model], model_dh_dintnl);
335  for (active_surface = active_surfaces_of_model.begin();
336  active_surface != active_surfaces_of_model.end();
337  ++active_surface)
338  dh_dintnl.push_back(model_dh_dintnl[*active_surface]);
339  }
340  }
341 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian().

◆ dhardPotential_dstress()

void MultiPlasticityRawComponentAssembler::dhardPotential_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  dh_dstress 
)
protectedvirtualinherited

The derivative of the active hardening potentials with respect to stress By assumption in the Userobjects, the h[a][alpha] is nonzero only for a = alpha, so we only calculate those here.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "dh_dstress"
[out]dh_dstressthe derivative. dh_dstress[a](i, j) = dh[a]/dstress(k, l)

Definition at line 289 of file MultiPlasticityRawComponentAssembler.C.

294 {
295  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
296  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
297 
298  dh_dstress.resize(0);
299  std::vector<unsigned int> active_surfaces_of_model;
300  std::vector<unsigned int>::iterator active_surface;
301  std::vector<RankTwoTensor> model_dh_dstress;
302  for (unsigned model = 0; model < _num_models; ++model)
303  {
304  activeModelSurfaces(model, active, active_surfaces_of_model);
305  if (active_surfaces_of_model.size() > 0)
306  {
307  _f[model]->dhardPotential_dstressV(stress, intnl[model], model_dh_dstress);
308  for (active_surface = active_surfaces_of_model.begin();
309  active_surface != active_surfaces_of_model.end();
310  ++active_surface)
311  dh_dstress.push_back(model_dh_dstress[*active_surface]);
312  }
313  }
314 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian().

◆ dof_included()

bool MultiPlasticityDebugger::dof_included ( unsigned int  dof,
const std::vector< bool > &  deactivated_due_to_ld 
)
private

Definition at line 349 of file MultiPlasticityDebugger.C.

351 {
352  if (dof < unsigned(6))
353  // these are the stress components
354  return true;
355  unsigned eff_dof = dof - 6;
356  if (eff_dof < _num_surfaces)
357  // these are the plastic multipliers, pm
358  return !deactivated_due_to_ld[eff_dof];
359  eff_dof -= _num_surfaces; // now we know the dof is an intnl parameter
360  std::vector<bool> active_surface(_num_surfaces);
361  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
362  active_surface[surface] = !deactivated_due_to_ld[surface];
363  return anyActiveSurfaces(eff_dof, active_surface);
364 }

Referenced by fdJacobian().

◆ dyieldFunction_dintnl()

void MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  df_dintnl 
)
protectedvirtualinherited

The derivative of active yield function(s) with respect to their internal parameters (the user objects assume there is exactly one internal param per yield function)

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "df_dintnl"
[out]df_dintnlthe derivatives. df_dstress[alpha] = dyieldFunction[alpha]/dintnl[alpha]

Definition at line 153 of file MultiPlasticityRawComponentAssembler.C.

157 {
158  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
159  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
160 
161  df_dintnl.resize(0);
162  std::vector<unsigned int> active_surfaces_of_model;
163  std::vector<unsigned int>::iterator active_surface;
164  std::vector<Real> model_df_dintnl;
165  for (unsigned model = 0; model < _num_models; ++model)
166  {
167  activeModelSurfaces(model, active, active_surfaces_of_model);
168  if (active_surfaces_of_model.size() > 0)
169  {
170  _f[model]->dyieldFunction_dintnlV(stress, intnl[model], model_df_dintnl);
171  for (active_surface = active_surfaces_of_model.begin();
172  active_surface != active_surfaces_of_model.end();
173  ++active_surface)
174  df_dintnl.push_back(model_df_dintnl[*active_surface]);
175  }
176  }
177 }

Referenced by MultiPlasticityLinearSystem::calculateJacobian(), checkDerivatives(), and ComputeMultiPlasticityStress::consistentTangentOperator().

◆ dyieldFunction_dstress()

void MultiPlasticityRawComponentAssembler::dyieldFunction_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  df_dstress 
)
protectedvirtualinherited

The derivative of the active yield function(s) with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
activeset of active constraints - only the active derivatives are put into "df_dstress"
[out]df_dstressthe derivative (or derivatives in the case of multisurface plasticity). df_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)

Definition at line 125 of file MultiPlasticityRawComponentAssembler.C.

130 {
131  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
132  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
133 
134  df_dstress.resize(0);
135  std::vector<unsigned int> active_surfaces_of_model;
136  std::vector<unsigned int>::iterator active_surface;
137  std::vector<RankTwoTensor> model_df_dstress;
138  for (unsigned model = 0; model < _num_models; ++model)
139  {
140  activeModelSurfaces(model, active, active_surfaces_of_model);
141  if (active_surfaces_of_model.size() > 0)
142  {
143  _f[model]->dyieldFunction_dstressV(stress, intnl[model], model_df_dstress);
144  for (active_surface = active_surfaces_of_model.begin();
145  active_surface != active_surfaces_of_model.end();
146  ++active_surface)
147  df_dstress.push_back(model_df_dstress[*active_surface]);
148  }
149  }
150 }

Referenced by ComputeMultiPlasticityStress::buildDumbOrder(), MultiPlasticityLinearSystem::calculateJacobian(), checkDerivatives(), ComputeMultiPlasticityStress::consistentTangentOperator(), and MultiPlasticityLinearSystem::eliminateLinearDependence().

◆ eliminateLinearDependence()

void MultiPlasticityLinearSystem::eliminateLinearDependence ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< Real > &  f,
const std::vector< RankTwoTensor > &  r,
const std::vector< bool > &  active,
std::vector< bool > &  deactivated_due_to_ld 
)
privatevirtualinherited

Performs a number of singular-value decompositions to check for linear-dependence of the active directions "r" If linear dependence is found, then deactivated_due_to_ld will contain 'true' entries where surfaces need to be deactivated_due_to_ld.

Parameters
stressthe current stress
intnlthe current values of internal parameters
fActive yield function values
rthe flow directions that for those yield functions that are active upon entry to this function
activetrue if active
[out]deactivated_due_to_ldYield functions deactivated due to linearly-dependent flow directions

Definition at line 107 of file MultiPlasticityLinearSystem.C.

113 {
114  deactivated_due_to_ld.resize(_num_surfaces, false);
115 
116  unsigned num_active = r.size();
117 
118  if (num_active <= 1)
119  return;
120 
121  std::vector<double> s;
122  int info = singularValuesOfR(r, s);
123  if (info != 0)
124  mooseError("In finding the SVD in the return-map algorithm, the PETSC LAPACK gesvd routine "
125  "returned with error code ",
126  info);
127 
128  // num_lin_dep are the number of linearly dependent
129  // "r vectors", if num_active <= 6
130  unsigned int num_lin_dep = 0;
131 
132  unsigned i = s.size();
133  while (i-- > 0)
134  if (s[i] < _svd_tol * s[0])
135  num_lin_dep++;
136  else
137  break;
138 
139  if (num_lin_dep == 0 && num_active <= 6)
140  return;
141 
142  // From here on, some flow directions are linearly dependent
143 
144  // Find the signed "distance" of the current (stress, internal) configuration
145  // from the yield surfaces. This distance will not be precise, but
146  // i want to preferentially deactivate yield surfaces that are close
147  // to the current stress point.
148  std::vector<RankTwoTensor> df_dstress;
149  dyieldFunction_dstress(stress, intnl, active, df_dstress);
150 
151  typedef std::pair<Real, unsigned> pair_for_sorting;
152  std::vector<pair_for_sorting> dist(num_active);
153  for (unsigned i = 0; i < num_active; ++i)
154  {
155  dist[i].first = f[i] / df_dstress[i].L2norm();
156  dist[i].second = i;
157  }
158  std::sort(dist.begin(), dist.end()); // sorted in ascending order
159 
160  // There is a potential problem when we have equal f[i], for it can give oscillations
161  bool equals_detected = false;
162  for (unsigned i = 0; i < num_active - 1; ++i)
163  if (std::abs(dist[i].first - dist[i + 1].first) < _min_f_tol)
164  {
165  equals_detected = true;
166  dist[i].first += _min_f_tol * (MooseRandom::rand() - 0.5);
167  }
168  if (equals_detected)
169  std::sort(dist.begin(), dist.end()); // sorted in ascending order
170 
171  std::vector<bool> scheduled_for_deactivation;
172  scheduled_for_deactivation.assign(num_active, false);
173 
174  // In the following loop we go through all the flow directions, from
175  // the one with the largest dist, to the one with the smallest dist,
176  // adding them one-by-one into r_tmp. Upon each addition we check
177  // for linear-dependence. if LD is found, we schedule the most
178  // recently added flow direction for deactivation, and pop it
179  // back off r_tmp
180  unsigned current_yf;
181  current_yf = dist[num_active - 1].second;
182  // the one with largest dist
183  std::vector<RankTwoTensor> r_tmp = {r[current_yf]};
184 
185  unsigned num_kept_active = 1;
186  for (unsigned yf_to_try = 2; yf_to_try <= num_active; ++yf_to_try)
187  {
188  current_yf = dist[num_active - yf_to_try].second;
189  if (num_active == 2) // shortcut to we don't have to singularValuesOfR
190  scheduled_for_deactivation[current_yf] = true;
191  else if (num_kept_active >= 6) // shortcut to we don't have to singularValuesOfR: there can
192  // never be > 6 linearly-independent r vectors
193  scheduled_for_deactivation[current_yf] = true;
194  else
195  {
196  r_tmp.push_back(r[current_yf]);
197  info = singularValuesOfR(r_tmp, s);
198  if (info != 0)
199  mooseError("In finding the SVD in the return-map algorithm, the PETSC LAPACK gesvd routine "
200  "returned with error code ",
201  info);
202  if (s[s.size() - 1] < _svd_tol * s[0])
203  {
204  scheduled_for_deactivation[current_yf] = true;
205  r_tmp.pop_back();
206  num_lin_dep--;
207  }
208  else
209  num_kept_active++;
210  if (num_lin_dep == 0 && num_active <= 6)
211  // have taken out all the vectors that were linearly dependent
212  // so no point continuing
213  break;
214  }
215  }
216 
217  unsigned int old_active_number = 0;
218  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
219  if (active[surface])
220  {
221  if (scheduled_for_deactivation[old_active_number])
222  deactivated_due_to_ld[surface] = true;
223  old_active_number++;
224  }
225 }

Referenced by MultiPlasticityLinearSystem::calculateRHS().

◆ fddflowPotential_dintnl()

void MultiPlasticityDebugger::fddflowPotential_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< RankTwoTensor > &  dr_dintnl 
)
privatevirtual

The finite-difference derivative of the flow potentials with respect to internal parameters.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
[out]dr_dintnlthe derivatives. dr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl[alpha]

Definition at line 607 of file MultiPlasticityDebugger.C.

610 {
611  dr_dintnl.resize(_num_surfaces);
612 
613  std::vector<bool> act;
614  act.assign(_num_surfaces, true);
615 
616  std::vector<RankTwoTensor> origr;
617  flowPotential(stress, intnl, act, origr);
618 
619  std::vector<Real> intnlep;
620  intnlep.resize(_num_models);
621  for (unsigned model = 0; model < _num_models; ++model)
622  intnlep[model] = intnl[model];
623  Real ep;
624  std::vector<RankTwoTensor> rep;
625  unsigned int model;
626  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
627  {
628  model = modelNumber(surface);
629  ep = _fspb_debug_intnl_change[model];
630  intnlep[model] += ep;
631  flowPotential(stress, intnlep, act, rep);
632  dr_dintnl[surface] = (rep[surface] - origr[surface]) / ep;
633  intnlep[model] -= ep;
634  }
635 }

Referenced by checkDerivatives().

◆ fddflowPotential_dstress()

void MultiPlasticityDebugger::fddflowPotential_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< RankFourTensor > &  dr_dstress 
)
privatevirtual

The finite-difference derivative of the flow potential(s) with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
[out]dr_dstressthe derivative. dr_dstress[alpha](i, j, k, l) = dr[alpha](i, j)/dstress(k, l)

Definition at line 578 of file MultiPlasticityDebugger.C.

581 {
582  dr_dstress.assign(_num_surfaces, RankFourTensor());
583 
584  std::vector<bool> act;
585  act.assign(_num_surfaces, true);
586 
587  Real ep = _fspb_debug_stress_change;
588  RankTwoTensor stressep;
589  std::vector<RankTwoTensor> rep, rep_minus;
590  for (unsigned i = 0; i < 3; ++i)
591  for (unsigned j = 0; j < 3; ++j)
592  {
593  stressep = stress;
594  // do a central difference
595  stressep(i, j) += ep / 2.0;
596  flowPotential(stressep, intnl, act, rep);
597  stressep(i, j) -= ep;
598  flowPotential(stressep, intnl, act, rep_minus);
599  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
600  for (unsigned k = 0; k < 3; ++k)
601  for (unsigned l = 0; l < 3; ++l)
602  dr_dstress[surface](k, l, i, j) = (rep[surface](k, l) - rep_minus[surface](k, l)) / ep;
603  }
604 }

Referenced by checkDerivatives().

◆ fddyieldFunction_dintnl()

void MultiPlasticityDebugger::fddyieldFunction_dintnl ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< Real > &  df_dintnl 
)
private

The finite-difference derivative of yield function(s) with respect to internal parameter(s)

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
[out]df_dintnlthe derivative (or derivatives in the case of multisurface plasticity). df_dintnl[alpha] = dyieldFunction[alpha]/dintnl[alpha]

Definition at line 547 of file MultiPlasticityDebugger.C.

550 {
551  df_dintnl.resize(_num_surfaces);
552 
553  std::vector<bool> act;
554  act.assign(_num_surfaces, true);
555 
556  std::vector<Real> origf;
557  yieldFunction(stress, intnl, act, origf);
558 
559  std::vector<Real> intnlep;
560  intnlep.resize(_num_models);
561  for (unsigned model = 0; model < _num_models; ++model)
562  intnlep[model] = intnl[model];
563  Real ep;
564  std::vector<Real> fep;
565  unsigned int model;
566  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
567  {
568  model = modelNumber(surface);
569  ep = _fspb_debug_intnl_change[model];
570  intnlep[model] += ep;
571  yieldFunction(stress, intnlep, act, fep);
572  df_dintnl[surface] = (fep[surface] - origf[surface]) / ep;
573  intnlep[model] -= ep;
574  }
575 }

Referenced by checkDerivatives().

◆ fddyieldFunction_dstress()

void MultiPlasticityDebugger::fddyieldFunction_dstress ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
std::vector< RankTwoTensor > &  df_dstress 
)
private

The finite-difference derivative of yield function(s) with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
[out]df_dstressthe derivative (or derivatives in the case of multisurface plasticity). df_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)

Definition at line 519 of file MultiPlasticityDebugger.C.

522 {
523  df_dstress.assign(_num_surfaces, RankTwoTensor());
524 
525  std::vector<bool> act;
526  act.assign(_num_surfaces, true);
527 
528  Real ep = _fspb_debug_stress_change;
529  RankTwoTensor stressep;
530  std::vector<Real> fep, fep_minus;
531  for (unsigned i = 0; i < 3; ++i)
532  for (unsigned j = 0; j < 3; ++j)
533  {
534  stressep = stress;
535  // do a central difference to attempt to capture discontinuities
536  // such as those encountered in tensile and Mohr-Coulomb
537  stressep(i, j) += ep / 2.0;
538  yieldFunction(stressep, intnl, act, fep);
539  stressep(i, j) -= ep;
540  yieldFunction(stressep, intnl, act, fep_minus);
541  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
542  df_dstress[surface](i, j) = (fep[surface] - fep_minus[surface]) / ep;
543  }
544 }

Referenced by checkDerivatives().

◆ fdJacobian()

void MultiPlasticityDebugger::fdJacobian ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankTwoTensor &  delta_dp,
const RankFourTensor &  E_inv,
bool  eliminate_ld,
std::vector< std::vector< Real >> &  jac 
)
privatevirtual

The Jacobian calculated using finite differences.

The output should be equal to calculateJacobian(...) if everything is coded correctly.

Parameters
stressthe stress at which to calculate the Jacobian
intnl_oldthe old values of internal variables (jacobian is inependent of these, but they are needed to do the finite-differencing cleanly)
intnlthe vector of internal parameters at which to calculate the Jacobian
pmthe plasticity multipliers at which to calculate the Jacobian
delta_dpplastic_strain - plastic_strain_old (Jacobian is independent of this, but it is needed to do the finite-differencing cleanly)
E_invinverse of the elasticity tensor
eliminate_ldonly calculate the Jacobian for the linearly independent constraints
[out]jacthe finite-difference Jacobian

Definition at line 221 of file MultiPlasticityDebugger.C.

229 {
230  std::vector<bool> active;
231  active.assign(_num_surfaces, true);
232 
233  std::vector<bool> deactivated_due_to_ld;
234  std::vector<bool> deactivated_due_to_ld_ep;
235 
236  std::vector<Real> orig_rhs;
237  calculateRHS(stress,
238  intnl_old,
239  intnl,
240  pm,
241  delta_dp,
242  orig_rhs,
243  active,
244  eliminate_ld,
245  deactivated_due_to_ld); // this calculates RHS, and also set deactivated_due_to_ld.
246  // The latter stays fixed for the rest of this routine
247 
248  unsigned int whole_system_size = 6 + _num_surfaces + _num_models;
249  unsigned int system_size =
250  orig_rhs.size(); // will be = whole_system_size if eliminate_ld = false, since all active=true
251  jac.resize(system_size);
252  for (unsigned row = 0; row < system_size; ++row)
253  jac[row].assign(system_size, 0);
254 
255  std::vector<Real> rhs_ep;
256  unsigned col = 0;
257 
258  RankTwoTensor stressep;
259  RankTwoTensor delta_dpep;
260  Real ep = _fspb_debug_stress_change;
261  for (unsigned i = 0; i < 3; ++i)
262  for (unsigned j = 0; j <= i; ++j)
263  {
264  stressep = stress;
265  stressep(i, j) += ep;
266  if (i != j)
267  stressep(j, i) += ep;
268  delta_dpep = delta_dp;
269  for (unsigned k = 0; k < 3; ++k)
270  for (unsigned l = 0; l < 3; ++l)
271  {
272  delta_dpep(k, l) -= E_inv(k, l, i, j) * ep;
273  if (i != j)
274  delta_dpep(k, l) -= E_inv(k, l, j, i) * ep;
275  }
276  active.assign(_num_surfaces, true);
277  calculateRHS(stressep,
278  intnl_old,
279  intnl,
280  pm,
281  delta_dpep,
282  rhs_ep,
283  active,
284  false,
285  deactivated_due_to_ld_ep);
286  unsigned row = 0;
287  for (unsigned dof = 0; dof < whole_system_size; ++dof)
288  if (dof_included(dof, deactivated_due_to_ld))
289  {
290  jac[row][col] =
291  -(rhs_ep[dof] - orig_rhs[row]) / ep; // remember jacobian = -d(rhs)/d(something)
292  row++;
293  }
294  col++; // all of the first 6 columns are dof_included since they're stresses
295  }
296 
297  std::vector<Real> pmep;
298  pmep.resize(_num_surfaces);
299  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
300  pmep[surface] = pm[surface];
301  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
302  {
303  if (!dof_included(6 + surface, deactivated_due_to_ld))
304  continue;
305  ep = _fspb_debug_pm_change[surface];
306  pmep[surface] += ep;
307  active.assign(_num_surfaces, true);
308  calculateRHS(
309  stress, intnl_old, intnl, pmep, delta_dp, rhs_ep, active, false, deactivated_due_to_ld_ep);
310  unsigned row = 0;
311  for (unsigned dof = 0; dof < whole_system_size; ++dof)
312  if (dof_included(dof, deactivated_due_to_ld))
313  {
314  jac[row][col] =
315  -(rhs_ep[dof] - orig_rhs[row]) / ep; // remember jacobian = -d(rhs)/d(something)
316  row++;
317  }
318  pmep[surface] -= ep;
319  col++;
320  }
321 
322  std::vector<Real> intnlep;
323  intnlep.resize(_num_models);
324  for (unsigned model = 0; model < _num_models; ++model)
325  intnlep[model] = intnl[model];
326  for (unsigned model = 0; model < _num_models; ++model)
327  {
328  if (!dof_included(6 + _num_surfaces + model, deactivated_due_to_ld))
329  continue;
330  ep = _fspb_debug_intnl_change[model];
331  intnlep[model] += ep;
332  active.assign(_num_surfaces, true);
333  calculateRHS(
334  stress, intnl_old, intnlep, pm, delta_dp, rhs_ep, active, false, deactivated_due_to_ld_ep);
335  unsigned row = 0;
336  for (unsigned dof = 0; dof < whole_system_size; ++dof)
337  if (dof_included(dof, deactivated_due_to_ld))
338  {
339  jac[row][col] =
340  -(rhs_ep[dof] - orig_rhs[row]) / ep; // remember jacobian = -d(rhs)/d(something)
341  row++;
342  }
343  intnlep[model] -= ep;
344  col++;
345  }
346 }

Referenced by checkJacobian(), and checkSolution().

◆ flowPotential()

void MultiPlasticityRawComponentAssembler::flowPotential ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< RankTwoTensor > &  r 
)
protectedvirtualinherited

The active flow potential(s) - one for each yield function.

Parameters
stressthe stress at which to calculate the flow potential
intnlvector of internal parameters
activeset of active constraints - only the active flow potentials are put into "r"
[out]rthe flow potential (flow potentials in the multi-surface case)

Definition at line 180 of file MultiPlasticityRawComponentAssembler.C.

184 {
185  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
186  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
187 
188  r.resize(0);
189  std::vector<unsigned int> active_surfaces_of_model;
190  std::vector<unsigned int>::iterator active_surface;
191  std::vector<RankTwoTensor> model_r;
192  for (unsigned model = 0; model < _num_models; ++model)
193  {
194  activeModelSurfaces(model, active, active_surfaces_of_model);
195  if (active_surfaces_of_model.size() > 0)
196  {
197  _f[model]->flowPotentialV(stress, intnl[model], model_r);
198  for (active_surface = active_surfaces_of_model.begin();
199  active_surface != active_surfaces_of_model.end();
200  ++active_surface)
201  r.push_back(model_r[*active_surface]);
202  }
203  }
204 }

Referenced by MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), ComputeMultiPlasticityStress::consistentTangentOperator(), fddflowPotential_dintnl(), and fddflowPotential_dstress().

◆ hardPotential()

void MultiPlasticityRawComponentAssembler::hardPotential ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  h 
)
protectedvirtualinherited

The active hardening potentials (one for each internal parameter and for each yield function) by assumption in the Userobjects, the h[a][alpha] is nonzero only if the surface alpha is part of model a, so we only calculate those here.

Parameters
stressthe stress at which to calculate the hardening potential
intnlvector of internal parameters
activeset of active constraints - only the active hardening potentials are put into "h"
[out]hthe hardening potentials. h[alpha] = hardening potential for yield fcn alpha (and, by the above assumption we know which hardening parameter, a, this belongs to)

Definition at line 262 of file MultiPlasticityRawComponentAssembler.C.

266 {
267  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
268  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
269 
270  h.resize(0);
271  std::vector<unsigned int> active_surfaces_of_model;
272  std::vector<unsigned int>::iterator active_surface;
273  std::vector<Real> model_h;
274  for (unsigned model = 0; model < _num_models; ++model)
275  {
276  activeModelSurfaces(model, active, active_surfaces_of_model);
277  if (active_surfaces_of_model.size() > 0)
278  {
279  _f[model]->hardPotentialV(stress, intnl[model], model_h);
280  for (active_surface = active_surfaces_of_model.begin();
281  active_surface != active_surfaces_of_model.end();
282  ++active_surface)
283  h.push_back(model_h[*active_surface]);
284  }
285  }
286 }

Referenced by MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), and ComputeMultiPlasticityStress::consistentTangentOperator().

◆ modelNumber()

unsigned int MultiPlasticityRawComponentAssembler::modelNumber ( unsigned int  surface)
protectedinherited

returns the model number, given the surface number

Definition at line 785 of file MultiPlasticityRawComponentAssembler.C.

786 {
787  return _model_given_surface[surface];
788 }

Referenced by ComputeMultiPlasticityStress::applyKuhnTucker(), MultiPlasticityLinearSystem::calculateJacobian(), ComputeMultiPlasticityStress::canAddConstraints(), ComputeMultiPlasticityStress::checkAdmissible(), ComputeMultiPlasticityStress::checkKuhnTucker(), ComputeMultiPlasticityStress::consistentTangentOperator(), fddflowPotential_dintnl(), fddyieldFunction_dintnl(), ComputeMultiPlasticityStress::residual2(), ComputeMultiPlasticityStress::returnMap(), and ComputeMultiPlasticityStress::singleStep().

◆ nrStep()

void MultiPlasticityLinearSystem::nrStep ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl_old,
const std::vector< Real > &  intnl,
const std::vector< Real > &  pm,
const RankFourTensor &  E_inv,
const RankTwoTensor &  delta_dp,
RankTwoTensor &  dstress,
std::vector< Real > &  dpm,
std::vector< Real > &  dintnl,
const std::vector< bool > &  active,
std::vector< bool > &  deactivated_due_to_ld 
)
protectedvirtualinherited

Performs one Newton-Raphson step.

The purpose here is to find the changes, dstress, dpm and dintnl according to the Newton-Raphson procedure

Parameters
stressCurrent value of stress
intnl_oldThe internal variables at the previous "time" step
intnlCurrent value of the internal variables
pmCurrent value of the plasticity multipliers (consistency parameters)
E_invinverse of the elasticity tensor
delta_dpCurrent value of the plastic-strain increment (ie plastic_strain - plastic_strain_old)
[out]dstressThe change in stress for a full Newton step
[out]dpmThe change in all plasticity multipliers for a full Newton step
[out]dintnlThe change in all internal variables for a full Newton step
activeThe active constraints
[out]deactivated_due_to_ldThe constraints deactivated due to linear-dependence of the flow directions

Definition at line 615 of file MultiPlasticityLinearSystem.C.

626 {
627  // Calculate RHS and Jacobian
628  std::vector<Real> rhs;
629  calculateRHS(stress, intnl_old, intnl, pm, delta_dp, rhs, active, true, deactivated_due_to_ld);
630 
631  std::vector<std::vector<Real>> jac;
632  calculateJacobian(stress, intnl, pm, E_inv, active, deactivated_due_to_ld, jac);
633 
634  // prepare for LAPACKgesv_ routine provided by PETSc
635  int system_size = rhs.size();
636 
637  std::vector<double> a(system_size * system_size);
638  // Fill in the a "matrix" by going down columns
639  unsigned ind = 0;
640  for (int col = 0; col < system_size; ++col)
641  for (int row = 0; row < system_size; ++row)
642  a[ind++] = jac[row][col];
643 
644  int nrhs = 1;
645  std::vector<int> ipiv(system_size);
646  int info;
647  LAPACKgesv_(&system_size, &nrhs, &a[0], &system_size, &ipiv[0], &rhs[0], &system_size, &info);
648 
649  if (info != 0)
650  mooseError("In solving the linear system in a Newton-Raphson process, the PETSC LAPACK gsev "
651  "routine returned with error code ",
652  info);
653 
654  // Extract the results back to dstress, dpm and dintnl
655  std::vector<bool> active_not_deact(_num_surfaces);
656  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
657  active_not_deact[surface] = (active[surface] && !deactivated_due_to_ld[surface]);
658 
659  unsigned int dim = 3;
660  ind = 0;
661 
662  for (unsigned i = 0; i < dim; ++i)
663  for (unsigned j = 0; j <= i; ++j)
664  dstress(i, j) = dstress(j, i) = rhs[ind++];
665  dpm.assign(_num_surfaces, 0);
666  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
667  if (active_not_deact[surface])
668  dpm[surface] = rhs[ind++];
669  dintnl.assign(_num_models, 0);
670  for (unsigned model = 0; model < _num_models; ++model)
671  if (anyActiveSurfaces(model, active_not_deact))
672  dintnl[model] = rhs[ind++];
673 
674  mooseAssert(static_cast<int>(ind) == system_size,
675  "Incorrect extracting of changes from NR solution in nrStep");
676 }

Referenced by checkSolution(), and ComputeMultiPlasticityStress::singleStep().

◆ outputAndCheckDebugParameters()

void MultiPlasticityDebugger::outputAndCheckDebugParameters ( )

Outputs the debug parameters: _fspb_debug_stress, _fspd_debug_pm, etc and checks that they are sized correctly.

Definition at line 55 of file MultiPlasticityDebugger.C.

56 {
57  Moose::err << "Debug Parameters are as follows\n";
58  Moose::err << "stress = \n";
59  _fspb_debug_stress.print();
60 
61  if (_fspb_debug_pm.size() != _num_surfaces || _fspb_debug_intnl.size() != _num_models ||
62  _fspb_debug_pm_change.size() != _num_surfaces ||
63  _fspb_debug_intnl_change.size() != _num_models)
64  mooseError("The debug parameters have the wrong size\n");
65 
66  Moose::err << "plastic multipliers =\n";
67  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
68  Moose::err << _fspb_debug_pm[surface] << "\n";
69 
70  Moose::err << "internal parameters =\n";
71  for (unsigned model = 0; model < _num_models; ++model)
72  Moose::err << _fspb_debug_intnl[model] << "\n";
73 
74  Moose::err << "finite-differencing parameter for stress-changes:\n"
75  << _fspb_debug_stress_change << "\n";
76  Moose::err << "finite-differencing parameter(s) for plastic-multiplier(s):\n";
77  for (unsigned surface = 0; surface < _num_surfaces; ++surface)
78  Moose::err << _fspb_debug_pm_change[surface] << "\n";
79  Moose::err << "finite-differencing parameter(s) for internal-parameter(s):\n";
80  for (unsigned model = 0; model < _num_models; ++model)
81  Moose::err << _fspb_debug_intnl_change[model] << "\n";
82 }

Referenced by checkDerivatives(), checkJacobian(), and checkSolution().

◆ returnMapAll()

bool MultiPlasticityRawComponentAssembler::returnMapAll ( const RankTwoTensor &  trial_stress,
const std::vector< Real > &  intnl_old,
const RankFourTensor &  E_ijkl,
Real  ep_plastic_tolerance,
RankTwoTensor &  stress,
std::vector< Real > &  intnl,
std::vector< Real > &  pm,
std::vector< Real > &  cumulative_pm,
RankTwoTensor &  delta_dp,
std::vector< Real > &  yf,
unsigned &  num_successful_plastic_returns,
unsigned &  custom_model 
)
protectedinherited

Performs a returnMap for each plastic model using their inbuilt returnMap functions.

Performs a returnMap for each plastic model.

This may be used to quickly ascertain whether a (trial_stress, intnl_old) configuration is admissible, or whether a single model's customized returnMap function can provide a solution to the return-map problem, or whether a full Newton-Raphson approach such as implemented in ComputeMultiPlasticityStress is needed.

There are three cases mentioned below: (A) The (trial_stress, intnl_old) configuration is admissible according to all plastic models (B) The (trial_stress, intnl_old) configuration is inadmissible to exactly one plastic model, and that model can successfully use its customized returnMap function to provide a returned (stress, intnl) configuration, and that configuration is admissible according to all plastic models (C) All other cases. This includes customized returnMap functions failing, or more than one plastic_model being inadmissible, etc

Parameters
trial_stressthe trial stress
intnl_oldthe old values of the internal parameters
E_ijklthe elasticity tensor
ep_plastic_tolerancethe tolerance on the plastic strain
[out]stressis set to trial_stress in case (A) or (C), and the returned value of stress in case (B).
[out]intnlis set to intnl_old in case (A) or (C), and the returned value of intnl in case (B)
[out]pmZero in case (A) or (C), otherwise the plastic multipliers needed to bring about the returnMap in case (B)
[in/out]cumulative_pm cumulative plastic multipliers, updated in case (B), otherwise left untouched
[out]delta_dpis unchanged in case (A) or (C), and is set to the change in plastic strain in case(B)
[out]yfwill contain the yield function values at (stress, intnl)
[out]num_successful_plastic_returnswill be 0 for (A) and (C), and 1 for (B)
Returns
true in case (A) and (B), and false in case (C)

If all models actually signal "elastic" by returning true from their returnMap, and by returning model_plastically_active=0, then yf will contain the yield function values num_successful_plastic_returns will be zero intnl = intnl_old delta_dp will be unchanged from its input value stress will be set to trial_stress pm will be zero cumulative_pm will be unchanged return value will be true num_successful_plastic_returns = 0

If only one model signals "plastically active" by returning true from its returnMap, and by returning model_plastically_active=1, then yf will contain the yield function values num_successful_plastic_returns will be one intnl will be set by the returnMap algorithm delta_dp will be set by the returnMap algorithm stress will be set by the returnMap algorithm pm will be nonzero for the single model, and zero for other models cumulative_pm will be updated return value will be true num_successful_plastic_returns = 1

If >1 model signals "plastically active" or if >=1 model's returnMap fails, then yf will contain the yield function values num_successful_plastic_returns will be set appropriately intnl = intnl_old delta_dp will be unchanged from its input value stress will be set to trial_stress pm will be zero cumulative_pm will be unchanged return value will be true if all returnMap functions returned true, otherwise it will be false num_successful_plastic_returns is set appropriately

Definition at line 599 of file MultiPlasticityRawComponentAssembler.C.

611 {
612  mooseAssert(intnl_old.size() == _num_models,
613  "returnMapAll: Incorrect size of internal parameters");
614  mooseAssert(intnl.size() == _num_models, "returnMapAll: Incorrect size of internal parameters");
615  mooseAssert(pm.size() == _num_surfaces, "returnMapAll: Incorrect size of pm");
616 
617  num_successful_plastic_returns = 0;
618  yf.resize(0);
619  pm.assign(_num_surfaces, 0.0);
620 
621  RankTwoTensor returned_stress; // each model will give a returned_stress. if only one model is
622  // plastically active, i set stress=returned_stress, so as to
623  // record this returned value
624  std::vector<Real> model_f;
625  RankTwoTensor model_delta_dp;
626  std::vector<Real> model_pm;
627  bool trial_stress_inadmissible;
628  bool successful_return = true;
629  unsigned the_single_plastic_model = 0;
630  bool using_custom_return_map = true;
631 
632  // run through all the plastic models, performing their
633  // returnMap algorithms.
634  // If one finds (trial_stress, intnl) inadmissible and
635  // successfully returns, break from the loop to evaluate
636  // all the others at that returned stress
637  for (unsigned model = 0; model < _num_models; ++model)
638  {
639  if (using_custom_return_map)
640  {
641  model_pm.assign(_f[model]->numberSurfaces(), 0.0);
642  bool model_returned = _f[model]->returnMap(trial_stress,
643  intnl_old[model],
644  E_ijkl,
645  ep_plastic_tolerance,
646  returned_stress,
647  intnl[model],
648  model_pm,
649  model_delta_dp,
650  model_f,
651  trial_stress_inadmissible);
652  if (!trial_stress_inadmissible)
653  {
654  // in the elastic zone: record the yield-function values (returned_stress, intnl, model_pm
655  // and model_delta_dp are undefined)
656  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces();
657  ++model_surface)
658  yf.push_back(model_f[model_surface]);
659  }
660  else if (trial_stress_inadmissible && !model_returned)
661  {
662  // in the plastic zone, and the customized returnMap failed
663  // for some reason (or wasn't implemented). The coder
664  // should have correctly returned model_f(trial_stress, intnl_old)
665  // so record them
666  // (returned_stress, intnl, model_pm and model_delta_dp are undefined)
667  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces();
668  ++model_surface)
669  yf.push_back(model_f[model_surface]);
670  // now there's almost zero point in using the custom
671  // returnMap functions
672  using_custom_return_map = false;
673  successful_return = false;
674  }
675  else
676  {
677  // in the plastic zone, and the customized returnMap
678  // succeeded.
679  // record the first returned_stress and delta_dp if everything is going OK
680  // as they could be the actual answer
681  if (trial_stress_inadmissible)
682  num_successful_plastic_returns++;
683  the_single_plastic_model = model;
684  stress = returned_stress;
685  // note that i break here, and don't push_back
686  // model_f to yf. So now yf contains only the values of
687  // model_f from previous models to the_single_plastic_model
688  // also i don't set delta_dp = model_delta_dp yet, because
689  // i might find problems later on
690  // also, don't increment cumulative_pm for the same reason
691 
692  break;
693  }
694  }
695  else
696  {
697  // not using custom returnMap functions because one
698  // has already failed and that one said trial_stress
699  // was inadmissible. So now calculate the yield functions
700  // at the trial stress
701  _f[model]->yieldFunctionV(trial_stress, intnl_old[model], model_f);
702  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
703  yf.push_back(model_f[model_surface]);
704  }
705  }
706 
707  if (num_successful_plastic_returns == 0)
708  {
709  // here either all the models were elastic (successful_return=true),
710  // or some were plastic and either the customized returnMap failed
711  // or wasn't implemented (successful_return=false).
712  // In either case, have to set the following:
713  stress = trial_stress;
714  for (unsigned model = 0; model < _num_models; ++model)
715  intnl[model] = intnl_old[model];
716  return successful_return;
717  }
718 
719  // Now we know that num_successful_plastic_returns == 1 and all the other
720  // models (with model number < the_single_plastic_model) must have been
721  // admissible at (trial_stress, intnl). However, all models might
722  // not be admissible at (trial_stress, intnl), so must check that
723  std::vector<Real> yf_at_returned_stress(0);
724  bool all_admissible = true;
725  for (unsigned model = 0; model < _num_models; ++model)
726  {
727  if (model == the_single_plastic_model)
728  {
729  // no need to spend time calculating the yield function: we know its admissible
730  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
731  yf_at_returned_stress.push_back(model_f[model_surface]);
732  continue;
733  }
734  _f[model]->yieldFunctionV(stress, intnl_old[model], model_f);
735  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
736  {
737  if (model_f[model_surface] > _f[model]->_f_tol)
738  // bummer, this model is not admissible at the returned_stress
739  all_admissible = false;
740  yf_at_returned_stress.push_back(model_f[model_surface]);
741  }
742  if (!all_admissible)
743  // no point in continuing computing yield functions
744  break;
745  }
746 
747  if (!all_admissible)
748  {
749  // we tried using the returned value of stress predicted by
750  // the_single_plastic_model, but it wasn't admissible according
751  // to other plastic models. We need to set:
752  stress = trial_stress;
753  for (unsigned model = 0; model < _num_models; ++model)
754  intnl[model] = intnl_old[model];
755  // and calculate the remainder of the yield functions at trial_stress
756  for (unsigned model = the_single_plastic_model; model < _num_models; ++model)
757  {
758  _f[model]->yieldFunctionV(trial_stress, intnl[model], model_f);
759  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
760  yf.push_back(model_f[model_surface]);
761  }
762  num_successful_plastic_returns = 0;
763  return false;
764  }
765 
766  // So the customized returnMap algorithm can provide a returned
767  // (stress, intnl) configuration, and that is admissible according
768  // to all plastic models
769  yf.resize(0);
770  for (unsigned surface = 0; surface < yf_at_returned_stress.size(); ++surface)
771  yf.push_back(yf_at_returned_stress[surface]);
772  delta_dp = model_delta_dp;
773  for (unsigned model_surface = 0; model_surface < _f[the_single_plastic_model]->numberSurfaces();
774  ++model_surface)
775  {
776  cumulative_pm[_surfaces_given_model[the_single_plastic_model][model_surface]] +=
777  model_pm[model_surface];
778  pm[_surfaces_given_model[the_single_plastic_model][model_surface]] = model_pm[model_surface];
779  }
780  custom_model = the_single_plastic_model;
781  return true;
782 }

Referenced by ComputeMultiPlasticityStress::quickStep().

◆ singularValuesOfR()

int MultiPlasticityLinearSystem::singularValuesOfR ( const std::vector< RankTwoTensor > &  r,
std::vector< Real > &  s 
)
privatevirtualinherited

Performs a singular-value decomposition of r and returns the singular values.

Example: If r has size 5 then the singular values of the following matrix are returned: ( r[0](0,0) r[0](0,1) r[0](0,2) r[0](1,1) r[0](1,2) r[0](2,2) ) ( r[1](0,0) r[1](0,1) r[1](0,2) r[1](1,1) r[1](1,2) r[1](2,2) ) a = ( r[2](0,0) r[2](0,1) r[2](0,2) r[2](1,1) r[2](1,2) r[2](2,2) ) ( r[3](0,0) r[3](0,1) r[3](0,2) r[3](1,1) r[3](1,2) r[3](2,2) ) ( r[4](0,0) r[4](0,1) r[4](0,2) r[4](1,1) r[4](1,2) r[4](2,2) )

Parameters
rThe flow directions
[out]sThe singular values
Returns
The return value from the PETSc LAPACK gesvd reoutine

Definition at line 47 of file MultiPlasticityLinearSystem.C.

49 {
50  int bm = r.size();
51  int bn = 6;
52 
53  s.resize(std::min(bm, bn));
54 
55  // prepare for gesvd or gesdd routine provided by PETSc
56  // Want to find the singular values of matrix
57  // ( r[0](0,0) r[0](0,1) r[0](0,2) r[0](1,1) r[0](1,2) r[0](2,2) )
58  // ( r[1](0,0) r[1](0,1) r[1](0,2) r[1](1,1) r[1](1,2) r[1](2,2) )
59  // a = ( r[2](0,0) r[2](0,1) r[2](0,2) r[2](1,1) r[2](1,2) r[2](2,2) )
60  // ( r[3](0,0) r[3](0,1) r[3](0,2) r[3](1,1) r[3](1,2) r[3](2,2) )
61  // ( r[4](0,0) r[4](0,1) r[4](0,2) r[4](1,1) r[4](1,2) r[4](2,2) )
62  // bm = 5
63 
64  std::vector<double> a(bm * 6);
65  // Fill in the a "matrix" by going down columns
66  unsigned ind = 0;
67  for (int col = 0; col < 3; ++col)
68  for (int row = 0; row < bm; ++row)
69  a[ind++] = r[row](0, col);
70  for (int col = 3; col < 5; ++col)
71  for (int row = 0; row < bm; ++row)
72  a[ind++] = r[row](1, col - 2);
73  for (int row = 0; row < bm; ++row)
74  a[ind++] = r[row](2, 2);
75 
76  // u and vt are dummy variables because they won't
77  // get referenced due to the "N" and "N" choices
78  int sizeu = 1;
79  std::vector<double> u(sizeu);
80  int sizevt = 1;
81  std::vector<double> vt(sizevt);
82 
83  int sizework = 16 * (bm + 6); // this is above the lowerbound specified in the LAPACK doco
84  std::vector<double> work(sizework);
85 
86  int info;
87 
88  LAPACKgesvd_("N",
89  "N",
90  &bm,
91  &bn,
92  &a[0],
93  &bm,
94  &s[0],
95  &u[0],
96  &sizeu,
97  &vt[0],
98  &sizevt,
99  &work[0],
100  &sizework,
101  &info);
102 
103  return info;
104 }

Referenced by MultiPlasticityLinearSystem::eliminateLinearDependence().

◆ validParams()

InputParameters MultiPlasticityDebugger::validParams ( )
static

Definition at line 17 of file MultiPlasticityDebugger.C.

18 {
19  InputParameters params = MultiPlasticityLinearSystem::validParams();
20  MooseEnum debug_fspb_type("none crash jacobian jacobian_and_linear_system", "none");
21  params.addParam<MooseEnum>("debug_fspb",
22  debug_fspb_type,
23  "Debug types for use by developers when creating new "
24  "plasticity models, not for general use. 2 = debug Jacobian "
25  "entries, 3 = check the entire Jacobian, and check Ax=b");
26  params.addParam<RealTensorValue>("debug_jac_at_stress",
27  RealTensorValue(),
28  "Debug Jacobian entries at this stress. For use by developers");
29  params.addParam<std::vector<Real>>("debug_jac_at_pm",
30  "Debug Jacobian entries at these plastic multipliers");
31  params.addParam<std::vector<Real>>("debug_jac_at_intnl",
32  "Debug Jacobian entries at these internal parameters");
33  params.addParam<Real>(
34  "debug_stress_change", 1.0, "Debug finite differencing parameter for the stress");
35  params.addParam<std::vector<Real>>(
36  "debug_pm_change", "Debug finite differencing parameters for the plastic multipliers");
37  params.addParam<std::vector<Real>>(
38  "debug_intnl_change", "Debug finite differencing parameters for the internal parameters");
39  return params;
40 }

Referenced by ComputeMultiPlasticityStress::validParams().

◆ yieldFunction()

void MultiPlasticityRawComponentAssembler::yieldFunction ( const RankTwoTensor &  stress,
const std::vector< Real > &  intnl,
const std::vector< bool > &  active,
std::vector< Real > &  f 
)
protectedvirtualinherited

The active yield function(s)

Parameters
stressthe stress at which to calculate the yield function
intnlvector of internal parameters
activeset of active constraints - only the active yield functions are put into "f"
[out]fthe yield function (or functions in the case of multisurface plasticity)

Definition at line 98 of file MultiPlasticityRawComponentAssembler.C.

102 {
103  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
104  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
105 
106  f.resize(0);
107  std::vector<unsigned int> active_surfaces_of_model;
108  std::vector<unsigned int>::iterator active_surface;
109  std::vector<Real> model_f;
110  for (unsigned model = 0; model < _num_models; ++model)
111  {
112  activeModelSurfaces(model, active, active_surfaces_of_model);
113  if (active_surfaces_of_model.size() > 0)
114  {
115  _f[model]->yieldFunctionV(stress, intnl[model], model_f);
116  for (active_surface = active_surfaces_of_model.begin();
117  active_surface != active_surfaces_of_model.end();
118  ++active_surface)
119  f.push_back(model_f[*active_surface]);
120  }
121  }
122 }

Referenced by ComputeMultiPlasticityStress::buildDumbOrder(), MultiPlasticityLinearSystem::calculateConstraints(), ComputeMultiPlasticityStress::checkAdmissible(), fddyieldFunction_dintnl(), fddyieldFunction_dstress(), and ComputeMultiPlasticityStress::returnMap().

Member Data Documentation

◆ _f

std::vector<const TensorMechanicsPlasticModel *> MultiPlasticityRawComponentAssembler::_f
protectedinherited

User objects that define the yield functions, flow potentials, etc.

Definition at line 75 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::activeModelSurfaces(), MultiPlasticityRawComponentAssembler::activeSurfaces(), MultiPlasticityRawComponentAssembler::anyActiveSurfaces(), ComputeMultiPlasticityStress::applyKuhnTucker(), MultiPlasticityRawComponentAssembler::buildActiveConstraints(), MultiPlasticityRawComponentAssembler::buildActiveConstraintsJoint(), MultiPlasticityRawComponentAssembler::buildActiveConstraintsRock(), ComputeMultiPlasticityStress::canAddConstraints(), ComputeMultiPlasticityStress::checkAdmissible(), ComputeMultiPlasticityStress::checkKuhnTucker(), MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), MultiPlasticityLinearSystem::MultiPlasticityLinearSystem(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), ComputeMultiPlasticityStress::quickStep(), ComputeMultiPlasticityStress::residual2(), ComputeMultiPlasticityStress::returnMap(), MultiPlasticityRawComponentAssembler::returnMapAll(), ComputeMultiPlasticityStress::singleStep(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ _fspb_debug

MooseEnum MultiPlasticityDebugger::_fspb_debug
protected

none - don't do any debugging crash - currently inactive jacobian - check the jacobian entries jacobian_and_linear_system - check entire jacobian and check that Ax=b

Definition at line 62 of file MultiPlasticityDebugger.h.

Referenced by ComputeMultiPlasticityStress::computeQpStress(), and ComputeMultiPlasticityStress::initQpStatefulProperties().

◆ _fspb_debug_intnl

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_intnl
protected

Debug the Jacobian entires at these internal parameters.

Definition at line 71 of file MultiPlasticityDebugger.h.

Referenced by checkDerivatives(), checkJacobian(), checkSolution(), outputAndCheckDebugParameters(), and ComputeMultiPlasticityStress::plasticStep().

◆ _fspb_debug_intnl_change

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_intnl_change
protected

Debug finite-differencing parameters for the internal parameters.

Definition at line 80 of file MultiPlasticityDebugger.h.

Referenced by fddflowPotential_dintnl(), fddyieldFunction_dintnl(), fdJacobian(), and outputAndCheckDebugParameters().

◆ _fspb_debug_pm

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_pm
protected

Debug the Jacobian entires at these plastic multipliers.

Definition at line 68 of file MultiPlasticityDebugger.h.

Referenced by checkJacobian(), checkSolution(), outputAndCheckDebugParameters(), and ComputeMultiPlasticityStress::plasticStep().

◆ _fspb_debug_pm_change

std::vector<Real> MultiPlasticityDebugger::_fspb_debug_pm_change
protected

Debug finite-differencing parameters for the plastic multipliers.

Definition at line 77 of file MultiPlasticityDebugger.h.

Referenced by fdJacobian(), and outputAndCheckDebugParameters().

◆ _fspb_debug_stress

RankTwoTensor MultiPlasticityDebugger::_fspb_debug_stress
protected

Debug the Jacobian entries at this stress.

Definition at line 65 of file MultiPlasticityDebugger.h.

Referenced by checkDerivatives(), checkJacobian(), checkSolution(), outputAndCheckDebugParameters(), and ComputeMultiPlasticityStress::plasticStep().

◆ _fspb_debug_stress_change

Real MultiPlasticityDebugger::_fspb_debug_stress_change
protected

Debug finite-differencing parameter for the stress.

Definition at line 74 of file MultiPlasticityDebugger.h.

Referenced by fddflowPotential_dstress(), fddyieldFunction_dstress(), fdJacobian(), and outputAndCheckDebugParameters().

◆ _min_f_tol

Real MultiPlasticityLinearSystem::_min_f_tol
protectedinherited

Minimum value of the _f_tol parameters for the Yield Function User Objects.

Definition at line 136 of file MultiPlasticityLinearSystem.h.

Referenced by MultiPlasticityLinearSystem::eliminateLinearDependence(), and MultiPlasticityLinearSystem::MultiPlasticityLinearSystem().

◆ _model_given_surface

std::vector<unsigned int> MultiPlasticityRawComponentAssembler::_model_given_surface
privateinherited

given a surface number, this returns the model number

Definition at line 295 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::modelNumber(), and MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _model_surface_given_surface

std::vector<unsigned int> MultiPlasticityRawComponentAssembler::_model_surface_given_surface
privateinherited

given a surface number, this returns the corresponding-model's internal surface number

Definition at line 298 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _num_models

unsigned int MultiPlasticityRawComponentAssembler::_num_models
protectedinherited

Number of plastic models for this material.

Definition at line 57 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints(), MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityLinearSystem::calculateRHS(), checkSolution(), MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), fddflowPotential_dintnl(), fddyieldFunction_dintnl(), fdJacobian(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), ComputeMultiPlasticityStress::initQpStatefulProperties(), MultiPlasticityLinearSystem::MultiPlasticityLinearSystem(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), MultiPlasticityLinearSystem::nrStep(), outputAndCheckDebugParameters(), ComputeMultiPlasticityStress::plasticStep(), ComputeMultiPlasticityStress::residual2(), ComputeMultiPlasticityStress::returnMap(), MultiPlasticityRawComponentAssembler::returnMapAll(), ComputeMultiPlasticityStress::singleStep(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ _num_surfaces

unsigned int MultiPlasticityRawComponentAssembler::_num_surfaces
protectedinherited

Number of surfaces within the plastic models.

For many situations this will be = _num_models since each model will contain just one surface. More generally it is >= _num_models. For instance, Mohr-Coulomb is a single model with 6 surfaces

Definition at line 66 of file MultiPlasticityRawComponentAssembler.h.

Referenced by ComputeMultiPlasticityStress::activeCombinationNumber(), ComputeMultiPlasticityStress::applyKuhnTucker(), MultiPlasticityRawComponentAssembler::buildActiveConstraints(), ComputeMultiPlasticityStress::buildDumbOrder(), MultiPlasticityLinearSystem::calculateConstraints(), MultiPlasticityLinearSystem::calculateJacobian(), MultiPlasticityLinearSystem::calculateRHS(), ComputeMultiPlasticityStress::canAddConstraints(), ComputeMultiPlasticityStress::canIncrementDumb(), ComputeMultiPlasticityStress::changeScheme(), ComputeMultiPlasticityStress::checkAdmissible(), checkDerivatives(), checkJacobian(), ComputeMultiPlasticityStress::checkKuhnTucker(), checkSolution(), ComputeMultiPlasticityStress::ComputeMultiPlasticityStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeMultiPlasticityStress::consistentTangentOperator(), MultiPlasticityRawComponentAssembler::dflowPotential_dintnl(), MultiPlasticityRawComponentAssembler::dflowPotential_dstress(), MultiPlasticityRawComponentAssembler::dhardPotential_dintnl(), MultiPlasticityRawComponentAssembler::dhardPotential_dstress(), dof_included(), MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl(), MultiPlasticityRawComponentAssembler::dyieldFunction_dstress(), MultiPlasticityLinearSystem::eliminateLinearDependence(), fddflowPotential_dintnl(), fddflowPotential_dstress(), fddyieldFunction_dintnl(), fddyieldFunction_dstress(), fdJacobian(), MultiPlasticityRawComponentAssembler::flowPotential(), MultiPlasticityRawComponentAssembler::hardPotential(), ComputeMultiPlasticityStress::incrementDumb(), ComputeMultiPlasticityStress::initQpStatefulProperties(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), MultiPlasticityLinearSystem::nrStep(), ComputeMultiPlasticityStress::numberActive(), outputAndCheckDebugParameters(), ComputeMultiPlasticityStress::plasticStep(), ComputeMultiPlasticityStress::reinstateLinearDependentConstraints(), ComputeMultiPlasticityStress::residual2(), ComputeMultiPlasticityStress::returnMap(), MultiPlasticityRawComponentAssembler::returnMapAll(), ComputeMultiPlasticityStress::singleStep(), and MultiPlasticityRawComponentAssembler::yieldFunction().

◆ _params

const InputParameters& MultiPlasticityRawComponentAssembler::_params
protectedinherited

Definition at line 54 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _specialIC

MooseEnum MultiPlasticityRawComponentAssembler::_specialIC
protectedinherited

Allows initial set of active constraints to be chosen optimally.

Definition at line 72 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::buildActiveConstraints(), and MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler().

◆ _surfaces_given_model

std::vector<std::vector<unsigned int> > MultiPlasticityRawComponentAssembler::_surfaces_given_model
protectedinherited

_surfaces_given_model[model_number] = vector of surface numbers for this model

Definition at line 69 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler::activeModelSurfaces(), MultiPlasticityRawComponentAssembler::activeSurfaces(), MultiPlasticityRawComponentAssembler::anyActiveSurfaces(), MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler(), ComputeMultiPlasticityStress::quickStep(), and MultiPlasticityRawComponentAssembler::returnMapAll().

◆ _svd_tol

Real MultiPlasticityLinearSystem::_svd_tol
protectedinherited

Tolerance on the minimum ratio of singular values before flow-directions are deemed linearly dependent.

Definition at line 133 of file MultiPlasticityLinearSystem.h.

Referenced by MultiPlasticityLinearSystem::eliminateLinearDependence().

The documentation for this class was generated from the following files:
MultiPlasticityRawComponentAssembler::dflowPotential_dintnl
virtual void dflowPotential_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dr_dintnl)
The derivative of the active flow potentials with respect to the active internal parameters The UserO...
Definition: MultiPlasticityRawComponentAssembler.C:235
MultiPlasticityRawComponentAssembler::modelNumber
unsigned int modelNumber(unsigned int surface)
returns the model number, given the surface number
Definition: MultiPlasticityRawComponentAssembler.C:785
MultiPlasticityLinearSystem::nrStep
virtual void nrStep(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const RankTwoTensor &delta_dp, RankTwoTensor &dstress, std::vector< Real > &dpm, std::vector< Real > &dintnl, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
Performs one Newton-Raphson step.
Definition: MultiPlasticityLinearSystem.C:615
MultiPlasticityRawComponentAssembler::dhardPotential_dintnl
virtual void dhardPotential_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &dh_dintnl)
The derivative of the active hardening potentials with respect to the active internal parameters.
Definition: MultiPlasticityRawComponentAssembler.C:317
MultiPlasticityRawComponentAssembler::_model_given_surface
std::vector< unsigned int > _model_given_surface
given a surface number, this returns the model number
Definition: MultiPlasticityRawComponentAssembler.h:295
MultiPlasticityRawComponentAssembler::flowPotential
virtual void flowPotential(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &r)
The active flow potential(s) - one for each yield function.
Definition: MultiPlasticityRawComponentAssembler.C:180
MultiPlasticityRawComponentAssembler::_specialIC
MooseEnum _specialIC
Allows initial set of active constraints to be chosen optimally.
Definition: MultiPlasticityRawComponentAssembler.h:72
MultiPlasticityDebugger::_fspb_debug_intnl
std::vector< Real > _fspb_debug_intnl
Debug the Jacobian entires at these internal parameters.
Definition: MultiPlasticityDebugger.h:71
MultiPlasticityDebugger::fdJacobian
virtual void fdJacobian(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, const RankFourTensor &E_inv, bool eliminate_ld, std::vector< std::vector< Real >> &jac)
The Jacobian calculated using finite differences.
Definition: MultiPlasticityDebugger.C:221
MultiPlasticityLinearSystem::calculateConstraints
virtual void calculateConstraints(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &f, std::vector< RankTwoTensor > &r, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active)
The constraints.
Definition: MultiPlasticityLinearSystem.C:228
MultiPlasticityRawComponentAssembler::_num_surfaces
unsigned int _num_surfaces
Number of surfaces within the plastic models.
Definition: MultiPlasticityRawComponentAssembler.h:66
MultiPlasticityLinearSystem::eliminateLinearDependence
virtual void eliminateLinearDependence(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &f, const std::vector< RankTwoTensor > &r, const std::vector< bool > &active, std::vector< bool > &deactivated_due_to_ld)
Performs a number of singular-value decompositions to check for linear-dependence of the active direc...
Definition: MultiPlasticityLinearSystem.C:107
MultiPlasticityRawComponentAssembler::hardPotential
virtual void hardPotential(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &h)
The active hardening potentials (one for each internal parameter and for each yield function) by assu...
Definition: MultiPlasticityRawComponentAssembler.C:262
MultiPlasticityDebugger::_fspb_debug
MooseEnum _fspb_debug
none - don't do any debugging crash - currently inactive jacobian - check the jacobian entries jacobi...
Definition: MultiPlasticityDebugger.h:62
MultiPlasticityDebugger::fddyieldFunction_dstress
void fddyieldFunction_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &df_dstress)
The finite-difference derivative of yield function(s) with respect to stress.
Definition: MultiPlasticityDebugger.C:519
MultiPlasticityLinearSystem::calculateJacobian
virtual void calculateJacobian(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankFourTensor &E_inv, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, std::vector< std::vector< Real >> &jac)
d(rhs)/d(dof)
Definition: MultiPlasticityLinearSystem.C:367
MultiPlasticityLinearSystem::calculateRHS
virtual void calculateRHS(const RankTwoTensor &stress, const std::vector< Real > &intnl_old, const std::vector< Real > &intnl, const std::vector< Real > &pm, const RankTwoTensor &delta_dp, std::vector< Real > &rhs, const std::vector< bool > &active, bool eliminate_ld, std::vector< bool > &deactivated_due_to_ld)
Calculate the RHS which is rhs = -(epp(0,0), epp(1,0), epp(1,1), epp(2,0), epp(2,1),...
Definition: MultiPlasticityLinearSystem.C:288
MultiPlasticityDebugger::fddflowPotential_dintnl
virtual void fddflowPotential_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankTwoTensor > &dr_dintnl)
The finite-difference derivative of the flow potentials with respect to internal parameters.
Definition: MultiPlasticityDebugger.C:607
MultiPlasticityDebugger::_fspb_debug_stress
RankTwoTensor _fspb_debug_stress
Debug the Jacobian entries at this stress.
Definition: MultiPlasticityDebugger.h:65
MultiPlasticityDebugger::dof_included
bool dof_included(unsigned int dof, const std::vector< bool > &deactivated_due_to_ld)
Definition: MultiPlasticityDebugger.C:349
MultiPlasticityRawComponentAssembler::_f
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
Definition: MultiPlasticityRawComponentAssembler.h:75
MultiPlasticityDebugger::_fspb_debug_stress_change
Real _fspb_debug_stress_change
Debug finite-differencing parameter for the stress.
Definition: MultiPlasticityDebugger.h:74
MultiPlasticityRawComponentAssembler::dyieldFunction_dintnl
virtual void dyieldFunction_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &df_dintnl)
The derivative of active yield function(s) with respect to their internal parameters (the user object...
Definition: MultiPlasticityRawComponentAssembler.C:153
MultiPlasticityRawComponentAssembler::dhardPotential_dstress
virtual void dhardPotential_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &dh_dstress)
The derivative of the active hardening potentials with respect to stress By assumption in the Userobj...
Definition: MultiPlasticityRawComponentAssembler.C:289
MultiPlasticityLinearSystem::validParams
static InputParameters validParams()
Definition: MultiPlasticityLinearSystem.C:22
MultiPlasticityDebugger::_fspb_debug_intnl_change
std::vector< Real > _fspb_debug_intnl_change
Debug finite-differencing parameters for the internal parameters.
Definition: MultiPlasticityDebugger.h:80
MultiPlasticityDebugger::fddflowPotential_dstress
virtual void fddflowPotential_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< RankFourTensor > &dr_dstress)
The finite-difference derivative of the flow potential(s) with respect to stress.
Definition: MultiPlasticityDebugger.C:578
MultiPlasticityRawComponentAssembler::dyieldFunction_dstress
virtual void dyieldFunction_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankTwoTensor > &df_dstress)
The derivative of the active yield function(s) with respect to stress.
Definition: MultiPlasticityRawComponentAssembler.C:125
MultiPlasticityDebugger::outputAndCheckDebugParameters
void outputAndCheckDebugParameters()
Outputs the debug parameters: _fspb_debug_stress, _fspd_debug_pm, etc and checks that they are sized ...
Definition: MultiPlasticityDebugger.C:55
MultiPlasticityRawComponentAssembler::activeModelSurfaces
void activeModelSurfaces(int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces_of_model)
Returns the internal surface number(s) of the active surfaces of the given model This may be of size=...
Definition: MultiPlasticityRawComponentAssembler.C:811
MultiPlasticityDebugger::_fspb_debug_pm
std::vector< Real > _fspb_debug_pm
Debug the Jacobian entires at these plastic multipliers.
Definition: MultiPlasticityDebugger.h:68
RankTwoTensor
RankTwoTensorTempl< Real > RankTwoTensor
Definition: ACGrGrElasticDrivingForce.h:17
MultiPlasticityRawComponentAssembler::_num_models
unsigned int _num_models
Number of plastic models for this material.
Definition: MultiPlasticityRawComponentAssembler.h:57
MultiPlasticityRawComponentAssembler::anyActiveSurfaces
bool anyActiveSurfaces(int model, const std::vector< bool > &active)
returns true if any internal surfaces of the given model are active according to 'active'
Definition: MultiPlasticityRawComponentAssembler.C:791
MultiPlasticityLinearSystem::MultiPlasticityLinearSystem
MultiPlasticityLinearSystem(const MooseObject *moose_object)
Definition: MultiPlasticityLinearSystem.C:34
RankFourTensorTempl
Definition: ACGrGrElasticDrivingForce.h:20
MultiPlasticityDebugger::fddyieldFunction_dintnl
void fddyieldFunction_dintnl(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &df_dintnl)
The finite-difference derivative of yield function(s) with respect to internal parameter(s)
Definition: MultiPlasticityDebugger.C:547
MultiPlasticityRawComponentAssembler::_params
const InputParameters & _params
Definition: MultiPlasticityRawComponentAssembler.h:54
MultiPlasticityRawComponentAssembler::buildActiveConstraintsJoint
void buildActiveConstraintsJoint(const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
"Joint" version Constructs a set of active constraints, given the yield functions,...
Definition: MultiPlasticityRawComponentAssembler.C:381
MultiPlasticityLinearSystem::_min_f_tol
Real _min_f_tol
Minimum value of the _f_tol parameters for the Yield Function User Objects.
Definition: MultiPlasticityLinearSystem.h:136
MultiPlasticityRawComponentAssembler::yieldFunction
virtual void yieldFunction(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< Real > &f)
The active yield function(s)
Definition: MultiPlasticityRawComponentAssembler.C:98
MultiPlasticityLinearSystem::singularValuesOfR
virtual int singularValuesOfR(const std::vector< RankTwoTensor > &r, std::vector< Real > &s)
Performs a singular-value decomposition of r and returns the singular values.
Definition: MultiPlasticityLinearSystem.C:47
MultiPlasticityDebugger::_fspb_debug_pm_change
std::vector< Real > _fspb_debug_pm_change
Debug finite-differencing parameters for the plastic multipliers.
Definition: MultiPlasticityDebugger.h:77
RankTwoTensorTempl< Real >
MultiPlasticityLinearSystem::_svd_tol
Real _svd_tol
Tolerance on the minimum ratio of singular values before flow-directions are deemed linearly dependen...
Definition: MultiPlasticityLinearSystem.h:133
MultiPlasticityRawComponentAssembler::activeSurfaces
void activeSurfaces(int model, const std::vector< bool > &active, std::vector< unsigned int > &active_surfaces)
Returns the external surface number(s) of the active surfaces of the given model This may be of size=...
Definition: MultiPlasticityRawComponentAssembler.C:800
MultiPlasticityRawComponentAssembler::dflowPotential_dstress
virtual void dflowPotential_dstress(const RankTwoTensor &stress, const std::vector< Real > &intnl, const std::vector< bool > &active, std::vector< RankFourTensor > &dr_dstress)
The derivative of the active flow potential(s) with respect to stress.
Definition: MultiPlasticityRawComponentAssembler.C:207
RankFourTensor
RankFourTensorTempl< Real > RankFourTensor
Definition: ACGrGrElasticDrivingForce.h:20
RankTwoScalarTools::L2norm
T L2norm(const RankTwoTensorTempl< T > &r2tensor)
Definition: RankTwoScalarTools.h:98
MultiPlasticityRawComponentAssembler::_surfaces_given_model
std::vector< std::vector< unsigned int > > _surfaces_given_model
_surfaces_given_model[model_number] = vector of surface numbers for this model
Definition: MultiPlasticityRawComponentAssembler.h:69
MultiPlasticityRawComponentAssembler::buildActiveConstraintsRock
void buildActiveConstraintsRock(const std::vector< Real > &f, const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &Eijkl, std::vector< bool > &act)
"Rock" version Constructs a set of active constraints, given the yield functions, f.
Definition: MultiPlasticityRawComponentAssembler.C:422
Generated on Sat Jan 25 2020 12:05:08 for www.mooseframework.org by   doxygen 1.8.16