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