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 "ReactionNetworkPhysicsBase.h"
      11             : #include <sstream>
      12             : 
      13             : InputParameters
      14          36 : ReactionNetworkPhysicsBase::validParams()
      15             : {
      16          36 :   InputParameters params = PhysicsBase::validParams();
      17          36 :   params.addClassDescription("Base class to create a reaction network Physics from.");
      18             : 
      19             :   // Variables parameters
      20          72 :   params.addRequiredParam<std::vector<VariableName>>("solver_variables", "Species to solve for");
      21          72 :   params.addRequiredParam<std::vector<AuxVariableName>>(
      22             :       "auxiliary_variables",
      23             :       "Additional species to output the concentration of, which do not require an additional "
      24             :       "equation to obtain");
      25          36 :   params.transferParam<MooseEnum>(MooseVariableBase::validParams(), "order", "variable_order");
      26          72 :   params.addParam<std::vector<Real>>("equation_scaling",
      27             :                                      "Scaling factor to apply to each equation");
      28             : 
      29             :   // Initial conditions
      30          72 :   params.addParam<std::vector<MooseFunctorName>>("initial_conditions",
      31             :                                                  "Initial conditions for the species to solve for");
      32             : 
      33             :   // Reaction network
      34          72 :   params.addParam<std::string>("reactions", "The list of chemical reactions");
      35             : 
      36          36 :   return params;
      37           0 : }
      38             : 
      39          36 : ReactionNetworkPhysicsBase::ReactionNetworkPhysicsBase(const InputParameters & parameters)
      40             :   : PhysicsBase(parameters),
      41          36 :     _solver_species(getParam<std::vector<VariableName>>("solver_variables")),
      42          36 :     _num_solver_species(_solver_species.size()),
      43          72 :     _aux_species(getParam<std::vector<AuxVariableName>>("auxiliary_variables")),
      44          36 :     _num_aux_species(_aux_species.size()),
      45          36 :     _reactions(
      46          72 :         ReactionNetworkUtils::parseReactionNetwork(getParam<std::string>("reactions"), _verbose)),
      47          36 :     _num_reactions(_reactions.size())
      48             : {
      49             :   // Keep track of variables
      50          90 :   for (const auto & var_name : _solver_species)
      51             :     saveSolverVariableName(var_name);
      52          99 :   for (const auto & var_name : _aux_species)
      53         126 :     saveAuxVariableName(var_name);
      54             : 
      55             :   // Parameter checking
      56          72 :   if (isParamValid("initial_conditions"))
      57          18 :     checkVectorParamsSameLength<VariableName, MooseFunctorName>("solver_variables",
      58             :                                                                 "initial_conditions");
      59          72 :   if (isParamValid("equation_scaling"))
      60          18 :     checkVectorParamsSameLength<VariableName, Real>("solver_variables", "equation_scaling");
      61             : 
      62          36 :   addRequiredPhysicsTask("add_ic");
      63          36 :   addRequiredPhysicsTask("add_variable");
      64          36 :   addRequiredPhysicsTask("add_aux_variable");
      65          36 :   addRequiredPhysicsTask("init_physics");
      66             : 
      67             :   // Parse the lines in the reaction
      68         108 :   std::stringstream reactions_param(getParam<std::string>("reactions"));
      69             :   std::string line;
      70          99 :   while (std::getline(reactions_param, line, '\n'))
      71          63 :     _reactions_input.push_back(line);
      72             :   mooseAssert(_num_reactions == _reactions_input.size(),
      73             :               "Should be the same size. Extra line break in the reaction network?");
      74          36 : }
      75             : 
      76             : void
      77          36 : ReactionNetworkPhysicsBase::addSolverVariables()
      78             : {
      79          90 :   for (const auto i : index_range(_solver_species))
      80             :   {
      81          54 :     const auto & var_name = _solver_species[i];
      82             :     // If the variable was added outside the Physics
      83          54 :     if (variableExists(var_name, /*error_if_aux*/ true))
      84             :     {
      85          45 :       reportPotentiallyMissedParameters(
      86             :           {"variable_order", "system_names", "equation_scaling"}, "MooseVariable", var_name);
      87          45 :       continue;
      88             :     }
      89             : 
      90           9 :     const std::string variable_type = "MooseVariable";
      91           9 :     InputParameters params = getFactory().getValidParams(variable_type);
      92          27 :     params.set<MooseEnum>("order") = getParam<MooseEnum>("variable_order");
      93          18 :     if (isParamValid("equation_scaling"))
      94             :       // all our variables have a single component
      95          18 :       params.set<std::vector<Real>>("scaling") = {
      96          36 :           getParam<std::vector<Real>>("equation_scaling")[i]};
      97           9 :     assignBlocks(params, _blocks);
      98           9 :     params.set<SolverSystemName>("solver_sys") = getSolverSystem(var_name);
      99             : 
     100           9 :     getProblem().addVariable(variable_type, var_name, params);
     101           9 :   }
     102          36 : }
     103             : 
     104             : void
     105          36 : ReactionNetworkPhysicsBase::addAuxiliaryVariables()
     106             : {
     107          99 :   for (const auto i : index_range(_aux_species))
     108             :   {
     109          63 :     const auto & var_name = _aux_species[i];
     110             :     // If the variable was added outside the Physics
     111         126 :     if (variableExists(var_name, /*error_if_aux*/ false))
     112             :     {
     113           0 :       reportPotentiallyMissedParameters({"variable_order"}, "MooseVariable", var_name);
     114           0 :       continue;
     115             :     }
     116             : 
     117          63 :     const std::string variable_type = "MooseVariable";
     118          63 :     InputParameters params = getFactory().getValidParams(variable_type);
     119         189 :     params.set<MooseEnum>("order") = getParam<MooseEnum>("variable_order");
     120          63 :     assignBlocks(params, _blocks);
     121          63 :     getProblem().addAuxVariable(variable_type, var_name, params);
     122          63 :   }
     123          36 : }
     124             : 
     125             : void
     126           0 : ReactionNetworkPhysicsBase::addPreconditioning()
     127             : {
     128             :   // TODO: identify fully independent or sequential groups in the reaction network
     129             :   // Solve using a segregated approach in the order of the sequence
     130           0 : }
     131             : 
     132             : void
     133           0 : ReactionNetworkPhysicsBase::addComponent(const ActionComponent & component)
     134             : {
     135           0 :   for (const auto & block : component.blocks())
     136           0 :     _blocks.push_back(block);
     137           0 : }
     138             : 
     139             : void
     140          36 : ReactionNetworkPhysicsBase::addInitialConditions()
     141             : {
     142          36 :   InputParameters params = getFactory().getValidParams("FunctorIC");
     143          36 :   assignBlocks(params, _blocks);
     144             :   std::vector<MooseFunctorName> empty_functor_vector;
     145             :   const auto & initial_functors =
     146          36 :       isParamValid("initial_conditions")
     147          54 :           ? getParam<std::vector<MooseFunctorName>>("initial_conditions")
     148             :           : empty_functor_vector;
     149             : 
     150          90 :   for (const auto i : index_range(_solver_species))
     151             :   {
     152          54 :     const auto & var_name = _solver_species[i];
     153             :     // always obey the user specification of initial conditions
     154             :     // Base class does not have a default but derived classes could set one
     155         162 :     if (isParamValid("initial_conditions") &&
     156           9 :         shouldCreateIC(var_name,
     157             :                        _blocks,
     158          72 :                        /*whether IC is a default*/ !isParamSetByUser("initial_conditions"),
     159             :                        /*error if already an IC*/ true))
     160             :     {
     161          18 :       params.set<VariableName>("variable") = var_name;
     162           9 :       params.set<MooseFunctorName>("functor") = initial_functors[i];
     163          36 :       getProblem().addInitialCondition("FunctorIC", prefix() + var_name + "_ic", params);
     164             :     }
     165             :   }
     166          36 : }