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 "AddCoupledSolidKinSpeciesAction.h"
      11             : #include "MooseUtils.h"
      12             : #include "FEProblem.h"
      13             : #include "Factory.h"
      14             : #include "MooseError.h"
      15             : #include "Parser.h"
      16             : #include <algorithm>
      17             : 
      18             : // Regular expression includes
      19             : #include "pcrecpp.h"
      20             : 
      21             : registerMooseAction("ChemicalReactionsApp", AddCoupledSolidKinSpeciesAction, "add_kernel");
      22             : 
      23             : registerMooseAction("ChemicalReactionsApp", AddCoupledSolidKinSpeciesAction, "add_aux_kernel");
      24             : 
      25             : InputParameters
      26          83 : AddCoupledSolidKinSpeciesAction::validParams()
      27             : {
      28          83 :   InputParameters params = Action::validParams();
      29         166 :   params.addRequiredParam<std::vector<NonlinearVariableName>>("primary_species",
      30             :                                                               "The list of primary species to add");
      31         166 :   params.addParam<std::vector<AuxVariableName>>(
      32             :       "secondary_species", "The list of solid kinetic species to be output as aux variables");
      33         166 :   params.addRequiredParam<std::string>("kin_reactions", "The list of solid kinetic reactions");
      34         166 :   params.addRequiredParam<std::vector<Real>>("log10_keq",
      35             :                                              "The list of equilibrium constants for all reactions");
      36         166 :   params.addRequiredParam<std::vector<Real>>(
      37             :       "specific_reactive_surface_area",
      38             :       "The list of specific reactive surface area for all minerals (m^2/L)");
      39         166 :   params.addRequiredParam<std::vector<Real>>(
      40             :       "kinetic_rate_constant", "The list of kinetic rate constant for all reactions (mol/m^2/s)");
      41         166 :   params.addRequiredParam<std::vector<Real>>(
      42             :       "activation_energy", "The list of activation energy values for all reactions (J/mol)");
      43         166 :   params.addParam<Real>("gas_constant", 8.314, "Gas constant. Default is 8.314 (J/mol/K)");
      44         166 :   params.addRequiredParam<std::vector<Real>>(
      45             :       "reference_temperature", "The list of reference temperatures for all reactions (K)");
      46         166 :   params.addRequiredCoupledVar("system_temperature",
      47             :                                "The system temperature for all reactions (K)");
      48          83 :   params.addClassDescription("Adds solid kinetic Kernels and AuxKernels for primary species");
      49          83 :   return params;
      50           0 : }
      51             : 
      52          83 : AddCoupledSolidKinSpeciesAction::AddCoupledSolidKinSpeciesAction(const InputParameters & params)
      53             :   : Action(params),
      54          83 :     _primary_species(getParam<std::vector<NonlinearVariableName>>("primary_species")),
      55         249 :     _secondary_species(getParam<std::vector<AuxVariableName>>("secondary_species")),
      56          83 :     _kinetic_species_involved(_primary_species.size()),
      57          83 :     _weights(_primary_species.size()),
      58         249 :     _input_reactions(getParam<std::string>("kin_reactions")),
      59         166 :     _logk(getParam<std::vector<Real>>("log10_keq")),
      60         166 :     _r_area(getParam<std::vector<Real>>("specific_reactive_surface_area")),
      61         166 :     _ref_kconst(getParam<std::vector<Real>>("kinetic_rate_constant")),
      62         166 :     _e_act(getParam<std::vector<Real>>("activation_energy")),
      63         166 :     _gas_const(getParam<Real>("gas_constant")),
      64         166 :     _ref_temp(getParam<std::vector<Real>>("reference_temperature")),
      65         332 :     _sys_temp(getParam<std::vector<VariableName>>("system_temperature"))
      66             : {
      67             :   // Note: as the reaction syntax has changed, check to see if the old syntax has
      68             :   // been used and throw an informative error. The number of = signs should be one
      69             :   // more than the number of commas, while the smallest number of spaces possible is 2
      70             :   bool old_syntax = false;
      71          83 :   if (std::count(_input_reactions.begin(), _input_reactions.end(), '=') !=
      72          83 :       std::count(_input_reactions.begin(), _input_reactions.end(), ',') + 1)
      73             :     old_syntax = true;
      74             : 
      75          83 :   if (std::count(_input_reactions.begin(), _input_reactions.end(), ' ') < 2)
      76             :     old_syntax = true;
      77             : 
      78          83 :   if (old_syntax)
      79           0 :     mooseError("Old solid kinetic reaction syntax present.\nReactions should now be comma "
      80             :                "separated, and must have spaces between species and +/-/= operators.\n"
      81             :                "See #9972 for details");
      82             : 
      83             :   // Parse the kinetic reactions
      84             :   pcrecpp::RE re_reactions("(.+?)" // A single reaction (any character until the comma delimiter)
      85             :                            "(?:,\\s*|$)" // comma or end of string
      86             :                            ,
      87          83 :                            pcrecpp::RE_Options().set_extended(true));
      88             : 
      89          83 :   pcrecpp::RE re_terms("(\\S+)");
      90             :   pcrecpp::RE re_coeff_and_species("(?: \\(? (.*?) \\)? )" // match the leading coefficent
      91             :                                    "([A-Za-z].*)"          // match the species
      92             :                                    ,
      93          83 :                                    pcrecpp::RE_Options().set_extended(true));
      94             : 
      95             :   pcrecpp::StringPiece input(_input_reactions);
      96             :   pcrecpp::StringPiece single_reaction, term;
      97             :   std::string single_reaction_str;
      98             : 
      99             :   // Parse reaction network to extract each individual reaction
     100         204 :   while (re_reactions.FindAndConsume(&input, &single_reaction_str))
     101         121 :     _reactions.push_back(single_reaction_str);
     102             : 
     103          83 :   _num_reactions = _reactions.size();
     104             : 
     105          83 :   if (_num_reactions == 0)
     106           0 :     mooseError("No solid kinetic reaction provided!");
     107             : 
     108             :   // Start parsing each reaction
     109         204 :   for (unsigned int i = 0; i < _num_reactions; ++i)
     110             :   {
     111         121 :     single_reaction = _reactions[i];
     112             : 
     113             :     // Capture all of the terms
     114             :     std::string species, coeff_str;
     115             :     Real coeff;
     116             :     int sign = 1;
     117             :     bool secondary = false;
     118             : 
     119             :     std::vector<Real> local_stos;
     120             :     std::vector<VariableName> local_species_list;
     121             : 
     122             :     // Find every single term in this reaction (species and operators)
     123         830 :     while (re_terms.FindAndConsume(&single_reaction, &term))
     124             :     {
     125             :       // Separating the stoichiometric coefficients from species
     126         709 :       if (re_coeff_and_species.PartialMatch(term, &coeff_str, &species))
     127             :       {
     128         415 :         if (coeff_str.length())
     129          95 :           coeff = std::stod(coeff_str);
     130             :         else
     131         320 :           coeff = 1.0;
     132             : 
     133         415 :         coeff *= sign;
     134             : 
     135         415 :         if (secondary)
     136         242 :           _solid_kinetic_species.push_back(species);
     137             :         else
     138             :         {
     139         294 :           local_stos.push_back(coeff);
     140         294 :           local_species_list.push_back(species);
     141             :         }
     142             :       }
     143             :       // Finding the operators and assign value of -1.0 to "-" sign
     144         294 :       else if (term == "+" || term == "=" || term == "-")
     145             :       {
     146         294 :         if (term == "-")
     147             :         {
     148             :           sign = -1;
     149          71 :           term = "+";
     150             :         }
     151             : 
     152         294 :         if (term == "=")
     153             :           secondary = true;
     154             :       }
     155             :       else
     156           0 :         mooseError("Error parsing term: ", term.as_string());
     157             :     }
     158             : 
     159         121 :     _stos.push_back(local_stos);
     160         121 :     _primary_species_involved.push_back(local_species_list);
     161         121 :   }
     162             : 
     163             :   // Start picking out primary species and coupled primary species and assigning
     164             :   // corresponding stoichiomentric coefficients
     165         339 :   for (unsigned int i = 0; i < _primary_species.size(); ++i)
     166         664 :     for (unsigned int j = 0; j < _num_reactions; ++j)
     167             :     {
     168        1380 :       for (unsigned int k = 0; k < _primary_species_involved[j].size(); ++k)
     169         972 :         if (_primary_species_involved[j][k] == _primary_species[i])
     170             :         {
     171         294 :           _weights[i].push_back(_stos[j][k]);
     172         294 :           _kinetic_species_involved[i].push_back(_solid_kinetic_species[j]);
     173             :         }
     174             :     }
     175             : 
     176             :   // Print out details of the solid kinetic reactions to the console
     177          83 :   _console << "Solid kinetic reactions:\n";
     178         204 :   for (unsigned int i = 0; i < _num_reactions; ++i)
     179         121 :     _console << "  Reaction " << i + 1 << ": " << _reactions[i] << "\n";
     180          83 :   _console << std::endl;
     181             : 
     182             :   // Check that all secondary species read from the reaction network have been added
     183             :   // as AuxVariables. Note: can't sort the _solid_kinetic_species vector as it throws
     184             :   // out the species and coefficient vectors so use std::is_permutation
     185          83 :   if (!std::is_permutation(
     186             :           _secondary_species.begin(), _secondary_species.end(), _solid_kinetic_species.begin()))
     187           0 :     mooseError("All solid kinetic species must be added as secondary species");
     188             : 
     189             :   // Check that the size of property vectors is equal to the number of reactions
     190          83 :   if (_logk.size() != _num_reactions)
     191           0 :     mooseError("The number of values entered for log10_keq is not equal to the number of solid "
     192             :                "kinetic reactions");
     193          83 :   if (_r_area.size() != _num_reactions)
     194           0 :     mooseError("The number of values entered for specific_reactive_surface_area is not equal to "
     195             :                "the number of solid kinetic reactions");
     196          83 :   if (_ref_kconst.size() != _num_reactions)
     197           0 :     mooseError("The number of values entered for kinetic_rate_constant is not equal to the number "
     198             :                "of solid kinetic reactions");
     199          83 :   if (_e_act.size() != _num_reactions)
     200           0 :     mooseError("The number of values entered for activation_energy is not equal to the number of "
     201             :                "solid kinetic reactions");
     202          83 :   if (_ref_temp.size() != _num_reactions)
     203           0 :     mooseError("The number of values entered for reference_temperature is not equal to the number "
     204             :                "of solid kinetic reactions");
     205          83 : }
     206             : 
     207             : void
     208         166 : AddCoupledSolidKinSpeciesAction::act()
     209             : {
     210         166 :   if (_current_task == "add_kernel")
     211             :   {
     212             :     // Add Kernels for each primary species
     213         339 :     for (unsigned int i = 0; i < _primary_species.size(); ++i)
     214             :     {
     215         512 :       InputParameters params_kin = _factory.getValidParams("CoupledBEKinetic");
     216         512 :       params_kin.set<NonlinearVariableName>("variable") = _primary_species[i];
     217         512 :       params_kin.set<std::vector<Real>>("weight") = _weights[i];
     218         256 :       params_kin.set<std::vector<VariableName>>("v") = _kinetic_species_involved[i];
     219         768 :       _problem->addKernel("CoupledBEKinetic", _primary_species[i] + "_" + "_kin", params_kin);
     220         256 :     }
     221             :   }
     222             : 
     223         166 :   if (_current_task == "add_aux_kernel")
     224             :   {
     225             :     // Add AuxKernels for each solid kinetic species
     226         204 :     for (unsigned int i = 0; i < _num_reactions; ++i)
     227             :     {
     228         121 :       InputParameters params_kin = _factory.getValidParams("KineticDisPreConcAux");
     229         242 :       params_kin.set<AuxVariableName>("variable") = _solid_kinetic_species[i];
     230         242 :       params_kin.defaultCoupledValue("log_k", _logk[i]);
     231         121 :       params_kin.set<Real>("r_area") = _r_area[i];
     232         121 :       params_kin.set<Real>("ref_kconst") = _ref_kconst[i];
     233         121 :       params_kin.set<Real>("e_act") = _e_act[i];
     234         121 :       params_kin.set<Real>("gas_const") = _gas_const;
     235         121 :       params_kin.set<Real>("ref_temp") = _ref_temp[i];
     236         242 :       params_kin.set<std::vector<VariableName>>("sys_temp") = _sys_temp;
     237         242 :       params_kin.set<std::vector<Real>>("sto_v") = _stos[i];
     238         121 :       params_kin.set<std::vector<VariableName>>("v") = _primary_species_involved[i];
     239         242 :       _problem->addAuxKernel(
     240         121 :           "KineticDisPreConcAux", "aux_" + _solid_kinetic_species[i], params_kin);
     241         121 :     }
     242             :   }
     243         166 : }