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

Rate-independent associative mean-cap tensile AND compressive failure with hardening/softening of the tensile and compressive strength. More...

#include <TensorMechanicsPlasticMeanCapTC.h>

Inheritance diagram for TensorMechanicsPlasticMeanCapTC:
[legend]

Public Member Functions

 TensorMechanicsPlasticMeanCapTC (const InputParameters &parameters)
 
virtual void activeConstraints (const std::vector< Real > &f, const RankTwoTensor &stress, Real intnl, const RankFourTensor &Eijkl, std::vector< bool > &act, RankTwoTensor &returned_stress) const override
 The active yield surfaces, given a vector of yield functions. More...
 
virtual bool useCustomReturnMap () const override
 Returns false. You will want to override this in your derived class if you write a custom returnMap function. More...
 
virtual bool useCustomCTO () const override
 Returns false. You will want to override this in your derived class if you write a custom consistent tangent operator function. More...
 
virtual bool returnMap (const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const override
 Performs a custom return-map. More...
 
virtual RankFourTensor consistentTangentOperator (const RankTwoTensor &trial_stress, Real intnl_old, const RankTwoTensor &stress, Real intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &cumulative_pm) const override
 Calculates a custom consistent tangent operator. More...
 
virtual std::string modelName () const override
 
void initialize ()
 
void execute ()
 
void finalize ()
 
virtual unsigned int numberSurfaces () const
 The number of yield surfaces for this plasticity model. More...
 
virtual void yieldFunctionV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &f) const
 Calculates the yield functions. More...
 
virtual void dyieldFunction_dstressV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &df_dstress) const
 The derivative of yield functions with respect to stress. More...
 
virtual void dyieldFunction_dintnlV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &df_dintnl) const
 The derivative of yield functions with respect to the internal parameter. More...
 
virtual void flowPotentialV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &r) const
 The flow potentials. More...
 
virtual void dflowPotential_dstressV (const RankTwoTensor &stress, Real intnl, std::vector< RankFourTensor > &dr_dstress) const
 The derivative of the flow potential with respect to stress. More...
 
virtual void dflowPotential_dintnlV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &dr_dintnl) const
 The derivative of the flow potential with respect to the internal parameter. More...
 
virtual void hardPotentialV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &h) const
 The hardening potential. More...
 
virtual void dhardPotential_dstressV (const RankTwoTensor &stress, Real intnl, std::vector< RankTwoTensor > &dh_dstress) const
 The derivative of the hardening potential with respect to stress. More...
 
virtual void dhardPotential_dintnlV (const RankTwoTensor &stress, Real intnl, std::vector< Real > &dh_dintnl) const
 The derivative of the hardening potential with respect to the internal parameter. More...
 
bool KuhnTuckerSingleSurface (Real yf, Real dpm, Real dpm_tol) const
 Returns true if the Kuhn-Tucker conditions for the single surface are satisfied. More...
 

Public Attributes

const Real _f_tol
 Tolerance on yield function. More...
 
const Real _ic_tol
 Tolerance on internal constraint. More...
 

Protected Member Functions

Real yieldFunction (const RankTwoTensor &stress, Real intnl) const override
 The following functions are what you should override when building single-plasticity models. More...
 
RankTwoTensor dyieldFunction_dstress (const RankTwoTensor &stress, Real intnl) const override
 The derivative of yield function with respect to stress. More...
 
Real dyieldFunction_dintnl (const RankTwoTensor &stress, Real intnl) const override
 The derivative of yield function with respect to the internal parameter. More...
 
RankTwoTensor flowPotential (const RankTwoTensor &stress, Real intnl) const override
 The flow potential. More...
 
RankFourTensor dflowPotential_dstress (const RankTwoTensor &stress, Real intnl) const override
 The derivative of the flow potential with respect to stress. More...
 
RankTwoTensor dflowPotential_dintnl (const RankTwoTensor &stress, Real intnl) const override
 The derivative of the flow potential with respect to the internal parameter. More...
 
Real hardPotential (const RankTwoTensor &stress, Real intnl) const override
 The hardening potential. More...
 
virtual RankTwoTensor dhardPotential_dstress (const RankTwoTensor &stress, Real intnl) const override
 The derivative of the hardening potential with respect to stress. More...
 
virtual Real dhardPotential_dintnl (const RankTwoTensor &stress, Real intnl) const override
 The derivative of the hardening potential with respect to the internal parameter. More...
 
RankTwoTensor df_dsig (const RankTwoTensor &stress, Real intnl) const
 Derivative of the yield function with respect to stress. More...
 
virtual Real tensile_strength (const Real internal_param) const
 tensile strength as a function of residual value, rate, and internal_param More...
 
virtual Real dtensile_strength (const Real internal_param) const
 d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param More...
 
virtual Real compressive_strength (const Real internal_param) const
 compressive strength as a function of residual value, rate, and internal_param More...
 
virtual Real dcompressive_strength (const Real internal_param) const
 d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param More...
 

Protected Attributes

const unsigned _max_iters
 max iters for custom return map loop More...
 
const bool _use_custom_returnMap
 Whether to use the custom return-map algorithm. More...
 
const bool _use_custom_cto
 Whether to use the custom consistent tangent operator algorithm. More...
 
const TensorMechanicsHardeningModel_strength
 the tensile strength More...
 
const TensorMechanicsHardeningModel_c_strength
 the compressive strength More...
 

Detailed Description

Rate-independent associative mean-cap tensile AND compressive failure with hardening/softening of the tensile and compressive strength.

The key point here is that the internal parameter is equal to the volumetric plastic strain. This means that upon tensile failure, the compressive strength can soften (using a TensorMechanicsHardening object) which physically means that a subsequent compressive stress can easily squash the material

Definition at line 30 of file TensorMechanicsPlasticMeanCapTC.h.

Constructor & Destructor Documentation

◆ TensorMechanicsPlasticMeanCapTC()

TensorMechanicsPlasticMeanCapTC::TensorMechanicsPlasticMeanCapTC ( const InputParameters &  parameters)

Definition at line 48 of file TensorMechanicsPlasticMeanCapTC.C.

49  : TensorMechanicsPlasticModel(parameters),
50  _max_iters(getParam<unsigned>("max_iterations")),
51  _use_custom_returnMap(getParam<bool>("use_custom_returnMap")),
52  _use_custom_cto(getParam<bool>("use_custom_cto")),
53  _strength(getUserObject<TensorMechanicsHardeningModel>("tensile_strength")),
54  _c_strength(getUserObject<TensorMechanicsHardeningModel>("compressive_strength"))
55 {
56  // cannot check the following for all values of the internal parameter, but this will catch most
57  // errors
58  if (_strength.value(0) <= _c_strength.value(0))
59  mooseError("MeanCapTC: tensile strength (which is usually positive) must not be less than "
60  "compressive strength (which is usually negative)");
61 }
TensorMechanicsPlasticModel(const InputParameters &parameters)
const TensorMechanicsHardeningModel & _c_strength
the compressive strength
const bool _use_custom_returnMap
Whether to use the custom return-map algorithm.
const unsigned _max_iters
max iters for custom return map loop
const bool _use_custom_cto
Whether to use the custom consistent tangent operator algorithm.
virtual Real value(Real intnl) const
const TensorMechanicsHardeningModel & _strength
the tensile strength

Member Function Documentation

◆ activeConstraints()

void TensorMechanicsPlasticMeanCapTC::activeConstraints ( const std::vector< Real > &  f,
const RankTwoTensor &  stress,
Real  intnl,
const RankFourTensor &  Eijkl,
std::vector< bool > &  act,
RankTwoTensor &  returned_stress 
) const
overridevirtual

The active yield surfaces, given a vector of yield functions.

This is used by FiniteStrainMultiPlasticity to determine the initial set of active constraints at the trial (stress, intnl) configuration. It is up to you (the coder) to determine how accurate you want the returned_stress to be. Currently it is only used by FiniteStrainMultiPlasticity to estimate a good starting value for the Newton-Rahson procedure, so currently it may not need to be super perfect.

Parameters
fvalues of the yield functions
stressstress tensor
intnlinternal parameter
Eijklelasticity tensor (stress = Eijkl*strain)
[out]actact[i] = true if the i_th yield function is active
[out]returned_stressApproximate value of the returned stress

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 246 of file TensorMechanicsPlasticMeanCapTC.C.

252 {
253  act.assign(1, false);
254 
255  if (f[0] <= _f_tol)
256  {
257  returned_stress = stress;
258  return;
259  }
260 
261  const Real tr = stress.trace();
262  const Real t_str = tensile_strength(intnl);
263  Real str;
264  Real dirn;
265  if (tr >= t_str)
266  {
267  str = t_str;
268  dirn = 1.0;
269  }
270  else
271  {
272  str = compressive_strength(intnl);
273  dirn = -1.0;
274  }
275 
276  RankTwoTensor n; // flow direction
277  for (unsigned i = 0; i < 3; ++i)
278  for (unsigned j = 0; j < 3; ++j)
279  for (unsigned k = 0; k < 3; ++k)
280  n(i, j) += dirn * Eijkl(i, j, k, k);
281 
282  // returned_stress = stress - gamma*n
283  // and taking the trace of this and using
284  // Tr(returned_stress) = str, gives
285  // gamma = (Tr(stress) - str)/Tr(n)
286  Real gamma = (stress.trace() - str) / n.trace();
287 
288  for (unsigned i = 0; i < 3; ++i)
289  for (unsigned j = 0; j < 3; ++j)
290  returned_stress(i, j) = stress(i, j) - gamma * n(i, j);
291 
292  act[0] = true;
293 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
const Real _f_tol
Tolerance on yield function.

◆ compressive_strength()

Real TensorMechanicsPlasticMeanCapTC::compressive_strength ( const Real  internal_param) const
protectedvirtual

compressive strength as a function of residual value, rate, and internal_param

Definition at line 234 of file TensorMechanicsPlasticMeanCapTC.C.

Referenced by activeConstraints(), df_dsig(), dflowPotential_dintnl(), dflowPotential_dstress(), dhardPotential_dintnl(), dhardPotential_dstress(), dyieldFunction_dintnl(), hardPotential(), returnMap(), and yieldFunction().

235 {
236  return _c_strength.value(internal_param);
237 }
const TensorMechanicsHardeningModel & _c_strength
the compressive strength
virtual Real value(Real intnl) const

◆ consistentTangentOperator()

RankFourTensor TensorMechanicsPlasticMeanCapTC::consistentTangentOperator ( const RankTwoTensor &  trial_stress,
Real  intnl_old,
const RankTwoTensor &  stress,
Real  intnl,
const RankFourTensor &  E_ijkl,
const std::vector< Real > &  cumulative_pm 
) const
overridevirtual

Calculates a custom consistent tangent operator.

You may choose to over-ride this in your derived TensorMechanicsPlasticXXXX class.

(Note, if you over-ride returnMap, you will probably want to override consistentTangentOpertor too, otherwise it will default to E_ijkl.)

Parameters
stress_oldtrial stress before returning
intnl_oldinternal parameter before returning
stresscurrent returned stress state
intnlinternal parameter
E_ijklelasticity tensor
cumulative_pmthe cumulative plastic multipliers
Returns
the consistent tangent operator: E_ijkl if not over-ridden

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 398 of file TensorMechanicsPlasticMeanCapTC.C.

405 {
406  if (!_use_custom_cto)
408  trial_stress, intnl_old, stress, intnl, E_ijkl, cumulative_pm);
409 
410  Real df_dq;
411  Real alpha;
412  if (trial_stress.trace() >= tensile_strength(intnl_old))
413  {
414  df_dq = -dtensile_strength(intnl);
415  alpha = 1.0;
416  }
417  else
418  {
419  df_dq = dcompressive_strength(intnl);
420  alpha = -1.0;
421  }
422 
423  RankTwoTensor elas;
424  for (unsigned int i = 0; i < 3; ++i)
425  for (unsigned int j = 0; j < 3; ++j)
426  for (unsigned int k = 0; k < 3; ++k)
427  elas(i, j) += E_ijkl(i, j, k, k);
428 
429  const Real hw = -df_dq + alpha * elas.trace();
430 
431  return E_ijkl - alpha / hw * elas.outerProduct(elas);
432 }
const bool _use_custom_cto
Whether to use the custom consistent tangent operator algorithm.
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
virtual RankFourTensor consistentTangentOperator(const RankTwoTensor &trial_stress, Real intnl_old, const RankTwoTensor &stress, Real intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &cumulative_pm) const
Calculates a custom consistent tangent operator.
virtual Real dtensile_strength(const Real internal_param) const
d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param ...
virtual Real dcompressive_strength(const Real internal_param) const
d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param ...

◆ dcompressive_strength()

Real TensorMechanicsPlasticMeanCapTC::dcompressive_strength ( const Real  internal_param) const
protectedvirtual

d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param

Definition at line 240 of file TensorMechanicsPlasticMeanCapTC.C.

Referenced by consistentTangentOperator(), dflowPotential_dintnl(), dhardPotential_dintnl(), dyieldFunction_dintnl(), and returnMap().

241 {
242  return _c_strength.derivative(internal_param);
243 }
const TensorMechanicsHardeningModel & _c_strength
the compressive strength
virtual Real derivative(Real intnl) const

◆ df_dsig()

RankTwoTensor TensorMechanicsPlasticMeanCapTC::df_dsig ( const RankTwoTensor &  stress,
Real  intnl 
) const
protected

Derivative of the yield function with respect to stress.

This is also the flow potential.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
Returns
the derivative

Definition at line 109 of file TensorMechanicsPlasticMeanCapTC.C.

Referenced by dyieldFunction_dstress(), and flowPotential().

110 {
111  const Real tr = stress.trace();
112  const Real t_str = tensile_strength(intnl);
113  if (tr >= t_str)
114  return stress.dtrace();
115 
116  const Real c_str = compressive_strength(intnl);
117  if (tr <= c_str)
118  return -stress.dtrace();
119 
120  return -std::cos(libMesh::pi * (tr - c_str) / (t_str - c_str)) * stress.dtrace();
121 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param

◆ dflowPotential_dintnl()

RankTwoTensor TensorMechanicsPlasticMeanCapTC::dflowPotential_dintnl ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The derivative of the flow potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
Returns
dr_dintnl(i, j) = dr(i, j)/dintnl

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 147 of file TensorMechanicsPlasticMeanCapTC.C.

149 {
150  const Real tr = stress.trace();
151  const Real t_str = tensile_strength(intnl);
152  if (tr >= t_str)
153  return RankTwoTensor();
154 
155  const Real c_str = compressive_strength(intnl);
156  if (tr <= c_str)
157  return RankTwoTensor();
158 
159  const Real dt = dtensile_strength(intnl);
160  const Real dc = dcompressive_strength(intnl);
161  return std::sin(libMesh::pi * (tr - c_str) / (t_str - c_str)) * stress.dtrace() * libMesh::pi /
162  Utility::pow<2>(t_str - c_str) * ((tr - t_str) * dc - (tr - c_str) * dt);
163 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
virtual Real dtensile_strength(const Real internal_param) const
d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param ...
virtual Real dcompressive_strength(const Real internal_param) const
d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param ...

◆ dflowPotential_dintnlV()

void TensorMechanicsPlasticModel::dflowPotential_dintnlV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< RankTwoTensor > &  dr_dintnl 
) const
virtualinherited

The derivative of the flow potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
[out]dr_dintnldr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl

Reimplemented in TensorMechanicsPlasticMohrCoulombMulti, and TensorMechanicsPlasticTensileMulti.

Definition at line 138 of file TensorMechanicsPlasticModel.C.

141 {
142  return dr_dintnl.assign(1, dflowPotential_dintnl(stress, intnl));
143 }
virtual RankTwoTensor dflowPotential_dintnl(const RankTwoTensor &stress, Real intnl) const
The derivative of the flow potential with respect to the internal parameter.

◆ dflowPotential_dstress()

RankFourTensor TensorMechanicsPlasticMeanCapTC::dflowPotential_dstress ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The derivative of the flow potential with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
Returns
dr_dstress(i, j, k, l) = dr(i, j)/dstress(k, l)

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 130 of file TensorMechanicsPlasticMeanCapTC.C.

132 {
133  const Real tr = stress.trace();
134  const Real t_str = tensile_strength(intnl);
135  if (tr >= t_str)
136  return RankFourTensor();
137 
138  const Real c_str = compressive_strength(intnl);
139  if (tr <= c_str)
140  return RankFourTensor();
141 
142  return libMesh::pi / (t_str - c_str) * std::sin(libMesh::pi * (tr - c_str) / (t_str - c_str)) *
143  stress.dtrace().outerProduct(stress.dtrace());
144 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param

◆ dflowPotential_dstressV()

void TensorMechanicsPlasticModel::dflowPotential_dstressV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< RankFourTensor > &  dr_dstress 
) const
virtualinherited

The derivative of the flow potential with respect to stress.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
[out]dr_dstressdr_dstress[alpha](i, j, k, l) = dr[alpha](i, j)/dstress(k, l)

Reimplemented in TensorMechanicsPlasticMohrCoulombMulti, and TensorMechanicsPlasticTensileMulti.

Definition at line 124 of file TensorMechanicsPlasticModel.C.

127 {
128  return dr_dstress.assign(1, dflowPotential_dstress(stress, intnl));
129 }
virtual RankFourTensor dflowPotential_dstress(const RankTwoTensor &stress, Real intnl) const
The derivative of the flow potential with respect to stress.

◆ dhardPotential_dintnl()

Real TensorMechanicsPlasticMeanCapTC::dhardPotential_dintnl ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The derivative of the hardening potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
Returns
the derivative

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 203 of file TensorMechanicsPlasticMeanCapTC.C.

205 {
206  const Real tr = stress.trace();
207  const Real t_str = tensile_strength(intnl);
208  if (tr >= t_str)
209  return 0.0;
210 
211  const Real c_str = compressive_strength(intnl);
212  if (tr <= c_str)
213  return 0.0;
214 
215  const Real dt = dtensile_strength(intnl);
216  const Real dc = dcompressive_strength(intnl);
217  return -std::sin(libMesh::pi * (tr - c_str) / (t_str - c_str)) * libMesh::pi /
218  Utility::pow<2>(t_str - c_str) * ((tr - t_str) * dc - (tr - c_str) * dt);
219 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
virtual Real dtensile_strength(const Real internal_param) const
d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param ...
virtual Real dcompressive_strength(const Real internal_param) const
d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param ...

◆ dhardPotential_dintnlV()

void TensorMechanicsPlasticModel::dhardPotential_dintnlV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< Real > &  dh_dintnl 
) const
virtualinherited

The derivative of the hardening potential with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
[out]dh_dintnldh_dintnl[alpha] = dh[alpha]/dintnl

Definition at line 179 of file TensorMechanicsPlasticModel.C.

182 {
183  dh_dintnl.resize(numberSurfaces(), dhardPotential_dintnl(stress, intnl));
184 }
virtual unsigned int numberSurfaces() const
The number of yield surfaces for this plasticity model.
virtual Real dhardPotential_dintnl(const RankTwoTensor &stress, Real intnl) const
The derivative of the hardening potential with respect to the internal parameter. ...

◆ dhardPotential_dstress()

RankTwoTensor TensorMechanicsPlasticMeanCapTC::dhardPotential_dstress ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The derivative of the hardening potential with respect to stress.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
Returns
dh_dstress(i, j) = dh/dstress(i, j)

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 186 of file TensorMechanicsPlasticMeanCapTC.C.

188 {
189  const Real tr = stress.trace();
190  const Real t_str = tensile_strength(intnl);
191  if (tr >= t_str)
192  return RankTwoTensor();
193 
194  const Real c_str = compressive_strength(intnl);
195  if (tr <= c_str)
196  return RankTwoTensor();
197 
198  return -std::sin(libMesh::pi * (tr - c_str) / (t_str - c_str)) * libMesh::pi / (t_str - c_str) *
199  stress.dtrace();
200 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param

◆ dhardPotential_dstressV()

void TensorMechanicsPlasticModel::dhardPotential_dstressV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< RankTwoTensor > &  dh_dstress 
) const
virtualinherited

The derivative of the hardening potential with respect to stress.

Parameters
stressthe stress at which to calculate the hardening potentials
intnlinternal parameter
[out]dh_dstressdh_dstress[alpha](i, j) = dh[alpha]/dstress(i, j)

Definition at line 165 of file TensorMechanicsPlasticModel.C.

168 {
169  dh_dstress.assign(numberSurfaces(), dhardPotential_dstress(stress, intnl));
170 }
virtual unsigned int numberSurfaces() const
The number of yield surfaces for this plasticity model.
virtual RankTwoTensor dhardPotential_dstress(const RankTwoTensor &stress, Real intnl) const
The derivative of the hardening potential with respect to stress.

◆ dtensile_strength()

Real TensorMechanicsPlasticMeanCapTC::dtensile_strength ( const Real  internal_param) const
protectedvirtual

d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param

Definition at line 228 of file TensorMechanicsPlasticMeanCapTC.C.

Referenced by consistentTangentOperator(), dflowPotential_dintnl(), dhardPotential_dintnl(), dyieldFunction_dintnl(), and returnMap().

229 {
230  return _strength.derivative(internal_param);
231 }
virtual Real derivative(Real intnl) const
const TensorMechanicsHardeningModel & _strength
the tensile strength

◆ dyieldFunction_dintnl()

Real TensorMechanicsPlasticMeanCapTC::dyieldFunction_dintnl ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The derivative of yield function with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
Returns
the derivative

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 89 of file TensorMechanicsPlasticMeanCapTC.C.

91 {
92  const Real tr = stress.trace();
93  const Real t_str = tensile_strength(intnl);
94  if (tr >= t_str)
95  return -dtensile_strength(intnl);
96 
97  const Real c_str = compressive_strength(intnl);
98  if (tr <= c_str)
99  return dcompressive_strength(intnl);
100 
101  const Real dt = dtensile_strength(intnl);
102  const Real dc = dcompressive_strength(intnl);
103  return (dc - dt) / libMesh::pi * std::sin(libMesh::pi * (tr - c_str) / (t_str - c_str)) +
104  1.0 / (t_str - c_str) * std::cos(libMesh::pi * (tr - c_str) / (t_str - c_str)) *
105  ((tr - c_str) * dt - (tr - t_str) * dc);
106 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
virtual Real dtensile_strength(const Real internal_param) const
d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param ...
virtual Real dcompressive_strength(const Real internal_param) const
d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param ...

◆ dyieldFunction_dintnlV()

void TensorMechanicsPlasticModel::dyieldFunction_dintnlV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< Real > &  df_dintnl 
) const
virtualinherited

The derivative of yield functions with respect to the internal parameter.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
[out]df_dintnldf_dintnl[alpha] = df[alpha]/dintnl

Reimplemented in TensorMechanicsPlasticMohrCoulombMulti, and TensorMechanicsPlasticTensileMulti.

Definition at line 97 of file TensorMechanicsPlasticModel.C.

100 {
101  return df_dintnl.assign(1, dyieldFunction_dintnl(stress, intnl));
102 }
virtual Real dyieldFunction_dintnl(const RankTwoTensor &stress, Real intnl) const
The derivative of yield function with respect to the internal parameter.

◆ dyieldFunction_dstress()

RankTwoTensor TensorMechanicsPlasticMeanCapTC::dyieldFunction_dstress ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The derivative of yield function with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
Returns
df_dstress(i, j) = dyieldFunction/dstress(i, j)

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 82 of file TensorMechanicsPlasticMeanCapTC.C.

84 {
85  return df_dsig(stress, intnl);
86 }
RankTwoTensor df_dsig(const RankTwoTensor &stress, Real intnl) const
Derivative of the yield function with respect to stress.

◆ dyieldFunction_dstressV()

void TensorMechanicsPlasticModel::dyieldFunction_dstressV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< RankTwoTensor > &  df_dstress 
) const
virtualinherited

The derivative of yield functions with respect to stress.

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
[out]df_dstressdf_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)

Reimplemented in TensorMechanicsPlasticMohrCoulombMulti, and TensorMechanicsPlasticTensileMulti.

Definition at line 83 of file TensorMechanicsPlasticModel.C.

86 {
87  df_dstress.assign(1, dyieldFunction_dstress(stress, intnl));
88 }
virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor &stress, Real intnl) const
The derivative of yield function with respect to stress.

◆ execute()

void TensorMechanicsPlasticModel::execute ( )
inherited

Definition at line 46 of file TensorMechanicsPlasticModel.C.

47 {
48 }

◆ finalize()

void TensorMechanicsPlasticModel::finalize ( )
inherited

Definition at line 51 of file TensorMechanicsPlasticModel.C.

52 {
53 }

◆ flowPotential()

RankTwoTensor TensorMechanicsPlasticMeanCapTC::flowPotential ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The flow potential.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
Returns
the flow potential

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 124 of file TensorMechanicsPlasticMeanCapTC.C.

125 {
126  return df_dsig(stress, intnl);
127 }
RankTwoTensor df_dsig(const RankTwoTensor &stress, Real intnl) const
Derivative of the yield function with respect to stress.

◆ flowPotentialV()

void TensorMechanicsPlasticModel::flowPotentialV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< RankTwoTensor > &  r 
) const
virtualinherited

The flow potentials.

Parameters
stressthe stress at which to calculate the flow potential
intnlinternal parameter
[out]rr[alpha] is the flow potential for the "alpha" yield function

Reimplemented in TensorMechanicsPlasticMohrCoulombMulti, and TensorMechanicsPlasticTensileMulti.

Definition at line 110 of file TensorMechanicsPlasticModel.C.

113 {
114  return r.assign(1, flowPotential(stress, intnl));
115 }
virtual RankTwoTensor flowPotential(const RankTwoTensor &stress, Real intnl) const
The flow potential.

◆ hardPotential()

Real TensorMechanicsPlasticMeanCapTC::hardPotential ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The hardening potential.

Note that it is -1 for stress.trace() > _strength, and +1 for stress.trace() < _c_strength. This implements the idea that tensile failure will cause a massive reduction in compressive strength

Parameters
stressthe stress at which to calculate the hardening potential
intnlinternal parameter
Returns
the hardening potential

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 166 of file TensorMechanicsPlasticMeanCapTC.C.

167 {
168  // This is the key for this whole class!
169  const Real tr = stress.trace();
170  const Real t_str = tensile_strength(intnl);
171 
172  if (tr >= t_str)
173  return -1.0; // this will serve to *increase* the internal parameter (so internal parameter will
174  // be a measure of volumetric strain)
175 
176  const Real c_str = compressive_strength(intnl);
177  if (tr <= c_str)
178  return 1.0; // this will serve to *decrease* the internal parameter (so internal parameter will
179  // be a measure of volumetric strain)
180 
181  return std::cos(libMesh::pi * (tr - c_str) /
182  (t_str - c_str)); // this interpolates C2 smoothly between 1 and -1
183 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param

◆ hardPotentialV()

void TensorMechanicsPlasticModel::hardPotentialV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< Real > &  h 
) const
virtualinherited

The hardening potential.

Parameters
stressthe stress at which to calculate the hardening potential
intnlinternal parameter
[out]hh[alpha] is the hardening potential for the "alpha" yield function

Definition at line 151 of file TensorMechanicsPlasticModel.C.

154 {
155  h.assign(numberSurfaces(), hardPotential(stress, intnl));
156 }
virtual Real hardPotential(const RankTwoTensor &stress, Real intnl) const
The hardening potential.
virtual unsigned int numberSurfaces() const
The number of yield surfaces for this plasticity model.

◆ initialize()

void TensorMechanicsPlasticModel::initialize ( )
inherited

Definition at line 41 of file TensorMechanicsPlasticModel.C.

42 {
43 }

◆ KuhnTuckerSingleSurface()

bool TensorMechanicsPlasticModel::KuhnTuckerSingleSurface ( Real  yf,
Real  dpm,
Real  dpm_tol 
) const
inherited

Returns true if the Kuhn-Tucker conditions for the single surface are satisfied.

Parameters
yfYield function value
dpmplastic multiplier
dpm_toltolerance on plastic multiplier: viz dpm>-dpm_tol means "dpm is non-negative"

Definition at line 247 of file TensorMechanicsPlasticModel.C.

Referenced by TensorMechanicsPlasticMohrCoulombMulti::KuhnTuckerOK(), TensorMechanicsPlasticTensileMulti::KuhnTuckerOK(), and TensorMechanicsPlasticModel::returnMap().

248 {
249  return (dpm == 0 && yf <= _f_tol) || (dpm > -dpm_tol && yf <= _f_tol && yf >= -_f_tol);
250 }
const Real _f_tol
Tolerance on yield function.

◆ modelName()

std::string TensorMechanicsPlasticMeanCapTC::modelName ( ) const
overridevirtual

Implements TensorMechanicsPlasticModel.

Definition at line 447 of file TensorMechanicsPlasticMeanCapTC.C.

448 {
449  return "MeanCapTC";
450 }

◆ numberSurfaces()

unsigned TensorMechanicsPlasticModel::numberSurfaces ( ) const
virtualinherited

◆ returnMap()

bool TensorMechanicsPlasticMeanCapTC::returnMap ( const RankTwoTensor &  trial_stress,
Real  intnl_old,
const RankFourTensor &  E_ijkl,
Real  ep_plastic_tolerance,
RankTwoTensor &  returned_stress,
Real &  returned_intnl,
std::vector< Real > &  dpm,
RankTwoTensor &  delta_dp,
std::vector< Real > &  yf,
bool &  trial_stress_inadmissible 
) const
overridevirtual

Performs a custom return-map.

You may choose to over-ride this in your derived TensorMechanicsPlasticXXXX class, and you may implement the return-map algorithm in any way that suits you. Eg, using a Newton-Raphson approach, or a radial-return, etc. This may also be used as a quick way of ascertaining whether (trial_stress, intnl_old) is in fact admissible.

For over-riding this function, please note the following.

(1) Denoting the return value of the function by "successful_return", the only possible output values should be: (A) trial_stress_inadmissible=false, successful_return=true. That is, (trial_stress, intnl_old) is in fact admissible (in the elastic domain). (B) trial_stress_inadmissible=true, successful_return=false. That is (trial_stress, intnl_old) is inadmissible (outside the yield surface), and you didn't return to the yield surface. (C) trial_stress_inadmissible=true, successful_return=true. That is (trial_stress, intnl_old) is inadmissible (outside the yield surface), but you did return to the yield surface. The default implementation only handles case (A) and (B): it does not attempt to do a return-map algorithm.

(2) you must correctly signal "successful_return" using the return value of this function. Don't assume the calling function will do Kuhn-Tucker checking and so forth!

(3) In cases (A) and (B) you needn't set returned_stress, returned_intnl, delta_dp, or dpm. This is for computational efficiency.

(4) In cases (A) and (B), you MUST place the yield function values at (trial_stress, intnl_old) into yf so the calling function can use this information optimally. You will have already calculated these yield function values, which can be quite expensive, and it's not very optimal for the calling function to have to re-calculate them.

(5) In case (C), you need to set: returned_stress (the returned value of stress) returned_intnl (the returned value of the internal variable) delta_dp (the change in plastic strain) dpm (the plastic multipliers needed to bring about the return) yf (yield function values at the returned configuration)

(Note, if you over-ride returnMap, you will probably want to override consistentTangentOpertor too, otherwise it will default to E_ijkl.)

Parameters
trial_stressThe trial stress
intnl_oldValue of the internal parameter
E_ijklElasticity tensor
ep_plastic_toleranceTolerance defined by the user for the plastic strain
[out]returned_stressIn case (C): lies on the yield surface after returning and produces the correct plastic strain (normality condition). Otherwise: not defined
[out]returned_intnlIn case (C): the value of the internal parameter after returning. Otherwise: not defined
[out]dpmIn case (C): the plastic multipliers needed to bring about the return. Otherwise: not defined
[out]delta_dpIn case (C): The change in plastic strain induced by the return process. Otherwise: not defined
[out]yfIn case (C): the yield function at (returned_stress, returned_intnl). Otherwise: the yield function at (trial_stress, intnl_old)
[out]trial_stress_inadmissibleShould be set to false if the trial_stress is admissible, and true if the trial_stress is inadmissible. This can be used by the calling prorgram
Returns
true if a successful return (or a return-map not needed), false if the trial_stress is inadmissible but the return process failed

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 296 of file TensorMechanicsPlasticMeanCapTC.C.

306 {
307  if (!(_use_custom_returnMap))
308  return TensorMechanicsPlasticModel::returnMap(trial_stress,
309  intnl_old,
310  E_ijkl,
311  ep_plastic_tolerance,
312  returned_stress,
313  returned_intnl,
314  dpm,
315  delta_dp,
316  yf,
317  trial_stress_inadmissible);
318 
319  yf.resize(1);
320 
321  Real yf_orig = yieldFunction(trial_stress, intnl_old);
322 
323  yf[0] = yf_orig;
324 
325  if (yf_orig < _f_tol)
326  {
327  // the trial_stress is admissible
328  trial_stress_inadmissible = false;
329  return true;
330  }
331 
332  trial_stress_inadmissible = true;
333 
334  // In the following we want to solve
335  // trial_stress - stress = E_ijkl * dpm * r ...... (1)
336  // and either
337  // stress.trace() = tensile_strength(intnl) ...... (2a)
338  // intnl = intnl_old + dpm ...... (3a)
339  // or
340  // stress.trace() = compressive_strength(intnl) ... (2b)
341  // intnl = intnl_old - dpm ...... (3b)
342  // The former (2a and 3a) are chosen if
343  // trial_stress.trace() > tensile_strength(intnl_old)
344  // while the latter (2b and 3b) are chosen if
345  // trial_stress.trace() < compressive_strength(intnl_old)
346  // The variables we want to solve for are stress, dpm
347  // and intnl. We do this using a Newton approach, starting
348  // with stress=trial_stress and intnl=intnl_old and dpm=0
349  const bool tensile_failure = (trial_stress.trace() >= tensile_strength(intnl_old));
350  const Real dirn = (tensile_failure ? 1.0 : -1.0);
351 
352  RankTwoTensor n; // flow direction, which is E_ijkl * r
353  for (unsigned i = 0; i < 3; ++i)
354  for (unsigned j = 0; j < 3; ++j)
355  for (unsigned k = 0; k < 3; ++k)
356  n(i, j) += dirn * E_ijkl(i, j, k, k);
357  const Real n_trace = n.trace();
358 
359  // Perform a Newton-Raphson to find dpm when
360  // residual = trial_stress.trace() - tensile_strength(intnl) - dpm * n.trace() [for
361  // tensile_failure=true]
362  // or
363  // residual = trial_stress.trace() - compressive_strength(intnl) - dpm * n.trace() [for
364  // tensile_failure=false]
365  Real trial_trace = trial_stress.trace();
366  Real residual;
367  Real jac;
368  dpm[0] = 0;
369  unsigned int iter = 0;
370  do
371  {
372  if (tensile_failure)
373  {
374  residual = trial_trace - tensile_strength(intnl_old + dpm[0]) - dpm[0] * n_trace;
375  jac = -dtensile_strength(intnl_old + dpm[0]) - n_trace;
376  }
377  else
378  {
379  residual = trial_trace - compressive_strength(intnl_old - dpm[0]) - dpm[0] * n_trace;
380  jac = -dcompressive_strength(intnl_old - dpm[0]) - n_trace;
381  }
382  dpm[0] += -residual / jac;
383  if (iter > _max_iters) // not converging
384  return false;
385  iter++;
386  } while (residual * residual > _f_tol * _f_tol);
387 
388  // set the returned values
389  yf[0] = 0;
390  returned_intnl = intnl_old + dirn * dpm[0];
391  returned_stress = trial_stress - dpm[0] * n;
392  delta_dp = dpm[0] * dirn * returned_stress.dtrace();
393 
394  return true;
395 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
const bool _use_custom_returnMap
Whether to use the custom return-map algorithm.
const unsigned _max_iters
max iters for custom return map loop
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
Real yieldFunction(const RankTwoTensor &stress, Real intnl) const override
The following functions are what you should override when building single-plasticity models...
const Real _f_tol
Tolerance on yield function.
virtual bool returnMap(const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const
Performs a custom return-map.
virtual Real dtensile_strength(const Real internal_param) const
d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param ...
virtual Real dcompressive_strength(const Real internal_param) const
d(compressive strength)/d(internal_param) as a function of residual value, rate, and internal_param ...

◆ tensile_strength()

Real TensorMechanicsPlasticMeanCapTC::tensile_strength ( const Real  internal_param) const
protectedvirtual

tensile strength as a function of residual value, rate, and internal_param

Definition at line 222 of file TensorMechanicsPlasticMeanCapTC.C.

Referenced by activeConstraints(), consistentTangentOperator(), df_dsig(), dflowPotential_dintnl(), dflowPotential_dstress(), dhardPotential_dintnl(), dhardPotential_dstress(), dyieldFunction_dintnl(), hardPotential(), returnMap(), and yieldFunction().

223 {
224  return _strength.value(internal_param);
225 }
virtual Real value(Real intnl) const
const TensorMechanicsHardeningModel & _strength
the tensile strength

◆ useCustomCTO()

bool TensorMechanicsPlasticMeanCapTC::useCustomCTO ( ) const
overridevirtual

Returns false. You will want to override this in your derived class if you write a custom consistent tangent operator function.

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 441 of file TensorMechanicsPlasticMeanCapTC.C.

442 {
443  return _use_custom_cto;
444 }
const bool _use_custom_cto
Whether to use the custom consistent tangent operator algorithm.

◆ useCustomReturnMap()

bool TensorMechanicsPlasticMeanCapTC::useCustomReturnMap ( ) const
overridevirtual

Returns false. You will want to override this in your derived class if you write a custom returnMap function.

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 435 of file TensorMechanicsPlasticMeanCapTC.C.

436 {
437  return _use_custom_returnMap;
438 }
const bool _use_custom_returnMap
Whether to use the custom return-map algorithm.

◆ yieldFunction()

Real TensorMechanicsPlasticMeanCapTC::yieldFunction ( const RankTwoTensor &  stress,
Real  intnl 
) const
overrideprotectedvirtual

The following functions are what you should override when building single-plasticity models.

The yield function

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
Returns
the yield function

Reimplemented from TensorMechanicsPlasticModel.

Definition at line 64 of file TensorMechanicsPlasticMeanCapTC.C.

Referenced by returnMap().

65 {
66  const Real tr = stress.trace();
67  const Real t_str = tensile_strength(intnl);
68  if (tr >= t_str)
69  return tr - t_str;
70 
71  const Real c_str = compressive_strength(intnl);
72  if (tr <= c_str)
73  return -(tr - c_str);
74  // the following is zero at tr = t_str, and at tr = c_str
75  // it also has derivative = 1 at tr = t_str, and derivative = -1 at tr = c_str
76  // it also has second derivative = 0, at these points.
77  // This makes the complete yield function C2 continuous.
78  return (c_str - t_str) / libMesh::pi * std::sin(libMesh::pi * (tr - c_str) / (t_str - c_str));
79 }
virtual Real compressive_strength(const Real internal_param) const
compressive strength as a function of residual value, rate, and internal_param
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param

◆ yieldFunctionV()

void TensorMechanicsPlasticModel::yieldFunctionV ( const RankTwoTensor &  stress,
Real  intnl,
std::vector< Real > &  f 
) const
virtualinherited

Calculates the yield functions.

Note that for single-surface plasticity you don't want to override this - override the private yieldFunction below

Parameters
stressthe stress at which to calculate the yield function
intnlinternal parameter
[out]fthe yield functions

Reimplemented in TensorMechanicsPlasticMohrCoulombMulti, and TensorMechanicsPlasticTensileMulti.

Definition at line 68 of file TensorMechanicsPlasticModel.C.

Referenced by TensorMechanicsPlasticModel::returnMap().

71 {
72  f.assign(1, yieldFunction(stress, intnl));
73 }
virtual Real yieldFunction(const RankTwoTensor &stress, Real intnl) const
The following functions are what you should override when building single-plasticity models...

Member Data Documentation

◆ _c_strength

const TensorMechanicsHardeningModel& TensorMechanicsPlasticMeanCapTC::_c_strength
protected

the compressive strength

Definition at line 81 of file TensorMechanicsPlasticMeanCapTC.h.

Referenced by compressive_strength(), dcompressive_strength(), and TensorMechanicsPlasticMeanCapTC().

◆ _f_tol

const Real TensorMechanicsPlasticModel::_f_tol
inherited

◆ _ic_tol

const Real TensorMechanicsPlasticModel::_ic_tol
inherited

Tolerance on internal constraint.

Definition at line 177 of file TensorMechanicsPlasticModel.h.

◆ _max_iters

const unsigned TensorMechanicsPlasticMeanCapTC::_max_iters
protected

max iters for custom return map loop

Definition at line 69 of file TensorMechanicsPlasticMeanCapTC.h.

Referenced by returnMap().

◆ _strength

const TensorMechanicsHardeningModel& TensorMechanicsPlasticMeanCapTC::_strength
protected

the tensile strength

Definition at line 78 of file TensorMechanicsPlasticMeanCapTC.h.

Referenced by dtensile_strength(), tensile_strength(), and TensorMechanicsPlasticMeanCapTC().

◆ _use_custom_cto

const bool TensorMechanicsPlasticMeanCapTC::_use_custom_cto
protected

Whether to use the custom consistent tangent operator algorithm.

Definition at line 75 of file TensorMechanicsPlasticMeanCapTC.h.

Referenced by consistentTangentOperator(), and useCustomCTO().

◆ _use_custom_returnMap

const bool TensorMechanicsPlasticMeanCapTC::_use_custom_returnMap
protected

Whether to use the custom return-map algorithm.

Definition at line 72 of file TensorMechanicsPlasticMeanCapTC.h.

Referenced by returnMap(), and useCustomReturnMap().


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