Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #pragma once
11 : /* MOOSE - Multiphysics Object Oriented Simulation Environment */
12 : /* */
13 : /* All contents are licensed under LGPL V2.1 */
14 : /* See LICENSE for full restrictions */
15 : /****************************************************************/
16 :
17 : #include "SolidMechanicsPlasticModel.h"
18 : #include "UserObjectInterface.h"
19 :
20 : /**
21 : * MultiPlasticityRawComponentAssembler holds and computes yield functions,
22 : * flow directions, etc, for use in FiniteStrainMultiPlasticity
23 : *
24 : * Here there are a number of numbering systems.
25 : *
26 : * There are _num_models of plastic models, read from the input file.
27 : * Each of these models can, in principal, contain multiple
28 : * 'internal surfaces'.
29 : * Models are numbered from 0 to _num_models - 1.
30 : *
31 : * There are _num_surfaces surfaces. This
32 : * = sum_{plastic_models} (numberSurfaces in model)
33 : * Evidently _num_surface >= _num_models
34 : * Surfaces are numbered from 0 to _num_surfaces - 1.
35 : *
36 : * The std::vectors _model_given_surface, _model_surface_given_surface
37 : * and _surfaces_given_model allow translation between these
38 : */
39 : class MultiPlasticityRawComponentAssembler : public UserObjectInterface
40 : {
41 : public:
42 : static InputParameters validParams();
43 :
44 : MultiPlasticityRawComponentAssembler(const MooseObject * moose_object);
45 :
46 6192 : virtual ~MultiPlasticityRawComponentAssembler() {}
47 :
48 : protected:
49 : const InputParameters & _params;
50 :
51 : /// Number of plastic models for this material
52 : unsigned int _num_models;
53 :
54 : /**
55 : * Number of surfaces within the plastic models.
56 : * For many situations this will be = _num_models
57 : * since each model will contain just one surface.
58 : * More generally it is >= _num_models. For instance,
59 : * Mohr-Coulomb is a single model with 6 surfaces
60 : */
61 : unsigned int _num_surfaces;
62 :
63 : /// _surfaces_given_model[model_number] = vector of surface numbers for this model
64 : std::vector<std::vector<unsigned int>> _surfaces_given_model;
65 :
66 : /// Allows initial set of active constraints to be chosen optimally
67 : MooseEnum _specialIC;
68 :
69 : /// User objects that define the yield functions, flow potentials, etc
70 : std::vector<const SolidMechanicsPlasticModel *> _f;
71 :
72 : /**
73 : * The active yield function(s)
74 : * @param stress the stress at which to calculate the yield function
75 : * @param intnl vector of internal parameters
76 : * @param active set of active constraints - only the active yield functions are put into "f"
77 : * @param[out] f the yield function (or functions in the case of multisurface plasticity)
78 : */
79 : virtual void yieldFunction(const RankTwoTensor & stress,
80 : const std::vector<Real> & intnl,
81 : const std::vector<bool> & active,
82 : std::vector<Real> & f);
83 :
84 : /**
85 : * The derivative of the active yield function(s) with respect to stress
86 : * @param stress the stress at which to calculate the yield function
87 : * @param intnl vector of internal parameters
88 : * @param active set of active constraints - only the active derivatives are put into "df_dstress"
89 : * @param[out] df_dstress the derivative (or derivatives in the case of multisurface plasticity).
90 : * df_dstress[alpha](i, j) = dyieldFunction[alpha]/dstress(i, j)
91 : */
92 : virtual void dyieldFunction_dstress(const RankTwoTensor & stress,
93 : const std::vector<Real> & intnl,
94 : const std::vector<bool> & active,
95 : std::vector<RankTwoTensor> & df_dstress);
96 :
97 : /**
98 : * The derivative of active yield function(s) with respect to their internal parameters (the user
99 : * objects assume there is exactly one internal param per yield function)
100 : * @param stress the stress at which to calculate the yield function
101 : * @param intnl vector of internal parameters
102 : * @param active set of active constraints - only the active derivatives are put into "df_dintnl"
103 : * @param[out] df_dintnl the derivatives. df_dstress[alpha] = dyieldFunction[alpha]/dintnl[alpha]
104 : */
105 : virtual void dyieldFunction_dintnl(const RankTwoTensor & stress,
106 : const std::vector<Real> & intnl,
107 : const std::vector<bool> & active,
108 : std::vector<Real> & df_dintnl);
109 :
110 : /**
111 : * The active flow potential(s) - one for each yield function
112 : * @param stress the stress at which to calculate the flow potential
113 : * @param intnl vector of internal parameters
114 : * @param active set of active constraints - only the active flow potentials are put into "r"
115 : * @param[out] r the flow potential (flow potentials in the multi-surface case)
116 : */
117 : virtual void flowPotential(const RankTwoTensor & stress,
118 : const std::vector<Real> & intnl,
119 : const std::vector<bool> & active,
120 : std::vector<RankTwoTensor> & r);
121 :
122 : /**
123 : * The derivative of the active flow potential(s) with respect to stress
124 : * @param stress the stress at which to calculate the flow potential
125 : * @param intnl vector of internal parameters
126 : * @param active set of active constraints - only the active derivatives are put into "dr_dstress"
127 : * @param[out] dr_dstress the derivative. dr_dstress[alpha](i, j, k, l) = dr[alpha](i,
128 : * j)/dstress(k, l)
129 : */
130 : virtual void dflowPotential_dstress(const RankTwoTensor & stress,
131 : const std::vector<Real> & intnl,
132 : const std::vector<bool> & active,
133 : std::vector<RankFourTensor> & dr_dstress);
134 :
135 : /**
136 : * The derivative of the active flow potentials with respect to the active internal parameters
137 : * The UserObjects explicitly assume that r[alpha] is only dependent on intnl[alpha]
138 : * @param stress the stress at which to calculate the flow potential
139 : * @param intnl vector of internal parameters
140 : * @param active set of active constraints - only the active derivatives are put into "dr_dintnl"
141 : * @param[out] dr_dintnl the derivatives. dr_dintnl[alpha](i, j) = dr[alpha](i, j)/dintnl[alpha]
142 : */
143 : virtual void dflowPotential_dintnl(const RankTwoTensor & stress,
144 : const std::vector<Real> & intnl,
145 : const std::vector<bool> & active,
146 : std::vector<RankTwoTensor> & dr_dintnl);
147 :
148 : /**
149 : * The active hardening potentials (one for each internal parameter and for each yield function)
150 : * by assumption in the Userobjects, the h[a][alpha] is nonzero only if the surface alpha is part
151 : * of model a, so we only calculate those here
152 : * @param stress the stress at which to calculate the hardening potential
153 : * @param intnl vector of internal parameters
154 : * @param active set of active constraints - only the active hardening potentials are put into "h"
155 : * @param[out] h the hardening potentials. h[alpha] = hardening potential for yield fcn alpha
156 : * (and, by the above assumption we know which hardening parameter, a, this belongs to)
157 : */
158 : virtual void hardPotential(const RankTwoTensor & stress,
159 : const std::vector<Real> & intnl,
160 : const std::vector<bool> & active,
161 : std::vector<Real> & h);
162 :
163 : /**
164 : * The derivative of the active hardening potentials with respect to stress
165 : * By assumption in the Userobjects, the h[a][alpha] is nonzero only for a = alpha, so we only
166 : * calculate those here
167 : * @param stress the stress at which to calculate the hardening potentials
168 : * @param intnl vector of internal parameters
169 : * @param active set of active constraints - only the active derivatives are put into "dh_dstress"
170 : * @param[out] dh_dstress the derivative. dh_dstress[a](i, j) = dh[a]/dstress(k, l)
171 : */
172 : virtual void dhardPotential_dstress(const RankTwoTensor & stress,
173 : const std::vector<Real> & intnl,
174 : const std::vector<bool> & active,
175 : std::vector<RankTwoTensor> & dh_dstress);
176 :
177 : /**
178 : * The derivative of the active hardening potentials with respect to the active internal
179 : * parameters
180 : * @param stress the stress at which to calculate the hardening potentials
181 : * @param intnl vector of internal parameters
182 : * @param active set of active constraints - only the active derivatives are put into "dh_dintnl"
183 : * @param[out] dh_dintnl the derivatives. dh_dintnl[a][alpha][b] = dh[a][alpha]/dintnl[b]. Note
184 : * that the userobjects assume that there is exactly one internal parameter per yield function, so
185 : * the derivative is only nonzero for a=alpha=b, so that is all we calculate
186 : */
187 : virtual void dhardPotential_dintnl(const RankTwoTensor & stress,
188 : const std::vector<Real> & intnl,
189 : const std::vector<bool> & active,
190 : std::vector<Real> & dh_dintnl);
191 :
192 : /**
193 : * Constructs a set of active constraints, given the yield functions, f.
194 : * This uses SolidMechanicsPlasticModel::activeConstraints to identify the active
195 : * constraints for each model.
196 : * @param f yield functions (should be _num_surfaces of these)
197 : * @param stress stress tensor
198 : * @param intnl internal parameters
199 : * @param Eijkl elasticity tensor (stress = Eijkl*strain)
200 : * @param[out] act the set of active constraints (will be resized to _num_surfaces)
201 : */
202 : virtual void buildActiveConstraints(const std::vector<Real> & f,
203 : const RankTwoTensor & stress,
204 : const std::vector<Real> & intnl,
205 : const RankFourTensor & Eijkl,
206 : std::vector<bool> & act);
207 :
208 : /// returns the model number, given the surface number
209 : unsigned int modelNumber(unsigned int surface);
210 :
211 : /// returns true if any internal surfaces of the given model are active according to 'active'
212 : bool anyActiveSurfaces(int model, const std::vector<bool> & active);
213 :
214 : /**
215 : * Returns the internal surface number(s) of the active surfaces of the given model
216 : * This may be of size=0 if there are no active surfaces of the given model
217 : * @param model the model number
218 : * @param active array with entries being 'true' if the surface is active
219 : * @param[out] active_surfaces_of_model the output
220 : */
221 : void activeModelSurfaces(int model,
222 : const std::vector<bool> & active,
223 : std::vector<unsigned int> & active_surfaces_of_model);
224 :
225 : /**
226 : * Returns the external surface number(s) of the active surfaces of the given model
227 : * This may be of size=0 if there are no active surfaces of the given model
228 : * @param model the model number
229 : * @param active array with entries being 'true' if the surface is active
230 : * @param[out] active_surfaces the output
231 : */
232 : void activeSurfaces(int model,
233 : const std::vector<bool> & active,
234 : std::vector<unsigned int> & active_surfaces);
235 :
236 : /**
237 : * Performs a returnMap for each plastic model using
238 : * their inbuilt returnMap functions. This may be used
239 : * to quickly ascertain whether a (trial_stress, intnl_old) configuration
240 : * is admissible, or whether a single model's customized returnMap
241 : * function can provide a solution to the return-map problem,
242 : * or whether a full Newton-Raphson approach such as implemented
243 : * in ComputeMultiPlasticityStress is needed.
244 : *
245 : * There are three cases mentioned below:
246 : * (A) The (trial_stress, intnl_old) configuration is admissible
247 : * according to all plastic models
248 : * (B) The (trial_stress, intnl_old) configuration is inadmissible
249 : * to exactly one plastic model, and that model can successfully
250 : * use its customized returnMap function to provide a returned
251 : * (stress, intnl) configuration, and that configuration is
252 : * admissible according to all plastic models
253 : * (C) All other cases. This includes customized returnMap
254 : * functions failing, or more than one plastic_model being
255 : * inadmissible, etc
256 : *
257 : * @param trial_stress the trial stress
258 : * @param intnl_old the old values of the internal parameters
259 : * @param E_ijkl the elasticity tensor
260 : * @param ep_plastic_tolerance the tolerance on the plastic strain
261 : * @param[out] stress is set to trial_stress in case (A) or (C), and the returned value of stress
262 : * in case (B).
263 : * @param[out] intnl is set to intnl_old in case (A) or (C), and the returned value of intnl in
264 : * case (B)
265 : * @param[out] pm Zero in case (A) or (C), otherwise the plastic multipliers needed to bring
266 : * about the returnMap in case (B)
267 : * @param[in/out] cumulative_pm cumulative plastic multipliers, updated in case (B), otherwise
268 : * left untouched
269 : * @param[out] delta_dp is unchanged in case (A) or (C), and is set to the change in plastic
270 : * strain in case(B)
271 : * @param[out] yf will contain the yield function values at (stress, intnl)
272 : * @param[out] num_successful_plastic_returns will be 0 for (A) and (C), and 1 for (B)
273 : * @return true in case (A) and (B), and false in case (C)
274 : */
275 : bool returnMapAll(const RankTwoTensor & trial_stress,
276 : const std::vector<Real> & intnl_old,
277 : const RankFourTensor & E_ijkl,
278 : Real ep_plastic_tolerance,
279 : RankTwoTensor & stress,
280 : std::vector<Real> & intnl,
281 : std::vector<Real> & pm,
282 : std::vector<Real> & cumulative_pm,
283 : RankTwoTensor & delta_dp,
284 : std::vector<Real> & yf,
285 : unsigned & num_successful_plastic_returns,
286 : unsigned & custom_model);
287 :
288 : private:
289 : /// given a surface number, this returns the model number
290 : std::vector<unsigned int> _model_given_surface;
291 :
292 : /// given a surface number, this returns the corresponding-model's internal surface number
293 : std::vector<unsigned int> _model_surface_given_surface;
294 :
295 : /**
296 : * "Rock" version
297 : * Constructs a set of active constraints, given the yield functions, f.
298 : * This uses SolidMechanicsPlasticModel::activeConstraints to identify the active
299 : * constraints for each model.
300 : * @param f yield functions (should be _num_surfaces of these)
301 : * @param stress stress tensor
302 : * @param intnl internal parameters
303 : * @param Eijkl elasticity tensor (stress = Eijkl*strain)
304 : * @param[out] act the set of active constraints (will be resized to _num_surfaces)
305 : */
306 : void buildActiveConstraintsRock(const std::vector<Real> & f,
307 : const RankTwoTensor & stress,
308 : const std::vector<Real> & intnl,
309 : const RankFourTensor & Eijkl,
310 : std::vector<bool> & act);
311 :
312 : /**
313 : * "Joint" version
314 : * Constructs a set of active constraints, given the yield functions, f.
315 : * This uses SolidMechanicsPlasticModel::activeConstraints to identify the active
316 : * constraints for each model.
317 : * @param f yield functions (should be _num_surfaces of these)
318 : * @param stress stress tensor
319 : * @param intnl internal parameters
320 : * @param Eijkl elasticity tensor (stress = Eijkl*strain)
321 : * @param[out] act the set of active constraints (will be resized to _num_surfaces)
322 : */
323 : void buildActiveConstraintsJoint(const std::vector<Real> & f,
324 : const RankTwoTensor & stress,
325 : const std::vector<Real> & intnl,
326 : const RankFourTensor & Eijkl,
327 : std::vector<bool> & act);
328 : };
|