LCOV - code coverage report
Current view: top level - src/materials - PorousFlowPorosityLinear.C (source / functions) Hit Total Coverage
Test: idaholab/moose porous_flow: #31405 (292dce) with base fef103 Lines: 93 94 98.9 %
Date: 2025-09-04 07:55:56 Functions: 5 5 100.0 %
Legend: Lines: hit not hit

          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 "PorousFlowPorosityLinear.h"
      11             : 
      12             : registerMooseObject("PorousFlowApp", PorousFlowPorosityLinear);
      13             : 
      14             : InputParameters
      15         455 : PorousFlowPorosityLinear::validParams()
      16             : {
      17         455 :   InputParameters params = PorousFlowPorosityBase::validParams();
      18         910 :   params.addParam<bool>("strain_at_nearest_qp",
      19         910 :                         false,
      20             :                         "When calculating nodal porosity that depends on strain, use the strain at "
      21             :                         "the nearest quadpoint.  This adds a small extra computational burden, and "
      22             :                         "is not necessary for simulations involving only linear lagrange elements. "
      23             :                         " If you set this to true, you will also want to set the same parameter to "
      24             :                         "true for related Kernels and Materials");
      25         910 :   params.addRequiredCoupledVar(
      26             :       "porosity_ref",
      27             :       "The porosity at reference effective porepressure, temperature and volumetric strain.  This "
      28             :       "must be a real number or a constant monomial variable (not a linear lagrange or other type "
      29             :       "of variable)");
      30         910 :   params.addCoupledVar(
      31             :       "P_ref",
      32             :       0.0,
      33             :       "The reference effective porepressure.  This should usually be a linear-lagrange variable");
      34         910 :   params.addCoupledVar(
      35             :       "T_ref",
      36             :       0.0,
      37             :       "The reference temperature.  This should usually be a linear-lagrange variable");
      38         910 :   params.addCoupledVar("epv_ref",
      39             :                        0.0,
      40             :                        "The reference volumetric strain.  This must be a real number or a constant "
      41             :                        "monomial variable (not a linear lagrange or other type of variable)");
      42         910 :   params.addParam<Real>("P_coeff", 0.0, "Effective porepressure coefficient");
      43         910 :   params.addParam<Real>("T_coeff", 0.0, "Temperature coefficient");
      44         910 :   params.addParam<Real>("epv_coeff", 0.0, "Volumetric-strain coefficient");
      45        1365 :   params.addRangeCheckedParam<Real>(
      46             :       "porosity_min",
      47         910 :       0.0,
      48             :       "porosity_min >= 0",
      49             :       "Minimum allowed value of the porosity: if the linear relationship gives "
      50             :       "values less than this value, then porosity is set to this value instead");
      51         910 :   params.addParam<Real>(
      52             :       "zero_modifier",
      53         910 :       1E-3,
      54             :       "If the linear relationship produces porosity < porosity_min, then porosity is set "
      55             :       "porosity_min.  This means the derivatives of it will be zero.  However, these zero "
      56             :       "derivatives often result in poor NR convergence, so the derivatives are set to "
      57             :       "_zero_modifier * (values that are relevant for min porosity) to hint to the NR that "
      58             :       "porosity is not always constant");
      59         910 :   params.addParamNamesToGroup("zero_modifier", "Advanced");
      60         455 :   params.addClassDescription(
      61             :       "This Material calculates the porosity in PorousFlow simulations using the relationship "
      62             :       "porosity_ref + P_coeff * (P - P_ref) + T_coeff * (T - T_ref) + epv_coeff * (epv - epv_ref), "
      63             :       "where P is the effective porepressure, T is the temperature and epv is the volumetric "
      64             :       "strain");
      65         455 :   return params;
      66           0 : }
      67             : 
      68         357 : PorousFlowPorosityLinear::PorousFlowPorosityLinear(const InputParameters & parameters)
      69             :   : PorousFlowPorosityBase(parameters),
      70         357 :     _strain_at_nearest_qp(getParam<bool>("strain_at_nearest_qp")),
      71         714 :     _porosity_min(getParam<Real>("porosity_min")),
      72         357 :     _phi_ref(coupledValue("porosity_ref")),
      73         615 :     _P_ref(_nodal_material ? coupledDofValues("P_ref") : coupledValue("P_ref")),
      74         615 :     _T_ref(_nodal_material ? coupledDofValues("T_ref") : coupledValue("T_ref")),
      75         357 :     _epv_ref(coupledValue("epv_ref")),
      76         714 :     _P_coeff(getParam<Real>("P_coeff")),
      77         714 :     _T_coeff(getParam<Real>("T_coeff")),
      78         714 :     _epv_coeff(getParam<Real>("epv_coeff")),
      79             : 
      80         507 :     _uses_volstrain(parameters.isParamSetByUser("epv_coeff") && _epv_coeff != 0.0),
      81         714 :     _vol_strain_qp(getOptionalMaterialProperty<Real>("PorousFlow_total_volumetric_strain_qp")),
      82         714 :     _dvol_strain_qp_dvar(getOptionalMaterialProperty<std::vector<RealGradient>>(
      83             :         "dPorousFlow_total_volumetric_strain_qp_dvar")),
      84             : 
      85         375 :     _uses_pf(parameters.isParamSetByUser("P_coeff") && _P_coeff != 0.0),
      86         714 :     _pf_nodal(getOptionalMaterialProperty<Real>("PorousFlow_effective_fluid_pressure_nodal")),
      87         714 :     _pf_qp(getOptionalMaterialProperty<Real>("PorousFlow_effective_fluid_pressure_qp")),
      88         714 :     _dpf_dvar_nodal(getOptionalMaterialProperty<std::vector<Real>>(
      89             :         "dPorousFlow_effective_fluid_pressure_nodal_dvar")),
      90         714 :     _dpf_dvar_qp(getOptionalMaterialProperty<std::vector<Real>>(
      91             :         "dPorousFlow_effective_fluid_pressure_qp_dvar")),
      92             : 
      93         507 :     _uses_T(parameters.isParamSetByUser("T_coeff") && _T_coeff != 0.0),
      94         714 :     _temperature_nodal(getOptionalMaterialProperty<Real>("PorousFlow_temperature_nodal")),
      95         714 :     _temperature_qp(getOptionalMaterialProperty<Real>("PorousFlow_temperature_qp")),
      96         357 :     _dtemperature_dvar_nodal(
      97         357 :         getOptionalMaterialProperty<std::vector<Real>>("dPorousFlow_temperature_nodal_dvar")),
      98         357 :     _dtemperature_dvar_qp(
      99         357 :         getOptionalMaterialProperty<std::vector<Real>>("dPorousFlow_temperature_qp_dvar")),
     100        1071 :     _zero_modifier(getParam<Real>("zero_modifier"))
     101             : {
     102         357 : }
     103             : 
     104             : void
     105         336 : PorousFlowPorosityLinear::initialSetup()
     106             : {
     107         336 :   _pf = _nodal_material ? _pf_nodal.get() : _pf_qp.get();
     108         336 :   _dpf_dvar = _nodal_material ? _dpf_dvar_nodal.get() : _dpf_dvar_qp.get();
     109         336 :   if (!(_pf && _dpf_dvar) && _uses_pf)
     110             :   {
     111           2 :     mooseError("PorousFlowPorosityLinear: When P_coeff is given you must have an effective fluid "
     112             :                "pressure Material");
     113             :   }
     114             : 
     115         334 :   _temperature = _nodal_material ? _temperature_nodal.get() : _temperature_qp.get();
     116         334 :   _dtemperature_dvar =
     117         334 :       _nodal_material ? _dtemperature_dvar_nodal.get() : _dtemperature_dvar_qp.get();
     118         334 :   if (!_temperature && _uses_T)
     119           2 :     mooseError(
     120             :         "PorousFlowPorosityLinear: When T_coeff is given you must have a temperature Material");
     121             : 
     122         332 :   if (!(_vol_strain_qp && _dvol_strain_qp_dvar) && _uses_volstrain)
     123           2 :     mooseError("PorousFlowPorosityLinear: When epv_coeff is given you must have a "
     124             :                "volumetric-strain Material");
     125             : 
     126         330 :   PorousFlowPorosityBase::initialSetup();
     127         330 : }
     128             : 
     129             : void
     130      179731 : PorousFlowPorosityLinear::initQpStatefulProperties()
     131             : {
     132      179731 :   _porosity[_qp] = _phi_ref[0];
     133      179731 :   (*_dporosity_dvar)[_qp].assign(_num_var, 0.0);
     134      179731 :   (*_dporosity_dgradvar)[_qp].assign(_num_var, 0.0);
     135      179731 : }
     136             : 
     137             : void
     138     1028640 : PorousFlowPorosityLinear::computeQpProperties()
     139             : {
     140             :   // note the [0] below: _phi_ref is a constant monomial and we use [0] regardless of
     141             :   // _nodal_material
     142     1028640 :   Real porosity = _phi_ref[0];
     143     1028640 :   if (_uses_pf)
     144     1028640 :     porosity += _P_coeff * ((*_pf)[_qp] - _P_ref[_qp]);
     145     1028640 :   if (_uses_T)
     146        2440 :     porosity += _T_coeff * ((*_temperature)[_qp] - _T_ref[_qp]);
     147     1028640 :   if (_uses_volstrain)
     148             :   {
     149             :     // Note that in the following _strain[_qp] is evaluated at q quadpoint
     150             :     // So _porosity_nodal[_qp], which should be the nodal value of porosity
     151             :     // actually uses the strain at a quadpoint.  This
     152             :     // is OK for LINEAR elements, as strain is constant over the element anyway.
     153             :     const unsigned qp_to_use =
     154        2440 :         (_nodal_material && (_bnd || _strain_at_nearest_qp) ? nearestQP(_qp) : _qp);
     155        2440 :     porosity += _epv_coeff * (_vol_strain_qp[qp_to_use] - _epv_ref[0]);
     156             :   }
     157             : 
     158     1028640 :   if (porosity < _porosity_min)
     159         205 :     _porosity[_qp] = _porosity_min;
     160             :   else
     161     1028435 :     _porosity[_qp] = porosity;
     162             : 
     163     1028640 :   (*_dporosity_dvar)[_qp].resize(_num_var);
     164     1028640 :   (*_dporosity_dgradvar)[_qp].resize(_num_var);
     165     3092760 :   for (unsigned int pvar = 0; pvar < _num_var; ++pvar)
     166             :   {
     167             :     Real deriv = 0.0;
     168     2064120 :     if (_uses_pf)
     169     2064120 :       deriv += _P_coeff * (*_dpf_dvar)[_qp][pvar];
     170     2064120 :     if (_uses_T)
     171       11720 :       deriv += _T_coeff * (*_dtemperature_dvar)[_qp][pvar];
     172     2064120 :     if (porosity < _porosity_min)
     173         590 :       (*_dporosity_dvar)[_qp][pvar] = _zero_modifier * deriv;
     174             :     else
     175     2063530 :       (*_dporosity_dvar)[_qp][pvar] = deriv;
     176             : 
     177             :     RealGradient deriv_grad(0.0, 0.0, 0.0);
     178     2064120 :     if (_uses_volstrain)
     179             :     {
     180             :       const unsigned qp_to_use =
     181       11720 :           (_nodal_material && (_bnd || _strain_at_nearest_qp) ? nearestQP(_qp) : _qp);
     182       11720 :       deriv_grad += _epv_coeff * _dvol_strain_qp_dvar[qp_to_use][pvar];
     183             :     }
     184     2064120 :     if (porosity < _porosity_min)
     185         590 :       (*_dporosity_dgradvar)[_qp][pvar] = _zero_modifier * deriv_grad;
     186             :     else
     187     2063530 :       (*_dporosity_dgradvar)[_qp][pvar] = deriv_grad;
     188             :   }
     189     1028640 : }

Generated by: LCOV version 1.14