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 "GeochemistryQuantityAux.h"
      11             : 
      12             : registerMooseObject("GeochemistryApp", GeochemistryQuantityAux);
      13             : 
      14             : InputParameters
      15      140617 : GeochemistryQuantityAux::validParams()
      16             : {
      17      140617 :   InputParameters params = AuxKernel::validParams();
      18      140617 :   params.addClassDescription(
      19             :       "Extracts information from the Reactor and records it in the AuxVariable");
      20      281234 :   params.addRequiredParam<std::string>("species", "Species name");
      21             :   MooseEnum quantity_choice("molal mg_per_kg free_mg free_cm3 neglog10a activity bulk_moles "
      22             :                             "surface_charge surface_potential temperature kinetic_moles "
      23             :                             "kinetic_additions moles_dumped transported_moles_in_original_basis",
      24      281234 :                             "molal");
      25      281234 :   params.addParam<MooseEnum>(
      26             :       "quantity",
      27             :       quantity_choice,
      28             :       "Type of quantity to output.  These are available for non-kinetic species: activity (which "
      29             :       "equals "
      30             :       "fugacity for gases); bulk moles (this will be zero if the species is not in the basis); "
      31             :       "neglog10a = -log10(activity); transported_moles_in_original_basis (will throw an error if "
      32             :       "species is not in original basis).  These are available only for non-kinetic non-minerals: "
      33             :       "molal = "
      34             :       "mol(species)/kg(solvent_water); mg_per_kg = mg(species)/kg(solvent_water).  These are "
      35             :       "available only for minerals: "
      36             :       "free_mg = free mg; free_cm3 = free cubic-centimeters; moles_dumped = moles dumped when "
      37             :       "special dump and flow-through modes are active.  These are available for minerals "
      38             :       "that host sorbing sites: surface_charge (C/m^2); surface_potential (V).  These are "
      39             :       "available for kinetic species: kinetic_moles; kinetic_additions (-dt * rate = mole "
      40             :       "increment in kinetic species for this timestep).  If "
      41             :       "quantity=temperature, then the 'species' is ignored and the AuxVariable will record the "
      42             :       "aqueous solution temperature in degC");
      43      281234 :   params.addRequiredParam<UserObjectName>("reactor",
      44             :                                           "The name of the Geochemistry*Reactor user object.");
      45      281234 :   params.declareControllable("species");
      46      140617 :   return params;
      47      140617 : }
      48             : 
      49       75777 : GeochemistryQuantityAux::GeochemistryQuantityAux(const InputParameters & parameters)
      50             :   : AuxKernel(parameters),
      51       75777 :     _species(getParam<std::string>("species")),
      52       75777 :     _reactor(getUserObject<GeochemistryReactorBase>("reactor")),
      53      151554 :     _quantity_choice(getParam<MooseEnum>("quantity").getEnum<QuantityChoiceEnum>()),
      54       75777 :     _surface_sorption_mineral_index(0)
      55             : {
      56             :   const ModelGeochemicalDatabase & mgd =
      57       75777 :       _reactor.getPertinentGeochemicalSystem().modelGeochemicalDatabase();
      58       75777 :   if (!(mgd.basis_species_index.count(_species) == 1 ||
      59       60038 :         mgd.eqm_species_index.count(_species) == 1 || mgd.kin_species_index.count(_species) == 1))
      60           2 :     paramError(
      61             :         "species",
      62           2 :         _species,
      63             :         " does not appear in the model's geochemical system either as a basis or equilibrium "
      64             :         "or kinetic species, but you requested an Aux involving it");
      65             : 
      66             :   bool is_mineral = false;
      67       75775 :   if (mgd.basis_species_index.count(_species) == 1)
      68       16349 :     is_mineral = mgd.basis_species_mineral[mgd.basis_species_index.at(_species)];
      69       59426 :   else if (mgd.eqm_species_index.count(_species) == 1)
      70       58818 :     is_mineral = mgd.eqm_species_mineral[mgd.eqm_species_index.at(_species)];
      71             :   else
      72         608 :     is_mineral = mgd.kin_species_mineral[mgd.kin_species_index.at(_species)];
      73       75775 :   bool is_kinetic = (mgd.kin_species_index.count(_species) == 1);
      74             : 
      75       75775 :   if (!is_mineral && (_quantity_choice == QuantityChoiceEnum::FREE_CM3 ||
      76       70825 :                       _quantity_choice == QuantityChoiceEnum::FREE_MG ||
      77             :                       _quantity_choice == QuantityChoiceEnum::SURFACE_CHARGE ||
      78             :                       _quantity_choice == QuantityChoiceEnum::SURFACE_POTENTIAL ||
      79             :                       _quantity_choice == QuantityChoiceEnum::MOLES_DUMPED))
      80          10 :     paramError(
      81             :         "quantity",
      82             :         "the free_mg, free_cm3, moles_dumped and surface-related quantities are only available for "
      83             :         "mineral species");
      84       75765 :   if ((is_mineral || is_kinetic) && (_quantity_choice == QuantityChoiceEnum::MOLAL ||
      85             :                                      _quantity_choice == QuantityChoiceEnum::MG_PER_KG))
      86           8 :     paramError("quantity",
      87             :                "the molal and mg_per_kg quantities are only available for "
      88             :                "non-kinetic, non-mineral species");
      89       75757 :   if (_quantity_choice == QuantityChoiceEnum::SURFACE_CHARGE ||
      90             :       _quantity_choice == QuantityChoiceEnum::SURFACE_POTENTIAL)
      91             :   {
      92         188 :     _surface_sorption_mineral_index = std::distance(
      93             :         mgd.surface_sorption_name.begin(),
      94         188 :         std::find(mgd.surface_sorption_name.begin(), mgd.surface_sorption_name.end(), _species));
      95         188 :     if (_surface_sorption_mineral_index >= mgd.surface_sorption_name.size())
      96           4 :       paramError("species",
      97             :                  "cannot record surface charge or surface potential for a species that is not "
      98             :                  "involved in surface sorption");
      99             :   }
     100       75753 :   if (!is_kinetic && (_quantity_choice == QuantityChoiceEnum::KINETIC_MOLES ||
     101             :                       _quantity_choice == QuantityChoiceEnum::KINETIC_ADDITIONS))
     102           4 :     paramError("quantity",
     103             :                "the kinetic_moles and kinetic_additions quantities are only available for kinetic "
     104             :                "species");
     105         604 :   if (is_kinetic && (_quantity_choice == QuantityChoiceEnum::NEGLOG10A ||
     106         604 :                      _quantity_choice == QuantityChoiceEnum::ACTIVITY ||
     107             :                      _quantity_choice == QuantityChoiceEnum::BULK_MOLES))
     108           6 :     paramError("species", "cannot record activity, neglog10a or bulk_moles for a kinetic species");
     109       75743 :   if (_quantity_choice == QuantityChoiceEnum::TRANSPORTED_MOLES_IN_ORIGINAL_BASIS)
     110          24 :     _reactor.getPertinentGeochemicalSystem().getIndexOfOriginalBasisSpecies(
     111             :         _species); // will throw error if species not in original basis
     112       75743 : }
     113             : 
     114             : Real
     115      640888 : GeochemistryQuantityAux::computeValue()
     116             : {
     117             :   Real ret = 0.0;
     118      640888 :   const GeochemicalSystem & egs = _reactor.getGeochemicalSystem(_current_node->id());
     119      640888 :   const ModelGeochemicalDatabase & mgd = egs.getModelGeochemicalDatabase();
     120      640888 :   const PertinentGeochemicalSystem & pgs = _reactor.getPertinentGeochemicalSystem();
     121      640888 :   if (_quantity_choice == QuantityChoiceEnum::SURFACE_CHARGE)
     122         408 :     ret = egs.getSurfaceCharge(_surface_sorption_mineral_index);
     123             :   else if (_quantity_choice == QuantityChoiceEnum::SURFACE_POTENTIAL)
     124         408 :     ret = egs.getSurfacePotential(_surface_sorption_mineral_index);
     125             :   else if (_quantity_choice == QuantityChoiceEnum::TEMPERATURE)
     126        5792 :     ret = egs.getTemperature();
     127             :   else if (_quantity_choice == QuantityChoiceEnum::MOLES_DUMPED)
     128         440 :     ret = _reactor.getMolesDumped(_current_node->id(), _species);
     129             :   else if (_quantity_choice == QuantityChoiceEnum::TRANSPORTED_MOLES_IN_ORIGINAL_BASIS)
     130         160 :     ret = egs.getTransportedBulkInOriginalBasis()(pgs.getIndexOfOriginalBasisSpecies(_species));
     131             :   else
     132             :   {
     133      633680 :     if (mgd.basis_species_index.count(_species) == 1)
     134             :     {
     135      146852 :       const unsigned basis_ind = mgd.basis_species_index.at(_species);
     136      146852 :       switch (_quantity_choice)
     137             :       {
     138       24932 :         case QuantityChoiceEnum::MG_PER_KG:
     139       24932 :           ret = egs.getSolventMassAndFreeMolalityAndMineralMoles()[basis_ind] *
     140             :                 mgd.basis_species_molecular_weight[basis_ind] * 1000.0;
     141       24932 :           break;
     142        6592 :         case QuantityChoiceEnum::FREE_MG:
     143        6592 :           ret = egs.getSolventMassAndFreeMolalityAndMineralMoles()[basis_ind] *
     144             :                 mgd.basis_species_molecular_weight[basis_ind] * 1000.0;
     145        6592 :           break;
     146        6592 :         case QuantityChoiceEnum::FREE_CM3:
     147        6592 :           ret = egs.getSolventMassAndFreeMolalityAndMineralMoles()[basis_ind] *
     148             :                 mgd.basis_species_molecular_volume[basis_ind];
     149        6592 :           break;
     150        3460 :         case QuantityChoiceEnum::NEGLOG10A:
     151        3460 :           ret = -std::log10(egs.getBasisActivity(basis_ind));
     152        3460 :           break;
     153       37286 :         case QuantityChoiceEnum::ACTIVITY:
     154       37286 :           ret = egs.getBasisActivity(basis_ind);
     155       37286 :           break;
     156       37286 :         case QuantityChoiceEnum::BULK_MOLES:
     157       37286 :           ret = egs.getBulkMolesOld()[basis_ind];
     158       37286 :           break;
     159       30704 :         default:
     160       30704 :           ret = egs.getSolventMassAndFreeMolalityAndMineralMoles()[basis_ind];
     161             :       }
     162             :     }
     163      486828 :     else if (mgd.eqm_species_index.count(_species) == 1)
     164             :     {
     165      479846 :       const unsigned eqm_ind = mgd.eqm_species_index.at(_species);
     166      479846 :       switch (_quantity_choice)
     167             :       {
     168      114470 :         case QuantityChoiceEnum::MG_PER_KG:
     169      114470 :           ret = egs.getEquilibriumMolality(eqm_ind) * mgd.eqm_species_molecular_weight[eqm_ind] *
     170             :                 1000.0;
     171      114470 :           break;
     172        5120 :         case QuantityChoiceEnum::FREE_CM3:
     173        5120 :           ret = egs.getEquilibriumMolality(eqm_ind) * mgd.eqm_species_molecular_volume[eqm_ind];
     174        5120 :           break;
     175        5120 :         case QuantityChoiceEnum::FREE_MG:
     176        5120 :           ret = egs.getEquilibriumMolality(eqm_ind) * mgd.eqm_species_molecular_weight[eqm_ind] *
     177             :                 1000.0;
     178        5120 :           break;
     179        1486 :         case QuantityChoiceEnum::NEGLOG10A:
     180        1486 :           ret = -std::log10(egs.getEquilibriumActivity(eqm_ind));
     181        1486 :           break;
     182      119590 :         case QuantityChoiceEnum::ACTIVITY:
     183      119590 :           ret = egs.getEquilibriumActivity(eqm_ind);
     184      119590 :           break;
     185             :         case QuantityChoiceEnum::BULK_MOLES:
     186             :           ret = 0.0;
     187             :           break;
     188      114470 :         default:
     189      114470 :           ret = egs.getEquilibriumMolality(eqm_ind);
     190             :       }
     191             :     }
     192             :     else
     193             :     {
     194        6982 :       const unsigned kin_ind = mgd.kin_species_index.at(_species);
     195        6982 :       switch (_quantity_choice)
     196             :       {
     197        1648 :         case QuantityChoiceEnum::FREE_CM3:
     198        1648 :           ret = egs.getKineticMoles(kin_ind) * mgd.kin_species_molecular_volume[kin_ind];
     199        1648 :           break;
     200        1658 :         case QuantityChoiceEnum::FREE_MG:
     201        1658 :           ret = egs.getKineticMoles(kin_ind) * mgd.kin_species_molecular_weight[kin_ind] * 1000.0;
     202        1658 :           break;
     203        1758 :         case QuantityChoiceEnum::KINETIC_ADDITIONS:
     204        1758 :           ret = _reactor.getMoleAdditions(_current_node->id())(kin_ind + egs.getNumInBasis());
     205        1758 :           break;
     206        1918 :         default:
     207        1918 :           ret = egs.getKineticMoles(kin_ind);
     208             :       }
     209             :     }
     210             :   }
     211      640888 :   return ret;
     212             : }