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 "PorousFlowPreDis.h"
      11             : 
      12             : registerMooseObject("PorousFlowApp", PorousFlowPreDis);
      13             : 
      14             : InputParameters
      15        1174 : PorousFlowPreDis::validParams()
      16             : {
      17        1174 :   InputParameters params = TimeKernel::validParams();
      18        2348 :   params.addRequiredParam<std::vector<Real>>(
      19             :       "mineral_density",
      20             :       "Density (kg(precipitate)/m^3(precipitate)) of each secondary species in the "
      21             :       "aqueous precipitation-dissolution reaction system");
      22        2348 :   params.addRequiredParam<UserObjectName>(
      23             :       "PorousFlowDictator", "The UserObject that holds the list of PorousFlow variable names.");
      24        2348 :   params.addRequiredParam<std::vector<Real>>("stoichiometry",
      25             :                                              "A vector of stoichiometric coefficients for the "
      26             :                                              "primary species that is the Variable of this Kernel: "
      27             :                                              "one for each precipitation-dissolution reaction "
      28             :                                              "(these are one columns of the 'reactions' matrix)");
      29        1174 :   params.addClassDescription("Precipitation-dissolution of chemical species");
      30        1174 :   return params;
      31           0 : }
      32             : 
      33         653 : PorousFlowPreDis::PorousFlowPreDis(const InputParameters & parameters)
      34             :   : TimeKernel(parameters),
      35         653 :     _mineral_density(getParam<std::vector<Real>>("mineral_density")),
      36         653 :     _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
      37         653 :     _aq_ph(_dictator.aqueousPhaseNumber()),
      38        1306 :     _porosity_old(getMaterialPropertyOld<Real>("PorousFlow_porosity_nodal")),
      39        1306 :     _saturation(getMaterialProperty<std::vector<Real>>("PorousFlow_saturation_nodal")),
      40         653 :     _dsaturation_dvar(
      41         653 :         getMaterialProperty<std::vector<std::vector<Real>>>("dPorousFlow_saturation_nodal_dvar")),
      42         653 :     _reaction_rate(
      43         653 :         getMaterialProperty<std::vector<Real>>("PorousFlow_mineral_reaction_rate_nodal")),
      44        1306 :     _dreaction_rate_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      45             :         "dPorousFlow_mineral_reaction_rate_nodal_dvar")),
      46        2612 :     _stoichiometry(getParam<std::vector<Real>>("stoichiometry"))
      47             : {
      48             :   /* Not needed due to PorousFlow_mineral_reaction_rate already checking this condition
      49             :   if (_dictator.numPhases() < 1)
      50             :     mooseError("PorousFlowPreDis: The number of fluid phases must not be zero");
      51             :   */
      52             : 
      53         653 :   if (_mineral_density.size() != _dictator.numAqueousKinetic())
      54           2 :     paramError("mineral_density",
      55             :                "The Dictator proclaims that the number of precipitation-dissolution secondary "
      56             :                "species in this simulation is ",
      57           2 :                _dictator.numAqueousKinetic(),
      58             :                " whereas you have provided ",
      59             :                _mineral_density.size(),
      60             :                ". The Dictator does not take such mistakes lightly");
      61             : 
      62         651 :   if (_stoichiometry.size() != _dictator.numAqueousKinetic())
      63           2 :     paramError("stoichiometry",
      64             :                "The Dictator proclaims that the number of precipitation-dissolution secondary "
      65             :                "species in this simulation is ",
      66           2 :                _dictator.numAqueousKinetic(),
      67             :                " whereas you have provided ",
      68             :                _stoichiometry.size(),
      69             :                ". The Dictator does not take such mistakes lightly");
      70         649 : }
      71             : 
      72             : Real
      73     1385432 : PorousFlowPreDis::computeQpResidual()
      74             : {
      75             :   /*
      76             :    *
      77             :    * Note the use of the OLD value of porosity here.
      78             :    * This strategy, which breaks the cyclic dependency between porosity
      79             :    * and mineral concentration, is used in
      80             :    * Kernel: PorousFlowPreDis
      81             :    * Material: PorousFlowPorosity
      82             :    * Material: PorousFlowAqueousPreDisChemistry
      83             :    * Material: PorousFlowAqueousPreDisMineral
      84             :    *
      85             :    */
      86             :   Real res = 0.0;
      87     2777344 :   for (unsigned r = 0; r < _dictator.numAqueousKinetic(); ++r)
      88     1391912 :     res += _stoichiometry[r] * _mineral_density[r] * _reaction_rate[_i][r];
      89     1385432 :   return _test[_i][_qp] * res * _porosity_old[_i] * _saturation[_i][_aq_ph];
      90             : }
      91             : 
      92             : Real
      93     2199496 : PorousFlowPreDis::computeQpJacobian()
      94             : {
      95             :   /// If the variable is not a PorousFlow variable (very unusual), the diag Jacobian terms are 0
      96     2199496 :   if (_dictator.notPorousFlowVariable(_var.number()))
      97             :     return 0.0;
      98     2199496 :   return computeQpJac(_dictator.porousFlowVariableNum(_var.number()));
      99             : }
     100             : 
     101             : Real
     102     5034600 : PorousFlowPreDis::computeQpOffDiagJacobian(unsigned int jvar)
     103             : {
     104             :   /// If the variable is not a PorousFlow variable, the OffDiag Jacobian terms are 0
     105     5034600 :   if (_dictator.notPorousFlowVariable(jvar))
     106             :     return 0.0;
     107     5034600 :   return computeQpJac(_dictator.porousFlowVariableNum(jvar));
     108             : }
     109             : 
     110             : Real
     111     7234096 : PorousFlowPreDis::computeQpJac(unsigned int pvar)
     112             : {
     113     7234096 :   if (_i != _j)
     114             :     return 0.0;
     115             : 
     116             :   Real res = 0.0;
     117             :   Real dres = 0.0;
     118     6027936 :   for (unsigned r = 0; r < _dictator.numAqueousKinetic(); ++r)
     119             :   {
     120     3015048 :     dres += _stoichiometry[r] * _mineral_density[r] * _dreaction_rate_dvar[_i][r][pvar];
     121     3015048 :     res += _stoichiometry[r] * _mineral_density[r] * _reaction_rate[_i][r];
     122             :   }
     123             : 
     124     3012888 :   return _test[_i][_qp] *
     125     3012888 :          (dres * _saturation[_i][_aq_ph] + res * _dsaturation_dvar[_i][_aq_ph][pvar]) *
     126     3012888 :          _porosity_old[_i];
     127             : }