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

MultiPlasticityRawComponentAssembler holds and computes yield functions, flow directions, etc, for use in FiniteStrainMultiPlasticity. More...

#include <MultiPlasticityRawComponentAssembler.h>

Inheritance diagram for MultiPlasticityRawComponentAssembler:
[legend]

Public Member Functions

 MultiPlasticityRawComponentAssembler (const MooseObject *moose_object)
 
virtual ~MultiPlasticityRawComponentAssembler ()
 

Protected Member Functions

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

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 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

MultiPlasticityRawComponentAssembler holds and computes yield functions, flow directions, etc, for use in FiniteStrainMultiPlasticity.

Here there are a number of numbering systems.

There are _num_models of plastic models, read from the input file. Each of these models can, in principal, contain multiple 'internal surfaces'. Models are numbered from 0 to _num_models - 1.

There are num_surfaces surfaces. This = sum{plastic_models} (numberSurfaces in model) Evidently _num_surface >= _num_models Surfaces are numbered from 0 to _num_surfaces - 1.

The std::vectors _model_given_surface, _model_surface_given_surface and _surfaces_given_model allow translation between these

Definition at line 44 of file MultiPlasticityRawComponentAssembler.h.

Constructor & Destructor Documentation

◆ MultiPlasticityRawComponentAssembler()

MultiPlasticityRawComponentAssembler::MultiPlasticityRawComponentAssembler ( const MooseObject *  moose_object)

Definition at line 37 of file MultiPlasticityRawComponentAssembler.C.

39  : UserObjectInterface(moose_object),
40  _params(moose_object->parameters()),
41  _num_models(_params.get<std::vector<UserObjectName>>("plastic_models").size()),
42  _num_surfaces(0),
43  _specialIC(_params.get<MooseEnum>("specialIC"))
44 {
45  _f.resize(_num_models);
46  for (unsigned model = 0; model < _num_models; ++model)
47  _f[model] = &getUserObjectByName<TensorMechanicsPlasticModel>(
48  _params.get<std::vector<UserObjectName>>("plastic_models")[model]);
49 
50  for (unsigned model = 0; model < _num_models; ++model)
51  _num_surfaces += _f[model]->numberSurfaces();
52 
55  unsigned int surface = 0;
56  for (unsigned model = 0; model < _num_models; ++model)
57  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
58  {
59  _model_given_surface[surface] = model;
60  _model_surface_given_surface[surface] = model_surface;
61  surface++;
62  }
63 
65  surface = 0;
66  for (unsigned model = 0; model < _num_models; ++model)
67  {
68  _surfaces_given_model[model].resize(_f[model]->numberSurfaces());
69  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
70  _surfaces_given_model[model][model_surface] = surface++;
71  }
72 
73  // check the plastic_models for specialIC
74  if (_specialIC == "rock")
75  {
76  if (_num_models != 2)
77  mooseError("Choosing specialIC=rock, you must have plasticity models of type 'TensileMulti "
78  "MohrCoulombMulti'\n");
79  if (!(_f[0]->modelName().compare("TensileMulti") == 0 &&
80  _f[1]->modelName().compare("MohrCoulombMulti") == 0))
81  mooseError("Choosing specialIC=rock, you must have plasticity models of type 'TensileMulti "
82  "MohrCoulombMulti'\n");
83  }
84  if (_specialIC == "joint")
85  {
86  if (_num_models != 2)
87  mooseError("Choosing specialIC=joint, you must have plasticity models of type "
88  "'WeakPlaneTensile WeakPlaneShear'\n");
89  if (!(_f[0]->modelName().compare("WeakPlaneTensile") == 0 &&
90  _f[1]->modelName().compare("WeakPlaneShear") == 0))
91  mooseError("Choosing specialIC=joint, you must have plasticity models of type "
92  "'WeakPlaneTensile WeakPlaneShear'\n");
93  }
94 }
std::vector< std::vector< unsigned int > > _surfaces_given_model
_surfaces_given_model[model_number] = vector of surface numbers for this model
std::vector< unsigned int > _model_surface_given_surface
given a surface number, this returns the corresponding-model&#39;s internal surface number ...
MooseEnum _specialIC
Allows initial set of active constraints to be chosen optimally.
std::vector< unsigned int > _model_given_surface
given a surface number, this returns the model number
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ ~MultiPlasticityRawComponentAssembler()

virtual MultiPlasticityRawComponentAssembler::~MultiPlasticityRawComponentAssembler ( )
inlinevirtual

Definition at line 49 of file MultiPlasticityRawComponentAssembler.h.

49 {}

Member Function Documentation

◆ activeModelSurfaces()

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

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 810 of file MultiPlasticityRawComponentAssembler.C.

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

814 {
815  active_surfaces_of_model.resize(0);
816  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
817  if (active[_surfaces_given_model[model][model_surface]])
818  active_surfaces_of_model.push_back(model_surface);
819 }
std::vector< std::vector< unsigned int > > _surfaces_given_model
_surfaces_given_model[model_number] = vector of surface numbers for this model
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.

◆ activeSurfaces()

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

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 799 of file MultiPlasticityRawComponentAssembler.C.

Referenced by MultiPlasticityLinearSystem::calculateConstraints().

802 {
803  active_surfaces.resize(0);
804  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
805  if (active[_surfaces_given_model[model][model_surface]])
806  active_surfaces.push_back(_surfaces_given_model[model][model_surface]);
807 }
std::vector< std::vector< unsigned int > > _surfaces_given_model
_surfaces_given_model[model_number] = vector of surface numbers for this model
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.

◆ anyActiveSurfaces()

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

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

Definition at line 790 of file MultiPlasticityRawComponentAssembler.C.

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

791 {
792  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
793  if (active[_surfaces_given_model[model][model_surface]])
794  return true;
795  return false;
796 }
std::vector< std::vector< unsigned int > > _surfaces_given_model
_surfaces_given_model[model_number] = vector of surface numbers for this model
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.

◆ buildActiveConstraints()

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

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 343 of file MultiPlasticityRawComponentAssembler.C.

Referenced by ComputeMultiPlasticityStress::returnMap().

348 {
349  mooseAssert(f.size() == _num_surfaces,
350  "buildActiveConstraints called with f.size = " << f.size() << " while there are "
351  << _num_surfaces << " surfaces");
352  mooseAssert(intnl.size() == _num_models,
353  "buildActiveConstraints called with intnl.size = "
354  << intnl.size() << " while there are " << _num_models << " models");
355 
356  if (_specialIC == "rock")
357  buildActiveConstraintsRock(f, stress, intnl, Eijkl, act);
358  else if (_specialIC == "joint")
359  buildActiveConstraintsJoint(f, stress, intnl, Eijkl, act);
360  else // no specialIC
361  {
362  act.resize(0);
363  unsigned ind = 0;
364  for (unsigned model = 0; model < _num_models; ++model)
365  {
366  std::vector<Real> model_f(0);
367  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
368  model_f.push_back(f[ind++]);
369  std::vector<bool> model_act;
370  RankTwoTensor returned_stress;
371  _f[model]->activeConstraints(
372  model_f, stress, intnl[model], Eijkl, model_act, returned_stress);
373  for (unsigned model_surface = 0; model_surface < _f[model]->numberSurfaces(); ++model_surface)
374  act.push_back(model_act[model_surface]);
375  }
376  }
377 }
MooseEnum _specialIC
Allows initial set of active constraints to be chosen optimally.
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...
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.
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ buildActiveConstraintsJoint()

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

"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 380 of file MultiPlasticityRawComponentAssembler.C.

Referenced by buildActiveConstraints().

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

◆ buildActiveConstraintsRock()

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

"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 421 of file MultiPlasticityRawComponentAssembler.C.

Referenced by buildActiveConstraints().

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

◆ dflowPotential_dintnl()

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

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 234 of file MultiPlasticityRawComponentAssembler.C.

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

238 {
239  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
240  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
241 
242  dr_dintnl.resize(0);
243  std::vector<unsigned int> active_surfaces_of_model;
244  std::vector<unsigned int>::iterator active_surface;
245  std::vector<RankTwoTensor> model_dr_dintnl;
246  for (unsigned model = 0; model < _num_models; ++model)
247  {
248  activeModelSurfaces(model, active, active_surfaces_of_model);
249  if (active_surfaces_of_model.size() > 0)
250  {
251  _f[model]->dflowPotential_dintnlV(stress, intnl[model], model_dr_dintnl);
252  for (active_surface = active_surfaces_of_model.begin();
253  active_surface != active_surfaces_of_model.end();
254  ++active_surface)
255  dr_dintnl.push_back(model_dr_dintnl[*active_surface]);
256  }
257  }
258 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ dflowPotential_dstress()

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

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 206 of file MultiPlasticityRawComponentAssembler.C.

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

211 {
212  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
213  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
214 
215  dr_dstress.resize(0);
216  std::vector<unsigned int> active_surfaces_of_model;
217  std::vector<unsigned int>::iterator active_surface;
218  std::vector<RankFourTensor> model_dr_dstress;
219  for (unsigned model = 0; model < _num_models; ++model)
220  {
221  activeModelSurfaces(model, active, active_surfaces_of_model);
222  if (active_surfaces_of_model.size() > 0)
223  {
224  _f[model]->dflowPotential_dstressV(stress, intnl[model], model_dr_dstress);
225  for (active_surface = active_surfaces_of_model.begin();
226  active_surface != active_surfaces_of_model.end();
227  ++active_surface)
228  dr_dstress.push_back(model_dr_dstress[*active_surface]);
229  }
230  }
231 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ dhardPotential_dintnl()

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

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 316 of file MultiPlasticityRawComponentAssembler.C.

Referenced by MultiPlasticityLinearSystem::calculateJacobian().

320 {
321  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
322  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
323 
324  dh_dintnl.resize(0);
325  std::vector<unsigned int> active_surfaces_of_model;
326  std::vector<unsigned int>::iterator active_surface;
327  std::vector<Real> model_dh_dintnl;
328  for (unsigned model = 0; model < _num_models; ++model)
329  {
330  activeModelSurfaces(model, active, active_surfaces_of_model);
331  if (active_surfaces_of_model.size() > 0)
332  {
333  _f[model]->dhardPotential_dintnlV(stress, intnl[model], model_dh_dintnl);
334  for (active_surface = active_surfaces_of_model.begin();
335  active_surface != active_surfaces_of_model.end();
336  ++active_surface)
337  dh_dintnl.push_back(model_dh_dintnl[*active_surface]);
338  }
339  }
340 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ dhardPotential_dstress()

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

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 288 of file MultiPlasticityRawComponentAssembler.C.

Referenced by MultiPlasticityLinearSystem::calculateJacobian().

293 {
294  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
295  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
296 
297  dh_dstress.resize(0);
298  std::vector<unsigned int> active_surfaces_of_model;
299  std::vector<unsigned int>::iterator active_surface;
300  std::vector<RankTwoTensor> model_dh_dstress;
301  for (unsigned model = 0; model < _num_models; ++model)
302  {
303  activeModelSurfaces(model, active, active_surfaces_of_model);
304  if (active_surfaces_of_model.size() > 0)
305  {
306  _f[model]->dhardPotential_dstressV(stress, intnl[model], model_dh_dstress);
307  for (active_surface = active_surfaces_of_model.begin();
308  active_surface != active_surfaces_of_model.end();
309  ++active_surface)
310  dh_dstress.push_back(model_dh_dstress[*active_surface]);
311  }
312  }
313 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ dyieldFunction_dintnl()

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

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 152 of file MultiPlasticityRawComponentAssembler.C.

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

156 {
157  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
158  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
159 
160  df_dintnl.resize(0);
161  std::vector<unsigned int> active_surfaces_of_model;
162  std::vector<unsigned int>::iterator active_surface;
163  std::vector<Real> model_df_dintnl;
164  for (unsigned model = 0; model < _num_models; ++model)
165  {
166  activeModelSurfaces(model, active, active_surfaces_of_model);
167  if (active_surfaces_of_model.size() > 0)
168  {
169  _f[model]->dyieldFunction_dintnlV(stress, intnl[model], model_df_dintnl);
170  for (active_surface = active_surfaces_of_model.begin();
171  active_surface != active_surfaces_of_model.end();
172  ++active_surface)
173  df_dintnl.push_back(model_df_dintnl[*active_surface]);
174  }
175  }
176 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ dyieldFunction_dstress()

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

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 124 of file MultiPlasticityRawComponentAssembler.C.

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

129 {
130  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
131  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
132 
133  df_dstress.resize(0);
134  std::vector<unsigned int> active_surfaces_of_model;
135  std::vector<unsigned int>::iterator active_surface;
136  std::vector<RankTwoTensor> model_df_dstress;
137  for (unsigned model = 0; model < _num_models; ++model)
138  {
139  activeModelSurfaces(model, active, active_surfaces_of_model);
140  if (active_surfaces_of_model.size() > 0)
141  {
142  _f[model]->dyieldFunction_dstressV(stress, intnl[model], model_df_dstress);
143  for (active_surface = active_surfaces_of_model.begin();
144  active_surface != active_surfaces_of_model.end();
145  ++active_surface)
146  df_dstress.push_back(model_df_dstress[*active_surface]);
147  }
148  }
149 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ flowPotential()

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

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 179 of file MultiPlasticityRawComponentAssembler.C.

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

183 {
184  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
185  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
186 
187  r.resize(0);
188  std::vector<unsigned int> active_surfaces_of_model;
189  std::vector<unsigned int>::iterator active_surface;
190  std::vector<RankTwoTensor> model_r;
191  for (unsigned model = 0; model < _num_models; ++model)
192  {
193  activeModelSurfaces(model, active, active_surfaces_of_model);
194  if (active_surfaces_of_model.size() > 0)
195  {
196  _f[model]->flowPotentialV(stress, intnl[model], model_r);
197  for (active_surface = active_surfaces_of_model.begin();
198  active_surface != active_surfaces_of_model.end();
199  ++active_surface)
200  r.push_back(model_r[*active_surface]);
201  }
202  }
203 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ hardPotential()

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

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 261 of file MultiPlasticityRawComponentAssembler.C.

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

265 {
266  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
267  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
268 
269  h.resize(0);
270  std::vector<unsigned int> active_surfaces_of_model;
271  std::vector<unsigned int>::iterator active_surface;
272  std::vector<Real> model_h;
273  for (unsigned model = 0; model < _num_models; ++model)
274  {
275  activeModelSurfaces(model, active, active_surfaces_of_model);
276  if (active_surfaces_of_model.size() > 0)
277  {
278  _f[model]->hardPotentialV(stress, intnl[model], model_h);
279  for (active_surface = active_surfaces_of_model.begin();
280  active_surface != active_surfaces_of_model.end();
281  ++active_surface)
282  h.push_back(model_h[*active_surface]);
283  }
284  }
285 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

◆ modelNumber()

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

◆ 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 
)
protected

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 598 of file MultiPlasticityRawComponentAssembler.C.

Referenced by ComputeMultiPlasticityStress::quickStep().

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

◆ yieldFunction()

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

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 97 of file MultiPlasticityRawComponentAssembler.C.

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

101 {
102  mooseAssert(intnl.size() == _num_models, "Incorrect size of internal parameters");
103  mooseAssert(active.size() == _num_surfaces, "Incorrect size of active");
104 
105  f.resize(0);
106  std::vector<unsigned int> active_surfaces_of_model;
107  std::vector<unsigned int>::iterator active_surface;
108  std::vector<Real> model_f;
109  for (unsigned model = 0; model < _num_models; ++model)
110  {
111  activeModelSurfaces(model, active, active_surfaces_of_model);
112  if (active_surfaces_of_model.size() > 0)
113  {
114  _f[model]->yieldFunctionV(stress, intnl[model], model_f);
115  for (active_surface = active_surfaces_of_model.begin();
116  active_surface != active_surfaces_of_model.end();
117  ++active_surface)
118  f.push_back(model_f[*active_surface]);
119  }
120  }
121 }
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=...
unsigned int _num_surfaces
Number of surfaces within the plastic models.
std::vector< const TensorMechanicsPlasticModel * > _f
User objects that define the yield functions, flow potentials, etc.
unsigned int _num_models
Number of plastic models for this material.

Member Data Documentation

◆ _f

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

◆ _model_given_surface

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

given a surface number, this returns the model number

Definition at line 293 of file MultiPlasticityRawComponentAssembler.h.

Referenced by modelNumber(), and MultiPlasticityRawComponentAssembler().

◆ _model_surface_given_surface

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

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

Definition at line 296 of file MultiPlasticityRawComponentAssembler.h.

Referenced by MultiPlasticityRawComponentAssembler().

◆ _num_models

unsigned int MultiPlasticityRawComponentAssembler::_num_models
protected

◆ _num_surfaces

unsigned int MultiPlasticityRawComponentAssembler::_num_surfaces
protected

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 64 of file MultiPlasticityRawComponentAssembler.h.

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

◆ _params

const InputParameters& MultiPlasticityRawComponentAssembler::_params
protected

◆ _specialIC

MooseEnum MultiPlasticityRawComponentAssembler::_specialIC
protected

Allows initial set of active constraints to be chosen optimally.

Definition at line 70 of file MultiPlasticityRawComponentAssembler.h.

Referenced by buildActiveConstraints(), and MultiPlasticityRawComponentAssembler().

◆ _surfaces_given_model

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

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

Definition at line 67 of file MultiPlasticityRawComponentAssembler.h.

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


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