www.mooseframework.org
ConservedAction.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 #include "ConservedAction.h"
9 // MOOSE includes
10 #include "Conversion.h"
11 #include "FEProblem.h"
12 #include "Factory.h"
13 #include "MooseObjectAction.h"
14 #include "MooseMesh.h"
15 #include "AddVariableAction.h"
16 
17 #include "libmesh/string_to_enum.h"
18 
19 template <>
20 InputParameters
22 {
23  InputParameters params = validParams<Action>();
24  params.addClassDescription(
25  "Set up the variable(s) and the kernels needed for a conserved phase field variable."
26  " Note that for a direct solve, the element family and order are overwritten with hermite "
27  "and third.");
28  MooseEnum solves("DIRECT REVERSE_SPLIT FORWARD_SPLIT");
29  params.addRequiredParam<MooseEnum>("solve_type", solves, "Split or direct solve?");
30  // Get MooseEnums for the possible order/family options for this variable
31  MooseEnum families(AddVariableAction::getNonlinearVariableFamilies());
32  MooseEnum orders(AddVariableAction::getNonlinearVariableOrders());
33  params.addParam<MooseEnum>("family",
34  families,
35  "Specifies the family of FE "
36  "shape functions to use for this variable");
37  params.addParam<MooseEnum>("order",
38  orders,
39  "Specifies the order of the FE "
40  "shape function to use for this variable");
41  params.addParam<Real>("scaling", 1.0, "Specifies a scaling factor to apply to this variable");
42  params.addParam<bool>("implicit", true, "Whether kernels are implicit or not");
43  params.addParam<bool>(
44  "use_displaced_mesh", false, "Whether to use displaced mesh in the kernels");
45  params.addParamNamesToGroup("scaling implicit use_displaced_mesh", "Advanced");
46  params.addRequiredParam<MaterialPropertyName>("mobility", "The mobility used with the kernel");
47  params.addParam<std::vector<VariableName>>("args",
48  "Vector of variable arguments this kernel depends on");
49  params.addRequiredParam<MaterialPropertyName>(
50  "free_energy", "Base name of the free energy function F defined in a free energy material");
51  params.addRequiredParam<MaterialPropertyName>("kappa", "The kappa used with the kernel");
52 
53  return params;
54 }
55 
56 ConservedAction::ConservedAction(const InputParameters & params)
57  : Action(params),
58  _solve_type(getParam<MooseEnum>("solve_type").getEnum<SolveType>()),
59  _var_name(name()),
60  _scaling(getParam<Real>("scaling"))
61 {
62  switch (_solve_type)
63  {
64  case SolveType::DIRECT:
65  _fe_type = FEType(Utility::string_to_enum<Order>("THIRD"),
66  Utility::string_to_enum<FEFamily>("HERMITE"));
67  if (!parameters().isParamSetByAddParam("order") &&
68  !parameters().isParamSetByAddParam("family"))
69  mooseWarning("Order and family autoset to third and hermite in ConservedAction");
70  break;
73  _fe_type = FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("order")),
74  Utility::string_to_enum<FEFamily>(getParam<MooseEnum>("family")));
75  // Set name of chemical potential variable
76  _chempot_name = "chem_pot_" + _var_name;
77  break;
78  default:
79  mooseError("Incorrect solve_type in ConservedAction");
80  }
81 }
82 
83 void
85 {
86  //
87  // Add variable(s)
88  //
89  if (_current_task == "add_variable")
90  {
91  // Create conserved variable _var_name
92  _problem->addVariable(_var_name, _fe_type, _scaling);
93 
94  // Create chemical potential variable for split form
95  switch (_solve_type)
96  {
97  case SolveType::DIRECT:
98  break;
101  _problem->addVariable(_chempot_name, _fe_type, _scaling);
102  }
103  }
104 
105  //
106  // Add Kernels
107  //
108  else if (_current_task == "add_kernel")
109  {
110  switch (_solve_type)
111  {
112  case SolveType::DIRECT:
113  // Add time derivative kernel
114  {
115  std::string kernel_type = "TimeDerivative";
116 
117  std::string kernel_name = _var_name + "_" + kernel_type;
118  InputParameters params = _factory.getValidParams(kernel_type);
119  params.set<NonlinearVariableName>("variable") = _var_name;
120  params.applyParameters(parameters());
121 
122  _problem->addKernel(kernel_type, kernel_name, params);
123  }
124 
125  // Add CahnHilliard kernel
126  {
127  std::string kernel_type = "CahnHilliard";
128 
129  std::string kernel_name = _var_name + "_" + kernel_type;
130  InputParameters params = _factory.getValidParams(kernel_type);
131  params.set<NonlinearVariableName>("variable") = _var_name;
132  params.set<MaterialPropertyName>("mob_name") = getParam<MaterialPropertyName>("mobility");
133  params.set<MaterialPropertyName>("f_name") =
134  getParam<MaterialPropertyName>("free_energy");
135  params.applyParameters(parameters());
136 
137  _problem->addKernel(kernel_type, kernel_name, params);
138  }
139 
140  // Add ACInterface kernel
141  {
142  std::string kernel_type = "CHInterface";
143 
144  std::string kernel_name = _var_name + "_" + kernel_type;
145  InputParameters params = _factory.getValidParams(kernel_type);
146  params.set<NonlinearVariableName>("variable") = _var_name;
147  params.set<MaterialPropertyName>("mob_name") = getParam<MaterialPropertyName>("mobility");
148  params.set<MaterialPropertyName>("kappa_name") = getParam<MaterialPropertyName>("kappa");
149  params.applyParameters(parameters());
150 
151  _problem->addKernel(kernel_type, kernel_name, params);
152  }
153  break;
154 
156  // Add time derivative kernel
157  {
158  std::string kernel_type = "CoupledTimeDerivative";
159 
160  std::string kernel_name = _var_name + "_" + kernel_type;
161  InputParameters params = _factory.getValidParams(kernel_type);
162  params.set<NonlinearVariableName>("variable") = _chempot_name;
163  params.set<std::vector<VariableName>>("v") = {_var_name};
164  params.applyParameters(parameters());
165 
166  _problem->addKernel(kernel_type, kernel_name, params);
167  }
168 
169  // Add SplitCHWRes kernel
170  {
171  std::string kernel_type = "SplitCHWRes";
172 
173  std::string kernel_name = _var_name + "_" + kernel_type;
174  InputParameters params = _factory.getValidParams(kernel_type);
175  params.set<NonlinearVariableName>("variable") = _chempot_name;
176  params.set<MaterialPropertyName>("mob_name") = getParam<MaterialPropertyName>("mobility");
177  params.applyParameters(parameters());
178 
179  _problem->addKernel(kernel_type, kernel_name, params);
180  }
181 
182  // Add SplitCHParsed kernel
183  {
184  std::string kernel_type = "SplitCHParsed";
185 
186  std::string kernel_name = _var_name + "_" + kernel_type;
187  InputParameters params = _factory.getValidParams(kernel_type);
188  params.set<NonlinearVariableName>("variable") = _var_name;
189  params.set<std::vector<VariableName>>("w") = {_chempot_name};
190  params.set<MaterialPropertyName>("f_name") =
191  getParam<MaterialPropertyName>("free_energy");
192  params.set<MaterialPropertyName>("kappa_name") = getParam<MaterialPropertyName>("kappa");
193  params.applyParameters(parameters());
194 
195  _problem->addKernel(kernel_type, kernel_name, params);
196  }
197  break;
198 
200  // Add time derivative kernel
201  {
202  std::string kernel_type = "TimeDerivative";
203 
204  std::string kernel_name = _var_name + "_" + kernel_type;
205  InputParameters params = _factory.getValidParams(kernel_type);
206  params.set<NonlinearVariableName>("variable") = _var_name;
207  params.applyParameters(parameters());
208 
209  _problem->addKernel(kernel_type, kernel_name, params);
210  }
211 
212  // Add MatDiffusion kernel for c residual
213  {
214  std::string kernel_type = "MatDiffusion";
215 
216  std::string kernel_name = _var_name + "_" + kernel_type;
217  InputParameters params = _factory.getValidParams(kernel_type);
218  params.set<NonlinearVariableName>("variable") = _var_name;
219  params.set<std::vector<VariableName>>("conc") = {_chempot_name};
220  params.set<MaterialPropertyName>("D_name") = getParam<MaterialPropertyName>("mobility");
221  params.applyParameters(parameters());
222 
223  _problem->addKernel(kernel_type, kernel_name, params);
224  }
225  // Add MatDiffusion kernel for chemical potential residual
226  {
227  std::string kernel_type = "MatDiffusion";
228 
229  std::string kernel_name = _chempot_name + "_" + kernel_type;
230  InputParameters params = _factory.getValidParams(kernel_type);
231  params.set<NonlinearVariableName>("variable") = _chempot_name;
232  params.set<std::vector<VariableName>>("conc") = {_var_name};
233  params.set<MaterialPropertyName>("D_name") = getParam<MaterialPropertyName>("kappa");
234  params.applyParameters(parameters());
235 
236  _problem->addKernel(kernel_type, kernel_name, params);
237  }
238 
239  // Add CoupledMaterialDerivative kernel
240  {
241  std::string kernel_type = "CoupledMaterialDerivative";
242 
243  std::string kernel_name = _chempot_name + "_" + kernel_type;
244  InputParameters params = _factory.getValidParams(kernel_type);
245  params.set<NonlinearVariableName>("variable") = _chempot_name;
246  params.set<std::vector<VariableName>>("v") = {_var_name};
247  params.set<MaterialPropertyName>("f_name") =
248  getParam<MaterialPropertyName>("free_energy");
249  params.applyParameters(parameters());
250 
251  _problem->addKernel(kernel_type, kernel_name, params);
252  }
253 
254  // Add CoefReaction kernel
255  {
256  std::string kernel_type = "CoefReaction";
257 
258  std::string kernel_name = _chempot_name + "_" + kernel_type;
259  InputParameters params = _factory.getValidParams(kernel_type);
260  params.set<NonlinearVariableName>("variable") = _chempot_name;
261  params.set<Real>("coefficient") = -1.0;
262  params.applyParameters(parameters());
263 
264  _problem->addKernel(kernel_type, kernel_name, params);
265  }
266  }
267  }
268 }
FEType _fe_type
FEType for the variable being created.
virtual void act() override
const Real _scaling
Scaling parameter.
std::string _chempot_name
Name of chemical potential variable for split solves.
SolveType
Type of solve.
const NonlinearVariableName _var_name
Name of the variable being created.
const SolveType _solve_type
Type of solve to use used in the action.
InputParameters validParams< ConservedAction >()
ConservedAction(const InputParameters &params)