Advection of heat via flux of component k in fluid phase alpha. More...

#include <PorousFlowHeatAdvection.h>

Inheritance diagram for PorousFlowHeatAdvection:

Public Member Functions
	PorousFlowHeatAdvection (const InputParameters &parameters)

Protected Types
enum	JacRes { JacRes::CALCULATE_RESIDUAL = 0, JacRes::CALCULATE_JACOBIAN = 1 }

enum	FallbackEnum { FallbackEnum::QUICK, FallbackEnum::HARMONIC }
	If full upwinding is failing due to nodes swapping between upwind and downwind in successive nonlinear iterations (within one timestep and element) a fallback scheme is used instead. More...

Protected Member Functions
virtual Real	mobility (unsigned nodenum, unsigned phase) const override
	The mobility of the fluid. More...

virtual Real	dmobility (unsigned nodenum, unsigned phase, unsigned pvar) const override
	The derivative of mobility with respect to PorousFlow variable pvar. More...

virtual void	timestepSetup () override

virtual Real	computeQpResidual () override

virtual void	computeResidual () override

virtual void	computeJacobian () override

virtual void	computeOffDiagJacobian (MooseVariableFEBase &jvar) override

virtual Real	darcyQp (unsigned int ph) const
	The Darcy part of the flux (this is the non-upwinded part) More...

virtual Real	darcyQpJacobian (unsigned int jvar, unsigned int ph) const
	Jacobian of the Darcy part of the flux. More...

void	computeResidualAndJacobian (JacRes res_or_jac, unsigned int jvar)
	Computation of the residual and Jacobian. More...

void	fullyUpwind (JacRes res_or_jac, unsigned int ph, unsigned int pvar)
	Calculate the residual or Jacobian using full upwinding. More...

void	quickUpwind (JacRes res_or_jac, unsigned int ph, unsigned int pvar)
	Calculate the residual or Jacobian using the nodal mobilities, but without conserving fluid mass. More...

void	harmonicMean (JacRes res_or_jac, unsigned int ph, unsigned int pvar)
	Calculate the residual or Jacobian by using the harmonic mean of the nodal mobilities for the entire element. More...

Protected Attributes
const MaterialProperty< std::vector< Real > > &	_enthalpy
	Enthalpy of each phase. More...

const MaterialProperty< std::vector< std::vector< Real > > > &	_denthalpy_dvar
	Derivative of the enthalpy wrt PorousFlow variables. More...

const MaterialProperty< std::vector< Real > > &	_relative_permeability
	Relative permeability of each phase. More...

const MaterialProperty< std::vector< std::vector< Real > > > &	_drelative_permeability_dvar
	Derivative of relative permeability of each phase wrt PorousFlow variables. More...

const MaterialProperty< RealTensorValue > &	_permeability
	Permeability of porous material. More...

const MaterialProperty< std::vector< RealTensorValue > > &	_dpermeability_dvar
	d(permeabiity)/d(PorousFlow variable) More...

const MaterialProperty< std::vector< std::vector< RealTensorValue > > > &	_dpermeability_dgradvar
	d(permeabiity)/d(grad(PorousFlow variable)) More...

const MaterialProperty< std::vector< Real > > &	_fluid_density_node
	Fluid density for each phase (at the node) More...

const MaterialProperty< std::vector< std::vector< Real > > > &	_dfluid_density_node_dvar
	Derivative of the fluid density for each phase wrt PorousFlow variables (at the node) More...

const MaterialProperty< std::vector< Real > > &	_fluid_density_qp
	Fluid density for each phase (at the qp) More...

const MaterialProperty< std::vector< std::vector< Real > > > &	_dfluid_density_qp_dvar
	Derivative of the fluid density for each phase wrt PorousFlow variables (at the qp) More...

const MaterialProperty< std::vector< Real > > &	_fluid_viscosity
	Viscosity of each component in each phase. More...

const MaterialProperty< std::vector< std::vector< Real > > > &	_dfluid_viscosity_dvar
	Derivative of the fluid viscosity for each phase wrt PorousFlow variables. More...

const MaterialProperty< std::vector< Real > > &	_pp
	Nodal pore pressure in each phase. More...

const MaterialProperty< std::vector< RealGradient > > &	_grad_p
	Gradient of the pore pressure in each phase. More...

const MaterialProperty< std::vector< std::vector< Real > > > &	_dgrad_p_dgrad_var
	Derivative of Grad porepressure in each phase wrt grad(PorousFlow variables) More...

const MaterialProperty< std::vector< std::vector< RealGradient > > > &	_dgrad_p_dvar
	Derivative of Grad porepressure in each phase wrt PorousFlow variables. More...

const PorousFlowDictator &	_dictator
	PorousFlowDictator UserObject. More...

const unsigned int	_num_phases
	The number of fluid phases. More...

const RealVectorValue	_gravity
	Gravity. Defaults to 9.81 m/s^2. More...

const bool	_perm_derivs
	Flag to check whether permeabiity derivatives are non-zero. More...

const unsigned	_full_upwind_threshold
	If the number of upwind-downwind swaps is less than this amount then full upwinding is used. More...

enum PorousFlowDarcyBase::FallbackEnum	_fallback_scheme

std::vector< std::vector< Real > >	_proto_flux
	The Darcy flux. More...

std::vector< std::vector< std::vector< Real > > >	_jacobian
	Derivative of _proto_flux with respect to nodal variables. More...

std::unordered_map< unsigned, std::vector< std::vector< unsigned > > >	_num_upwinds
	Number of nonlinear iterations (in this timestep and this element) that a node is an upwind node for a given fluid phase. More...

std::unordered_map< unsigned, std::vector< std::vector< unsigned > > >	_num_downwinds
	Number of nonlinear iterations (in this timestep and this element) that a node is an downwind node for a given fluid phase. More...

Detailed Description

Advection of heat via flux of component k in fluid phase alpha.

A fully-updwinded version is implemented, where the mobility of the upstream nodes is used.

Definition at line 24 of file PorousFlowHeatAdvection.h.

Member Enumeration Documentation

◆ FallbackEnum

enum PorousFlowDarcyBase::FallbackEnum

strongprotectedinherited

If full upwinding is failing due to nodes swapping between upwind and downwind in successive nonlinear iterations (within one timestep and element) a fallback scheme is used instead.

QUICK: use the nodal mobilities, as in full upwinding, but do not attempt to conserve fluid mass. HARMONIC: assign the harmonic mean of the mobilities to the entire element.

Enumerator
QUICK
HARMONIC

Definition at line 152 of file PorousFlowDarcyBase.h.

152 { QUICK, HARMONIC } _fallback_scheme;

◆ JacRes

enum PorousFlowDarcyBase::JacRes

strongprotectedinherited

Enumerator
CALCULATE_RESIDUAL
CALCULATE_JACOBIAN

Definition at line 63 of file PorousFlowDarcyBase.h.

   {
     CALCULATE_RESIDUAL = 0,
     CALCULATE_JACOBIAN = 1
   };

Constructor & Destructor Documentation

◆ PorousFlowHeatAdvection()

PorousFlowHeatAdvection::PorousFlowHeatAdvection ( const InputParameters & parameters )

Definition at line 23 of file PorousFlowHeatAdvection.C.

   : PorousFlowDarcyBase(parameters),
     _enthalpy(getMaterialProperty<std::vector<Real>>("PorousFlow_fluid_phase_enthalpy_nodal")),
     _denthalpy_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
         "dPorousFlow_fluid_phase_enthalpy_nodal_dvar")),
     _relative_permeability(
         getMaterialProperty<std::vector<Real>>("PorousFlow_relative_permeability_nodal")),
     _drelative_permeability_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
         "dPorousFlow_relative_permeability_nodal_dvar"))
 {
 }

Member Function Documentation

◆ computeJacobian()

void PorousFlowDarcyBase::computeJacobian ( )

overrideprotectedvirtualinherited

Definition at line 141 of file PorousFlowDarcyBase.C.

 {
   computeResidualAndJacobian(JacRes::CALCULATE_JACOBIAN, _var.number());
 }

◆ computeOffDiagJacobian()

void PorousFlowDarcyBase::computeOffDiagJacobian ( MooseVariableFEBase & jvar )

overrideprotectedvirtualinherited

Definition at line 147 of file PorousFlowDarcyBase.C.

 {
   computeResidualAndJacobian(JacRes::CALCULATE_JACOBIAN, jvar.number());
 }

◆ computeQpResidual()

Real PorousFlowDarcyBase::computeQpResidual ( )

overrideprotectedvirtualinherited

Definition at line 128 of file PorousFlowDarcyBase.C.

 {
   mooseError("PorousFlowDarcyBase: computeQpResidual called");
   return 0.0;
 }

◆ computeResidual()

void PorousFlowDarcyBase::computeResidual ( )

overrideprotectedvirtualinherited

Definition at line 135 of file PorousFlowDarcyBase.C.

 {
   computeResidualAndJacobian(JacRes::CALCULATE_RESIDUAL, 0.0);
 }

◆ computeResidualAndJacobian()

void PorousFlowDarcyBase::computeResidualAndJacobian	(	JacRes	res_or_jac,
		unsigned int	jvar
	)

protectedinherited

Computation of the residual and Jacobian.

If res_or_jac=CALCULATE_JACOBIAN then the residual gets calculated anyway (becuase we have to know which nodes are upwind and which are downwind) except we don't actually need to store the residual in _assembly.residualBlock and we don't need to save in save_in.

If res_or_jac=CALCULATE_RESIDUAL then we don't have to do all the Jacobian calculations.

Parameters

res_or_jac	Whether to calculate the residual or the jacobian.
jvar	PorousFlow variable to take derivatives wrt to

Definition at line 153 of file PorousFlowDarcyBase.C.

 {
   if ((res_or_jac == JacRes::CALCULATE_JACOBIAN) && _dictator.notPorousFlowVariable(jvar))
     return;
  
   // The PorousFlow variable index corresponding to the variable number jvar
   const unsigned int pvar =
       ((res_or_jac == JacRes::CALCULATE_JACOBIAN) ? _dictator.porousFlowVariableNum(jvar) : 0);
  
   DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
   if ((ke.n() == 0) && (res_or_jac == JacRes::CALCULATE_JACOBIAN)) // this removes a problem
                                                                    // encountered in the initial
                                                                    // timestep when
                                                                    // use_displaced_mesh=true
     return;
  
   // The number of nodes in the element
   const unsigned int num_nodes = _test.size();
  
   // Compute the residual and jacobian without the mobility terms. Even if we are computing the
   // Jacobian we still need this in order to see which nodes are upwind and which are downwind.
   for (unsigned ph = 0; ph < _num_phases; ++ph)
   {
     _proto_flux[ph].assign(num_nodes, 0);
     for (_qp = 0; _qp < _qrule->n_points(); _qp++)
     {
       for (_i = 0; _i < num_nodes; ++_i)
         _proto_flux[ph][_i] += _JxW[_qp] * _coord[_qp] * darcyQp(ph);
     }
   }
  
   // for this element, record whether each node is "upwind" or "downwind" (or neither)
   const unsigned elem = _current_elem->id();
   if (_num_upwinds.find(elem) == _num_upwinds.end())
   {
     _num_upwinds[elem] = std::vector<std::vector<unsigned>>(_num_phases);
     _num_downwinds[elem] = std::vector<std::vector<unsigned>>(_num_phases);
     for (unsigned ph = 0; ph < _num_phases; ++ph)
     {
       _num_upwinds[elem][ph].assign(num_nodes, 0);
       _num_downwinds[elem][ph].assign(num_nodes, 0);
     }
   }
   // record the information once per nonlinear iteration
   if (res_or_jac == JacRes::CALCULATE_JACOBIAN && jvar == _var.number())
   {
     for (unsigned ph = 0; ph < _num_phases; ++ph)
     {
       for (unsigned nod = 0; nod < num_nodes; ++nod)
       {
         if (_proto_flux[ph][nod] > 0)
           _num_upwinds[elem][ph][nod]++;
         else if (_proto_flux[ph][nod] < 0)
           _num_downwinds[elem][ph][nod]++;
       }
     }
   }
  
   // based on _num_upwinds and _num_downwinds, calculate the maximum number
   // of upwind-downwind swaps that have been encountered in this timestep
   // for this element
   std::vector<unsigned> max_swaps(_num_phases, 0);
   for (unsigned ph = 0; ph < _num_phases; ++ph)
   {
     for (unsigned nod = 0; nod < num_nodes; ++nod)
       max_swaps[ph] = std::max(
           max_swaps[ph], std::min(_num_upwinds[elem][ph][nod], _num_downwinds[elem][ph][nod]));
   }
  
   // size the _jacobian correctly and calculate it for the case residual = _proto_flux
   if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
   {
     for (unsigned ph = 0; ph < _num_phases; ++ph)
     {
       _jacobian[ph].resize(ke.m());
       for (_i = 0; _i < _test.size(); _i++)
       {
         _jacobian[ph][_i].assign(ke.n(), 0.0);
         for (_j = 0; _j < _phi.size(); _j++)
           for (_qp = 0; _qp < _qrule->n_points(); _qp++)
             _jacobian[ph][_i][_j] += _JxW[_qp] * _coord[_qp] * darcyQpJacobian(jvar, ph);
       }
     }
   }
  
   // Loop over all the phases, computing the mass flux, which
   // gets placed into _proto_flux, and the derivative of this
   // which gets placed into _jacobian
   for (unsigned int ph = 0; ph < _num_phases; ++ph)
   {
     if (max_swaps[ph] < _full_upwind_threshold)
       fullyUpwind(res_or_jac, ph, pvar);
     else
     {
       switch (_fallback_scheme)
       {
         case FallbackEnum::QUICK:
           quickUpwind(res_or_jac, ph, pvar);
           break;
         case FallbackEnum::HARMONIC:
           harmonicMean(res_or_jac, ph, pvar);
           break;
       }
     }
   }
  
   // Add results to the Residual or Jacobian
   if (res_or_jac == JacRes::CALCULATE_RESIDUAL)
   {
     DenseVector<Number> & re = _assembly.residualBlock(_var.number());
  
     _local_re.resize(re.size());
     _local_re.zero();
     for (_i = 0; _i < _test.size(); _i++)
       for (unsigned int ph = 0; ph < _num_phases; ++ph)
         _local_re(_i) += _proto_flux[ph][_i];
  
     re += _local_re;
  
     if (_has_save_in)
     {
       Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
       for (unsigned int i = 0; i < _save_in.size(); i++)
         _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
     }
   }
  
   if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
   {
     _local_ke.resize(ke.m(), ke.n());
     _local_ke.zero();
  
     for (_i = 0; _i < _test.size(); _i++)
       for (_j = 0; _j < _phi.size(); _j++)
         for (unsigned int ph = 0; ph < _num_phases; ++ph)
           _local_ke(_i, _j) += _jacobian[ph][_i][_j];
  
     ke += _local_ke;
  
     if (_has_diag_save_in && jvar == _var.number())
     {
       unsigned int rows = ke.m();
       DenseVector<Number> diag(rows);
       for (unsigned int i = 0; i < rows; i++)
         diag(i) = _local_ke(i, i);
  
       Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
       for (unsigned int i = 0; i < _diag_save_in.size(); i++)
         _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
     }
   }
 }

Referenced by PorousFlowDarcyBase::computeJacobian(), PorousFlowDarcyBase::computeOffDiagJacobian(), and PorousFlowDarcyBase::computeResidual().

◆ darcyQp()

Real PorousFlowDarcyBase::darcyQp ( unsigned int ph ) const

protectedvirtualinherited

The Darcy part of the flux (this is the non-upwinded part)

Definition at line 95 of file PorousFlowDarcyBase.C.

 {
   return _grad_test[_i][_qp] *
          (_permeability[_qp] * (_grad_p[_qp][ph] - _fluid_density_qp[_qp][ph] * _gravity));
 }

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ darcyQpJacobian()

Real PorousFlowDarcyBase::darcyQpJacobian	(	unsigned int	jvar,
		unsigned int	ph
	)		const

protectedvirtualinherited

Jacobian of the Darcy part of the flux.

Definition at line 102 of file PorousFlowDarcyBase.C.

 {
   if (_dictator.notPorousFlowVariable(jvar))
     return 0.0;
  
   const unsigned int pvar = _dictator.porousFlowVariableNum(jvar);
  
   RealVectorValue deriv =
       _permeability[_qp] * (_grad_phi[_j][_qp] * _dgrad_p_dgrad_var[_qp][ph][pvar] -
                             _phi[_j][_qp] * _dfluid_density_qp_dvar[_qp][ph][pvar] * _gravity);
  
   deriv += _permeability[_qp] * (_dgrad_p_dvar[_qp][ph][pvar] * _phi[_j][_qp]);
  
   if (_perm_derivs)
   {
     deriv += _dpermeability_dvar[_qp][pvar] * _phi[_j][_qp] *
              (_grad_p[_qp][ph] - _fluid_density_qp[_qp][ph] * _gravity);
     for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
       deriv += _dpermeability_dgradvar[_qp][i][pvar] * _grad_phi[_j][_qp](i) *
                (_grad_p[_qp][ph] - _fluid_density_qp[_qp][ph] * _gravity);
   }
  
   return _grad_test[_i][_qp] * deriv;
 }

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ dmobility()

Real PorousFlowHeatAdvection::dmobility	(	unsigned	nodenum,
		unsigned	phase,
		unsigned	pvar
	)		const

overrideprotectedvirtual

The derivative of mobility with respect to PorousFlow variable pvar.

Parameters

nodenum	The node-number to evaluate the mobility for
phase	the fluid phase number
pvar	the PorousFlow variable pvar

Reimplemented from PorousFlowDarcyBase.

Definition at line 43 of file PorousFlowHeatAdvection.C.

 {
   Real dm = _denthalpy_dvar[nodenum][phase][pvar] * _fluid_density_node[nodenum][phase] *
             _relative_permeability[nodenum][phase] / _fluid_viscosity[nodenum][phase];
   dm += _enthalpy[nodenum][phase] * _dfluid_density_node_dvar[nodenum][phase][pvar] *
         _relative_permeability[nodenum][phase] / _fluid_viscosity[nodenum][phase];
   dm += _enthalpy[nodenum][phase] * _fluid_density_node[nodenum][phase] *
         _drelative_permeability_dvar[nodenum][phase][pvar] / _fluid_viscosity[nodenum][phase];
   dm -= _enthalpy[nodenum][phase] * _fluid_density_node[nodenum][phase] *
         _relative_permeability[nodenum][phase] * _dfluid_viscosity_dvar[nodenum][phase][pvar] /
         std::pow(_fluid_viscosity[nodenum][phase], 2);
   return dm;
 }

◆ fullyUpwind()

void PorousFlowDarcyBase::fullyUpwind	(	JacRes	res_or_jac,
		unsigned int	ph,
		unsigned int	pvar
	)

protectedinherited

Calculate the residual or Jacobian using full upwinding.

Parameters

res_or_jac	whether to compute the residual or jacobian
ph	fluid phase number
pvar	differentiate wrt to this PorousFlow variable (when computing the jacobian)

Perform the full upwinding by multiplying the residuals at the upstream nodes by their mobilities. Mobility is different for each phase, and in each situation: mobility = density / viscosity for single-component Darcy flow mobility = mass_fraction * density * relative_perm / viscosity for multi-component, multiphase flow mobility = enthalpy * density * relative_perm / viscosity for heat convection

The residual for the kernel is the sum over Darcy fluxes for each phase. The Darcy flux for a particular phase is R_i = int{mobility*flux_no_mob} = int{mobility*grad(pot)*permeability*grad(test_i)} for node i. where int is the integral over the element. However, in fully-upwind, the first step is to take the mobility outside the integral, which was done in the _proto_flux calculation above.

NOTE: Physically _proto_flux[i][ph] is a measure of fluid of phase ph flowing out of node i. If we had left in mobility, it would be exactly the component mass flux flowing out of node i.

This leads to the definition of upwinding:

If _proto_flux(i)[ph] is positive then we use mobility_i. That is we use the upwind value of mobility.

The final subtle thing is we must also conserve fluid mass: the total component mass flowing out of node i must be the sum of the masses flowing into the other nodes.

note the -= means the result is positive

Definition at line 307 of file PorousFlowDarcyBase.C.

 {
   // The number of nodes in the element
   const unsigned int num_nodes = _test.size();
  
   Real mob;
   Real dmob;
   // Define variables used to ensure mass conservation
   Real total_mass_out = 0.0;
   Real total_in = 0.0;
  
   // The following holds derivatives of these
   std::vector<Real> dtotal_mass_out;
   std::vector<Real> dtotal_in;
   if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
   {
     dtotal_mass_out.assign(num_nodes, 0.0);
     dtotal_in.assign(num_nodes, 0.0);
   }
  
   // Perform the upwinding using the mobility
   std::vector<bool> upwind_node(num_nodes);
   for (unsigned int n = 0; n < num_nodes; ++n)
   {
     if (_proto_flux[ph][n] >= 0.0) // upstream node
     {
       upwind_node[n] = true;
       // The mobility at the upstream node
       mob = mobility(n, ph);
       if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
       {
         // The derivative of the mobility wrt the PorousFlow variable
         dmob = dmobility(n, ph, pvar);
  
         for (_j = 0; _j < _phi.size(); _j++)
           _jacobian[ph][n][_j] *= mob;
  
         if (_test.size() == _phi.size())
           /* mobility at node=n depends only on the variables at node=n, by construction.  For
            * linear-lagrange variables, this means that Jacobian entries involving the derivative
            * of mobility will only be nonzero for derivatives wrt variables at node=n.  Hence the
            * [n][n] in the line below.  However, for other variable types (eg constant monomials)
            * I cannot tell what variable number contributes to the derivative.  However, in all
            * cases I can possibly imagine, the derivative is zero anyway, since in the full
            * upwinding scheme, mobility shouldn't depend on these other sorts of variables.
            */
           _jacobian[ph][n][n] += dmob * _proto_flux[ph][n];
  
         for (_j = 0; _j < _phi.size(); _j++)
           dtotal_mass_out[_j] += _jacobian[ph][n][_j];
       }
       _proto_flux[ph][n] *= mob;
       total_mass_out += _proto_flux[ph][n];
     }
     else
     {
       upwind_node[n] = false;
       total_in -= _proto_flux[ph][n]; 
       if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
         for (_j = 0; _j < _phi.size(); _j++)
           dtotal_in[_j] -= _jacobian[ph][n][_j];
     }
   }
  
   // Conserve mass over all phases by proportioning the total_mass_out mass to the inflow nodes,
   // weighted by their proto_flux values
   for (unsigned int n = 0; n < num_nodes; ++n)
   {
     if (!upwind_node[n]) // downstream node
     {
       if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
         for (_j = 0; _j < _phi.size(); _j++)
         {
           _jacobian[ph][n][_j] *= total_mass_out / total_in;
           _jacobian[ph][n][_j] +=
               _proto_flux[ph][n] * (dtotal_mass_out[_j] / total_in -
                                     dtotal_in[_j] * total_mass_out / total_in / total_in);
         }
       _proto_flux[ph][n] *= total_mass_out / total_in;
     }
   }
 }

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ harmonicMean()

void PorousFlowDarcyBase::harmonicMean	(	JacRes	res_or_jac,
		unsigned int	ph,
		unsigned int	pvar
	)

protectedinherited

Calculate the residual or Jacobian by using the harmonic mean of the nodal mobilities for the entire element.

Parameters

res_or_jac	whether to compute the residual or jacobian
ph	fluid phase number
pvar	differentiate wrt to this PorousFlow variable (when computing the jacobian)

Definition at line 457 of file PorousFlowDarcyBase.C.

 {
   // The number of nodes in the element
   const unsigned int num_nodes = _test.size();
  
   std::vector<Real> mob(num_nodes);
   unsigned num_zero = 0;
   unsigned zero_mobility_node = std::numeric_limits<unsigned>::max();
   Real harmonic_mob = 0;
   for (unsigned n = 0; n < num_nodes; ++n)
   {
     mob[n] = mobility(n, ph);
     if (mob[n] == 0.0)
     {
       zero_mobility_node = n;
       num_zero++;
     }
     else
       harmonic_mob += 1.0 / mob[n];
   }
   if (num_zero > 0)
     harmonic_mob = 0.0;
   else
     harmonic_mob = (1.0 * num_nodes) / harmonic_mob;
  
   // d(harmonic_mob)/d(PorousFlow variable at node n)
   std::vector<Real> dharmonic_mob(num_nodes, 0.0);
   if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
   {
     const Real harm2 = std::pow(harmonic_mob, 2) / (1.0 * num_nodes);
     if (num_zero == 0)
       for (unsigned n = 0; n < num_nodes; ++n)
         dharmonic_mob[n] = dmobility(n, ph, pvar) * harm2 / std::pow(mob[n], 2);
     else if (num_zero == 1)
       dharmonic_mob[zero_mobility_node] =
           num_nodes * dmobility(zero_mobility_node, ph, pvar); // other derivs are zero
     // if num_zero > 1 then all dharmonic_mob = 0.0
   }
  
   if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
     for (unsigned n = 0; n < num_nodes; ++n)
       for (_j = 0; _j < _phi.size(); _j++)
       {
         _jacobian[ph][n][_j] *= harmonic_mob;
         if (_test.size() == _phi.size())
           _jacobian[ph][n][_j] += dharmonic_mob[_j] * _proto_flux[ph][n];
       }
  
   if (res_or_jac == JacRes::CALCULATE_RESIDUAL)
     for (unsigned n = 0; n < num_nodes; ++n)
       _proto_flux[ph][n] *= harmonic_mob;
 }

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ mobility()

Real PorousFlowHeatAdvection::mobility	(	unsigned	nodenum,
		unsigned	phase
	)		const

overrideprotectedvirtual

The mobility of the fluid.

For multi-component Darcy flow this is mass_fraction * fluid_density * relative_permeability / fluid_viscosity

Parameters

nodenum	The node-number to evaluate the mobility for
phase	the fluid phase number

Reimplemented from PorousFlowDarcyBase.

Definition at line 36 of file PorousFlowHeatAdvection.C.

 {
   return _enthalpy[nodenum][phase] * _fluid_density_node[nodenum][phase] *
          _relative_permeability[nodenum][phase] / _fluid_viscosity[nodenum][phase];
 }

◆ quickUpwind()

void PorousFlowDarcyBase::quickUpwind	(	JacRes	res_or_jac,
		unsigned int	ph,
		unsigned int	pvar
	)

protectedinherited

Calculate the residual or Jacobian using the nodal mobilities, but without conserving fluid mass.

Parameters

res_or_jac	whether to compute the residual or jacobian
ph	fluid phase number
pvar	differentiate wrt to this PorousFlow variable (when computing the jacobian)

Definition at line 420 of file PorousFlowDarcyBase.C.

 {
   // The number of nodes in the element
   const unsigned int num_nodes = _test.size();
  
   Real mob;
   Real dmob;
  
   // Use the raw nodal mobility
   for (unsigned int n = 0; n < num_nodes; ++n)
   {
     // The mobility at the node
     mob = mobility(n, ph);
     if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
     {
       // The derivative of the mobility wrt the PorousFlow variable
       dmob = dmobility(n, ph, pvar);
  
       for (_j = 0; _j < _phi.size(); _j++)
         _jacobian[ph][n][_j] *= mob;
  
       if (_test.size() == _phi.size())
         /* mobility at node=n depends only on the variables at node=n, by construction.  For
          * linear-lagrange variables, this means that Jacobian entries involving the derivative
          * of mobility will only be nonzero for derivatives wrt variables at node=n.  Hence the
          * [n][n] in the line below.  However, for other variable types (eg constant monomials)
          * I cannot tell what variable number contributes to the derivative.  However, in all
          * cases I can possibly imagine, the derivative is zero anyway, since in the full
          * upwinding scheme, mobility shouldn't depend on these other sorts of variables.
          */
         _jacobian[ph][n][n] += dmob * _proto_flux[ph][n];
     }
     _proto_flux[ph][n] *= mob;
   }
 }

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ timestepSetup()

void PorousFlowDarcyBase::timestepSetup ( )

overrideprotectedvirtualinherited

Definition at line 87 of file PorousFlowDarcyBase.C.

 {
   Kernel::timestepSetup();
   _num_upwinds = std::unordered_map<unsigned, std::vector<std::vector<unsigned>>>();
   _num_downwinds = std::unordered_map<unsigned, std::vector<std::vector<unsigned>>>();
 }

Member Data Documentation

◆ _denthalpy_dvar

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowHeatAdvection::_denthalpy_dvar

protected

Derivative of the enthalpy wrt PorousFlow variables.

Definition at line 37 of file PorousFlowHeatAdvection.h.

Referenced by dmobility().

◆ _dfluid_density_node_dvar

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowDarcyBase::_dfluid_density_node_dvar

protectedinherited

Derivative of the fluid density for each phase wrt PorousFlow variables (at the node)

Definition at line 99 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowAdvectiveFlux::dmobility(), dmobility(), and PorousFlowDarcyBase::dmobility().

◆ _dfluid_density_qp_dvar

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowDarcyBase::_dfluid_density_qp_dvar

protectedinherited

Derivative of the fluid density for each phase wrt PorousFlow variables (at the qp)

Definition at line 105 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQpJacobian().

◆ _dfluid_viscosity_dvar

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowDarcyBase::_dfluid_viscosity_dvar

protectedinherited

Derivative of the fluid viscosity for each phase wrt PorousFlow variables.

Definition at line 111 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowAdvectiveFlux::dmobility(), dmobility(), and PorousFlowDarcyBase::dmobility().

◆ _dgrad_p_dgrad_var

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowDarcyBase::_dgrad_p_dgrad_var

protectedinherited

Derivative of Grad porepressure in each phase wrt grad(PorousFlow variables)

Definition at line 120 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQpJacobian().

◆ _dgrad_p_dvar

const MaterialProperty<std::vector<std::vector<RealGradient> > >& PorousFlowDarcyBase::_dgrad_p_dvar

protectedinherited

Derivative of Grad porepressure in each phase wrt PorousFlow variables.

Definition at line 123 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQpJacobian().

◆ _dictator

const PorousFlowDictator& PorousFlowDarcyBase::_dictator

protectedinherited

PorousFlowDictator UserObject.

Definition at line 126 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian(), and PorousFlowDarcyBase::darcyQpJacobian().

◆ _dpermeability_dgradvar

const MaterialProperty<std::vector<std::vector<RealTensorValue> > >& PorousFlowDarcyBase::_dpermeability_dgradvar

protectedinherited

d(permeabiity)/d(grad(PorousFlow variable))

Definition at line 93 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQpJacobian().

◆ _dpermeability_dvar

const MaterialProperty<std::vector<RealTensorValue> >& PorousFlowDarcyBase::_dpermeability_dvar

protectedinherited

d(permeabiity)/d(PorousFlow variable)

Definition at line 90 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQpJacobian().

◆ _drelative_permeability_dvar

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowHeatAdvection::_drelative_permeability_dvar

protected

Derivative of relative permeability of each phase wrt PorousFlow variables.

Definition at line 43 of file PorousFlowHeatAdvection.h.

Referenced by dmobility().

◆ _enthalpy

const MaterialProperty<std::vector<Real> >& PorousFlowHeatAdvection::_enthalpy

protected

Enthalpy of each phase.

Definition at line 34 of file PorousFlowHeatAdvection.h.

Referenced by dmobility(), and mobility().

◆ _fallback_scheme

enum PorousFlowDarcyBase::FallbackEnum PorousFlowDarcyBase::_fallback_scheme

protectedinherited

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ _fluid_density_node

const MaterialProperty<std::vector<Real> >& PorousFlowDarcyBase::_fluid_density_node

protectedinherited

Fluid density for each phase (at the node)

Definition at line 96 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowAdvectiveFlux::dmobility(), dmobility(), PorousFlowDarcyBase::dmobility(), PorousFlowAdvectiveFlux::mobility(), mobility(), and PorousFlowDarcyBase::mobility().

◆ _fluid_density_qp

const MaterialProperty<std::vector<Real> >& PorousFlowDarcyBase::_fluid_density_qp

protectedinherited

Fluid density for each phase (at the qp)

Definition at line 102 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQp(), and PorousFlowDarcyBase::darcyQpJacobian().

◆ _fluid_viscosity

const MaterialProperty<std::vector<Real> >& PorousFlowDarcyBase::_fluid_viscosity

protectedinherited

Viscosity of each component in each phase.

Definition at line 108 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowAdvectiveFlux::dmobility(), dmobility(), PorousFlowDarcyBase::dmobility(), PorousFlowAdvectiveFlux::mobility(), mobility(), and PorousFlowDarcyBase::mobility().

◆ _full_upwind_threshold

const unsigned PorousFlowDarcyBase::_full_upwind_threshold

protectedinherited

If the number of upwind-downwind swaps is less than this amount then full upwinding is used.

Otherwise the fallback scheme is employed

Definition at line 141 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ _grad_p

const MaterialProperty<std::vector<RealGradient> >& PorousFlowDarcyBase::_grad_p

protectedinherited

Gradient of the pore pressure in each phase.

Definition at line 117 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQp(), and PorousFlowDarcyBase::darcyQpJacobian().

◆ _gravity

const RealVectorValue PorousFlowDarcyBase::_gravity

protectedinherited

Gravity. Defaults to 9.81 m/s^2.

Definition at line 132 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQp(), and PorousFlowDarcyBase::darcyQpJacobian().

◆ _jacobian

std::vector<std::vector<std::vector<Real> > > PorousFlowDarcyBase::_jacobian

protectedinherited

Derivative of _proto_flux with respect to nodal variables.

This gets modified with the mobility and its derivative during computation of MOOSE's Jacobian

Definition at line 167 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian(), PorousFlowDarcyBase::fullyUpwind(), PorousFlowDarcyBase::harmonicMean(), and PorousFlowDarcyBase::quickUpwind().

◆ _num_downwinds

std::unordered_map<unsigned, std::vector<std::vector<unsigned> > > PorousFlowDarcyBase::_num_downwinds

protectedinherited

Number of nonlinear iterations (in this timestep and this element) that a node is an downwind node for a given fluid phase.

_num_upwinds[element_number][phase][node_number_in_element]

Definition at line 181 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian(), and PorousFlowDarcyBase::timestepSetup().

◆ _num_phases

const unsigned int PorousFlowDarcyBase::_num_phases

protectedinherited

The number of fluid phases.

Definition at line 129 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian().

◆ _num_upwinds

std::unordered_map<unsigned, std::vector<std::vector<unsigned> > > PorousFlowDarcyBase::_num_upwinds

protectedinherited

Number of nonlinear iterations (in this timestep and this element) that a node is an upwind node for a given fluid phase.

_num_upwinds[element_number][phase][node_number_in_element]

Definition at line 174 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian(), and PorousFlowDarcyBase::timestepSetup().

◆ _perm_derivs

const bool PorousFlowDarcyBase::_perm_derivs

protectedinherited

Flag to check whether permeabiity derivatives are non-zero.

Definition at line 135 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQpJacobian().

◆ _permeability

const MaterialProperty<RealTensorValue>& PorousFlowDarcyBase::_permeability

protectedinherited

Permeability of porous material.

Definition at line 87 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::darcyQp(), and PorousFlowDarcyBase::darcyQpJacobian().

◆ _pp

const MaterialProperty<std::vector<Real> >& PorousFlowDarcyBase::_pp

protectedinherited

Nodal pore pressure in each phase.

Definition at line 114 of file PorousFlowDarcyBase.h.

◆ _proto_flux

std::vector<std::vector<Real> > PorousFlowDarcyBase::_proto_flux

protectedinherited

The Darcy flux.

When multiplied by the mobility, this is the mass flux of fluid moving out of a node. This multiplication occurs during the formation of the residual and Jacobian. _proto_flux[num_phases][num_nodes]

Definition at line 160 of file PorousFlowDarcyBase.h.

Referenced by PorousFlowDarcyBase::computeResidualAndJacobian(), PorousFlowDarcyBase::fullyUpwind(), PorousFlowDarcyBase::harmonicMean(), and PorousFlowDarcyBase::quickUpwind().

◆ _relative_permeability

const MaterialProperty<std::vector<Real> >& PorousFlowHeatAdvection::_relative_permeability

protected

Relative permeability of each phase.

Definition at line 40 of file PorousFlowHeatAdvection.h.

Referenced by dmobility(), and mobility().

The documentation for this class was generated from the following files:

porous_flow/include/kernels/PorousFlowHeatAdvection.h
porous_flow/src/kernels/PorousFlowHeatAdvection.C

Public Member Functions

Protected Types

Protected Member Functions

Protected Attributes

Detailed Description

Member Enumeration Documentation

◆ FallbackEnum

◆ JacRes

Constructor & Destructor Documentation

◆ PorousFlowHeatAdvection()

Member Function Documentation

◆ computeJacobian()

◆ computeOffDiagJacobian()

◆ computeQpResidual()

◆ computeResidual()

◆ computeResidualAndJacobian()

◆ darcyQp()

◆ darcyQpJacobian()

◆ dmobility()

◆ fullyUpwind()

◆ harmonicMean()

◆ mobility()

◆ quickUpwind()

◆ timestepSetup()

Member Data Documentation

◆ _denthalpy_dvar

◆ _dfluid_density_node_dvar

◆ _dfluid_density_qp_dvar

◆ _dfluid_viscosity_dvar

◆ _dgrad_p_dgrad_var

◆ _dgrad_p_dvar

◆ _dictator

◆ _dpermeability_dgradvar

◆ _dpermeability_dvar

◆ _drelative_permeability_dvar

◆ _enthalpy

◆ _fallback_scheme

◆ _fluid_density_node

◆ _fluid_density_qp

◆ _fluid_viscosity

◆ _full_upwind_threshold

◆ _grad_p

◆ _gravity

◆ _jacobian

◆ _num_downwinds

◆ _num_phases

◆ _num_upwinds

◆ _perm_derivs

◆ _permeability

◆ _pp

◆ _proto_flux

◆ _relative_permeability