www.mooseframework.org
DomainIntegralAction.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 // MOOSE includes
9 #include "DomainIntegralAction.h"
10 #include "Factory.h"
11 #include "FEProblem.h"
12 #include "Parser.h"
13 #include "CrackFrontDefinition.h"
14 #include "MooseMesh.h"
15 
16 #include "libmesh/string_to_enum.h"
17 
18 template <>
19 InputParameters
21 {
22  InputParameters params = validParams<Action>();
24  MultiMooseEnum integral_vec("JIntegral InteractionIntegralKI InteractionIntegralKII "
25  "InteractionIntegralKIII InteractionIntegralT");
26  params.addRequiredParam<MultiMooseEnum>("integrals",
27  integral_vec,
28  "Domain integrals to calculate. Choices are: " +
29  integral_vec.getRawNames());
30  params.addParam<std::vector<BoundaryName>>(
31  "boundary", "The list of boundary IDs from the mesh where this boundary condition applies");
32  params.addParam<std::vector<Point>>("crack_front_points", "Set of points to define crack front");
33  params.addParam<std::string>(
34  "order", "FIRST", "Specifies the order of the FE shape function to use for q AuxVariables");
35  params.addParam<std::string>(
36  "family", "LAGRANGE", "Specifies the family of FE shape functions to use for q AuxVariables");
37  params.addParam<std::vector<Real>>("radius_inner", "Inner radius for volume integral domain");
38  params.addParam<std::vector<Real>>("radius_outer", "Outer radius for volume integral domain");
39  params.addParam<unsigned int>("ring_first",
40  "The first ring of elements for volume integral domain");
41  params.addParam<unsigned int>("ring_last",
42  "The last ring of elements for volume integral domain");
43  params.addParam<std::vector<VariableName>>(
44  "output_variable", "Variable values to be reported along the crack front");
45  params.addParam<bool>(
46  "convert_J_to_K", false, "Convert J-integral to stress intensity factor K.");
47  params.addParam<Real>("poissons_ratio", "Poisson's ratio");
48  params.addParam<Real>("youngs_modulus", "Young's modulus");
49  params.addParam<std::vector<SubdomainName>>("block", "The block ids where integrals are defined");
50 
51  params.addParam<std::vector<VariableName>>(
52  "displacements",
53  "The displacements appropriate for the simulation geometry and coordinate system");
54  params.addParam<VariableName>("disp_x", "The x displacement");
55  params.addParam<VariableName>("disp_y", "The y displacement");
56  params.addParam<VariableName>("disp_z", "The z displacement");
57  params.addParam<VariableName>("temp", "", "The temperature");
58  MooseEnum position_type("Angle Distance", "Distance");
59  params.addParam<MooseEnum>(
60  "position_type",
61  position_type,
62  "The method used to calculate position along crack front. Options are: " +
63  position_type.getRawNames());
64  MooseEnum q_function_type("Geometry Topology", "Geometry");
65  params.addParam<MooseEnum>("q_function_type",
66  q_function_type,
67  "The method used to define the integration domain. Options are: " +
68  q_function_type.getRawNames());
69  params.addParam<bool>(
70  "equivalent_k",
71  false,
72  "Calculate an equivalent K from KI, KII and KIII, assuming self-similar crack growth.");
73  params.addParam<bool>("output_q", true, "Output q");
74  params.addParam<bool>("solid_mechanics",
75  false,
76  "Set to true if the solid_mechanics system is "
77  "used. This option is only needed for "
78  "interaction integrals.");
79  params.addParam<std::vector<MaterialPropertyName>>(
80  "eigenstrain_names", "List of eigenstrains applied in the strain calculation");
81  return params;
82 }
83 
84 DomainIntegralAction::DomainIntegralAction(const InputParameters & params)
85  : Action(params),
86  _boundary_names(getParam<std::vector<BoundaryName>>("boundary")),
87  _closed_loop(getParam<bool>("closed_loop")),
88  _use_crack_front_points_provider(false),
89  _order(getParam<std::string>("order")),
90  _family(getParam<std::string>("family")),
91  _direction_method_moose_enum(getParam<MooseEnum>("crack_direction_method")),
92  _end_direction_method_moose_enum(getParam<MooseEnum>("crack_end_direction_method")),
93  _have_crack_direction_vector(false),
94  _have_crack_direction_vector_end_1(false),
95  _have_crack_direction_vector_end_2(false),
96  _treat_as_2d(getParam<bool>("2d")),
97  _axis_2d(getParam<unsigned int>("axis_2d")),
98  _convert_J_to_K(false),
99  _has_symmetry_plane(isParamValid("symmetry_plane")),
100  _symmetry_plane(_has_symmetry_plane ? getParam<unsigned int>("symmetry_plane")
101  : std::numeric_limits<unsigned int>::max()),
102  _position_type(getParam<MooseEnum>("position_type")),
103  _q_function_type(getParam<MooseEnum>("q_function_type")),
104  _get_equivalent_k(getParam<bool>("equivalent_k")),
105  _use_displaced_mesh(false),
106  _output_q(getParam<bool>("output_q")),
107  _solid_mechanics(getParam<bool>("solid_mechanics"))
108 {
109  if (_q_function_type == GEOMETRY)
110  {
111  if (isParamValid("radius_inner") && isParamValid("radius_outer"))
112  {
113  _radius_inner = getParam<std::vector<Real>>("radius_inner");
114  _radius_outer = getParam<std::vector<Real>>("radius_outer");
115  }
116  else
117  mooseError("DomainIntegral error: must set radius_inner and radius_outer.");
118  for (unsigned int i = 0; i < _radius_inner.size(); ++i)
119  _ring_vec.push_back(i + 1);
120  }
121  else if (_q_function_type == TOPOLOGY)
122  {
123  if (isParamValid("ring_first") && isParamValid("ring_last"))
124  {
125  _ring_first = getParam<unsigned int>("ring_first");
126  _ring_last = getParam<unsigned int>("ring_last");
127  }
128  else
129  mooseError(
130  "DomainIntegral error: must set ring_first and ring_last if q_function_type = Topology.");
131  for (unsigned int i = _ring_first; i <= _ring_last; ++i)
132  _ring_vec.push_back(i);
133  }
134  else
135  mooseError("DomainIntegral error: invalid q_function_type.");
136 
137  if (isParamValid("crack_front_points"))
138  _crack_front_points = getParam<std::vector<Point>>("crack_front_points");
139 
140  if (isParamValid("crack_front_points_provider"))
141  {
142  if (!isParamValid("number_points_from_provider"))
143  mooseError("DomainIntegral error: when crack_front_points_provider is used, the "
144  "number_points_from_provider must be "
145  "provided.");
147  _crack_front_points_provider = getParam<UserObjectName>("crack_front_points_provider");
148  }
149  else if (isParamValid("number_points_from_provider"))
150  mooseError("DomainIntegral error: number_points_from_provider is provided but "
151  "crack_front_points_provider cannot "
152  "be found.");
153  if (isParamValid("crack_direction_vector"))
154  {
155  _crack_direction_vector = getParam<RealVectorValue>("crack_direction_vector");
157  }
158  if (isParamValid("crack_direction_vector_end_1"))
159  {
160  _crack_direction_vector_end_1 = getParam<RealVectorValue>("crack_direction_vector_end_1");
162  }
163  if (isParamValid("crack_direction_vector_end_2"))
164  {
165  _crack_direction_vector_end_2 = getParam<RealVectorValue>("crack_direction_vector_end_2");
167  }
168  if (isParamValid("crack_mouth_boundary"))
169  _crack_mouth_boundary_names = getParam<std::vector<BoundaryName>>("crack_mouth_boundary");
170  if (isParamValid("intersecting_boundary"))
171  _intersecting_boundary_names = getParam<std::vector<BoundaryName>>("intersecting_boundary");
172  if (_radius_inner.size() != _radius_outer.size())
173  mooseError("Number of entries in 'radius_inner' and 'radius_outer' must match.");
174 
175  bool youngs_modulus_set(false);
176  bool poissons_ratio_set(false);
177  MultiMooseEnum integral_moose_enums = getParam<MultiMooseEnum>("integrals");
178  if (integral_moose_enums.size() == 0)
179  mooseError("Must specify at least one domain integral to perform.");
180  for (unsigned int i = 0; i < integral_moose_enums.size(); ++i)
181  {
182  if (integral_moose_enums[i] != "JIntegral")
183  {
184  // Check that parameters required for interaction integrals are defined
185  if (isParamValid("displacements"))
186  {
187  _displacements = getParam<std::vector<VariableName>>("displacements");
188 
189  if (_displacements.size() < 2)
190  mooseError(
191  "DomainIntegral error: The size of the displacements vector should atleast be 2.");
192  }
193  else
194  {
195  if (isParamValid("disp_x") || isParamValid("disp_y") || isParamValid("disp_z"))
196  mooseDeprecated("DomainIntegral Warning: disp_x, disp_y and disp_z are deprecated. "
197  "Please specify displacements using the `dispalcements` parameter.");
198 
199  if (!isParamValid("disp_x") || !isParamValid("disp_y"))
200  mooseError(
201  "DomainIntegral error: Specify displacements using the `displacements` parameter.");
202  else
203  {
204  _displacements.clear();
205  _displacements.push_back(getParam<VariableName>("disp_x"));
206  _displacements.push_back(getParam<VariableName>("disp_y"));
207  if (isParamValid("disp_z"))
208  _displacements.push_back(getParam<VariableName>("disp_z"));
209  }
210  }
211 
212  if (!(isParamValid("poissons_ratio")) || !(isParamValid("youngs_modulus")))
213  mooseError(
214  "DomainIntegral error: must set Poisson's ratio and Young's modulus for integral: ",
215  integral_moose_enums[i]);
216 
217  if (!(isParamValid("block")))
218  mooseError("DomainIntegral error: must set block ID or name for integral: ",
219  integral_moose_enums[i]);
220 
221  _poissons_ratio = getParam<Real>("poissons_ratio");
222  poissons_ratio_set = true;
223  _youngs_modulus = getParam<Real>("youngs_modulus");
224  youngs_modulus_set = true;
225  _blocks = getParam<std::vector<SubdomainName>>("block");
226  }
227 
228  _integrals.insert(INTEGRAL(int(integral_moose_enums.get(i))));
229  }
230 
231  if (isParamValid("temp"))
232  _temp = getParam<VariableName>("temp");
233 
234  if (_temp != "" && !isParamValid("eigenstrain_names") && !_solid_mechanics)
235  mooseError(
236  "DomainIntegral error: must provide `eigenstrain_names` when temperature is coupled.");
237 
239  _integrals.count(INTERACTION_INTEGRAL_KII) == 0 ||
241  mooseError("DomainIntegral error: must calculate KI, KII and KIII to get equivalent K.");
242 
243  if (isParamValid("output_variable"))
244  {
245  _output_variables = getParam<std::vector<VariableName>>("output_variable");
246  if (_crack_front_points.size() > 0)
247  mooseError("'output_variables' not yet supported with 'crack_front_points'");
248  }
249 
250  if (isParamValid("convert_J_to_K"))
251  _convert_J_to_K = getParam<bool>("convert_J_to_K");
252  if (_convert_J_to_K)
253  {
254  if (!isParamValid("youngs_modulus") || !isParamValid("poissons_ratio"))
255  mooseError("DomainIntegral error: must set Young's modulus and Poisson's ratio for "
256  "J-integral if convert_J_to_K = true.");
257  if (!youngs_modulus_set)
258  _youngs_modulus = getParam<Real>("youngs_modulus");
259  if (!poissons_ratio_set)
260  _poissons_ratio = getParam<Real>("poissons_ratio");
261  }
262 }
263 
265 
266 void
268 {
269  const std::string uo_name("crackFrontDefinition");
270  const std::string ak_base_name("q");
271  const std::string av_base_name("q");
272  const unsigned int num_crack_front_points = calcNumCrackFrontPoints();
273  const std::string aux_stress_base_name("aux_stress");
274  const std::string aux_grad_disp_base_name("aux_grad_disp");
275 
276  if (_current_task == "add_user_object")
277  {
278  const std::string uo_type_name("CrackFrontDefinition");
279 
280  InputParameters params = _factory.getValidParams(uo_type_name);
281  params.set<MultiMooseEnum>("execute_on") = "initial timestep_end";
282  params.set<MooseEnum>("crack_direction_method") = _direction_method_moose_enum;
283  params.set<MooseEnum>("crack_end_direction_method") = _end_direction_method_moose_enum;
285  params.set<RealVectorValue>("crack_direction_vector") = _crack_direction_vector;
287  params.set<RealVectorValue>("crack_direction_vector_end_1") = _crack_direction_vector_end_1;
289  params.set<RealVectorValue>("crack_direction_vector_end_2") = _crack_direction_vector_end_2;
290  if (_crack_mouth_boundary_names.size() != 0)
291  params.set<std::vector<BoundaryName>>("crack_mouth_boundary") = _crack_mouth_boundary_names;
292  if (_intersecting_boundary_names.size() != 0)
293  params.set<std::vector<BoundaryName>>("intersecting_boundary") = _intersecting_boundary_names;
294  params.set<bool>("2d") = _treat_as_2d;
295  params.set<unsigned int>("axis_2d") = _axis_2d;
297  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
298  if (_boundary_names.size() != 0)
299  params.set<std::vector<BoundaryName>>("boundary") = _boundary_names;
300  if (_crack_front_points.size() != 0)
301  params.set<std::vector<Point>>("crack_front_points") = _crack_front_points;
303  params.applyParameters(parameters(),
304  {"crack_front_points_provider, number_points_from_provider"});
305  if (_closed_loop)
306  params.set<bool>("closed_loop") = _closed_loop;
307  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
308  if (_integrals.count(INTERACTION_INTEGRAL_T) != 0)
309  {
310  params.set<VariableName>("disp_x") = _displacements[0];
311  params.set<VariableName>("disp_y") = _displacements[1];
312  if (_displacements.size() == 3)
313  params.set<VariableName>("disp_z") = _displacements[2];
314  params.set<bool>("t_stress") = true;
315  }
316 
317  unsigned int nrings = 0;
318  if (_q_function_type == TOPOLOGY)
319  {
320  params.set<bool>("q_function_rings") = true;
321  params.set<unsigned int>("last_ring") = _ring_last;
322  params.set<unsigned int>("first_ring") = _ring_first;
323  nrings = _ring_last - _ring_first + 1;
324  }
325  else if (_q_function_type == GEOMETRY)
326  {
327  params.set<std::vector<Real>>("j_integral_radius_inner") = _radius_inner;
328  params.set<std::vector<Real>>("j_integral_radius_outer") = _radius_outer;
329  nrings = _ring_vec.size();
330  }
331 
332  params.set<unsigned int>("nrings") = nrings;
333  params.set<MooseEnum>("q_function_type") = _q_function_type;
334 
335  _problem->addUserObject(uo_type_name, uo_name, params);
336  }
337  else if (_current_task == "add_aux_variable" && _output_q)
338  {
339  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
340  {
341  if (_treat_as_2d)
342  {
343  std::ostringstream av_name_stream;
344  av_name_stream << av_base_name << "_" << _ring_vec[ring_index];
345  _problem->addAuxVariable(av_name_stream.str(),
346  FEType(Utility::string_to_enum<Order>(_order),
347  Utility::string_to_enum<FEFamily>(_family)));
348  }
349  else
350  {
351  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
352  {
353  std::ostringstream av_name_stream;
354  av_name_stream << av_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
355  _problem->addAuxVariable(av_name_stream.str(),
356  FEType(Utility::string_to_enum<Order>(_order),
357  Utility::string_to_enum<FEFamily>(_family)));
358  }
359  }
360  }
361  }
362 
363  else if (_current_task == "add_aux_kernel" && _output_q)
364  {
365  std::string ak_type_name;
366  unsigned int nrings = 0;
367  if (_q_function_type == GEOMETRY)
368  {
369  ak_type_name = "DomainIntegralQFunction";
370  nrings = _ring_vec.size();
371  }
372  else if (_q_function_type == TOPOLOGY)
373  {
374  ak_type_name = "DomainIntegralTopologicalQFunction";
375  nrings = _ring_last - _ring_first + 1;
376  }
377 
378  InputParameters params = _factory.getValidParams(ak_type_name);
379  params.set<MultiMooseEnum>("execute_on") = "initial timestep_end";
380  params.set<UserObjectName>("crack_front_definition") = uo_name;
381  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
382 
383  for (unsigned int ring_index = 0; ring_index < nrings; ++ring_index)
384  {
385  if (_q_function_type == GEOMETRY)
386  {
387  params.set<Real>("j_integral_radius_inner") = _radius_inner[ring_index];
388  params.set<Real>("j_integral_radius_outer") = _radius_outer[ring_index];
389  }
390  else if (_q_function_type == TOPOLOGY)
391  {
392  params.set<unsigned int>("ring_index") = _ring_first + ring_index;
393  }
394 
395  if (_treat_as_2d)
396  {
397  std::ostringstream ak_name_stream;
398  ak_name_stream << ak_base_name << "_" << _ring_vec[ring_index];
399  std::ostringstream av_name_stream;
400  av_name_stream << av_base_name << "_" << _ring_vec[ring_index];
401  params.set<AuxVariableName>("variable") = av_name_stream.str();
402  _problem->addAuxKernel(ak_type_name, ak_name_stream.str(), params);
403  }
404  else
405  {
406  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
407  {
408  std::ostringstream ak_name_stream;
409  ak_name_stream << ak_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
410  std::ostringstream av_name_stream;
411  av_name_stream << av_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
412  params.set<AuxVariableName>("variable") = av_name_stream.str();
413  params.set<unsigned int>("crack_front_point_index") = cfp_index;
414  _problem->addAuxKernel(ak_type_name, ak_name_stream.str(), params);
415  }
416  }
417  }
418  }
419 
420  else if (_current_task == "add_postprocessor")
421  {
422  if (_integrals.count(J_INTEGRAL) != 0)
423  {
424  std::string pp_base_name;
425  if (_convert_J_to_K)
426  pp_base_name = "K";
427  else
428  pp_base_name = "J";
429  const std::string pp_type_name("JIntegral");
430  InputParameters params = _factory.getValidParams(pp_type_name);
431  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
432  params.set<UserObjectName>("crack_front_definition") = uo_name;
433  params.set<bool>("convert_J_to_K") = _convert_J_to_K;
434  if (_convert_J_to_K)
435  {
436  params.set<Real>("youngs_modulus") = _youngs_modulus;
437  params.set<Real>("poissons_ratio") = _poissons_ratio;
438  }
440  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
441  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
442  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
443  {
444  params.set<unsigned int>("ring_index") = _ring_vec[ring_index];
445  if (_q_function_type == TOPOLOGY)
446  params.set<unsigned int>("ring_first") = _ring_first;
447  params.set<MooseEnum>("q_function_type") = _q_function_type;
448 
449  if (_treat_as_2d)
450  {
451  std::ostringstream av_name_stream;
452  av_name_stream << av_base_name << "_" << _ring_vec[ring_index];
453  std::ostringstream pp_name_stream;
454  pp_name_stream << pp_base_name << "_" << _ring_vec[ring_index];
455  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
456  }
457  else
458  {
459  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
460  {
461  std::ostringstream av_name_stream;
462  av_name_stream << av_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
463  std::ostringstream pp_name_stream;
464  pp_name_stream << pp_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
465  params.set<unsigned int>("crack_front_point_index") = cfp_index;
466  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
467  }
468  }
469  }
470  }
471  if (_integrals.count(INTERACTION_INTEGRAL_KI) != 0 ||
472  _integrals.count(INTERACTION_INTEGRAL_KII) != 0 ||
475  {
476 
479  mooseError("In DomainIntegral, symmetry_plane option cannot be used with mode-II or "
480  "mode-III interaction integral");
481 
482  const std::string pp_base_name("II");
483  std::string pp_type_name("InteractionIntegral");
484 
485  if (_solid_mechanics)
486  pp_type_name = "InteractionIntegralSM";
487 
488  InputParameters params = _factory.getValidParams(pp_type_name);
489  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
490  params.set<UserObjectName>("crack_front_definition") = uo_name;
491  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
493  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
494  params.set<Real>("poissons_ratio") = _poissons_ratio;
495  params.set<Real>("youngs_modulus") = _youngs_modulus;
496  params.set<std::vector<VariableName>>("displacements") = _displacements;
497  if (_temp != "")
498  params.set<std::vector<VariableName>>("temp") = {_temp};
500  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
501 
502  for (std::set<INTEGRAL>::iterator sit = _integrals.begin(); sit != _integrals.end(); ++sit)
503  {
504  std::string pp_base_name;
505  std::string aux_mode_name;
506  switch (*sit)
507  {
508  case J_INTEGRAL:
509  continue;
510 
512  pp_base_name = "II_KI";
513  aux_mode_name = "_I_";
514  params.set<Real>("K_factor") =
515  0.5 * _youngs_modulus / (1.0 - std::pow(_poissons_ratio, 2.0));
516  params.set<MooseEnum>("sif_mode") = "KI";
517  break;
518 
520  pp_base_name = "II_KII";
521  aux_mode_name = "_II_";
522  params.set<Real>("K_factor") =
523  0.5 * _youngs_modulus / (1.0 - std::pow(_poissons_ratio, 2.0));
524  params.set<MooseEnum>("sif_mode") = "KII";
525  break;
526 
528  pp_base_name = "II_KIII";
529  aux_mode_name = "_III_";
530  params.set<Real>("K_factor") = 0.5 * _youngs_modulus / (1.0 + _poissons_ratio);
531  params.set<MooseEnum>("sif_mode") = "KIII";
532  break;
533 
535  pp_base_name = "II_T";
536  aux_mode_name = "_T_";
537  params.set<Real>("K_factor") = _youngs_modulus / (1 - std::pow(_poissons_ratio, 2));
538  params.set<MooseEnum>("sif_mode") = "T";
539  break;
540  }
541 
542  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
543  {
544  params.set<unsigned int>("ring_index") = _ring_vec[ring_index];
545  params.set<unsigned int>("ring_first") = _ring_first;
546  params.set<MooseEnum>("q_function_type") = _q_function_type;
547 
548  if (_treat_as_2d)
549  {
550  std::ostringstream av_name_stream;
551  av_name_stream << av_base_name << "_" << _ring_vec[ring_index];
552  std::ostringstream pp_name_stream;
553  pp_name_stream << pp_base_name << "_" << _ring_vec[ring_index];
554 
555  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
556  }
557  else
558  {
559  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
560  {
561  std::ostringstream av_name_stream;
562  av_name_stream << av_base_name << "_" << cfp_index + 1 << "_"
563  << _ring_vec[ring_index];
564  std::ostringstream pp_name_stream;
565  pp_name_stream << pp_base_name << "_" << cfp_index + 1 << "_"
566  << _ring_vec[ring_index];
567  std::ostringstream cfn_index_stream;
568  cfn_index_stream << cfp_index + 1;
569 
570  params.set<unsigned int>("crack_front_point_index") = cfp_index;
571  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
572  }
573  }
574  }
575  }
576  }
577  for (unsigned int i = 0; i < _output_variables.size(); ++i)
578  {
579  const std::string ov_base_name(_output_variables[i]);
580  const std::string pp_type_name("CrackFrontData");
581  InputParameters params = _factory.getValidParams(pp_type_name);
582  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
583  params.set<UserObjectName>("crack_front_definition") = uo_name;
584  if (_treat_as_2d)
585  {
586  std::ostringstream pp_name_stream;
587  pp_name_stream << ov_base_name << "_crack";
588  params.set<VariableName>("variable") = _output_variables[i];
589  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
590  }
591  else
592  {
593  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
594  {
595  std::ostringstream pp_name_stream;
596  pp_name_stream << ov_base_name << "_crack_" << cfp_index + 1;
597  params.set<VariableName>("variable") = _output_variables[i];
598  params.set<unsigned int>("crack_front_point_index") = cfp_index;
599  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
600  }
601  }
602  }
603  if (_get_equivalent_k)
604  {
605  std::string pp_base_name("Keq");
606  const std::string pp_type_name("MixedModeEquivalentK");
607  InputParameters params = _factory.getValidParams(pp_type_name);
608  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
609  params.set<Real>("poissons_ratio") = _poissons_ratio;
610  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
611  {
612  if (_treat_as_2d)
613  {
614  std::ostringstream ki_name_stream;
615  ki_name_stream << "II_KI_" << _ring_vec[ring_index];
616  std::ostringstream kii_name_stream;
617  kii_name_stream << "II_KII_" << _ring_vec[ring_index];
618  std::ostringstream kiii_name_stream;
619  kiii_name_stream << "II_KIII_" << _ring_vec[ring_index];
620  params.set<PostprocessorName>("KI_name") = ki_name_stream.str();
621  params.set<PostprocessorName>("KII_name") = kii_name_stream.str();
622  params.set<PostprocessorName>("KIII_name") = kiii_name_stream.str();
623  std::ostringstream pp_name_stream;
624  pp_name_stream << pp_base_name << "_" << _ring_vec[ring_index];
625  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
626  }
627  else
628  {
629  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
630  {
631  std::ostringstream ki_name_stream;
632  ki_name_stream << "II_KI_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
633  std::ostringstream kii_name_stream;
634  kii_name_stream << "II_KII_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
635  std::ostringstream kiii_name_stream;
636  kiii_name_stream << "II_KIII_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
637  params.set<PostprocessorName>("KI_name") = ki_name_stream.str();
638  params.set<PostprocessorName>("KII_name") = kii_name_stream.str();
639  params.set<PostprocessorName>("KIII_name") = kiii_name_stream.str();
640  std::ostringstream pp_name_stream;
641  pp_name_stream << pp_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
642  params.set<unsigned int>("crack_front_point_index") = cfp_index;
643  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
644  }
645  }
646  }
647  }
648  }
649 
650  else if (_current_task == "add_vector_postprocessor")
651  {
652  if (!_treat_as_2d)
653  {
654  for (std::set<INTEGRAL>::iterator sit = _integrals.begin(); sit != _integrals.end(); ++sit)
655  {
656  std::string pp_base_name;
657  switch (*sit)
658  {
659  case J_INTEGRAL:
660  if (_convert_J_to_K)
661  pp_base_name = "K";
662  else
663  pp_base_name = "J";
664  break;
666  pp_base_name = "II_KI";
667  break;
669  pp_base_name = "II_KII";
670  break;
672  pp_base_name = "II_KIII";
673  break;
675  pp_base_name = "II_T";
676  break;
677  }
678  const std::string vpp_type_name("CrackDataSampler");
679  InputParameters params = _factory.getValidParams(vpp_type_name);
680  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
681  params.set<UserObjectName>("crack_front_definition") = uo_name;
682  params.set<MooseEnum>("sort_by") = "id";
683  params.set<MooseEnum>("position_type") = _position_type;
684  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
685  {
686  std::vector<PostprocessorName> postprocessor_names;
687  std::ostringstream vpp_name_stream;
688  vpp_name_stream << pp_base_name << "_" << _ring_vec[ring_index];
689  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
690  {
691  std::ostringstream pp_name_stream;
692  pp_name_stream << pp_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
693  postprocessor_names.push_back(pp_name_stream.str());
694  }
695  params.set<std::vector<PostprocessorName>>("postprocessors") = postprocessor_names;
696  _problem->addVectorPostprocessor(vpp_type_name, vpp_name_stream.str(), params);
697  }
698  }
699 
700  for (unsigned int i = 0; i < _output_variables.size(); ++i)
701  {
702  const std::string vpp_type_name("VectorOfPostprocessors");
703  InputParameters params = _factory.getValidParams(vpp_type_name);
704  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
705  std::ostringstream vpp_name_stream;
706  vpp_name_stream << _output_variables[i] << "_crack";
707  std::vector<PostprocessorName> postprocessor_names;
708  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
709  {
710  std::ostringstream pp_name_stream;
711  pp_name_stream << vpp_name_stream.str() << "_" << cfp_index + 1;
712  postprocessor_names.push_back(pp_name_stream.str());
713  }
714  params.set<std::vector<PostprocessorName>>("postprocessors") = postprocessor_names;
715  _problem->addVectorPostprocessor(vpp_type_name, vpp_name_stream.str(), params);
716  }
717  }
719  {
720  std::string pp_base_name("Keq");
721  const std::string vpp_type_name("CrackDataSampler");
722  InputParameters params = _factory.getValidParams(vpp_type_name);
723  params.set<MultiMooseEnum>("execute_on") = "timestep_end";
724  params.set<UserObjectName>("crack_front_definition") = uo_name;
725  params.set<MooseEnum>("sort_by") = "id";
726  params.set<MooseEnum>("position_type") = _position_type;
727  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
728  {
729  std::vector<PostprocessorName> postprocessor_names;
730  std::ostringstream vpp_name_stream;
731  vpp_name_stream << pp_base_name << "_" << _ring_vec[ring_index];
732  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
733  {
734  std::ostringstream pp_name_stream;
735  pp_name_stream << pp_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
736  postprocessor_names.push_back(pp_name_stream.str());
737  }
738  params.set<std::vector<PostprocessorName>>("postprocessors") = postprocessor_names;
739  _problem->addVectorPostprocessor(vpp_type_name, vpp_name_stream.str(), params);
740  }
741  }
742  }
743 
744  else if (_current_task == "add_material")
745  {
746  if (_temp != "" && !_solid_mechanics)
747  {
748  std::string mater_name;
749  const std::string mater_type_name("ThermalFractureIntegral");
750  if (isParamValid("blocks"))
751  {
752  _blocks = getParam<std::vector<SubdomainName>>("blocks");
753  mater_name = "ThermalFractureIntegral" + _blocks[0];
754  }
755  else
756  mater_name = "ThermalFractureIntegral";
757 
758  InputParameters params = _factory.getValidParams(mater_type_name);
759  params.set<std::vector<MaterialPropertyName>>("eigenstrain_names") =
760  getParam<std::vector<MaterialPropertyName>>("eigenstrain_names");
761  params.set<std::vector<VariableName>>("temperature") = {_temp};
762 
763  _problem->addMaterial(mater_type_name, mater_name, params);
764  }
765  }
766 }
767 
768 unsigned int
770 {
771  unsigned int num_points = 0;
772  if (_boundary_names.size() != 0)
773  {
774  std::vector<BoundaryID> bids = _mesh->getBoundaryIDs(_boundary_names, true);
775  std::set<unsigned int> nodes;
776 
777  ConstBndNodeRange & bnd_nodes = *_mesh->getBoundaryNodeRange();
778  for (ConstBndNodeRange::const_iterator nd = bnd_nodes.begin(); nd != bnd_nodes.end(); ++nd)
779  {
780  const BndNode * bnode = *nd;
781  BoundaryID boundary_id = bnode->_bnd_id;
782 
783  for (unsigned int ibid = 0; ibid < bids.size(); ++ibid)
784  {
785  if (boundary_id == bids[ibid])
786  {
787  nodes.insert(bnode->_node->id());
788  break;
789  }
790  }
791  }
792  num_points = nodes.size();
793  }
794  else if (_crack_front_points.size() != 0)
795  num_points = _crack_front_points.size();
797  num_points = getParam<unsigned int>("number_points_from_provider");
798  else
799  mooseError("Must define either 'boundary' or 'crack_front_points'");
800  return num_points;
801 }
std::vector< VariableName > _output_variables
const std::vector< BoundaryName > & _boundary_names
UserObjectName _crack_front_points_provider
const std::string _order
RealVectorValue _crack_direction_vector_end_1
void addCrackFrontDefinitionParams(InputParameters &params)
std::vector< BoundaryName > _intersecting_boundary_names
unsigned int calcNumCrackFrontPoints()
InputParameters validParams< DomainIntegralAction >()
RealVectorValue _crack_direction_vector
std::vector< Real > _radius_inner
DomainIntegralAction(const InputParameters &params)
const std::string _family
MooseEnum _direction_method_moose_enum
std::vector< VariableName > _displacements
std::vector< BoundaryName > _crack_mouth_boundary_names
std::vector< SubdomainName > _blocks
MooseEnum _end_direction_method_moose_enum
std::vector< unsigned int > _ring_vec
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
std::vector< Real > _radius_outer
std::vector< Point > _crack_front_points
RealVectorValue _crack_direction_vector_end_2
std::set< INTEGRAL > _integrals