Skip to content

Commit 33defc3

Browse files
WilkAndyWilkAndy
committed
MultiParameterPlasticityStressUpdate and all its derived classes, namely MohrCoulombStressUpdate, TensileStressUpdate and TwoParameterPlasticityStressUpdate, and derived classes from those, such as CappedDruckerPragerStressUpdate and CappedDruckerPragerCosseratStressUpdate, now have no variables that are heap-allocated in functions that are relevant to explicit time stepping (ie, functions relevant to Jacobian-only calculations have largely been untouched). This will hopefully speed up explicit simulations where it has been observed that heap alloc/dealloc takes considerable time. Refs #32178 .
Most of the strategy boils down to making a few mutable scratch class variables to take the place of the old heap-allocated variables. Mutable variables have been used so that const could be retained, which i think makes most sense in this situation. Most of the changes occur within MultiParameterPlasticityStressUpdate, and do not impact derived classes. However, the signatures of two important functions - dstress_param_dstress and d2stress_param_dstress - have been changed. The derived classes within solid_mechanics have all had their signatures updated, but external Apps such as Blackbear should be updated to take advantage of all the speedups associated with no heap alloc/dealloc. When that is completed, the old versions may be deleted from MultiParameterPlasticityStressUpdate, and I've included notes on how to do that.
1 parent 8af10d8 commit 33defc3

10 files changed

+392
-144
lines changed

‎modules/solid_mechanics/include/materials/CappedMohrCoulombCosseratStressUpdate.h‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -74,4 +74,7 @@ class CappedMohrCoulombCosseratStressUpdate : public CappedMohrCoulombStressUpda
7474
bool compute_full_tangent_operator,
7575
const std::vector<std::vector<Real>> & dvar_dtrial,
7676
RankFourTensor & cto) override;
77+
78+
private:
79+
std::vector<Real> _eigvals_scratch; // used as temporary storage to avoid heap allocation
7780
};

‎modules/solid_mechanics/include/materials/CappedMohrCoulombStressUpdate.h‎

Lines changed: 29 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -65,9 +65,11 @@ class CappedMohrCoulombStressUpdate : public MultiParameterPlasticityStressUpdat
6565
void computeStressParams(const RankTwoTensor & stress,
6666
std::vector<Real> & stress_params) const override;
6767

68-
std::vector<RankTwoTensor> dstress_param_dstress(const RankTwoTensor & stress) const override;
68+
void dstressparam_dstress(const RankTwoTensor & stress,
69+
std::vector<RankTwoTensor> & dsp) const override;
6970

70-
std::vector<RankFourTensor> d2stress_param_dstress(const RankTwoTensor & stress) const override;
71+
void d2stressparam_dstress(const RankTwoTensor & stress,
72+
std::vector<RankFourTensor> & d2sp) const override;
7173

7274
virtual void setStressAfterReturnV(const RankTwoTensor & stress_trial,
7375
const std::vector<Real> & stress_params,
@@ -119,4 +121,29 @@ class CappedMohrCoulombStressUpdate : public MultiParameterPlasticityStressUpdat
119121
bool compute_full_tangent_operator,
120122
const std::vector<std::vector<Real>> & dvar_dtrial,
121123
RankFourTensor & cto) override;
124+
125+
private:
126+
/**
127+
* this is d(stress_param[:])/d(stress) which is calculated by dstressparam_dstress. Using this
128+
* mutable means repeated allocation/deallocation is unnecessary.
129+
*/
130+
mutable std::vector<RankTwoTensor> _dsp_trial_scratch;
131+
132+
/**
133+
* eigenvalues of the stress, used in dstressparam_dstress and preReturnMapV. Using this mutable
134+
* means repeated allocation/deallocation is unnecessary.
135+
*/
136+
mutable std::vector<Real> _eigvals_scratch;
137+
138+
/**
139+
* derivative of ga_shear w.r.t. the stress parameters, used in setIntnlDerivativesV. Using
140+
* this mutable means repeated allocation/deallocation is unnecessary.
141+
*/
142+
mutable std::vector<Real> _dga_shear_scratch;
143+
144+
/**
145+
* scratch vector used in setIntnlDerivativesV. Using
146+
* this mutable means repeated allocation/deallocation is unnecessary.
147+
*/
148+
mutable std::vector<Real> _dshear_correction_scratch;
122149
};

‎modules/solid_mechanics/include/materials/MultiParameterPlasticityStressUpdate.h‎

Lines changed: 112 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -121,7 +121,7 @@ class MultiParameterPlasticityStressUpdate : public StressUpdateBase
121121
return TangentCalculationMethod::FULL;
122122
}
123123

124-
/// Internal dimensionality of tensors (currently this is 3 throughout tensor_mechanics)
124+
/// Internal dimensionality of tensors (currently this is 3 throughout solid mechanics)
125125
constexpr static unsigned _tensor_dimensionality = 3;
126126

127127
/// Number of stress parameters
@@ -228,6 +228,19 @@ class MultiParameterPlasticityStressUpdate : public StressUpdateBase
228228
{
229229
}
230230

231+
// Zero everything without resizing/reallocating
232+
void reset()
233+
{
234+
f = 0.0;
235+
std::fill(df.begin(), df.end(), 0.0);
236+
std::fill(df_di.begin(), df_di.end(), 0.0);
237+
std::fill(dg.begin(), dg.end(), 0.0);
238+
for (auto & row : d2g)
239+
std::fill(row.begin(), row.end(), 0.0);
240+
for (auto & row : d2g_di)
241+
std::fill(row.begin(), row.end(), 0.0);
242+
}
243+
231244
// this may be involved in the smoothing of a group of yield functions
232245
bool operator<(const yieldAndFlow & fd) const { return f < fd.f; }
233246
};
@@ -328,7 +341,7 @@ class MultiParameterPlasticityStressUpdate : public StressUpdateBase
328341
const std::vector<Real> & intnl_ok,
329342
std::vector<Real> & intnl,
330343
std::vector<Real> & rhs,
331-
Real & linesearch_needed) const;
344+
Real & linesearch_needed);
332345

333346
/**
334347
* Performs a Newton-Raphson step to attempt to zero rhs
@@ -566,13 +579,12 @@ class MultiParameterPlasticityStressUpdate : public StressUpdateBase
566579
* @param inelastic_strain_increment[out] The inelastic strain increment resulting from this
567580
* return-map
568581
*/
569-
virtual void
570-
setInelasticStrainIncrementAfterReturn(const RankTwoTensor & stress_trial,
571-
Real gaE,
572-
const yieldAndFlow & smoothed_q,
573-
const RankFourTensor & elasticity_tensor,
574-
const RankTwoTensor & returned_stress,
575-
RankTwoTensor & inelastic_strain_increment) const;
582+
virtual void setInelasticStrainIncrementAfterReturn(const RankTwoTensor & stress_trial,
583+
Real gaE,
584+
const yieldAndFlow & smoothed_q,
585+
const RankFourTensor & elasticity_tensor,
586+
const RankTwoTensor & returned_stress,
587+
RankTwoTensor & inelastic_strain_increment);
576588

577589
/**
578590
* Calculates the consistent tangent operator.
@@ -608,18 +620,39 @@ class MultiParameterPlasticityStressUpdate : public StressUpdateBase
608620

609621
/**
610622
* d(stress_param[i])/d(stress) at given stress
623+
* TODO: remove this function when Blackbear and others have implemented the void version
611624
* @param stress stress tensor
612625
* @return d(stress_param[:])/d(stress)
613626
*/
614-
virtual std::vector<RankTwoTensor> dstress_param_dstress(const RankTwoTensor & stress) const = 0;
627+
virtual std::vector<RankTwoTensor> dstress_param_dstress(const RankTwoTensor & stress) const;
628+
629+
/**
630+
* d(stress_param[i])/d(stress) at given stress. When overriding, ensure to assert that dsp has
631+
* size _num_sp
632+
* TODO: when Blackbear and others have this implemented, make this pure virtual
633+
* @param stress[in] stress tensor
634+
* @param dsp[out] d(stress_param[:])/d(stress)
635+
*/
636+
virtual void dstressparam_dstress(const RankTwoTensor & stress,
637+
std::vector<RankTwoTensor> & dsp) const;
615638

616639
/**
617640
* d2(stress_param[i])/d(stress)/d(stress) at given stress
641+
* TODO: remove this function when Blackbear and others have implemented the void version
618642
* @param stress stress tensor
619643
* @return d2(stress_param[:])/d(stress)/d(stress)
620644
*/
621-
virtual std::vector<RankFourTensor>
622-
d2stress_param_dstress(const RankTwoTensor & stress) const = 0;
645+
virtual std::vector<RankFourTensor> d2stress_param_dstress(const RankTwoTensor & stress) const;
646+
647+
/**
648+
* d2(stress_param[i])/d(stress)/d(stress) at given stress. When overriding, ensure to assert
649+
* that d2sp has size _num_sp
650+
* TODO: when Blackbear and others have this implemented, make this pure virtual
651+
* @param stress[in] stress tensor
652+
* @param d2sp[out] d2(stress_param[:])/d(stress)/d(stress)
653+
*/
654+
virtual void d2stressparam_dstress(const RankTwoTensor & stress,
655+
std::vector<RankFourTensor> & d2sp) const;
623656

624657
/// Sets _Eij and _En and _Cij
625658
virtual void setEffectiveElasticity(const RankFourTensor & Eijkl) = 0;
@@ -732,6 +765,73 @@ class MultiParameterPlasticityStressUpdate : public StressUpdateBase
732765
*/
733766
std::vector<Real> _current_intnl;
734767

768+
/**
769+
* Current values of the yield function, derivatives, etc during updateState.
770+
* During Newton-Raphson convergence, this is potentially re-calculated multiple times within the
771+
* loops of updateState and by the lineSearch function that is called within those loops.
772+
*/
773+
yieldAndFlow _smoothed_q;
774+
775+
/**
776+
* d(stress_param[:])/d(stress) used in dstressparam_dstress
777+
* to avoid repeatedly allocating/deallocating this vector
778+
*/
779+
mutable std::vector<RankTwoTensor> _dsp_scratch;
780+
781+
/**
782+
* d(stress_param[:])/d(stress_trial) used in dstressparam_dstress
783+
* to avoid repeatedly allocating/deallocating this vector
784+
*/
785+
mutable std::vector<RankTwoTensor> _dsp_trial_scratch;
786+
787+
/**
788+
* d2(stress_param[:])/d(stress)/d(stress) used in d2stressparam_dstress
789+
* to avoid repeatedly allocating/deallocating this vector
790+
*/
791+
mutable std::vector<RankFourTensor> _d2sp_scratch;
792+
793+
/**
794+
* values of yield functions in yieldF,
795+
* to avoid repeatedly allocating/deallocating this vector
796+
*/
797+
mutable std::vector<Real> _yfs_scratch;
798+
799+
/**
800+
* container to hold old values of old stress_params in lineSearch
801+
* to avoid repeatedly allocating/deallocating this vector
802+
*/
803+
mutable std::vector<Real> _sp_params_old_scratch;
804+
805+
/**
806+
* container to hold initial values of rhs in lineSearch
807+
* to avoid repeatedly allocating/deallocating this vector
808+
*/
809+
mutable std::vector<Real> _delta_nr_params_scratch;
810+
811+
/**
812+
* derivative of internal parameters with respect to stress parameters,
813+
* to avoid repeatedly allocating/deallocating this vector
814+
*/
815+
mutable std::vector<std::vector<Real>> _dintnl_scratch;
816+
817+
/**
818+
* jacobian in the NR process,
819+
* to avoid repeatedly allocating/deallocating this vector
820+
*/
821+
mutable std::vector<double> _jac_scratch;
822+
823+
/**
824+
* container to hold LAPACK ipiv integers,
825+
* to avoid repeatedly allocating/deallocating this vector
826+
*/
827+
mutable std::vector<PetscBLASInt> _ipiv_scratch;
828+
829+
/**
830+
* all the yield function information, for use in smoothAllQuantities,
831+
* to avoid repeatedly allocating/deallocating this vector
832+
*/
833+
mutable std::vector<yieldAndFlow> _all_q_scratch;
834+
735835
private:
736836
/**
737837
* The type of smoother function

‎modules/solid_mechanics/include/materials/TensileStressUpdate.h‎

Lines changed: 17 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -47,9 +47,11 @@ class TensileStressUpdate : public MultiParameterPlasticityStressUpdate
4747
void computeStressParams(const RankTwoTensor & stress,
4848
std::vector<Real> & stress_params) const override;
4949

50-
std::vector<RankTwoTensor> dstress_param_dstress(const RankTwoTensor & stress) const override;
50+
void dstressparam_dstress(const RankTwoTensor & stress,
51+
std::vector<RankTwoTensor> & dsp) const override;
5152

52-
std::vector<RankFourTensor> d2stress_param_dstress(const RankTwoTensor & stress) const override;
53+
void d2stressparam_dstress(const RankTwoTensor & stress,
54+
std::vector<RankFourTensor> & d2sp) const override;
5355

5456
virtual void setStressAfterReturnV(const RankTwoTensor & stress_trial,
5557
const std::vector<Real> & stress_params,
@@ -101,4 +103,17 @@ class TensileStressUpdate : public MultiParameterPlasticityStressUpdate
101103
bool compute_full_tangent_operator,
102104
const std::vector<std::vector<Real>> & dvar_dtrial,
103105
RankFourTensor & cto) override;
106+
107+
private:
108+
/**
109+
* this is d(stress_param[:])/d(stress) which is calculated by dstressparam_dstress. Using this
110+
* mutable means repeated allocation/deallocation is unnecessary.
111+
*/
112+
mutable std::vector<RankTwoTensor> _dsp_trial_scratch;
113+
114+
/**
115+
* eigenvalues of the stress, used in dstressparam_dstress and preReturnMapV. Using this mutable
116+
* means repeated allocation/deallocation is unnecessary.
117+
*/
118+
mutable std::vector<Real> _eigvals_scratch;
104119
};

‎modules/solid_mechanics/include/materials/TwoParameterPlasticityStressUpdate.h‎

Lines changed: 10 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -248,13 +248,12 @@ class TwoParameterPlasticityStressUpdate : public MultiParameterPlasticityStress
248248
const RankFourTensor & Eijkl,
249249
RankTwoTensor & stress) const override;
250250

251-
void
252-
setInelasticStrainIncrementAfterReturn(const RankTwoTensor & stress_trial,
253-
Real gaE,
254-
const yieldAndFlow & smoothed_q,
255-
const RankFourTensor & elasticity_tensor,
256-
const RankTwoTensor & returned_stress,
257-
RankTwoTensor & inelastic_strain_increment) const override;
251+
void setInelasticStrainIncrementAfterReturn(const RankTwoTensor & stress_trial,
252+
Real gaE,
253+
const yieldAndFlow & smoothed_q,
254+
const RankFourTensor & elasticity_tensor,
255+
const RankTwoTensor & returned_stress,
256+
RankTwoTensor & inelastic_strain_increment) override;
258257

259258
/**
260259
* Calculates the consistent tangent operator.
@@ -298,10 +297,10 @@ class TwoParameterPlasticityStressUpdate : public MultiParameterPlasticityStress
298297
const std::vector<std::vector<Real>> & dvar_dtrial,
299298
RankFourTensor & cto) override;
300299

301-
virtual std::vector<RankTwoTensor>
302-
dstress_param_dstress(const RankTwoTensor & stress) const override;
303-
virtual std::vector<RankFourTensor>
304-
d2stress_param_dstress(const RankTwoTensor & stress) const override;
300+
virtual void dstressparam_dstress(const RankTwoTensor & stress,
301+
std::vector<RankTwoTensor> & dsp) const override;
302+
virtual void d2stressparam_dstress(const RankTwoTensor & stress,
303+
std::vector<RankFourTensor> & d2sp) const override;
305304
/**
306305
* d(p)/d(stress)
307306
* Derived classes must override this

‎modules/solid_mechanics/src/materials/CappedMohrCoulombCosseratStressUpdate.C‎

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -36,7 +36,8 @@ CappedMohrCoulombCosseratStressUpdate::CappedMohrCoulombCosseratStressUpdate(
3636
_host_young(getParam<Real>("host_youngs_modulus")),
3737
_host_poisson(getParam<Real>("host_poissons_ratio")),
3838
_host_E0011(_host_young * _host_poisson / (1.0 + _host_poisson) / (1.0 - 2.0 * _host_poisson)),
39-
_host_E0000(_host_E0011 + _host_young / (1.0 + _host_poisson))
39+
_host_E0000(_host_E0011 + _host_young / (1.0 + _host_poisson)),
40+
_eigvals_scratch(_tensor_dimensionality)
4041
{
4142
}
4243

@@ -48,8 +49,10 @@ CappedMohrCoulombCosseratStressUpdate::preReturnMapV(
4849
const std::vector<Real> & /*yf*/,
4950
const RankFourTensor & /*Eijkl*/)
5051
{
51-
std::vector<Real> eigvals;
52-
stress_trial.symmetricEigenvaluesEigenvectors(eigvals, _eigvecs);
52+
mooseAssert(
53+
_eigvals_scratch.size() == _tensor_dimensionality,
54+
"_eigvals_scratch incorrectly sized in CappedMohCoulombCosseratStressUpdate::preReturnMapV");
55+
stress_trial.symmetricEigenvaluesEigenvectors(_eigvals_scratch, _eigvecs);
5356
_poissons_ratio = _host_poisson;
5457
}
5558

0 commit comments

Comments
 (0)