www.mooseframework.org
CHPFCRFFSplitVariablesAction.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 /****************************************************************/
8 #include "Factory.h"
9 #include "FEProblem.h"
10 #include "Conversion.h"
11 #include "AddVariableAction.h"
12 
13 #include "libmesh/string_to_enum.h"
14 
15 template <>
16 InputParameters
18 {
19  InputParameters params = validParams<Action>();
20  MooseEnum familyEnum = AddVariableAction::getNonlinearVariableFamilies();
21  params.addParam<MooseEnum>(
22  "family",
23  familyEnum,
24  "Specifies the family of FE shape functions to use for the L variables");
25  MooseEnum orderEnum = AddVariableAction::getNonlinearVariableOrders();
26  params.addParam<MooseEnum>(
27  "order",
28  orderEnum,
29  "Specifies the order of the FE shape function to use for the L variables");
30  params.addParam<Real>("scaling", 1.0, "Specifies a scaling factor to apply to the L variables");
31  params.addRequiredParam<unsigned int>(
32  "num_L", "specifies the number of complex L variables will be solved for");
33  params.addRequiredParam<std::string>("L_name_base", "Base name for the complex L variables");
34  params.addRequiredParam<std::vector<FileName>>("sub_filenames",
35  "This is the filename of the sub.i file");
36  params.addRequiredParam<AuxVariableName>("n_name", "Name of atomic density variable");
37 
38  return params;
39 }
40 
42  : Action(params),
43  _num_L(getParam<unsigned int>("num_L")),
44  _L_name_base(getParam<std::string>("L_name_base")),
45  _sub_filenames(getParam<std::vector<FileName>>("sub_filenames")),
46  _n_name(getParam<AuxVariableName>("n_name"))
47 {
48 }
49 
50 void
52 {
53  MultiMooseEnum execute_options(SetupInterface::getExecuteOptions());
54  execute_options = "timestep_begin";
55 
56  // Setup MultiApp
57  InputParameters poly_params = _factory.getValidParams("TransientMultiApp");
58  poly_params.set<MooseEnum>("app_type") = "PhaseFieldApp";
59  poly_params.set<MultiMooseEnum>("execute_on") = execute_options;
60  poly_params.set<std::vector<FileName>>("input_files") = _sub_filenames;
61  poly_params.set<unsigned int>("max_procs_per_app") = 1;
62  poly_params.set<std::vector<Point>>("positions") = {Point()};
63  _problem->addMultiApp("TransientMultiApp", "HHEquationSolver", poly_params);
64 
65  poly_params = _factory.getValidParams("MultiAppNearestNodeTransfer");
66  poly_params.set<MooseEnum>("direction") = "to_multiapp";
67  poly_params.set<MultiMooseEnum>("execute_on") = execute_options;
68  poly_params.set<AuxVariableName>("variable") = _n_name;
69  poly_params.set<VariableName>("source_variable") = _n_name;
70  poly_params.set<MultiAppName>("multi_app") = "HHEquationSolver";
71  _problem->addTransfer("MultiAppNearestNodeTransfer", _n_name + "_trans", poly_params);
72 
73  // Loop through the number of L variables
74  for (unsigned int l = 0; l < _num_L; ++l)
75  {
76  // Create L base name
77  std::string L_name = _L_name_base + Moose::stringify(l);
78 
79  // Create real L variable
80  std::string real_name = L_name + "_real";
81 
82  _problem->addAuxVariable(
83  real_name,
84  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("order")),
85  Utility::string_to_enum<FEFamily>(getParam<MooseEnum>("family"))));
86 
87  poly_params = _factory.getValidParams("MultiAppNearestNodeTransfer");
88  poly_params.set<MooseEnum>("direction") = "from_multiapp";
89  poly_params.set<AuxVariableName>("variable") = real_name;
90  poly_params.set<VariableName>("source_variable") = real_name;
91  poly_params.set<MultiAppName>("multi_app") = "HHEquationSolver";
92  _problem->addTransfer("MultiAppNearestNodeTransfer", real_name + "_trans", poly_params);
93 
94  if (l > 0)
95  {
96  // Create imaginary L variable IF l > 0
97  std::string imag_name = L_name + "_imag";
98 
99  _problem->addAuxVariable(
100  imag_name,
101  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("order")),
102  Utility::string_to_enum<FEFamily>(getParam<MooseEnum>("family"))));
103 
104  poly_params = _factory.getValidParams("MultiAppNearestNodeTransfer");
105  poly_params.set<MooseEnum>("direction") = "from_multiapp";
106  poly_params.set<AuxVariableName>("variable") = imag_name;
107  poly_params.set<VariableName>("source_variable") = imag_name;
108  poly_params.set<MultiAppName>("multi_app") = "HHEquationSolver";
109  _problem->addTransfer("MultiAppNearestNodeTransfer", imag_name + "_trans", poly_params);
110  }
111  }
112 }
const std::vector< FileName > _sub_filenames
CHPFCRFFSplitVariablesAction(const InputParameters &params)
InputParameters validParams< CHPFCRFFSplitVariablesAction >()