www.mooseframework.org
PorousFlowEnergyTimeDerivative.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 
11 
12 #include "MooseVariable.h"
13 
15 
16 template <>
17 InputParameters
19 {
20  InputParameters params = validParams<TimeKernel>();
21  params.addParam<bool>("strain_at_nearest_qp",
22  false,
23  "When calculating nodal porosity that depends on strain, use the strain at "
24  "the nearest quadpoint. This adds a small extra computational burden, and "
25  "is not necessary for simulations involving only linear lagrange elements. "
26  " If you set this to true, you will also want to set the same parameter to "
27  "true for related Kernels and Materials");
28  params.addRequiredParam<UserObjectName>(
29  "PorousFlowDictator", "The UserObject that holds the list of PorousFlow variable names.");
30  params.addClassDescription("Derivative of heat-energy-density wrt time");
31  return params;
32 }
33 
35  : TimeKernel(parameters),
36  _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
37  _var_is_porflow_var(_dictator.isPorousFlowVariable(_var.number())),
38  _num_phases(_dictator.numPhases()),
39  _fluid_present(_num_phases > 0),
40  _strain_at_nearest_qp(getParam<bool>("strain_at_nearest_qp")),
41  _porosity(getMaterialProperty<Real>("PorousFlow_porosity_nodal")),
42  _porosity_old(getMaterialPropertyOld<Real>("PorousFlow_porosity_nodal")),
43  _dporosity_dvar(getMaterialProperty<std::vector<Real>>("dPorousFlow_porosity_nodal_dvar")),
44  _dporosity_dgradvar(
45  getMaterialProperty<std::vector<RealGradient>>("dPorousFlow_porosity_nodal_dgradvar")),
46  _nearest_qp(_strain_at_nearest_qp
47  ? &getMaterialProperty<unsigned int>("PorousFlow_nearestqp_nodal")
48  : nullptr),
49  _rock_energy_nodal(getMaterialProperty<Real>("PorousFlow_matrix_internal_energy_nodal")),
50  _rock_energy_nodal_old(getMaterialPropertyOld<Real>("PorousFlow_matrix_internal_energy_nodal")),
51  _drock_energy_nodal_dvar(
52  getMaterialProperty<std::vector<Real>>("dPorousFlow_matrix_internal_energy_nodal_dvar")),
53  _fluid_density(_fluid_present ? &getMaterialProperty<std::vector<Real>>(
54  "PorousFlow_fluid_phase_density_nodal")
55  : nullptr),
56  _fluid_density_old(_fluid_present ? &getMaterialPropertyOld<std::vector<Real>>(
57  "PorousFlow_fluid_phase_density_nodal")
58  : nullptr),
59  _dfluid_density_dvar(_fluid_present ? &getMaterialProperty<std::vector<std::vector<Real>>>(
60  "dPorousFlow_fluid_phase_density_nodal_dvar")
61  : nullptr),
62  _fluid_saturation_nodal(
63  _fluid_present ? &getMaterialProperty<std::vector<Real>>("PorousFlow_saturation_nodal")
64  : nullptr),
65  _fluid_saturation_nodal_old(
66  _fluid_present ? &getMaterialPropertyOld<std::vector<Real>>("PorousFlow_saturation_nodal")
67  : nullptr),
68  _dfluid_saturation_nodal_dvar(_fluid_present
69  ? &getMaterialProperty<std::vector<std::vector<Real>>>(
70  "dPorousFlow_saturation_nodal_dvar")
71  : nullptr),
72  _energy_nodal(_fluid_present ? &getMaterialProperty<std::vector<Real>>(
73  "PorousFlow_fluid_phase_internal_energy_nodal")
74  : nullptr),
75  _energy_nodal_old(_fluid_present ? &getMaterialPropertyOld<std::vector<Real>>(
76  "PorousFlow_fluid_phase_internal_energy_nodal")
77  : nullptr),
78  _denergy_nodal_dvar(_fluid_present ? &getMaterialProperty<std::vector<std::vector<Real>>>(
79  "dPorousFlow_fluid_phase_internal_energy_nodal_dvar")
80  : nullptr)
81 {
82 }
83 
84 Real
86 {
87  Real energy = (1.0 - _porosity[_i]) * _rock_energy_nodal[_i];
88  Real energy_old = (1.0 - _porosity_old[_i]) * _rock_energy_nodal_old[_i];
89  for (unsigned ph = 0; ph < _num_phases; ++ph)
90  {
91  energy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
92  (*_energy_nodal)[_i][ph] * _porosity[_i];
93  energy_old += (*_fluid_density_old)[_i][ph] * (*_fluid_saturation_nodal_old)[_i][ph] *
94  (*_energy_nodal_old)[_i][ph] * _porosity_old[_i];
95  }
96 
97  return _test[_i][_qp] * (energy - energy_old) / _dt;
98 }
99 
100 Real
102 {
103  // If the variable is not a PorousFlow variable (very unusual), the diag Jacobian terms are 0
104  if (!_var_is_porflow_var)
105  return 0.0;
106  return computeQpJac(_dictator.porousFlowVariableNum(_var.number()));
107 }
108 
109 Real
111 {
112  // If the variable is not a PorousFlow variable, the OffDiag Jacobian terms are 0
114  return 0.0;
116 }
117 
118 Real
120 {
121  const unsigned nearest_qp = (_strain_at_nearest_qp ? (*_nearest_qp)[_i] : _i);
122 
123  // porosity is dependent on variables that are lumped to the nodes,
124  // but it can depend on the gradient
125  // of variables, which are NOT lumped to the nodes, hence:
126  Real denergy = -_dporosity_dgradvar[_i][pvar] * _grad_phi[_j][_i] * _rock_energy_nodal[_i];
127  for (unsigned ph = 0; ph < _num_phases; ++ph)
128  denergy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
129  (*_energy_nodal)[_i][ph] * _dporosity_dgradvar[_i][pvar] * _grad_phi[_j][nearest_qp];
130 
131  if (_i != _j)
132  return _test[_i][_qp] * denergy / _dt;
133 
134  // As the fluid energy is lumped to the nodes, only non-zero terms are for _i==_j
135  denergy += -_dporosity_dvar[_i][pvar] * _rock_energy_nodal[_i];
136  denergy += (1.0 - _porosity[_i]) * _drock_energy_nodal_dvar[_i][pvar];
137  for (unsigned ph = 0; ph < _num_phases; ++ph)
138  {
139  denergy += (*_dfluid_density_dvar)[_i][ph][pvar] * (*_fluid_saturation_nodal)[_i][ph] *
140  (*_energy_nodal)[_i][ph] * _porosity[_i];
141  denergy += (*_fluid_density)[_i][ph] * (*_dfluid_saturation_nodal_dvar)[_i][ph][pvar] *
142  (*_energy_nodal)[_i][ph] * _porosity[_i];
143  denergy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
144  (*_denergy_nodal_dvar)[_i][ph][pvar] * _porosity[_i];
145  denergy += (*_fluid_density)[_i][ph] * (*_fluid_saturation_nodal)[_i][ph] *
146  (*_energy_nodal)[_i][ph] * _dporosity_dvar[_i][pvar];
147  }
148  return _test[_i][_qp] * denergy / _dt;
149 }
Real computeQpJac(unsigned int pvar) const
Derivative of residual with respect to PorousFlow variable number pvar This is used by both computeQp...
VectorValue< Real > RealGradient
bool notPorousFlowVariable(unsigned int moose_var_num) const
Returns true if moose_var_num is not a porous flow variabe.
const MaterialProperty< std::vector< Real > > *const _fluid_density
Nodal fluid density.
const MaterialProperty< std::vector< Real > > & _drock_energy_nodal_dvar
d(nodal rock energy density)/d(PorousFlow variable)
InputParameters validParams< PorousFlowEnergyTimeDerivative >()
const PorousFlowDictator & _dictator
PorousFlowDictator UserObject.
PorousFlowEnergyTimeDerivative(const InputParameters &parameters)
const MaterialProperty< Real > & _porosity
Porosity at the nodes, but it can depend on grad(variables) which are actually evaluated at the qps...
const MaterialProperty< std::vector< RealGradient > > & _dporosity_dgradvar
d(porosity)/d(grad PorousFlow variable) - remember these derivatives will be wrt grad(vars) at qps ...
virtual Real computeQpOffDiagJacobian(unsigned int jvar) override
const MaterialProperty< Real > & _rock_energy_nodal
Nodal rock energy density.
const unsigned int _num_phases
Number of fluid phases.
Kernel = (heat_energy - heat_energy_old)/dt It is lumped to the nodes.
const MaterialProperty< Real > & _porosity_old
Old value of porosity.
This holds maps between the nonlinear variables used in a PorousFlow simulation and the variable numb...
unsigned int porousFlowVariableNum(unsigned int moose_var_num) const
The PorousFlow variable number.
registerMooseObject("PorousFlowApp", PorousFlowEnergyTimeDerivative)
const MaterialProperty< Real > & _rock_energy_nodal_old
Old value of nodal rock energy density.
const MaterialProperty< std::vector< Real > > & _dporosity_dvar
d(porosity)/d(PorousFlow variable) - these derivatives will be wrt variables at the nodes ...
const bool _var_is_porflow_var
Whether the Variable for this Kernel is a PorousFlow variable according to the Dictator.
const bool _strain_at_nearest_qp
Whether the porosity uses the volumetric strain at the closest quadpoint.