www.mooseframework.org
Public Member Functions | Private Attributes | List of all members
Q2PAction Class Reference

#include <Q2PAction.h>

Inheritance diagram for Q2PAction:
[legend]

Public Member Functions

 Q2PAction (const InputParameters &params)
 
virtual void act ()
 

Private Attributes

VariableName _pp_var
 
VariableName _sat_var
 
UserObjectName _water_density
 
UserObjectName _water_relperm
 
UserObjectName _water_relperm_for_diffusivity
 
Real _water_viscosity
 
UserObjectName _gas_density
 
UserObjectName _gas_relperm
 
Real _gas_viscosity
 
Real _diffusivity
 
std::vector< OutputName > _output_nodal_masses_to
 
std::vector< OutputName > _output_total_masses_to
 
bool _save_gas_flux_in_Q2PGasFluxResidual
 
bool _save_water_flux_in_Q2PWaterFluxResidual
 
bool _save_gas_Jacobian_in_Q2PGasJacobian
 
bool _save_water_Jacobian_in_Q2PWaterJacobian
 
bool _nodal_masses_not_outputted
 
bool _total_masses_not_outputted
 
bool _no_mass_calculations
 

Detailed Description

Definition at line 20 of file Q2PAction.h.

Constructor & Destructor Documentation

◆ Q2PAction()

Q2PAction::Q2PAction ( const InputParameters &  params)

Definition at line 89 of file Q2PAction.C.

90  : Action(params),
91  _pp_var(getParam<NonlinearVariableName>("porepressure")),
92  _sat_var(getParam<NonlinearVariableName>("saturation")),
93  _water_density(getParam<UserObjectName>("water_density")),
94  _water_relperm(getParam<UserObjectName>("water_relperm")),
95  _water_relperm_for_diffusivity(isParamValid("water_relperm_for_diffusivity")
96  ? getParam<UserObjectName>("water_relperm_for_diffusivity")
97  : getParam<UserObjectName>("water_relperm")),
98  _water_viscosity(getParam<Real>("water_viscosity")),
99  _gas_density(getParam<UserObjectName>("gas_density")),
100  _gas_relperm(getParam<UserObjectName>("gas_relperm")),
101  _gas_viscosity(getParam<Real>("gas_viscosity")),
102  _diffusivity(getParam<Real>("diffusivity")),
103  _output_nodal_masses_to(getParam<std::vector<OutputName>>("output_nodal_masses_to")),
104  _output_total_masses_to(getParam<std::vector<OutputName>>("output_total_masses_to")),
105  _save_gas_flux_in_Q2PGasFluxResidual(getParam<bool>("save_gas_flux_in_Q2PGasFluxResidual")),
107  getParam<bool>("save_water_flux_in_Q2PWaterFluxResidual")),
108  _save_gas_Jacobian_in_Q2PGasJacobian(getParam<bool>("save_gas_Jacobian_in_Q2PGasJacobian")),
110  getParam<bool>("save_water_Jacobian_in_Q2PWaterJacobian"))
111 {
113  if (_output_nodal_masses_to.size() == 0)
115 
117  if (_output_total_masses_to.size() == 0)
119 
121 }
Real _diffusivity
Definition: Q2PAction.h:37
Real _gas_viscosity
Definition: Q2PAction.h:36
UserObjectName _water_density
Definition: Q2PAction.h:30
bool _total_masses_not_outputted
Definition: Q2PAction.h:45
std::vector< OutputName > _output_total_masses_to
Definition: Q2PAction.h:39
UserObjectName _water_relperm
Definition: Q2PAction.h:31
std::vector< OutputName > _output_nodal_masses_to
Definition: Q2PAction.h:38
bool _save_gas_Jacobian_in_Q2PGasJacobian
Definition: Q2PAction.h:42
UserObjectName _water_relperm_for_diffusivity
Definition: Q2PAction.h:32
bool _save_water_Jacobian_in_Q2PWaterJacobian
Definition: Q2PAction.h:43
bool _save_water_flux_in_Q2PWaterFluxResidual
Definition: Q2PAction.h:41
Real _water_viscosity
Definition: Q2PAction.h:33
bool _save_gas_flux_in_Q2PGasFluxResidual
Definition: Q2PAction.h:40
UserObjectName _gas_density
Definition: Q2PAction.h:34
VariableName _pp_var
Definition: Q2PAction.h:28
VariableName _sat_var
Definition: Q2PAction.h:29
bool _nodal_masses_not_outputted
Definition: Q2PAction.h:44
bool _no_mass_calculations
Definition: Q2PAction.h:46
UserObjectName _gas_relperm
Definition: Q2PAction.h:35

Member Function Documentation

◆ act()

void Q2PAction::act ( )
virtual

Definition at line 124 of file Q2PAction.C.

125 {
126  // add the kernels
127  if (_current_task == "add_kernel")
128  {
129  std::string kernel_name;
130  std::string kernel_type;
131  InputParameters params = _factory.getValidParams("Q2PNodalMass");
132 
133  kernel_name = "Q2P_nodal_water_mass";
134  kernel_type = "Q2PNodalMass";
135  params = _factory.getValidParams(kernel_type);
136  params.set<NonlinearVariableName>("variable") = _sat_var;
137  params.set<std::vector<VariableName>>("other_var") = {_pp_var};
138  params.set<bool>("var_is_porepressure") = false;
140  params.set<std::vector<AuxVariableName>>("save_in") = {"Q2P_nodal_water_mass_divided_by_dt"};
141  params.set<UserObjectName>("fluid_density") = _water_density;
142  _problem->addKernel(kernel_type, kernel_name, params);
143 
144  kernel_name = "Q2P_nodal_gas_mass";
145  kernel_type = "Q2PNodalMass";
146  params = _factory.getValidParams(kernel_type);
147  params.set<NonlinearVariableName>("variable") = _pp_var;
148  params.set<std::vector<VariableName>>("other_var") = {_sat_var};
149  params.set<bool>("var_is_porepressure") = true;
151  params.set<std::vector<AuxVariableName>>("save_in") = {"Q2P_nodal_gas_mass_divided_by_dt"};
152  params.set<UserObjectName>("fluid_density") = _gas_density;
153  _problem->addKernel(kernel_type, kernel_name, params);
154 
155  kernel_name = "Q2P_nodal_water_mass_old";
156  kernel_type = "Q2PNegativeNodalMassOld";
157  params = _factory.getValidParams(kernel_type);
158  params.set<NonlinearVariableName>("variable") = _sat_var;
159  params.set<std::vector<VariableName>>("other_var") = {_pp_var};
160  params.set<bool>("var_is_porepressure") = false;
161  params.set<UserObjectName>("fluid_density") = _water_density;
162  _problem->addKernel(kernel_type, kernel_name, params);
163 
164  kernel_name = "Q2P_nodal_gas_mass_old";
165  kernel_type = "Q2PNegativeNodalMassOld";
166  params = _factory.getValidParams(kernel_type);
167  params.set<NonlinearVariableName>("variable") = _pp_var;
168  params.set<std::vector<VariableName>>("other_var") = {_sat_var};
169  params.set<bool>("var_is_porepressure") = true;
170  params.set<UserObjectName>("fluid_density") = _gas_density;
171  _problem->addKernel(kernel_type, kernel_name, params);
172 
173  kernel_name = "Q2P_water_flux";
174  kernel_type = "Q2PSaturationFlux";
175  params = _factory.getValidParams(kernel_type);
176  params.set<NonlinearVariableName>("variable") = _sat_var;
177  params.set<std::vector<VariableName>>("porepressure_variable") = {_pp_var};
178  params.set<UserObjectName>("fluid_density") = _water_density;
179  params.set<UserObjectName>("fluid_relperm") = _water_relperm;
180  params.set<Real>("fluid_viscosity") = _water_viscosity;
182  params.set<std::vector<AuxVariableName>>("save_in") = {"Q2PWaterFluxResidual"};
184  params.set<std::vector<AuxVariableName>>("diag_save_in") = {"Q2PWaterJacobian"};
185  _problem->addKernel(kernel_type, kernel_name, params);
186 
187  kernel_name = "Q2P_gas_flux";
188  kernel_type = "Q2PPorepressureFlux";
189  params = _factory.getValidParams(kernel_type);
190  params.set<NonlinearVariableName>("variable") = _pp_var;
191  params.set<std::vector<VariableName>>("saturation_variable") = {_sat_var};
192  params.set<UserObjectName>("fluid_density") = _gas_density;
193  params.set<UserObjectName>("fluid_relperm") = _gas_relperm;
194  params.set<Real>("fluid_viscosity") = _gas_viscosity;
196  params.set<std::vector<AuxVariableName>>("save_in") = {"Q2PGasFluxResidual"};
198  params.set<std::vector<AuxVariableName>>("diag_save_in") = {"Q2PGasJacobian"};
199  _problem->addKernel(kernel_type, kernel_name, params);
200 
201  kernel_name = "Q2P_liquid_diffusion";
202  kernel_type = "Q2PSaturationDiffusion";
203  params = _factory.getValidParams(kernel_type);
204  params.set<NonlinearVariableName>("variable") = _sat_var;
205  params.set<std::vector<VariableName>>("porepressure_variable") = {_pp_var};
206  params.set<UserObjectName>("fluid_density") = _water_density;
207  params.set<UserObjectName>("fluid_relperm") = _water_relperm_for_diffusivity;
208  params.set<Real>("fluid_viscosity") = _water_viscosity;
209  params.set<Real>("diffusivity") = _diffusivity;
210  _problem->addKernel(kernel_type, kernel_name, params);
211  }
212 
213  if (_current_task == "add_aux_variable")
214  {
216  {
217  // user wants nodal masses or total masses
218  _problem->addAuxVariable("Q2P_nodal_water_mass_divided_by_dt",
219  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("ORDER")),
220  Utility::string_to_enum<FEFamily>("LAGRANGE")));
221  _problem->addAuxVariable("Q2P_nodal_gas_mass_divided_by_dt",
222  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("ORDER")),
223  Utility::string_to_enum<FEFamily>("LAGRANGE")));
224  }
226  _problem->addAuxVariable("Q2PGasFluxResidual",
227  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("ORDER")),
228  Utility::string_to_enum<FEFamily>("LAGRANGE")));
230  _problem->addAuxVariable("Q2PWaterFluxResidual",
231  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("ORDER")),
232  Utility::string_to_enum<FEFamily>("LAGRANGE")));
234  _problem->addAuxVariable("Q2PGasJacobian",
235  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("ORDER")),
236  Utility::string_to_enum<FEFamily>("LAGRANGE")));
238  _problem->addAuxVariable("Q2PWaterJacobian",
239  FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("ORDER")),
240  Utility::string_to_enum<FEFamily>("LAGRANGE")));
241  }
242 
243  if (_current_task == "add_function" && _output_total_masses_to.size() > 0)
244  {
245  // user wants total masses, so need to build Functions to do this
246  InputParameters params = _factory.getValidParams("ParsedFunction");
247 
248  params.set<std::string>("value") = "a*b";
249 
250  std::vector<std::string> vars;
251  vars.push_back("a");
252  vars.push_back("b");
253  params.set<std::vector<std::string>>("vars") = vars;
254 
255  std::vector<std::string> vals_water;
256  vals_water.push_back("Q2P_mass_water_divided_by_dt");
257  vals_water.push_back("Q2P_dt");
258  params.set<std::vector<std::string>>("vals") = vals_water;
259  _problem->addFunction("ParsedFunction", "Q2P_water_mass_fcn", params);
260 
261  std::vector<std::string> vals_gas;
262  vals_gas.push_back("Q2P_mass_gas_divided_by_dt");
263  vals_gas.push_back("Q2P_dt");
264  params.set<std::vector<std::string>>("vals") = vals_gas;
265  _problem->addFunction("ParsedFunction", "Q2P_gas_mass_fcn", params);
266  }
267 
268  if (_current_task == "add_postprocessor" && _output_total_masses_to.size() > 0)
269  {
270  // user wants total masses, so need to build Postprocessors to do this
271 
272  InputParameters params = _factory.getValidParams("TimestepSize");
273 
274  params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_BEGIN;
275  params.set<std::vector<OutputName>>("outputs") = {"none"};
276  _problem->addPostprocessor("TimestepSize", "Q2P_dt", params);
277 
278  params = _factory.getValidParams("NodalSum");
279  params.set<std::vector<OutputName>>("outputs") = {"none"};
280  params.set<std::vector<VariableName>>("variable") = {"Q2P_nodal_water_mass_divided_by_dt"};
281  _problem->addPostprocessor("NodalSum", "Q2P_mass_water_divided_by_dt", params);
282 
283  params = _factory.getValidParams("FunctionValuePostprocessor");
284  params.set<FunctionName>("function") = "Q2P_water_mass_fcn";
285  params.set<std::vector<OutputName>>("outputs") = _output_total_masses_to;
286  _problem->addPostprocessor("FunctionValuePostprocessor", "mass_water", params);
287 
288  params = _factory.getValidParams("NodalSum");
289  params.set<std::vector<OutputName>>("outputs") = {"none"};
290  params.set<std::vector<VariableName>>("variable") = {"Q2P_nodal_gas_mass_divided_by_dt"};
291  _problem->addPostprocessor("NodalSum", "Q2P_mass_gas_divided_by_dt", params);
292 
293  params = _factory.getValidParams("FunctionValuePostprocessor");
294  params.set<FunctionName>("function") = "Q2P_gas_mass_fcn";
295  params.set<std::vector<OutputName>>("outputs") = _output_total_masses_to;
296  _problem->addPostprocessor("FunctionValuePostprocessor", "mass_gas", params);
297  }
298 }
Real _diffusivity
Definition: Q2PAction.h:37
Real _gas_viscosity
Definition: Q2PAction.h:36
UserObjectName _water_density
Definition: Q2PAction.h:30
std::vector< OutputName > _output_total_masses_to
Definition: Q2PAction.h:39
UserObjectName _water_relperm
Definition: Q2PAction.h:31
bool _save_gas_Jacobian_in_Q2PGasJacobian
Definition: Q2PAction.h:42
UserObjectName _water_relperm_for_diffusivity
Definition: Q2PAction.h:32
bool _save_water_Jacobian_in_Q2PWaterJacobian
Definition: Q2PAction.h:43
bool _save_water_flux_in_Q2PWaterFluxResidual
Definition: Q2PAction.h:41
Real _water_viscosity
Definition: Q2PAction.h:33
bool _save_gas_flux_in_Q2PGasFluxResidual
Definition: Q2PAction.h:40
UserObjectName _gas_density
Definition: Q2PAction.h:34
VariableName _pp_var
Definition: Q2PAction.h:28
VariableName _sat_var
Definition: Q2PAction.h:29
bool _no_mass_calculations
Definition: Q2PAction.h:46
UserObjectName _gas_relperm
Definition: Q2PAction.h:35

Member Data Documentation

◆ _diffusivity

Real Q2PAction::_diffusivity
private

Definition at line 37 of file Q2PAction.h.

Referenced by act().

◆ _gas_density

UserObjectName Q2PAction::_gas_density
private

Definition at line 34 of file Q2PAction.h.

Referenced by act().

◆ _gas_relperm

UserObjectName Q2PAction::_gas_relperm
private

Definition at line 35 of file Q2PAction.h.

Referenced by act().

◆ _gas_viscosity

Real Q2PAction::_gas_viscosity
private

Definition at line 36 of file Q2PAction.h.

Referenced by act().

◆ _no_mass_calculations

bool Q2PAction::_no_mass_calculations
private

Definition at line 46 of file Q2PAction.h.

Referenced by act(), and Q2PAction().

◆ _nodal_masses_not_outputted

bool Q2PAction::_nodal_masses_not_outputted
private

Definition at line 44 of file Q2PAction.h.

Referenced by Q2PAction().

◆ _output_nodal_masses_to

std::vector<OutputName> Q2PAction::_output_nodal_masses_to
private

Definition at line 38 of file Q2PAction.h.

Referenced by Q2PAction().

◆ _output_total_masses_to

std::vector<OutputName> Q2PAction::_output_total_masses_to
private

Definition at line 39 of file Q2PAction.h.

Referenced by act(), and Q2PAction().

◆ _pp_var

VariableName Q2PAction::_pp_var
private

Definition at line 28 of file Q2PAction.h.

Referenced by act().

◆ _sat_var

VariableName Q2PAction::_sat_var
private

Definition at line 29 of file Q2PAction.h.

Referenced by act().

◆ _save_gas_flux_in_Q2PGasFluxResidual

bool Q2PAction::_save_gas_flux_in_Q2PGasFluxResidual
private

Definition at line 40 of file Q2PAction.h.

Referenced by act().

◆ _save_gas_Jacobian_in_Q2PGasJacobian

bool Q2PAction::_save_gas_Jacobian_in_Q2PGasJacobian
private

Definition at line 42 of file Q2PAction.h.

Referenced by act().

◆ _save_water_flux_in_Q2PWaterFluxResidual

bool Q2PAction::_save_water_flux_in_Q2PWaterFluxResidual
private

Definition at line 41 of file Q2PAction.h.

Referenced by act().

◆ _save_water_Jacobian_in_Q2PWaterJacobian

bool Q2PAction::_save_water_Jacobian_in_Q2PWaterJacobian
private

Definition at line 43 of file Q2PAction.h.

Referenced by act().

◆ _total_masses_not_outputted

bool Q2PAction::_total_masses_not_outputted
private

Definition at line 45 of file Q2PAction.h.

Referenced by Q2PAction().

◆ _water_density

UserObjectName Q2PAction::_water_density
private

Definition at line 30 of file Q2PAction.h.

Referenced by act().

◆ _water_relperm

UserObjectName Q2PAction::_water_relperm
private

Definition at line 31 of file Q2PAction.h.

Referenced by act().

◆ _water_relperm_for_diffusivity

UserObjectName Q2PAction::_water_relperm_for_diffusivity
private

Definition at line 32 of file Q2PAction.h.

Referenced by act().

◆ _water_viscosity

Real Q2PAction::_water_viscosity
private

Definition at line 33 of file Q2PAction.h.

Referenced by act().


The documentation for this class was generated from the following files: