Line data Source code
1 : /*************************************************/
2 : /* DO NOT MODIFY THIS HEADER */
3 : /* */
4 : /* MASTODON */
5 : /* */
6 : /* (c) 2015 Battelle Energy Alliance, LLC */
7 : /* ALL RIGHTS RESERVED */
8 : /* */
9 : /* Prepared by Battelle Energy Alliance, LLC */
10 : /* With the U. S. Department of Energy */
11 : /* */
12 : /* See COPYRIGHT for full restrictions */
13 : /*************************************************/
14 :
15 : // This code was implemented in colloboration with Ozgun Numanoglu
16 : // ([email protected]) and Omar Baltaji ([email protected]) from UIUC.
17 : #include "ISoilAction.h"
18 :
19 : #include "ColumnMajorMatrix.h"
20 : #include "Conversion.h"
21 : #include "FEProblem.h"
22 : #include "Factory.h"
23 : #include "Parser.h"
24 :
25 : registerMooseAction("MastodonApp", ISoilAction, "add_material");
26 :
27 : InputParameters
28 74 : ISoilAction::validParams()
29 : {
30 74 : InputParameters params = Action::validParams();
31 74 : params.addClassDescription("Set up I-Soil material model.");
32 148 : params.addRequiredParam<std::vector<VariableName>>(
33 : "displacements", "The vector of displacement variables in the problem.");
34 148 : params.addRequiredParam<std::vector<SubdomainName>>(
35 : "block", "The blocks where this material model is applied.");
36 148 : params.addRequiredParam<std::vector<VariableName>>(
37 : "layer_variable", "The auxvariable providing the soil layer identification.");
38 148 : params.addRequiredParam<std::vector<unsigned int>>(
39 : "layer_ids",
40 : "Vector of layer ids that map one-to-one to the rest of the "
41 : "soil layer parameters provided as input.");
42 148 : params.addRequiredParam<std::vector<Real>>("poissons_ratio",
43 : "Poissons's ratio for the soil layers. The "
44 : "size of the vector should be same as the size of "
45 : "layer_ids.");
46 148 : params.addParam<Real>("b_exp",
47 148 : 0.0,
48 : "The exponential factors for pressure "
49 : "dependent stiffness for all the soil "
50 : "layers.");
51 148 : params.addParam<std::vector<Real>>("p_ref",
52 : {},
53 : "The reference pressure at which "
54 : "the parameters are defined for "
55 : "each soil layer. If 'soil_type = "
56 : "darendeli', then the reference "
57 : "pressure must be input in kilopascals.");
58 148 : params.addParam<Real>("a0",
59 148 : 1.0,
60 : "The first coefficient for pressure dependent yield strength "
61 : "calculation for all the soil layers. If a0 = 1, a1 = 0 and "
62 : "a2=0 for one soil layer, then the yield strength of that "
63 : "layer is independent of pressure.");
64 148 : params.addParam<Real>("a1",
65 148 : 0.0,
66 : "The second coefficient for pressure dependent yield "
67 : "strength calculation for all the soil layers. If a0 = "
68 : "1, a1 = 0, a2 = 0 for one soil layer, then the yield "
69 : "strength of that layer is independent of pressure.");
70 148 : params.addParam<Real>("a2",
71 148 : 0.0,
72 : "The third coefficient for pressure dependent yield "
73 : "strength calculation for all the soil layers. If a0 = "
74 : "1, a1=0 and a2=0 for one soil layer, then the yield "
75 : "strength of that layer is independent of pressure.");
76 148 : params.addParam<Real>("tension_pressure_cut_off",
77 148 : -1.0,
78 : "The tension cut-off for all the soil layers. If the "
79 : "pressure becomes lower than this value, then the "
80 : "stiffness of the soil reduces to zero. A negative "
81 : "pressure indicates tension. The default "
82 : "value is -1.0 for all the soil layers.");
83 148 : params.addParam<bool>("pressure_dependency",
84 148 : false,
85 : "Set to true to turn on pressure dependent stiffness "
86 : "and yield strength calculation.");
87 148 : params.addParam<bool>("implicit", true, "Set to false to use the explicit solver.");
88 148 : params.addParam<std::vector<FunctionName>>(
89 : "initial_soil_stress",
90 : {},
91 : "The function values for the initial stress distribution. 9 function "
92 : "names have to be provided corresponding to stress_xx, stress_xy, "
93 : "stress_xz, stress_yx, stress_yy, stress_yz, stress_zx, stress_zy, "
94 : "stress_zz. Each function can be a function of space.");
95 : // params for specific backbone types
96 148 : MooseEnum soil_type("user_defined darendeli gqh thin_layer");
97 148 : params.addRequiredParam<MooseEnum>(
98 : "soil_type",
99 : soil_type,
100 : "This parameter determines the type of backbone curve used. Use 'user_defined' "
101 : "for a user defined backbone curve provided in a data file, 'darendeli' "
102 : "if the backbone curve is to be determined using Darandeli equations, 'gqh' "
103 : "if the backbone curve is determined using the GQ/H approach and 'thin_layer' if the soil is "
104 : "being used to simulate a thin-layer friction interface.");
105 : // params required for user_defined backbone curve: soil_type = 'user_defined'
106 148 : params.addParam<std::vector<FileName>>(
107 : "backbone_curve_files",
108 : "The vector of file names of the files containing "
109 : "stress-strain backbone curves for the different soil layers. The "
110 : "size of the vector should be same as the size of layer_ids. All files "
111 : "should contain the same number of stress-strain points. Headers are not "
112 : "expected and it is assumed that the first column corresponds to strain values "
113 : "and the second column corresponds to the stress values. Additionally, two "
114 : "segments of a backbone curve cannot have the same slope.");
115 : // params required for soil_type = 'darendeli' and 'GQ/H'
116 148 : params.addRequiredParam<std::vector<Real>>("initial_shear_modulus",
117 : "The initial shear modulus of the soil layers.");
118 148 : params.addParam<unsigned int>("number_of_points",
119 : "The total number of data points in which the "
120 : "backbone curve needs to be split for all soil "
121 : "layers (required for Darandeli or GQ/H type backbone curves).");
122 : // params required for Darandeli backbone curve: soil_type = 'darendeli'
123 148 : params.addParam<std::vector<Real>>("over_consolidation_ratio",
124 : "The over consolidation ratio of the soil "
125 : "layers. Required for Darandeli backbone curve.");
126 148 : params.addParam<std::vector<Real>>(
127 : "plasticity_index",
128 : "The plasticity index of the soil layers. Required for Darandeli backbone curve.");
129 : // params required for GQ/H backbone curve: soil_type = "gqh"
130 148 : params.addParam<std::vector<Real>>("theta_1",
131 : "The curve fit coefficient for "
132 : "GQ/H model"
133 : "for each soil layer.");
134 148 : params.addParam<std::vector<Real>>("theta_2",
135 : "The curve fit coefficient for "
136 : "GQ/H model"
137 : "for each soil layer.");
138 148 : params.addParam<std::vector<Real>>("theta_3",
139 : "The curve fit coefficient for "
140 : "GQ/H model"
141 : "for each soil layer.");
142 148 : params.addParam<std::vector<Real>>("theta_4",
143 : "The curve fit coefficient for "
144 : "GQ/H model"
145 : "for each soil layer.");
146 148 : params.addParam<std::vector<Real>>("theta_5",
147 : "The curve fit coefficient for "
148 : "GQ/H model"
149 : "for each soil layer.");
150 148 : params.addParam<std::vector<Real>>("taumax",
151 : "The ultimate shear strength of "
152 : "the soil layers. Required for "
153 : "GQ/H model");
154 : // params required for thin_layer contact model soil_type = "thin_layer"
155 148 : params.addParam<std::vector<Real>>("friction_coefficient",
156 : "Friction coefficients of the thin layers.");
157 148 : params.addParam<std::vector<Real>>("hardening_ratio",
158 : "Post-yield hardening ratios of the layers. "
159 : "Defaults to 0.01, which corresponds to 1 percent hardening.");
160 :
161 : // flag to use finite strain calculator
162 148 : params.addParam<bool>("finite_strain",
163 148 : false,
164 : "Set to true to use finite "
165 : "strain calculator instead of "
166 : "incremental small strain.");
167 :
168 148 : params.addRequiredParam<std::vector<Real>>(
169 : "density",
170 : "Vector of density values that map one-to-one with the number "
171 : "'layer_ids' parameter.");
172 148 : params.addParam<bool>("use_automatic_differentiation",
173 148 : false,
174 : "Flag to use automatic differentiation (AD) objects when possible");
175 74 : return params;
176 74 : }
177 :
178 74 : ISoilAction::ISoilAction(const InputParameters & params) : Action(params) {}
179 :
180 : void
181 74 : ISoilAction::act()
182 : {
183 148 : const bool use_ad = getParam<bool>("use_automatic_differentiation");
184 74 : std::string ad_prepend = "";
185 74 : if (use_ad)
186 : ad_prepend = "AD";
187 :
188 148 : std::vector<SubdomainName> block = getParam<std::vector<SubdomainName>>("block");
189 148 : std::vector<VariableName> layer_variable = getParam<std::vector<VariableName>>("layer_variable");
190 148 : std::vector<unsigned int> layer_ids = getParam<std::vector<unsigned int>>("layer_ids");
191 222 : std::vector<Real> poissons_ratio = getParam<std::vector<Real>>("poissons_ratio");
192 74 : if (poissons_ratio.size() != layer_ids.size())
193 0 : mooseError("Error in" + name() + ". Poisson's ratio should be of the same size as layer_ids.");
194 222 : std::vector<Real> density = getParam<std::vector<Real>>("density");
195 74 : if (density.size() != layer_ids.size())
196 0 : mooseError("Error in" + name() + ". Density should be of the same size as layer_ids.");
197 148 : MooseEnum soil_type = getParam<MooseEnum>("soil_type");
198 : // Stress calculation
199 74 : auto params = _factory.getValidParams(ad_prepend + "ComputeISoilStress");
200 74 : params.set<std::vector<unsigned int>>("layer_ids") = layer_ids;
201 74 : params.set<std::vector<VariableName>>("layer_variable") = layer_variable;
202 74 : params.set<std::vector<SubdomainName>>("block") = block;
203 74 : params.set<std::vector<Real>>("poissons_ratio") = poissons_ratio;
204 148 : params.set<Real>("b_exp") = getParam<Real>("b_exp");
205 222 : params.set<std::vector<Real>>("p_ref") = getParam<std::vector<Real>>("p_ref");
206 148 : params.set<Real>("a0") = getParam<Real>("a0");
207 148 : params.set<Real>("a1") = getParam<Real>("a1");
208 148 : params.set<Real>("a2") = getParam<Real>("a2");
209 148 : params.set<Real>("tension_pressure_cut_off") = getParam<Real>("tension_pressure_cut_off");
210 148 : params.set<bool>("pressure_dependency") = getParam<bool>("pressure_dependency");
211 148 : params.set<bool>("implicit") = getParam<bool>("implicit");
212 74 : params.set<bool>("wave_speed_calculation") = true;
213 148 : params.set<std::vector<FunctionName>>("initial_soil_stress") =
214 148 : getParam<std::vector<FunctionName>>("initial_soil_stress");
215 74 : params.set<MooseEnum>("soil_type") = soil_type;
216 : // setting soiltype specific parameters
217 74 : if (soil_type == "user_defined")
218 46 : params.set<std::vector<FileName>>("backbone_curve_files") =
219 69 : getParam<std::vector<FileName>>("backbone_curve_files");
220 148 : params.set<std::vector<Real>>("initial_shear_modulus") =
221 148 : getParam<std::vector<Real>>("initial_shear_modulus");
222 74 : if (soil_type == "darendeli")
223 : {
224 36 : params.set<std::vector<Real>>("initial_shear_modulus") =
225 36 : getParam<std::vector<Real>>("initial_shear_modulus");
226 36 : params.set<unsigned int>("number_of_points") = getParam<unsigned int>("number_of_points");
227 36 : params.set<std::vector<Real>>("over_consolidation_ratio") =
228 36 : getParam<std::vector<Real>>("over_consolidation_ratio");
229 36 : params.set<std::vector<Real>>("plasticity_index") =
230 54 : getParam<std::vector<Real>>("plasticity_index");
231 : }
232 74 : if (soil_type == "gqh")
233 : {
234 54 : params.set<std::vector<Real>>("initial_shear_modulus") =
235 54 : getParam<std::vector<Real>>("initial_shear_modulus");
236 54 : params.set<unsigned int>("number_of_points") = getParam<unsigned int>("number_of_points");
237 81 : params.set<std::vector<Real>>("theta_1") = getParam<std::vector<Real>>("theta_1");
238 81 : params.set<std::vector<Real>>("theta_2") = getParam<std::vector<Real>>("theta_2");
239 81 : params.set<std::vector<Real>>("theta_3") = getParam<std::vector<Real>>("theta_3");
240 81 : params.set<std::vector<Real>>("theta_4") = getParam<std::vector<Real>>("theta_4");
241 81 : params.set<std::vector<Real>>("theta_5") = getParam<std::vector<Real>>("theta_5");
242 81 : params.set<std::vector<Real>>("taumax") = getParam<std::vector<Real>>("taumax");
243 : }
244 74 : if (soil_type == "thin_layer")
245 : {
246 12 : params.set<std::vector<Real>>("initial_shear_modulus") =
247 12 : getParam<std::vector<Real>>("initial_shear_modulus");
248 12 : params.set<std::vector<Real>>("friction_coefficient") =
249 12 : getParam<std::vector<Real>>("friction_coefficient");
250 12 : params.set<std::vector<Real>>("hardening_ratio") =
251 18 : getParam<std::vector<Real>>("hardening_ratio");
252 : }
253 74 : if (use_ad)
254 : {
255 36 : _problem->addMaterial(ad_prepend + "ComputeISoilStress", name() + "_stress" + block[0], params);
256 18 : _problem->haveADObjects(true);
257 : }
258 : else
259 110 : _problem->addMaterial("ComputeISoilStress", "stress_" + block[0], params);
260 :
261 : // strain calculation
262 144 : std::vector<VariableName> displacements = getParam<std::vector<VariableName>>("displacements");
263 144 : bool finite_strain = getParam<bool>("finite_strain");
264 72 : if (!finite_strain && soil_type != "thin_layer")
265 : {
266 : // create small incremental strain block
267 132 : params = _factory.getValidParams(ad_prepend + "ComputeIncrementalSmallStrain");
268 66 : std::string unique_strain_name = "strain_" + block[0];
269 66 : params.set<std::vector<SubdomainName>>("block") = block;
270 66 : params.set<std::vector<VariableName>>("displacements") = displacements;
271 : // params.set<bool>("stateful_displacements") = true; // deprecated
272 66 : if (use_ad)
273 : {
274 30 : _problem->addMaterial(
275 30 : ad_prepend + "ComputeIncrementalSmallStrain", name() + "_strain" + block[0], params);
276 15 : _problem->haveADObjects(true);
277 : }
278 : else
279 102 : _problem->addMaterial("ComputeIncrementalSmallStrain", unique_strain_name, params);
280 : }
281 : else
282 : {
283 6 : if (soil_type == "thin_layer")
284 : // create finite strain block
285 12 : params = _factory.getValidParams(ad_prepend + "ComputeFiniteStrain");
286 6 : std::string unique_strain_name = "strain_" + block[0];
287 6 : params.set<std::vector<SubdomainName>>("block") = block;
288 6 : params.set<std::vector<VariableName>>("displacements") = displacements;
289 : // params.set<bool>("stateful_displacements") = true; // deprecated
290 6 : if (use_ad)
291 : {
292 6 : _problem->addMaterial(
293 6 : ad_prepend + "ComputeFiniteStrain", name() + "_strain" + block[0], params);
294 3 : _problem->haveADObjects(true);
295 : }
296 : else
297 6 : _problem->addMaterial("ComputeFiniteStrain", unique_strain_name, params);
298 : }
299 :
300 144 : params = _factory.getValidParams("ComputeIsotropicElasticityTensorSoil");
301 72 : std::vector<Real> elastic_mod(layer_ids.size());
302 144 : std::vector<Real> shear_mod = getParam<std::vector<Real>>("initial_shear_modulus");
303 318 : for (std::size_t i = 0; i < layer_ids.size(); i++)
304 246 : elastic_mod[i] = 2 * shear_mod[i] * (1 + poissons_ratio[i]);
305 72 : std::string unique_elasticity_name = "elasticity_" + block[0];
306 72 : params.set<std::vector<SubdomainName>>("block") = block;
307 72 : params.set<std::vector<Real>>("elastic_modulus") = elastic_mod;
308 72 : params.set<std::vector<Real>>("poissons_ratio") = poissons_ratio;
309 72 : params.set<std::vector<Real>>("density") = density;
310 72 : params.set<bool>("wave_speed_calculation") = false;
311 72 : params.set<std::vector<unsigned int>>("layer_ids") = layer_ids;
312 72 : params.set<std::vector<VariableName>>("layer_variable") = layer_variable;
313 72 : _problem->addMaterial(
314 144 : ad_prepend + "ComputeIsotropicElasticityTensorSoil", unique_elasticity_name, params);
315 144 : }
|
