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             : //  "Monomial" form of relative permeability
      11             : //
      12             : #include "RichardsRelPermMonomial.h"
      13             : 
      14             : registerMooseObject("RichardsApp", RichardsRelPermMonomial);
      15             : 
      16             : InputParameters
      17          26 : RichardsRelPermMonomial::validParams()
      18             : {
      19          26 :   InputParameters params = RichardsRelPerm::validParams();
      20          52 :   params.addRequiredRangeCheckedParam<Real>(
      21             :       "simm",
      22             :       "simm >= 0 & simm < 1",
      23             :       "Immobile saturation.  Must be between 0 and 1.   Define s = "
      24             :       "(seff - simm)/(1 - simm).  Then relperm = s^n");
      25          52 :   params.addRequiredRangeCheckedParam<Real>(
      26             :       "n",
      27             :       "n >= 0",
      28             :       "Exponent.  Must be >= 0.   Define s = (seff - simm)/(1 - simm).  Then relperm = s^n");
      29          52 :   params.addParam<Real>(
      30          52 :       "zero_to_the_zero", 0.0, "If n=0, this is the value of relative permeability for s<=simm");
      31          26 :   params.addClassDescription("Monomial form of relative permeability.  Define s = (seff - simm)/(1 "
      32             :                              "- simm).  Then relperm = s^n if s<simm, otherwise relperm=1");
      33          26 :   return params;
      34           0 : }
      35             : 
      36           9 : RichardsRelPermMonomial::RichardsRelPermMonomial(const InputParameters & parameters)
      37             :   : RichardsRelPerm(parameters),
      38           9 :     _simm(getParam<Real>("simm")),
      39          18 :     _n(getParam<Real>("n")),
      40          27 :     _zero_to_the_zero(getParam<Real>("zero_to_the_zero"))
      41             : {
      42           9 : }
      43             : 
      44             : Real
      45     1598740 : RichardsRelPermMonomial::relperm(Real seff) const
      46             : {
      47     1598740 :   if (seff >= 1.0)
      48             :     return 1.0;
      49             : 
      50     1598740 :   if (_n == 0 && seff <= _simm)
      51          60 :     return _zero_to_the_zero;
      52             : 
      53     1598680 :   if (seff <= _simm)
      54             :     return 0.0;
      55             : 
      56     1598680 :   Real s_internal = (seff - _simm) / (1.0 - _simm);
      57     1598680 :   Real krel = std::pow(s_internal, _n);
      58             : 
      59             :   // bound, just in case
      60     1598680 :   if (krel < 0)
      61             :   {
      62             :     krel = 0;
      63             :   }
      64     1598680 :   if (krel > 1)
      65             :   {
      66             :     krel = 1;
      67             :   }
      68             :   return krel;
      69             : }
      70             : 
      71             : Real
      72     2532268 : RichardsRelPermMonomial::drelperm(Real seff) const
      73             : {
      74     2532268 :   if (seff >= 1.0)
      75             :     return 0.0;
      76             : 
      77     2532268 :   if (seff <= _simm)
      78             :     return 0.0;
      79             : 
      80     2532208 :   if (_n == 0)
      81             :     return 0.0;
      82             : 
      83     2531662 :   Real s_internal = (seff - _simm) / (1.0 - _simm);
      84     2531662 :   Real krelp = _n * std::pow(s_internal, _n - 1);
      85     2531662 :   return krelp / (1.0 - _simm);
      86             : }
      87             : 
      88             : Real
      89      934740 : RichardsRelPermMonomial::d2relperm(Real seff) const
      90             : {
      91      934740 :   if (seff >= 1.0)
      92             :     return 0.0;
      93             : 
      94      934740 :   if (seff <= _simm)
      95             :     return 0.0;
      96             : 
      97      934680 :   if (_n == 0)
      98             :     return 0.0;
      99             : 
     100      934134 :   Real s_internal = (seff - _simm) / (1.0 - _simm);
     101      934134 :   Real krelpp = _n * (_n - 1) * std::pow(s_internal, _n - 2);
     102      934134 :   return krelpp / std::pow(1.0 - _simm, 2);
     103             : }