Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "AddGeochemistrySolverAction.h"
11 : #include "GeochemicalModelDefinition.h"
12 : #include "GeochemistryReactorBase.h"
13 : #include "NearestNodeNumber.h"
14 : #include "GeochemistryConsoleOutput.h"
15 : #include "BlockRestrictable.h"
16 : #include "BoundaryRestrictable.h"
17 :
18 : registerMooseAction("GeochemistryApp", AddGeochemistrySolverAction, "add_output");
19 : registerMooseAction("GeochemistryApp", AddGeochemistrySolverAction, "add_user_object");
20 : registerMooseAction("GeochemistryApp",
21 : AddGeochemistrySolverAction,
22 : "add_geochemistry_molality_aux");
23 :
24 : InputParameters
25 675 : AddGeochemistrySolverAction::validParams()
26 : {
27 675 : InputParameters params = Action::validParams();
28 1350 : params.addParam<UserObjectName>(
29 : "geochemistry_reactor_name",
30 : "geochemistry_reactor",
31 : "The name that will be given to the GeochemistryReactor UserObject built by this action");
32 1350 : params.addParam<bool>("include_moose_solve",
33 1350 : false,
34 : "Include a usual MOOSE solve involving Variables and Kernels. In pure "
35 : "reaction systems (without transport) include_moose_solve = false is "
36 : "appropriate, but with transport 'true' must be used");
37 :
38 1350 : params.addRequiredParam<UserObjectName>("model_definition",
39 : "The name of the GeochemicalModelDefinition user object "
40 : "(you must create this UserObject yourself)");
41 675 : params += GeochemistryReactorBase::sharedParams();
42 :
43 675 : params += BlockRestrictable::validParams();
44 675 : params += BoundaryRestrictable::validParams();
45 :
46 2025 : params.addRangeCheckedParam<Real>(
47 : "stoichiometry_tolerance",
48 1350 : 1E-6,
49 : "stoichiometry_tolerance >= 0.0",
50 : "Swapping involves inverting matrices via a singular value decomposition. During this "
51 : "process: (1) if abs(singular value) < stoi_tol * L1norm(singular values), then the "
52 : "matrix is deemed singular (so the basis swap is deemed invalid); (2) if abs(any "
53 : "stoichiometric coefficient) < stoi_tol then it is set to zero.");
54 :
55 : // following are exclusively for the GeochemistryConsoleOutput
56 675 : params += GeochemistryConsoleOutput::sharedParams();
57 675 : ExecFlagEnum exec_enum = MooseUtils::getDefaultExecFlagEnum();
58 2700 : exec_enum = {EXEC_INITIAL, EXEC_FINAL};
59 1350 : params.addParam<ExecFlagEnum>(
60 : "execute_console_output_on", exec_enum, "When to execute the geochemistry console output");
61 675 : params.addParam<Point>("point",
62 675 : Point(0.0, 0.0, 0.0),
63 : "The geochemistry console output will be regarding the aqueous "
64 : "solution at node that is closest to this point");
65 :
66 : // following are the Aux possibilities
67 1350 : params.addParam<bool>(
68 : "add_aux_solvent_kg",
69 1350 : true,
70 : "Add AuxVariable, called kg_solvent_H2O, that records the kg of solvent water");
71 1350 : params.addParam<bool>(
72 1350 : "add_aux_pH", true, "Add AuxVariable, called pH, that records the pH of solvent water");
73 1350 : params.addParam<bool>(
74 : "add_aux_molal",
75 1350 : true,
76 : "Add AuxVariables measured in molal units (ie mol(species)/kg(solvent_water)). These are "
77 : "named molal_name, where 'name' is the species name. AuxVariables are added for all species "
78 : "except minerals");
79 1350 : params.addParam<bool>("add_aux_mg_per_kg",
80 1350 : true,
81 : "Add AuxVariables measured in mg(species)/kg(solvent_water). These are "
82 : "named mg_per_kg_name, where 'name' is the species name. AuxVariables are "
83 : "added for all species except minerals");
84 1350 : params.addParam<bool>("add_aux_free_mg",
85 1350 : true,
86 : "Add AuxVariables for all minerals measured in free mg. These are named "
87 : "free_mg_name, where 'name' is the species name");
88 1350 : params.addParam<bool>("add_aux_free_cm3",
89 1350 : true,
90 : "Add AuxVariables for all minerals measured in free cm^3. These are named "
91 : "free_cm3_name, where 'name' is the species name");
92 1350 : params.addParam<bool>(
93 : "add_aux_activity",
94 1350 : true,
95 : "Add AuxVariables that record the activity for all species (for gas species this equals the "
96 : "gas fugacity). These are called activity_name where 'name' is the species name.");
97 1350 : params.addParam<bool>(
98 : "add_aux_bulk_moles",
99 1350 : true,
100 : "Add AuxVariables that record the number of bulk-composition moles for all species. Note "
101 : "that these will be zero for any species not currently in the basis. These are called "
102 : "bulk_moles_name where 'name' is the species name.");
103 1350 : params.addParam<bool>("add_aux_surface_charge",
104 1350 : true,
105 : "Add AuxVariables, measured in C/m^2, corresponding to specific surface "
106 : "charge for each mineral involved in surface sorption. These are "
107 : "surface_charge_name, where 'name' is the mineral name");
108 1350 : params.addParam<bool>("add_aux_surface_potential",
109 1350 : true,
110 : "Add AuxVariables, measured in V, corresponding to surface potential "
111 : "for each mineral involved in surface sorption. These are "
112 : "surface_potential_name, where 'name' is the mineral name");
113 1350 : params.addParam<bool>("add_aux_temperature",
114 1350 : true,
115 : "Add AuxVariable, called solution_temperature, that records the "
116 : "temperature of the aqueous solution in degC");
117 1350 : params.addParam<bool>(
118 : "add_aux_kinetic_moles",
119 1350 : true,
120 : "Add AuxVariables that record the number of moles for all kinetic species. These are called "
121 : "moles_name where 'name' is the species name.");
122 1350 : params.addParam<bool>("add_aux_kinetic_additions",
123 1350 : true,
124 : "Add AuxVariables that record the rate-of-change (-reaction_rate * dt) for "
125 : "all kinetic species. These are called "
126 : "mol_change_name where 'name' is the species name.");
127 1350 : params.addClassDescription("Base class for an Action that sets up a reaction solver. This class "
128 : "adds a GeochemistryConsoleOutput and AuxVariables corresponding to "
129 : "molalities, etc. Derived classes will create the solver.");
130 :
131 675 : return params;
132 675 : }
133 :
134 675 : AddGeochemistrySolverAction::AddGeochemistrySolverAction(const InputParameters & params)
135 675 : : Action(params)
136 : {
137 675 : }
138 :
139 : void
140 4698 : AddGeochemistrySolverAction::act()
141 : {
142 6048 : if (_current_task == "add_user_object" && isParamValid("execute_console_output_on"))
143 : {
144 459 : const std::string class_name = "NearestNodeNumberUO";
145 459 : auto params = _factory.getValidParams(class_name);
146 918 : params.set<Point>("point") = getParam<Point>("point");
147 918 : if (isParamValid("block"))
148 0 : params.set<std::vector<SubdomainName>>("block") =
149 0 : getParam<std::vector<SubdomainName>>("block");
150 918 : if (isParamValid("boundary"))
151 0 : params.set<std::vector<BoundaryName>>("boundary") =
152 0 : getParam<std::vector<BoundaryName>>("boundary");
153 459 : params.set<ExecFlagEnum>("execute_on") = EXEC_INITIAL; // NOTE: adaptivity not active yet
154 459 : _problem->addUserObject(class_name, "geochemistry_nearest_node_number", params);
155 459 : }
156 5541 : else if (_current_task == "add_output" && isParamValid("execute_console_output_on"))
157 : {
158 435 : const std::string class_name = "GeochemistryConsoleOutput";
159 435 : auto params = _factory.getValidParams(class_name);
160 870 : params.set<UserObjectName>("geochemistry_reactor") =
161 435 : getParam<UserObjectName>("geochemistry_reactor_name");
162 870 : params.set<unsigned>("precision") = getParam<unsigned>("precision");
163 870 : params.set<Real>("mol_cutoff") = getParam<Real>("mol_cutoff");
164 870 : params.set<Real>("stoichiometry_tolerance") = getParam<Real>("stoichiometry_tolerance");
165 870 : params.set<bool>("solver_info") = getParam<bool>("solver_info");
166 870 : params.set<UserObjectName>("nearest_node_number_UO") = "geochemistry_nearest_node_number";
167 1305 : params.set<ExecFlagEnum>("execute_on") = getParam<ExecFlagEnum>("execute_console_output_on");
168 435 : _problem->addOutput(class_name, "geochemistry_console_output", params);
169 435 : }
170 3804 : else if (_current_task == "add_geochemistry_molality_aux")
171 : {
172 651 : const ModelGeochemicalDatabase & mgd = _problem
173 651 : ->getUserObject<GeochemicalModelDefinition>(
174 651 : getParam<UserObjectName>("model_definition"))
175 651 : .getDatabase();
176 : // add temperature aux, if requested
177 1302 : if (getParam<bool>("add_aux_temperature"))
178 1298 : addAuxSpecies("solution_temperature", "H2O", "temperature");
179 : // add water, if requested
180 1302 : if (getParam<bool>("add_aux_solvent_kg"))
181 1298 : addAuxSpecies("kg_solvent_H2O", "H2O", "molal");
182 1302 : if (getParam<bool>("add_aux_activity"))
183 1298 : addAuxSpecies("activity_H2O", "H2O", "activity");
184 1302 : if (getParam<bool>("add_aux_bulk_moles"))
185 1298 : addAuxSpecies("bulk_moles_H2O", "H2O", "bulk_moles");
186 : // add pH, if requested
187 1302 : if (getParam<bool>("add_aux_pH"))
188 1230 : addAuxSpecies("pH", "H+", "neglog10a");
189 :
190 : // add the remaining ones
191 651 : const unsigned num_basis = mgd.basis_species_name.size();
192 4218 : for (unsigned i = 1; i < num_basis; ++i)
193 : {
194 10701 : if (getParam<bool>("add_aux_molal") && !mgd.basis_species_mineral[i])
195 7126 : addAuxSpecies("molal_" + mgd.basis_species_name[i], mgd.basis_species_name[i], "molal");
196 10701 : if (getParam<bool>("add_aux_mg_per_kg") && !mgd.basis_species_mineral[i])
197 10689 : addAuxSpecies(
198 7126 : "mg_per_kg_" + mgd.basis_species_name[i], mgd.basis_species_name[i], "mg_per_kg");
199 10701 : if (getParam<bool>("add_aux_free_cm3") && mgd.basis_species_mineral[i])
200 0 : addAuxSpecies(
201 0 : "free_cm3_" + mgd.basis_species_name[i], mgd.basis_species_name[i], "free_cm3");
202 10701 : if (getParam<bool>("add_aux_free_mg") && mgd.basis_species_mineral[i])
203 0 : addAuxSpecies("free_mg_" + mgd.basis_species_name[i], mgd.basis_species_name[i], "free_mg");
204 7134 : if (getParam<bool>("add_aux_activity"))
205 10689 : addAuxSpecies(
206 7126 : "activity_" + mgd.basis_species_name[i], mgd.basis_species_name[i], "activity");
207 7134 : if (getParam<bool>("add_aux_bulk_moles"))
208 10689 : addAuxSpecies(
209 7126 : "bulk_moles_" + mgd.basis_species_name[i], mgd.basis_species_name[i], "bulk_moles");
210 : }
211 651 : const unsigned num_eqm = mgd.eqm_species_name.size();
212 17128 : for (unsigned j = 0; j < num_eqm; ++j)
213 : {
214 49431 : if (getParam<bool>("add_aux_molal") && !mgd.eqm_species_mineral[j])
215 30550 : addAuxSpecies("molal_" + mgd.eqm_species_name[j], mgd.eqm_species_name[j], "molal");
216 49431 : if (getParam<bool>("add_aux_mg_per_kg") && !mgd.eqm_species_mineral[j])
217 30550 : addAuxSpecies("mg_per_kg_" + mgd.eqm_species_name[j], mgd.eqm_species_name[j], "mg_per_kg");
218 49431 : if (getParam<bool>("add_aux_free_cm3") && mgd.eqm_species_mineral[j])
219 2396 : addAuxSpecies("free_cm3_" + mgd.eqm_species_name[j], mgd.eqm_species_name[j], "free_cm3");
220 49431 : if (getParam<bool>("add_aux_free_mg") && mgd.eqm_species_mineral[j])
221 2396 : addAuxSpecies("free_mg_" + mgd.eqm_species_name[j], mgd.eqm_species_name[j], "free_mg");
222 32954 : if (getParam<bool>("add_aux_activity"))
223 32946 : addAuxSpecies("activity_" + mgd.eqm_species_name[j], mgd.eqm_species_name[j], "activity");
224 32954 : if (getParam<bool>("add_aux_bulk_moles"))
225 49419 : addAuxSpecies(
226 32946 : "bulk_moles_" + mgd.eqm_species_name[j], mgd.eqm_species_name[j], "bulk_moles");
227 : }
228 : // add the kinetic aux variables
229 651 : const unsigned num_kin = mgd.kin_species_name.size();
230 777 : for (unsigned k = 0; k < num_kin; ++k)
231 : {
232 378 : if (getParam<bool>("add_aux_free_cm3") && mgd.kin_species_mineral[k])
233 216 : addAuxSpecies("free_cm3_" + mgd.kin_species_name[k], mgd.kin_species_name[k], "free_cm3");
234 378 : if (getParam<bool>("add_aux_free_mg") && mgd.kin_species_mineral[k])
235 216 : addAuxSpecies("free_mg_" + mgd.kin_species_name[k], mgd.kin_species_name[k], "free_mg");
236 252 : if (getParam<bool>("add_aux_kinetic_moles"))
237 252 : addAuxSpecies("moles_" + mgd.kin_species_name[k], mgd.kin_species_name[k], "kinetic_moles");
238 252 : if (getParam<bool>("add_aux_kinetic_additions"))
239 378 : addAuxSpecies(
240 252 : "mol_change_" + mgd.kin_species_name[k], mgd.kin_species_name[k], "kinetic_additions");
241 : }
242 :
243 : // add surface stuff
244 736 : for (const auto & mineral : mgd.surface_sorption_name)
245 : {
246 170 : if (getParam<bool>("add_aux_surface_charge"))
247 170 : addAuxSpecies("surface_charge_" + mineral, mineral, "surface_charge");
248 170 : if (getParam<bool>("add_aux_surface_potential"))
249 170 : addAuxSpecies("surface_potential_" + mineral, mineral, "surface_potential");
250 : }
251 : }
252 4698 : }
253 :
254 : void
255 83993 : AddGeochemistrySolverAction::addAuxSpecies(const std::string & var_name,
256 : const std::string & species_name,
257 : const std::string & quantity)
258 : {
259 : // add AuxVariable
260 83993 : auto var_params = _factory.getValidParams("MooseVariable");
261 83993 : _problem->addAuxVariable("MooseVariable", var_name, var_params);
262 : // add AuxKernel
263 83993 : const std::string class_name = "GeochemistryQuantityAux";
264 83993 : auto params = _factory.getValidParams(class_name);
265 83993 : params.set<std::string>("species") = species_name;
266 83993 : params.set<MooseEnum>("quantity") = quantity;
267 251979 : params.set<UserObjectName>("reactor") = getParam<UserObjectName>("geochemistry_reactor_name");
268 167986 : params.set<AuxVariableName>("variable") = var_name;
269 83993 : params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_END;
270 83993 : _problem->addAuxKernel(class_name, var_name, params);
271 167986 : }
|
