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 "PorousFlowHeatVolumetricExpansion.h"
      11             : 
      12             : #include "MooseVariable.h"
      13             : 
      14             : registerMooseObject("PorousFlowApp", PorousFlowHeatVolumetricExpansion);
      15             : 
      16             : InputParameters
      17         430 : PorousFlowHeatVolumetricExpansion::validParams()
      18             : {
      19         430 :   InputParameters params = TimeKernel::validParams();
      20         860 :   params.addParam<bool>("strain_at_nearest_qp",
      21         860 :                         false,
      22             :                         "When calculating nodal porosity that depends on strain, use the strain at "
      23             :                         "the nearest quadpoint.  This adds a small extra computational burden, and "
      24             :                         "is not necessary for simulations involving only linear lagrange elements. "
      25             :                         " If you set this to true, you will also want to set the same parameter to "
      26             :                         "true for related Kernels and Materials");
      27         860 :   params.addRequiredParam<UserObjectName>(
      28             :       "PorousFlowDictator", "The UserObject that holds the list of PorousFlow variable names.");
      29         430 :   params.set<bool>("use_displaced_mesh") = false;
      30         430 :   params.suppressParameter<bool>("use_displaced_mesh");
      31         430 :   params.addClassDescription("Energy-density*rate_of_solid_volumetric_expansion.  The "
      32             :                              "energy-density is lumped to the nodes");
      33         430 :   return params;
      34           0 : }
      35             : 
      36         232 : PorousFlowHeatVolumetricExpansion::PorousFlowHeatVolumetricExpansion(
      37         232 :     const InputParameters & parameters)
      38             :   : TimeKernel(parameters),
      39         232 :     _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
      40         232 :     _var_is_porflow_var(_dictator.isPorousFlowVariable(_var.number())),
      41         232 :     _num_phases(_dictator.numPhases()),
      42         232 :     _fluid_present(_num_phases > 0),
      43         464 :     _strain_at_nearest_qp(getParam<bool>("strain_at_nearest_qp")),
      44         464 :     _porosity(getMaterialProperty<Real>("PorousFlow_porosity_nodal")),
      45         464 :     _dporosity_dvar(getMaterialProperty<std::vector<Real>>("dPorousFlow_porosity_nodal_dvar")),
      46         232 :     _dporosity_dgradvar(
      47         232 :         getMaterialProperty<std::vector<RealGradient>>("dPorousFlow_porosity_nodal_dgradvar")),
      48         464 :     _nearest_qp(_strain_at_nearest_qp
      49         243 :                     ? &getMaterialProperty<unsigned int>("PorousFlow_nearestqp_nodal")
      50             :                     : nullptr),
      51         464 :     _rock_energy_nodal(getMaterialProperty<Real>("PorousFlow_matrix_internal_energy_nodal")),
      52         232 :     _drock_energy_nodal_dvar(
      53         232 :         getMaterialProperty<std::vector<Real>>("dPorousFlow_matrix_internal_energy_nodal_dvar")),
      54         464 :     _fluid_density(_fluid_present ? &getMaterialProperty<std::vector<Real>>(
      55             :                                         "PorousFlow_fluid_phase_density_nodal")
      56             :                                   : nullptr),
      57         464 :     _dfluid_density_dvar(_fluid_present ? &getMaterialProperty<std::vector<std::vector<Real>>>(
      58             :                                               "dPorousFlow_fluid_phase_density_nodal_dvar")
      59             :                                         : nullptr),
      60         232 :     _fluid_saturation_nodal(
      61         464 :         _fluid_present ? &getMaterialProperty<std::vector<Real>>("PorousFlow_saturation_nodal")
      62             :                        : nullptr),
      63         464 :     _dfluid_saturation_nodal_dvar(_fluid_present
      64         464 :                                       ? &getMaterialProperty<std::vector<std::vector<Real>>>(
      65             :                                             "dPorousFlow_saturation_nodal_dvar")
      66             :                                       : nullptr),
      67         464 :     _energy_nodal(_fluid_present ? &getMaterialProperty<std::vector<Real>>(
      68             :                                        "PorousFlow_fluid_phase_internal_energy_nodal")
      69             :                                  : nullptr),
      70         464 :     _denergy_nodal_dvar(_fluid_present ? &getMaterialProperty<std::vector<std::vector<Real>>>(
      71             :                                              "dPorousFlow_fluid_phase_internal_energy_nodal_dvar")
      72             :                                        : nullptr),
      73         464 :     _strain_rate_qp(getMaterialProperty<Real>("PorousFlow_volumetric_strain_rate_qp")),
      74         464 :     _dstrain_rate_qp_dvar(getMaterialProperty<std::vector<RealGradient>>(
      75         232 :         "dPorousFlow_volumetric_strain_rate_qp_dvar"))
      76             : {
      77         232 : }
      78             : 
      79             : Real
      80      275296 : PorousFlowHeatVolumetricExpansion::computeQpResidual()
      81             : {
      82      275296 :   Real energy = (1.0 - _porosity[_i]) * _rock_energy_nodal[_i];
      83      568352 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      84      293056 :     energy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
      85      293056 :               (*_energy_nodal)[_i][ph] * _porosity[_i];
      86             : 
      87      275296 :   return _test[_i][_qp] * energy * _strain_rate_qp[_qp];
      88             : }
      89             : 
      90             : Real
      91      981824 : PorousFlowHeatVolumetricExpansion::computeQpJacobian()
      92             : {
      93      981824 :   return computedEnergyQpJac(_var.number()) + computedVolQpJac(_var.number());
      94             : }
      95             : 
      96             : Real
      97     3392960 : PorousFlowHeatVolumetricExpansion::computeQpOffDiagJacobian(unsigned int jvar)
      98             : {
      99     3392960 :   return computedEnergyQpJac(jvar) + computedVolQpJac(jvar);
     100             : }
     101             : 
     102             : Real
     103     4374784 : PorousFlowHeatVolumetricExpansion::computedVolQpJac(unsigned int jvar)
     104             : {
     105     4374784 :   if (_dictator.notPorousFlowVariable(jvar))
     106             :     return 0.0;
     107             : 
     108     4374784 :   Real energy = (1.0 - _porosity[_i]) * _rock_energy_nodal[_i];
     109     8979968 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
     110     4605184 :     energy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
     111     4605184 :               (*_energy_nodal)[_i][ph] * _porosity[_i];
     112             : 
     113     4374784 :   const unsigned int pvar = _dictator.porousFlowVariableNum(jvar);
     114     4374784 :   Real dvol = _dstrain_rate_qp_dvar[_qp][pvar] * _grad_phi[_j][_qp];
     115             : 
     116     4374784 :   return _test[_i][_qp] * energy * dvol;
     117             : }
     118             : Real
     119     4374784 : PorousFlowHeatVolumetricExpansion::computedEnergyQpJac(unsigned int jvar)
     120             : {
     121     4374784 :   if (_dictator.notPorousFlowVariable(jvar))
     122             :     return 0.0;
     123             : 
     124     4374784 :   const unsigned int pvar = _dictator.porousFlowVariableNum(jvar);
     125     4374784 :   const unsigned nearest_qp = (_strain_at_nearest_qp ? (*_nearest_qp)[_i] : _i);
     126             : 
     127     4374784 :   Real denergy = -_dporosity_dgradvar[_i][pvar] * _grad_phi[_j][_i] * _rock_energy_nodal[_i];
     128     8979968 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
     129     4605184 :     denergy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
     130     4605184 :                (*_energy_nodal)[_i][ph] * _dporosity_dgradvar[_i][pvar] * _grad_phi[_j][nearest_qp];
     131             : 
     132     4374784 :   if (_i != _j)
     133     3563328 :     return _test[_i][_qp] * denergy * _strain_rate_qp[_qp];
     134             : 
     135      811456 :   denergy += _drock_energy_nodal_dvar[_i][pvar] * (1.0 - _porosity[_i]);
     136      811456 :   denergy -= _rock_energy_nodal[_i] * _dporosity_dvar[_i][pvar];
     137     1680512 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
     138             :   {
     139      869056 :     denergy += (*_dfluid_density_dvar)[_i][ph][pvar] * (*_fluid_saturation_nodal)[_i][ph] *
     140      869056 :                (*_energy_nodal)[_i][ph] * _porosity[_i];
     141      869056 :     denergy += (*_fluid_density)[_i][ph] * (*_dfluid_saturation_nodal_dvar)[_i][ph][pvar] *
     142      869056 :                (*_energy_nodal)[_i][ph] * _porosity[_i];
     143      869056 :     denergy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
     144      869056 :                (*_denergy_nodal_dvar)[_i][ph][pvar] * _porosity[_i];
     145      869056 :     denergy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
     146      869056 :                (*_energy_nodal)[_i][ph] * _dporosity_dvar[_i][pvar];
     147             :   }
     148             : 
     149      811456 :   return _test[_i][_qp] * denergy * _strain_rate_qp[_qp];
     150             : }