www.mooseframework.org
Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | List of all members
PorousFlowAdvectiveFluxCalculatorBase Class Referenceabstract

Base class to compute the advective flux of fluid in PorousFlow situations using the Kuzmin-Turek FEM-TVD multidimensional stabilization scheme. More...

#include <PorousFlowAdvectiveFluxCalculatorBase.h>

Inheritance diagram for PorousFlowAdvectiveFluxCalculatorBase:
[legend]

Public Member Functions

 PorousFlowAdvectiveFluxCalculatorBase (const InputParameters &parameters)
 
const std::map< dof_id_type, std::vector< Real > > & getdFluxOut_dvars (unsigned node_id) const
 Returns d(flux_out)/d(porous_flow_variables. More...
 
virtual void meshChanged () override
 
Real getFluxOut (dof_id_type node_i) const
 Returns the flux out of lobal node id. More...
 
const std::map< dof_id_type, Real > & getdFluxOutdu (dof_id_type node_i) const
 Returns r where r[j] = d(flux out of global node i)/du(global node j) used in Jacobian computations. More...
 
const std::vector< std::vector< Real > > & getdFluxOutdKjk (dof_id_type node_i) const
 Returns r where r[j][k] = d(flux out of global node i)/dK[connected node j][connected node k] used in Jacobian computations. More...
 
unsigned getValence (dof_id_type node_i) const
 Returns the valence of the global node i Valence is the number of times the node is encountered in a loop over elements (that have appropriate subdomain_id, if the user has employed the "blocks=" parameter) seen by this processor (including ghosted elements) More...
 

Protected Types

enum  FluxLimiterTypeEnum {
  FluxLimiterTypeEnum::MinMod, FluxLimiterTypeEnum::VanLeer, FluxLimiterTypeEnum::MC, FluxLimiterTypeEnum::superbee,
  FluxLimiterTypeEnum::None
}
 Determines Flux Limiter type (Page 135 of Kuzmin and Turek) "None" means that limitFlux=0 always, which implies zero antidiffusion will be added. More...
 
enum  PQPlusMinusEnum { PQPlusMinusEnum::PPlus, PQPlusMinusEnum::PMinus, PQPlusMinusEnum::QPlus, PQPlusMinusEnum::QMinus }
 Signals to the PQPlusMinus method what should be computed. More...
 

Protected Member Functions

virtual void timestepSetup () override
 
virtual void initialize () override
 
virtual void execute () override
 
virtual void finalize () override
 
virtual void threadJoin (const UserObject &uo) override
 
virtual Real computeVelocity (unsigned i, unsigned j, unsigned qp) const override
 Computes the transfer velocity between current node i and current node j at the current qp in the current element (_current_elem). More...
 
virtual void executeOnElement (dof_id_type global_i, dof_id_type global_j, unsigned local_i, unsigned local_j, unsigned qp) override
 This is called by multiple times in execute() in a double loop over _current_elem's nodes (local_i and local_j) nested in a loop over each of _current_elem's quadpoints (qp). More...
 
virtual Real computedU_dvar (unsigned i, unsigned pvar) const =0
 Compute d(u)/d(porous_flow_variable) More...
 
const std::map< dof_id_type, std::vector< Real > > & getdK_dvar (dof_id_type node_i, dof_id_type node_j) const
 Returns, r, where r[global node k][a] = d(K[node_i][node_j])/d(porous_flow_variable[global node k][porous_flow_variable a]) More...
 
virtual Real computeU (unsigned i) const =0
 Computes the value of u at the local node id of the current element (_current_elem) More...
 
void limitFlux (Real a, Real b, Real &limited, Real &dlimited_db) const
 flux limiter, L, on Page 135 of Kuzmin and Turek More...
 
Real rPlus (dof_id_type sequential_i, std::vector< Real > &dlimited_du, std::vector< Real > &dlimited_dk) const
 Returns the value of R_{i}^{+}, Eqn (49) of KT. More...
 
Real rMinus (dof_id_type sequential_i, std::vector< Real > &dlimited_du, std::vector< Real > &dlimited_dk) const
 Returns the value of R_{i}^{-}, Eqn (49) of KT. More...
 
Real PQPlusMinus (dof_id_type sequential_i, const PQPlusMinusEnum pq_plus_minus, std::vector< Real > &derivs, std::vector< Real > &dpq_dk) const
 Returns the value of P_{i}^{+}, P_{i}^{-}, Q_{i}^{+} or Q_{i}^{-} (depending on pq_plus_minus) which are defined in Eqns (47) and (48) of KT. More...
 
void zeroedConnection (std::map< dof_id_type, Real > &the_map, dof_id_type node_i) const
 Clears the_map, then, using _kij, constructs the_map so that the_map[node_id] = 0.0 for all node_id connected with node_i. More...
 

Protected Attributes

const PorousFlowDictator_dictator
 PorousFlowDictator UserObject. More...
 
const unsigned _num_vars
 Number of PorousFlow variables. More...
 
const RealVectorValue _gravity
 Gravity. More...
 
const unsigned int _phase
 The phase. 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_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< 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 FEType _fe_type
 FEType to use. More...
 
const VariablePhiValue & _phi
 Kuzmin-Turek shape function. More...
 
const VariablePhiGradient & _grad_phi
 grad(Kuzmin-Turek shape function) More...
 
std::vector< std::vector< Real > > _du_dvar
 _du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a]) More...
 
std::vector< bool > _du_dvar_computed_by_thread
 Whether _du_dvar has been computed by the local thread. More...
 
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
 _dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][porous_flow_variable a]) Here j is the j^th connection to sequential node sequential_i, and node k must be connected to node i More...
 
std::vector< std::map< dof_id_type, std::vector< Real > > > _dflux_out_dvars
 _dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous_flow_variable pvar at global node j) More...
 
bool _resizing_needed
 whether _kij, etc, need to be sized appropriately (and valence recomputed) at the start of the timestep More...
 
enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum _flux_limiter_type
 
std::vector< std::vector< Real > > _kij
 Kuzmin-Turek K_ij matrix. More...
 
std::vector< Real > _flux_out
 _flux_out[i] = flux of "heat" from sequential node i More...
 
std::vector< std::map< dof_id_type, Real > > _dflux_out_du
 _dflux_out_du[i][j] = d(flux_out[i])/d(u[j]). More...
 
std::vector< std::vector< std::vector< Real > > > _dflux_out_dKjk
 _dflux_out_dKjk[sequential_i][j][k] = d(flux_out[sequential_i])/d(K[j][k]). More...
 
std::vector< unsigned > _valence
 _valence[i] = number of times, in a loop over elements seen by this processor (viz, including ghost elements) and are part of the block-restricted blocks of this UserObject, that the sequential node i is encountered More...
 
std::vector< Real > _u_nodal
 _u_nodal[i] = value of _u at sequential node number i More...
 
std::vector< bool > _u_nodal_computed_by_thread
 _u_nodal_computed_by_thread(i) = true if _u_nodal[i] has been computed in execute() by the thread on this processor More...
 
PorousFlowConnectedNodes _connections
 Holds the sequential and global nodal IDs, and info regarding mesh connections between them. More...
 
std::size_t _number_of_nodes
 Number of nodes held by the _connections object. More...
 

Detailed Description

Base class to compute the advective flux of fluid in PorousFlow situations using the Kuzmin-Turek FEM-TVD multidimensional stabilization scheme.

The velocity is U * (-permeability * (grad(P) - density * gravity)) with derived classes defining U.

Much of this base class is spent with defining derivatives of the velocity with respect to the PorousFlow variables. As a first-time reader, if you ignore these derivatives, you will find this class is quite simple.

Definition at line 33 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Member Enumeration Documentation

◆ FluxLimiterTypeEnum

enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum
strongprotectedinherited

Determines Flux Limiter type (Page 135 of Kuzmin and Turek) "None" means that limitFlux=0 always, which implies zero antidiffusion will be added.

Enumerator
MinMod 
VanLeer 
MC 
superbee 
None 

Definition at line 154 of file AdvectiveFluxCalculatorBase.h.

154 { MinMod, VanLeer, MC, superbee, None } _flux_limiter_type;
enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum _flux_limiter_type

◆ PQPlusMinusEnum

enum AdvectiveFluxCalculatorBase::PQPlusMinusEnum
strongprotectedinherited

Signals to the PQPlusMinus method what should be computed.

Enumerator
PPlus 
PMinus 
QPlus 
QMinus 

Definition at line 195 of file AdvectiveFluxCalculatorBase.h.

196  {
197  PPlus,
198  PMinus,
199  QPlus,
200  QMinus
201  };

Constructor & Destructor Documentation

◆ PorousFlowAdvectiveFluxCalculatorBase()

PorousFlowAdvectiveFluxCalculatorBase::PorousFlowAdvectiveFluxCalculatorBase ( const InputParameters &  parameters)

Definition at line 42 of file PorousFlowAdvectiveFluxCalculatorBase.C.

44  : AdvectiveFluxCalculatorBase(parameters),
45  _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
47  _gravity(getParam<RealVectorValue>("gravity")),
48  _phase(getParam<unsigned int>("phase")),
49  _permeability(getMaterialProperty<RealTensorValue>("PorousFlow_permeability_qp")),
51  getMaterialProperty<std::vector<RealTensorValue>>("dPorousFlow_permeability_qp_dvar")),
52  _dpermeability_dgradvar(getMaterialProperty<std::vector<std::vector<RealTensorValue>>>(
53  "dPorousFlow_permeability_qp_dgradvar")),
54  _fluid_density_qp(getMaterialProperty<std::vector<Real>>("PorousFlow_fluid_phase_density_qp")),
55  _dfluid_density_qp_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
56  "dPorousFlow_fluid_phase_density_qp_dvar")),
57  _grad_p(getMaterialProperty<std::vector<RealGradient>>("PorousFlow_grad_porepressure_qp")),
58  _dgrad_p_dgrad_var(getMaterialProperty<std::vector<std::vector<Real>>>(
59  "dPorousFlow_grad_porepressure_qp_dgradvar")),
60  _dgrad_p_dvar(getMaterialProperty<std::vector<std::vector<RealGradient>>>(
61  "dPorousFlow_grad_porepressure_qp_dvar")),
62  _fe_type(isParamValid("fe_family") && isParamValid("fe_order")
63  ? FEType(Utility::string_to_enum<Order>(getParam<MooseEnum>("fe_order")),
64  Utility::string_to_enum<FEFamily>(getParam<MooseEnum>("fe_family")))
65  : _dictator.feType()),
66  _phi(_assembly.fePhi<Real>(_fe_type)),
67  _grad_phi(_assembly.feGradPhi<Real>(_fe_type)),
68  _du_dvar(),
70  _dkij_dvar(),
72 
73 {
74  if (_phase >= _dictator.numPhases())
75  paramError("phase",
76  "Phase number entered is greater than the number of phases specified in the "
77  "Dictator. Remember that indexing starts at 0");
78 
79  if (isParamValid("fe_family") && !isParamValid("fe_order"))
80  paramError("fe_order", "If you specify fe_family you must also specify fe_order");
81  if (isParamValid("fe_order") && !isParamValid("fe_family"))
82  paramError("fe_family", "If you specify fe_order you must also specify fe_family");
83  if (!_dictator.consistentFEType() && !isParamValid("fe_family"))
84  paramError("fe_family",
85  "The PorousFlowDictator cannot determine the appropriate FE type to use because "
86  "your porous_flow_vars are of different types. You must specify the appropriate "
87  "fe_family and fe_order to use.");
88 }
const MaterialProperty< RealTensorValue > & _permeability
Permeability of porous material.
const VariablePhiGradient & _grad_phi
grad(Kuzmin-Turek shape function)
const MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _dpermeability_dgradvar
d(permeabiity)/d(grad(PorousFlow variable))
const VariablePhiValue & _phi
Kuzmin-Turek shape function.
const MaterialProperty< std::vector< RealTensorValue > > & _dpermeability_dvar
d(permeabiity)/d(PorousFlow variable)
const MaterialProperty< std::vector< RealGradient > > & _grad_p
Gradient of the pore pressure in each phase.
std::vector< bool > _du_dvar_computed_by_thread
Whether _du_dvar has been computed by the local thread.
const MaterialProperty< std::vector< Real > > & _fluid_density_qp
Fluid density for each phase (at the qp)
std::vector< std::map< dof_id_type, std::vector< Real > > > _dflux_out_dvars
_dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous...
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) ...
std::vector< std::vector< Real > > _du_dvar
_du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a]) ...
unsigned int numPhases() const
The number of fluid phases.
unsigned int numVariables() const
The number of PorousFlow variables.
const MaterialProperty< std::vector< std::vector< Real > > > & _dgrad_p_dgrad_var
Derivative of Grad porepressure in each phase wrt grad(PorousFlow variables)
AdvectiveFluxCalculatorBase(const InputParameters &parameters)
const unsigned _num_vars
Number of PorousFlow variables.
const PorousFlowDictator & _dictator
PorousFlowDictator UserObject.
const MaterialProperty< std::vector< std::vector< RealGradient > > > & _dgrad_p_dvar
Derivative of Grad porepressure in each phase wrt PorousFlow variables.
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][por...
bool consistentFEType() const
Whether the porous_flow_vars all have the same FEType or if no porous_flow_vars were provided...

Member Function Documentation

◆ computedU_dvar()

virtual Real PorousFlowAdvectiveFluxCalculatorBase::computedU_dvar ( unsigned  i,
unsigned  pvar 
) const
protectedpure virtual

◆ computeU()

virtual Real AdvectiveFluxCalculatorBase::computeU ( unsigned  i) const
protectedpure virtualinherited

◆ computeVelocity()

Real PorousFlowAdvectiveFluxCalculatorBase::computeVelocity ( unsigned  i,
unsigned  j,
unsigned  qp 
) const
overrideprotectedvirtual

Computes the transfer velocity between current node i and current node j at the current qp in the current element (_current_elem).

For instance, (_grad_phi[i][qp] * _velocity) * _phi[j][qp];

Parameters
inode number in the current element
jnode number in the current element
qpquadpoint number in the current element

Implements AdvectiveFluxCalculatorBase.

Definition at line 91 of file PorousFlowAdvectiveFluxCalculatorBase.C.

92 {
93  // The following is but one choice for PorousFlow situations
94  // If you change this, you will probably have to change
95  // - the derivative in executeOnElement
96  // - computeU
97  // - computedU_dvar
98  return -_grad_phi[i][qp] *
100  _phi[j][qp];
101 }
const MaterialProperty< RealTensorValue > & _permeability
Permeability of porous material.
const VariablePhiGradient & _grad_phi
grad(Kuzmin-Turek shape function)
const VariablePhiValue & _phi
Kuzmin-Turek shape function.
const MaterialProperty< std::vector< RealGradient > > & _grad_p
Gradient of the pore pressure in each phase.
const MaterialProperty< std::vector< Real > > & _fluid_density_qp
Fluid density for each phase (at the qp)

◆ execute()

void PorousFlowAdvectiveFluxCalculatorBase::execute ( )
overrideprotectedvirtual

Reimplemented from AdvectiveFluxCalculatorBase.

Definition at line 182 of file PorousFlowAdvectiveFluxCalculatorBase.C.

183 {
185 
186  // compute d(U)/d(porous_flow_variables) for nodes in _current_elem and for this
187  // execution thread. In threadJoin all these computations get gathered
188  // using _du_dvar_computed_by_thread
189  for (unsigned i = 0; i < _current_elem->n_nodes(); ++i)
190  {
191  const dof_id_type global_i = _current_elem->node_id(i);
192  const dof_id_type sequential_i = _connections.sequentialID(global_i);
193  if (_du_dvar_computed_by_thread[sequential_i])
194  continue;
195  for (unsigned pvar = 0; pvar < _num_vars; ++pvar)
196  _du_dvar[sequential_i][pvar] = computedU_dvar(i, pvar);
197  _du_dvar_computed_by_thread[sequential_i] = true;
198  }
199 }
virtual Real computedU_dvar(unsigned i, unsigned pvar) const =0
Compute d(u)/d(porous_flow_variable)
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
std::vector< bool > _du_dvar_computed_by_thread
Whether _du_dvar has been computed by the local thread.
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
std::vector< std::vector< Real > > _du_dvar
_du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a]) ...
const unsigned _num_vars
Number of PorousFlow variables.

◆ executeOnElement()

void PorousFlowAdvectiveFluxCalculatorBase::executeOnElement ( dof_id_type  global_i,
dof_id_type  global_j,
unsigned  local_i,
unsigned  local_j,
unsigned  qp 
)
overrideprotectedvirtual

This is called by multiple times in execute() in a double loop over _current_elem's nodes (local_i and local_j) nested in a loop over each of _current_elem's quadpoints (qp).

It is used to compute _kij and its derivatives

Parameters
global_iglobal node id corresponding to the local node local_i
global_jglobal node id corresponding to the local node local_j
local_ilocal node number of the _current_elem
local_jlocal node number of the _current_elem
qpquadpoint number of the _current_elem

Reimplemented from AdvectiveFluxCalculatorBase.

Definition at line 104 of file PorousFlowAdvectiveFluxCalculatorBase.C.

106 {
107  AdvectiveFluxCalculatorBase::executeOnElement(global_i, global_j, local_i, local_j, qp);
108  const dof_id_type sequential_i = _connections.sequentialID(global_i);
109  const unsigned j = _connections.indexOfGlobalConnection(global_i, global_j);
110 
111  // compute d(Kij)/d(porous_flow_variables)
112  for (unsigned local_k = 0; local_k < _current_elem->n_nodes(); ++local_k)
113  {
114  const dof_id_type global_k = _current_elem->node_id(local_k);
115  for (unsigned pvar = 0; pvar < _num_vars; ++pvar)
116  {
117  RealVectorValue deriv = _dpermeability_dvar[qp][pvar] * _phi[local_k][qp] *
119  for (unsigned i = 0; i < LIBMESH_DIM; ++i)
120  deriv += _dpermeability_dgradvar[qp][i][pvar] * _grad_phi[local_k][qp](i) *
122  deriv += _permeability[qp] *
123  (_grad_phi[local_k][qp] * _dgrad_p_dgrad_var[qp][_phase][pvar] -
124  _phi[local_k][qp] * _dfluid_density_qp_dvar[qp][_phase][pvar] * _gravity);
125  deriv += _permeability[qp] * (_dgrad_p_dvar[qp][_phase][pvar] * _phi[local_k][qp]);
126  _dkij_dvar[sequential_i][j][global_k][pvar] +=
127  _JxW[qp] * _coord[qp] * (-_grad_phi[local_i][qp] * deriv * _phi[local_j][qp]);
128  }
129  }
130 }
const MaterialProperty< RealTensorValue > & _permeability
Permeability of porous material.
const VariablePhiGradient & _grad_phi
grad(Kuzmin-Turek shape function)
const MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _dpermeability_dgradvar
d(permeabiity)/d(grad(PorousFlow variable))
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
unsigned indexOfGlobalConnection(dof_id_type global_node_ID_from, dof_id_type global_node_ID_to) const
Return the index of global_node_ID_to in the globalConnectionsToGlobalID(global_node_ID_from) vector...
const VariablePhiValue & _phi
Kuzmin-Turek shape function.
virtual void executeOnElement(dof_id_type global_i, dof_id_type global_j, unsigned local_i, unsigned local_j, unsigned qp)
This is called by multiple times in execute() in a double loop over _current_elem&#39;s nodes (local_i an...
const MaterialProperty< std::vector< RealTensorValue > > & _dpermeability_dvar
d(permeabiity)/d(PorousFlow variable)
const MaterialProperty< std::vector< RealGradient > > & _grad_p
Gradient of the pore pressure in each phase.
const MaterialProperty< std::vector< Real > > & _fluid_density_qp
Fluid density for each phase (at the qp)
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
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) ...
const MaterialProperty< std::vector< std::vector< Real > > > & _dgrad_p_dgrad_var
Derivative of Grad porepressure in each phase wrt grad(PorousFlow variables)
const unsigned _num_vars
Number of PorousFlow variables.
const MaterialProperty< std::vector< std::vector< RealGradient > > > & _dgrad_p_dvar
Derivative of Grad porepressure in each phase wrt PorousFlow variables.
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][por...

◆ finalize()

void PorousFlowAdvectiveFluxCalculatorBase::finalize ( )
overrideprotectedvirtual

_dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous_flow_variable pvar at global node j)

Reimplemented from AdvectiveFluxCalculatorBase.

Definition at line 243 of file PorousFlowAdvectiveFluxCalculatorBase.C.

244 {
246 
247  // compute d(flux_out)/d(porous flow variable)
249  for (const auto & node_i : _connections.globalIDs())
250  {
251  const dof_id_type sequential_i = _connections.sequentialID(node_i);
252  _dflux_out_dvars[sequential_i].clear();
253 
254  const std::map<dof_id_type, Real> dflux_out_du =
256  for (const auto & node_du : dflux_out_du)
257  {
258  const dof_id_type j = node_du.first;
259  const Real dflux_out_du_j = node_du.second;
260  _dflux_out_dvars[sequential_i][j] = _du_dvar[_connections.sequentialID(j)];
261  for (unsigned pvar = 0; pvar < _num_vars; ++pvar)
262  _dflux_out_dvars[sequential_i][j][pvar] *= dflux_out_du_j;
263  }
264 
265  // _dflux_out_dvars is now sized correctly, because getdFluxOutdu(i) contains all nodes
266  // connected to i and all nodes connected to nodes connected to i. The
267  // getdFluxOutdKij contains no extra nodes, so just += the dflux/dK terms
268  const std::vector<std::vector<Real>> dflux_out_dKjk =
270  const std::vector<dof_id_type> con_i = _connections.globalConnectionsToGlobalID(node_i);
271  for (std::size_t index_j = 0; index_j < con_i.size(); ++index_j)
272  {
273  const dof_id_type node_j = con_i[index_j];
274  const std::vector<dof_id_type> con_j = _connections.globalConnectionsToGlobalID(node_j);
275  for (std::size_t index_k = 0; index_k < con_j.size(); ++index_k)
276  {
277  const dof_id_type node_k = con_j[index_k];
278  const Real dflux_out_dK_jk = dflux_out_dKjk[index_j][index_k];
279  const std::map<dof_id_type, std::vector<Real>> dkj_dvarl = getdK_dvar(node_j, node_k);
280  for (const auto & nodel_deriv : dkj_dvarl)
281  {
282  const dof_id_type l = nodel_deriv.first;
283  for (unsigned pvar = 0; pvar < _num_vars; ++pvar)
284  _dflux_out_dvars[sequential_i][l][pvar] += dflux_out_dK_jk * nodel_deriv.second[pvar];
285  }
286  }
287  }
288  }
289 }
const std::vector< dof_id_type > & globalIDs() const
Vector of all global node IDs (node numbers in the mesh)
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
std::vector< std::map< dof_id_type, std::vector< Real > > > _dflux_out_dvars
_dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
const std::map< dof_id_type, Real > & getdFluxOutdu(dof_id_type node_i) const
Returns r where r[j] = d(flux out of global node i)/du(global node j) used in Jacobian computations...
std::vector< std::vector< Real > > _du_dvar
_du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a]) ...
const std::vector< std::vector< Real > > & getdFluxOutdKjk(dof_id_type node_i) const
Returns r where r[j][k] = d(flux out of global node i)/dK[connected node j][connected node k] used in...
const unsigned _num_vars
Number of PorousFlow variables.
const std::vector< dof_id_type > & globalConnectionsToGlobalID(dof_id_type global_node_ID) const
Return all the nodes (global node IDs) connected to the given global node ID.
const std::map< dof_id_type, std::vector< Real > > & getdK_dvar(dof_id_type node_i, dof_id_type node_j) const
Returns, r, where r[global node k][a] = d(K[node_i][node_j])/d(porous_flow_variable[global node k][po...

◆ getdFluxOut_dvars()

const std::map< dof_id_type, std::vector< Real > > & PorousFlowAdvectiveFluxCalculatorBase::getdFluxOut_dvars ( unsigned  node_id) const

Returns d(flux_out)/d(porous_flow_variables.

Parameters
[in]node_idglobal node id
Returns
deriv[j][pvar] = d(flux_out[node_id])/d(porous_flow_variable pvar at global node j)

Definition at line 237 of file PorousFlowAdvectiveFluxCalculatorBase.C.

Referenced by PorousFlowFluxLimitedTVDAdvection::computeJacobian().

238 {
239  return _dflux_out_dvars[_connections.sequentialID(node_id)];
240 }
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
std::vector< std::map< dof_id_type, std::vector< Real > > > _dflux_out_dvars
_dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...

◆ getdFluxOutdKjk()

const std::vector< std::vector< Real > > & AdvectiveFluxCalculatorBase::getdFluxOutdKjk ( dof_id_type  node_i) const
inherited

Returns r where r[j][k] = d(flux out of global node i)/dK[connected node j][connected node k] used in Jacobian computations.

Parameters
node_iglobal id of node
Returns
the derivatives (after applying the KT procedure)

Definition at line 737 of file AdvectiveFluxCalculatorBase.C.

Referenced by finalize().

738 {
739  return _dflux_out_dKjk[_connections.sequentialID(node_i)];
740 }
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
std::vector< std::vector< std::vector< Real > > > _dflux_out_dKjk
_dflux_out_dKjk[sequential_i][j][k] = d(flux_out[sequential_i])/d(K[j][k]).
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...

◆ getdFluxOutdu()

const std::map< dof_id_type, Real > & AdvectiveFluxCalculatorBase::getdFluxOutdu ( dof_id_type  node_i) const
inherited

Returns r where r[j] = d(flux out of global node i)/du(global node j) used in Jacobian computations.

Parameters
node_iglobal id of node
Returns
the derivatives (after applying the KT procedure)

Definition at line 731 of file AdvectiveFluxCalculatorBase.C.

Referenced by FluxLimitedTVDAdvection::computeJacobian(), and finalize().

732 {
733  return _dflux_out_du[_connections.sequentialID(node_i)];
734 }
std::vector< std::map< dof_id_type, Real > > _dflux_out_du
_dflux_out_du[i][j] = d(flux_out[i])/d(u[j]).
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...

◆ getdK_dvar()

const std::map< dof_id_type, std::vector< Real > > & PorousFlowAdvectiveFluxCalculatorBase::getdK_dvar ( dof_id_type  node_i,
dof_id_type  node_j 
) const
protected

Returns, r, where r[global node k][a] = d(K[node_i][node_j])/d(porous_flow_variable[global node k][porous_flow_variable a])

Parameters
node_iglobal node id param node_j global node id

Definition at line 229 of file PorousFlowAdvectiveFluxCalculatorBase.C.

Referenced by finalize().

230 {
231  const dof_id_type sequential_i = _connections.sequentialID(node_i);
232  const unsigned j = _connections.indexOfGlobalConnection(node_i, node_j);
233  return _dkij_dvar[sequential_i][j];
234 }
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
unsigned indexOfGlobalConnection(dof_id_type global_node_ID_from, dof_id_type global_node_ID_to) const
Return the index of global_node_ID_to in the globalConnectionsToGlobalID(global_node_ID_from) vector...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][por...

◆ getFluxOut()

Real AdvectiveFluxCalculatorBase::getFluxOut ( dof_id_type  node_i) const
inherited

Returns the flux out of lobal node id.

Parameters
node_iid of node
Returns
advective flux out of node after applying the KT procedure

Definition at line 743 of file AdvectiveFluxCalculatorBase.C.

Referenced by FluxLimitedTVDAdvection::computeResidual(), and PorousFlowFluxLimitedTVDAdvection::computeResidual().

744 {
745  return _flux_out[_connections.sequentialID(node_i)];
746 }
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
std::vector< Real > _flux_out
_flux_out[i] = flux of "heat" from sequential node i
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...

◆ getValence()

unsigned AdvectiveFluxCalculatorBase::getValence ( dof_id_type  node_i) const
inherited

Returns the valence of the global node i Valence is the number of times the node is encountered in a loop over elements (that have appropriate subdomain_id, if the user has employed the "blocks=" parameter) seen by this processor (including ghosted elements)

Parameters
node_igloal id of i^th node
Returns
valence of the node

Definition at line 749 of file AdvectiveFluxCalculatorBase.C.

Referenced by FluxLimitedTVDAdvection::computeJacobian(), PorousFlowFluxLimitedTVDAdvection::computeJacobian(), FluxLimitedTVDAdvection::computeResidual(), and PorousFlowFluxLimitedTVDAdvection::computeResidual().

750 {
751  return _valence[_connections.sequentialID(node_i)];
752 }
dof_id_type sequentialID(dof_id_type global_node_ID) const
Return the sequential node ID corresponding to the global node ID This is guaranteed to lie in the ra...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
std::vector< unsigned > _valence
_valence[i] = number of times, in a loop over elements seen by this processor (viz, including ghost elements) and are part of the block-restricted blocks of this UserObject, that the sequential node i is encountered

◆ initialize()

void PorousFlowAdvectiveFluxCalculatorBase::initialize ( )
overrideprotectedvirtual

Reimplemented from AdvectiveFluxCalculatorBase.

Definition at line 165 of file PorousFlowAdvectiveFluxCalculatorBase.C.

166 {
168  const std::size_t num_nodes = _connections.numNodes();
169  _du_dvar_computed_by_thread.assign(num_nodes, false);
170  for (dof_id_type sequential_i = 0; sequential_i < num_nodes; ++sequential_i)
171  {
172  const std::vector<dof_id_type> con_i =
174  const std::size_t num_con_i = con_i.size();
175  for (unsigned j = 0; j < num_con_i; ++j)
176  for (const auto & global_neighbor_to_i : con_i)
177  _dkij_dvar[sequential_i][j][global_neighbor_to_i] = std::vector<Real>(_num_vars, 0.0);
178  }
179 }
const std::vector< dof_id_type > & globalConnectionsToSequentialID(dof_id_type sequential_node_ID) const
Return all the nodes (global node IDs) connected to the given sequential node ID. ...
std::vector< bool > _du_dvar_computed_by_thread
Whether _du_dvar has been computed by the local thread.
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
std::size_t numNodes() const
number of nodes known by this class
const unsigned _num_vars
Number of PorousFlow variables.
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][por...

◆ limitFlux()

void AdvectiveFluxCalculatorBase::limitFlux ( Real  a,
Real  b,
Real &  limited,
Real &  dlimited_db 
) const
protectedinherited

flux limiter, L, on Page 135 of Kuzmin and Turek

Parameters
aKT's "a" parameter
bKT's "b" parameter
limited[out]The value of the flux limiter, L
dlimited_db[out]The derivative dL/db

Definition at line 631 of file AdvectiveFluxCalculatorBase.C.

Referenced by AdvectiveFluxCalculatorBase::rMinus(), and AdvectiveFluxCalculatorBase::rPlus().

632 {
633  limited = 0.0;
634  dlimited_db = 0.0;
636  return;
637 
638  if ((a >= 0.0 && b <= 0.0) || (a <= 0.0 && b >= 0.0))
639  return;
640  const Real s = (a > 0.0 ? 1.0 : -1.0);
641 
642  const Real lal = std::abs(a);
643  const Real lbl = std::abs(b);
644  const Real dlbl = (b >= 0.0 ? 1.0 : -1.0); // d(lbl)/db
645  switch (_flux_limiter_type)
646  {
648  {
649  if (lal <= lbl)
650  {
651  limited = s * lal;
652  dlimited_db = 0.0;
653  }
654  else
655  {
656  limited = s * lbl;
657  dlimited_db = s * dlbl;
658  }
659  return;
660  }
662  {
663  limited = s * 2 * lal * lbl / (lal + lbl);
664  dlimited_db = s * 2 * lal * (dlbl / (lal + lbl) - lbl * dlbl / std::pow(lal + lbl, 2));
665  return;
666  }
668  {
669  const Real av = 0.5 * std::abs(a + b);
670  if (2 * lal <= av && lal <= lbl)
671  {
672  // 2 * lal is the smallest
673  limited = s * 2.0 * lal;
674  dlimited_db = 0.0;
675  }
676  else if (2 * lbl <= av && lbl <= lal)
677  {
678  // 2 * lbl is the smallest
679  limited = s * 2.0 * lbl;
680  dlimited_db = s * 2.0 * dlbl;
681  }
682  else
683  {
684  // av is the smallest
685  limited = s * av;
686  // if (a>0 and b>0) then d(av)/db = 0.5 = 0.5 * dlbl
687  // if (a<0 and b<0) then d(av)/db = -0.5 = 0.5 * dlbl
688  // if a and b have different sign then limited=0, above
689  dlimited_db = s * 0.5 * dlbl;
690  }
691  return;
692  }
694  {
695  const Real term1 = std::min(2.0 * lal, lbl);
696  const Real term2 = std::min(lal, 2.0 * lbl);
697  if (term1 >= term2)
698  {
699  if (2.0 * lal <= lbl)
700  {
701  limited = s * 2 * lal;
702  dlimited_db = 0.0;
703  }
704  else
705  {
706  limited = s * lbl;
707  dlimited_db = s * dlbl;
708  }
709  }
710  else
711  {
712  if (lal <= 2.0 * lbl)
713  {
714  limited = s * lal;
715  dlimited_db = 0.0;
716  }
717  else
718  {
719  limited = s * 2.0 * lbl;
720  dlimited_db = s * 2.0 * dlbl;
721  }
722  }
723  return;
724  }
725  default:
726  return;
727  }
728 }
enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum _flux_limiter_type
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ meshChanged()

void AdvectiveFluxCalculatorBase::meshChanged ( )
overridevirtualinherited

Definition at line 149 of file AdvectiveFluxCalculatorBase.C.

150 {
151  ElementUserObject::meshChanged();
152 
153  // Signal that _kij, _valence, etc need to be rebuilt
154  _resizing_needed = true;
155 }
bool _resizing_needed
whether _kij, etc, need to be sized appropriately (and valence recomputed) at the start of the timest...

◆ PQPlusMinus()

Real AdvectiveFluxCalculatorBase::PQPlusMinus ( dof_id_type  sequential_i,
const PQPlusMinusEnum  pq_plus_minus,
std::vector< Real > &  derivs,
std::vector< Real > &  dpq_dk 
) const
protectedinherited

Returns the value of P_{i}^{+}, P_{i}^{-}, Q_{i}^{+} or Q_{i}^{-} (depending on pq_plus_minus) which are defined in Eqns (47) and (48) of KT.

Parameters
sequential_isequential nodal ID
pq_plus_minusindicates whether P_{i}^{+}, P_{i}^{-}, Q_{i}^{+} or Q_{i}^{-} should be returned
derivs[out]derivs[j] = d(result)/d(u[sequential_j]). Here sequential_j is the j^th connection to sequential_i
dpq_dk[out]dpq_dk[j] = d(result)/d(K[node_i][j]). Here j indexes a connection to sequential_i. Recall that d(result)/d(K[l][m]) are zero unless l=sequential_i

Definition at line 764 of file AdvectiveFluxCalculatorBase.C.

Referenced by AdvectiveFluxCalculatorBase::rMinus(), and AdvectiveFluxCalculatorBase::rPlus().

768 {
769  // Find the value of u at sequential_i
770  const Real u_i = _u_nodal[sequential_i];
771 
772  // Connections to sequential_i
773  const std::vector<dof_id_type> con_i =
775  const std::size_t num_con = con_i.size();
776  // The neighbor number of sequential_i to sequential_i
777  const unsigned i_index_i = _connections.indexOfSequentialConnection(sequential_i, sequential_i);
778 
779  // Initialize the results
780  Real result = 0.0;
781  derivs.assign(num_con, 0.0);
782  dpqdk.assign(num_con, 0.0);
783 
784  // Sum over all nodes connected with node_i.
785  for (std::size_t j = 0; j < num_con; ++j)
786  {
787  const dof_id_type sequential_j = con_i[j];
788  if (sequential_j == sequential_i)
789  continue;
790  const Real kentry = _kij[sequential_i][j];
791 
792  // Find the value of u at node_j
793  const Real u_j = _u_nodal[sequential_j];
794  const Real ujminusi = u_j - u_i;
795 
796  // Evaluate the i-j contribution to the result
797  switch (pq_plus_minus)
798  {
800  {
801  if (ujminusi < 0.0 && kentry < 0.0)
802  {
803  result += kentry * ujminusi;
804  derivs[j] += kentry;
805  derivs[i_index_i] -= kentry;
806  dpqdk[j] += ujminusi;
807  }
808  break;
809  }
811  {
812  if (ujminusi > 0.0 && kentry < 0.0)
813  {
814  result += kentry * ujminusi;
815  derivs[j] += kentry;
816  derivs[i_index_i] -= kentry;
817  dpqdk[j] += ujminusi;
818  }
819  break;
820  }
822  {
823  if (ujminusi > 0.0 && kentry > 0.0)
824  {
825  result += kentry * ujminusi;
826  derivs[j] += kentry;
827  derivs[i_index_i] -= kentry;
828  dpqdk[j] += ujminusi;
829  }
830  break;
831  }
833  {
834  if (ujminusi < 0.0 && kentry > 0.0)
835  {
836  result += kentry * ujminusi;
837  derivs[j] += kentry;
838  derivs[i_index_i] -= kentry;
839  dpqdk[j] += ujminusi;
840  }
841  break;
842  }
843  }
844  }
845 
846  return result;
847 }
const std::vector< dof_id_type > & sequentialConnectionsToSequentialID(dof_id_type sequential_node_ID) const
Return all the nodes (sequential node IDs) connected to the given sequential node ID All elements of ...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
std::vector< Real > _u_nodal
_u_nodal[i] = value of _u at sequential node number i
std::vector< std::vector< Real > > _kij
Kuzmin-Turek K_ij matrix.
unsigned indexOfSequentialConnection(dof_id_type sequential_node_ID_from, dof_id_type sequential_node_ID_to) const
Return the index of sequential_node_ID_to in the sequentialConnectionsToSequentialID(sequential_node_...

◆ rMinus()

Real AdvectiveFluxCalculatorBase::rMinus ( dof_id_type  sequential_i,
std::vector< Real > &  dlimited_du,
std::vector< Real > &  dlimited_dk 
) const
protectedinherited

Returns the value of R_{i}^{-}, Eqn (49) of KT.

Parameters
sequential_iSequential nodal ID
dlimited_du[out]dlimited_du[j] = d(R_{sequential_i}^{-})/du[sequential_j]. Here sequential_j is the j^th connection to sequential_i
dlimited_dk[out]dlimited_dk[j] = d(R_{sequential_i}^{-})/d(K[sequential_i][j]). Note Derivatives w.r.t. K[l][m] with l!=i are zero

Definition at line 595 of file AdvectiveFluxCalculatorBase.C.

Referenced by AdvectiveFluxCalculatorBase::finalize().

598 {
599  const std::size_t num_con = _connections.sequentialConnectionsToSequentialID(sequential_i).size();
600  dlimited_du.assign(num_con, 0.0);
601  dlimited_dk.assign(num_con, 0.0);
603  return 0.0;
604  std::vector<Real> dp_du;
605  std::vector<Real> dp_dk;
606  const Real p = PQPlusMinus(sequential_i, PQPlusMinusEnum::PMinus, dp_du, dp_dk);
607  if (p == 0.0)
608  // Comment after Eqn (49): if P=0 then there's no antidiffusion, so no need to remove it
609  return 1.0;
610  std::vector<Real> dq_du;
611  std::vector<Real> dq_dk;
612  const Real q = PQPlusMinus(sequential_i, PQPlusMinusEnum::QMinus, dq_du, dq_dk);
613 
614  const Real r = q / p;
615  Real limited;
616  Real dlimited_dr;
617  limitFlux(1.0, r, limited, dlimited_dr);
618 
619  const Real p2 = std::pow(p, 2);
620  for (std::size_t j = 0; j < num_con; ++j)
621  {
622  const Real dr_du = dq_du[j] / p - q * dp_du[j] / p2;
623  const Real dr_dk = dq_dk[j] / p - q * dp_dk[j] / p2;
624  dlimited_du[j] = dlimited_dr * dr_du;
625  dlimited_dk[j] = dlimited_dr * dr_dk;
626  }
627  return limited;
628 }
void limitFlux(Real a, Real b, Real &limited, Real &dlimited_db) const
flux limiter, L, on Page 135 of Kuzmin and Turek
const std::vector< dof_id_type > & sequentialConnectionsToSequentialID(dof_id_type sequential_node_ID) const
Return all the nodes (sequential node IDs) connected to the given sequential node ID All elements of ...
Real PQPlusMinus(dof_id_type sequential_i, const PQPlusMinusEnum pq_plus_minus, std::vector< Real > &derivs, std::vector< Real > &dpq_dk) const
Returns the value of P_{i}^{+}, P_{i}^{-}, Q_{i}^{+} or Q_{i}^{-} (depending on pq_plus_minus) which ...
enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum _flux_limiter_type
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ rPlus()

Real AdvectiveFluxCalculatorBase::rPlus ( dof_id_type  sequential_i,
std::vector< Real > &  dlimited_du,
std::vector< Real > &  dlimited_dk 
) const
protectedinherited

Returns the value of R_{i}^{+}, Eqn (49) of KT.

Parameters
node_inodal id
dlimited_du[out]dlimited_du[j] = d(R_{sequential_i}^{+})/du[sequential_j]. Here sequential_j is the j^th connection to sequential_i
dlimited_dk[out]dlimited_dk[j] = d(R_{sequential_i}^{+})/d(K[sequential_i][j]). Note Derivatives w.r.t. K[l][m] with l!=i are zero

Definition at line 559 of file AdvectiveFluxCalculatorBase.C.

Referenced by AdvectiveFluxCalculatorBase::finalize().

562 {
563  const std::size_t num_con = _connections.sequentialConnectionsToSequentialID(sequential_i).size();
564  dlimited_du.assign(num_con, 0.0);
565  dlimited_dk.assign(num_con, 0.0);
567  return 0.0;
568  std::vector<Real> dp_du;
569  std::vector<Real> dp_dk;
570  const Real p = PQPlusMinus(sequential_i, PQPlusMinusEnum::PPlus, dp_du, dp_dk);
571  if (p == 0.0)
572  // Comment after Eqn (49): if P=0 then there's no antidiffusion, so no need to remove it
573  return 1.0;
574  std::vector<Real> dq_du;
575  std::vector<Real> dq_dk;
576  const Real q = PQPlusMinus(sequential_i, PQPlusMinusEnum::QPlus, dq_du, dq_dk);
577 
578  const Real r = q / p;
579  Real limited;
580  Real dlimited_dr;
581  limitFlux(1.0, r, limited, dlimited_dr);
582 
583  const Real p2 = std::pow(p, 2);
584  for (std::size_t j = 0; j < num_con; ++j)
585  {
586  const Real dr_du = dq_du[j] / p - q * dp_du[j] / p2;
587  const Real dr_dk = dq_dk[j] / p - q * dp_dk[j] / p2;
588  dlimited_du[j] = dlimited_dr * dr_du;
589  dlimited_dk[j] = dlimited_dr * dr_dk;
590  }
591  return limited;
592 }
void limitFlux(Real a, Real b, Real &limited, Real &dlimited_db) const
flux limiter, L, on Page 135 of Kuzmin and Turek
const std::vector< dof_id_type > & sequentialConnectionsToSequentialID(dof_id_type sequential_node_ID) const
Return all the nodes (sequential node IDs) connected to the given sequential node ID All elements of ...
Real PQPlusMinus(dof_id_type sequential_i, const PQPlusMinusEnum pq_plus_minus, std::vector< Real > &derivs, std::vector< Real > &dpq_dk) const
Returns the value of P_{i}^{+}, P_{i}^{-}, Q_{i}^{+} or Q_{i}^{-} (depending on pq_plus_minus) which ...
enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum _flux_limiter_type
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ threadJoin()

void PorousFlowAdvectiveFluxCalculatorBase::threadJoin ( const UserObject &  uo)
overrideprotectedvirtual

Reimplemented from AdvectiveFluxCalculatorBase.

Definition at line 202 of file PorousFlowAdvectiveFluxCalculatorBase.C.

203 {
206  static_cast<const PorousFlowAdvectiveFluxCalculatorBase &>(uo);
207  // add the values of _dkij_dvar computed by different threads
208  const std::size_t num_nodes = _connections.numNodes();
209  for (dof_id_type sequential_i = 0; sequential_i < num_nodes; ++sequential_i)
210  {
211  const std::vector<dof_id_type> con_i =
213  const std::size_t num_con_i = con_i.size();
214  for (unsigned j = 0; j < num_con_i; ++j)
215  for (const auto & global_derivs : pfafc._dkij_dvar[sequential_i][j])
216  for (unsigned pvar = 0; pvar < _num_vars; ++pvar)
217  _dkij_dvar[sequential_i][j][global_derivs.first][pvar] += global_derivs.second[pvar];
218  }
219 
220  // gather the values of _du_dvar computed by different threads
221  for (dof_id_type sequential_i = 0; sequential_i < _number_of_nodes; ++sequential_i)
222  if (!_du_dvar_computed_by_thread[sequential_i] &&
223  pfafc._du_dvar_computed_by_thread[sequential_i])
224  _du_dvar[sequential_i] = pfafc._du_dvar[sequential_i];
225 }
std::size_t _number_of_nodes
Number of nodes held by the _connections object.
Base class to compute the advective flux of fluid in PorousFlow situations using the Kuzmin-Turek FEM...
const std::vector< dof_id_type > & globalConnectionsToSequentialID(dof_id_type sequential_node_ID) const
Return all the nodes (global node IDs) connected to the given sequential node ID. ...
std::vector< bool > _du_dvar_computed_by_thread
Whether _du_dvar has been computed by the local thread.
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
std::vector< std::vector< Real > > _du_dvar
_du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a]) ...
std::size_t numNodes() const
number of nodes known by this class
const unsigned _num_vars
Number of PorousFlow variables.
virtual void threadJoin(const UserObject &uo) override
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][por...

◆ timestepSetup()

void PorousFlowAdvectiveFluxCalculatorBase::timestepSetup ( )
overrideprotectedvirtual

Reimplemented from AdvectiveFluxCalculatorBase.

Definition at line 133 of file PorousFlowAdvectiveFluxCalculatorBase.C.

134 {
135  const bool resizing_was_needed =
136  _resizing_needed; // _resizing_needed gets set to false at the end of
137  // AdvectiveFluxCalculatorBase::timestepSetup()
139 
140  // clear and appropriately size all the derivative info
141  // d(U)/d(porous_flow_variables) and
142  // d(Kij)/d(porous_flow_variables) and
143  // d(flux_out)/d(porous_flow_variables)
144  if (resizing_was_needed)
145  {
146  const std::size_t num_nodes = _connections.numNodes();
147  _du_dvar.assign(num_nodes, std::vector<Real>(_num_vars, 0.0));
148  _du_dvar_computed_by_thread.assign(num_nodes, false);
149  _dflux_out_dvars.assign(num_nodes, std::map<dof_id_type, std::vector<Real>>());
150  _dkij_dvar.resize(num_nodes);
151  for (dof_id_type sequential_i = 0; sequential_i < num_nodes; ++sequential_i)
152  {
153  const std::vector<dof_id_type> con_i =
155  const std::size_t num_con_i = con_i.size();
156  _dkij_dvar[sequential_i].assign(num_con_i, std::map<dof_id_type, std::vector<Real>>());
157  for (unsigned j = 0; j < num_con_i; ++j)
158  for (const auto & global_neighbor_to_i : con_i)
159  _dkij_dvar[sequential_i][j][global_neighbor_to_i] = std::vector<Real>(_num_vars, 0.0);
160  }
161  }
162 }
const std::vector< dof_id_type > & globalConnectionsToSequentialID(dof_id_type sequential_node_ID) const
Return all the nodes (global node IDs) connected to the given sequential node ID. ...
std::vector< bool > _du_dvar_computed_by_thread
Whether _du_dvar has been computed by the local thread.
std::vector< std::map< dof_id_type, std::vector< Real > > > _dflux_out_dvars
_dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous...
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
bool _resizing_needed
whether _kij, etc, need to be sized appropriately (and valence recomputed) at the start of the timest...
std::vector< std::vector< Real > > _du_dvar
_du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a]) ...
std::size_t numNodes() const
number of nodes known by this class
const unsigned _num_vars
Number of PorousFlow variables.
std::vector< std::vector< std::map< dof_id_type, std::vector< Real > > > > _dkij_dvar
_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][por...

◆ zeroedConnection()

void AdvectiveFluxCalculatorBase::zeroedConnection ( std::map< dof_id_type, Real > &  the_map,
dof_id_type  node_i 
) const
protectedinherited

Clears the_map, then, using _kij, constructs the_map so that the_map[node_id] = 0.0 for all node_id connected with node_i.

Parameters
[out]the_mapthe map to be zeroed appropriately
[in]node_inodal id

Definition at line 755 of file AdvectiveFluxCalculatorBase.C.

Referenced by AdvectiveFluxCalculatorBase::timestepSetup().

757 {
758  the_map.clear();
759  for (const auto & node_j : _connections.globalConnectionsToGlobalID(node_i))
760  the_map[node_j] = 0.0;
761 }
PorousFlowConnectedNodes _connections
Holds the sequential and global nodal IDs, and info regarding mesh connections between them...
const std::vector< dof_id_type > & globalConnectionsToGlobalID(dof_id_type global_node_ID) const
Return all the nodes (global node IDs) connected to the given global node ID.

Member Data Documentation

◆ _connections

PorousFlowConnectedNodes AdvectiveFluxCalculatorBase::_connections
protectedinherited

◆ _dfluid_density_qp_dvar

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowAdvectiveFluxCalculatorBase::_dfluid_density_qp_dvar
protected

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

Definition at line 100 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by executeOnElement().

◆ _dflux_out_dKjk

std::vector<std::vector<std::vector<Real> > > AdvectiveFluxCalculatorBase::_dflux_out_dKjk
protectedinherited

_dflux_out_dKjk[sequential_i][j][k] = d(flux_out[sequential_i])/d(K[j][k]).

Here sequential_i is a sequential node number according to the _connections object, and j represents the j^th node connected to i according to the _connections object, and k represents the k^th node connected to j according to the _connections object. Here j must be connected to i (this does include (the sequential version of j) == i), and k must be connected to j (this does include (the sequential version of k) = i and (the sequential version of k) == (sequential version of j)))

Definition at line 175 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::finalize(), AdvectiveFluxCalculatorBase::getdFluxOutdKjk(), and AdvectiveFluxCalculatorBase::timestepSetup().

◆ _dflux_out_du

std::vector<std::map<dof_id_type, Real> > AdvectiveFluxCalculatorBase::_dflux_out_du
protectedinherited

_dflux_out_du[i][j] = d(flux_out[i])/d(u[j]).

Here i is a sequential node number according to the _connections object, and j (global ID) must be connected to i, or to a node that is connected to i.

Definition at line 168 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::finalize(), AdvectiveFluxCalculatorBase::getdFluxOutdu(), and AdvectiveFluxCalculatorBase::timestepSetup().

◆ _dflux_out_dvars

std::vector<std::map<dof_id_type, std::vector<Real> > > PorousFlowAdvectiveFluxCalculatorBase::_dflux_out_dvars
protected

_dflux_out_dvars[sequential_i][global_j][pvar] = d(flux_out[global version of sequential_i])/d(porous_flow_variable pvar at global node j)

Definition at line 134 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by finalize(), getdFluxOut_dvars(), and timestepSetup().

◆ _dgrad_p_dgrad_var

const MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowAdvectiveFluxCalculatorBase::_dgrad_p_dgrad_var
protected

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

Definition at line 106 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by executeOnElement().

◆ _dgrad_p_dvar

const MaterialProperty<std::vector<std::vector<RealGradient> > >& PorousFlowAdvectiveFluxCalculatorBase::_dgrad_p_dvar
protected

Derivative of Grad porepressure in each phase wrt PorousFlow variables.

Definition at line 109 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by executeOnElement().

◆ _dictator

const PorousFlowDictator& PorousFlowAdvectiveFluxCalculatorBase::_dictator
protected

◆ _dkij_dvar

std::vector<std::vector<std::map<dof_id_type, std::vector<Real> > > > PorousFlowAdvectiveFluxCalculatorBase::_dkij_dvar
protected

_dkij_dvar[sequential_i][j][global_k][a] = d(K[sequential_i][j])/d(porous_flow_variable[global_k][porous_flow_variable a]) Here j is the j^th connection to sequential node sequential_i, and node k must be connected to node i

Definition at line 131 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by executeOnElement(), getdK_dvar(), initialize(), threadJoin(), and timestepSetup().

◆ _dpermeability_dgradvar

const MaterialProperty<std::vector<std::vector<RealTensorValue> > >& PorousFlowAdvectiveFluxCalculatorBase::_dpermeability_dgradvar
protected

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

Definition at line 94 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by executeOnElement().

◆ _dpermeability_dvar

const MaterialProperty<std::vector<RealTensorValue> >& PorousFlowAdvectiveFluxCalculatorBase::_dpermeability_dvar
protected

d(permeabiity)/d(PorousFlow variable)

Definition at line 91 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by executeOnElement().

◆ _du_dvar

std::vector<std::vector<Real> > PorousFlowAdvectiveFluxCalculatorBase::_du_dvar
protected

_du_dvar[sequential_i][a] = d(u[global version of sequential node i])/d(porous_flow_variable[a])

Definition at line 121 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by execute(), finalize(), threadJoin(), and timestepSetup().

◆ _du_dvar_computed_by_thread

std::vector<bool> PorousFlowAdvectiveFluxCalculatorBase::_du_dvar_computed_by_thread
protected

Whether _du_dvar has been computed by the local thread.

Definition at line 124 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by execute(), initialize(), threadJoin(), and timestepSetup().

◆ _fe_type

const FEType PorousFlowAdvectiveFluxCalculatorBase::_fe_type
protected

FEType to use.

Definition at line 112 of file PorousFlowAdvectiveFluxCalculatorBase.h.

◆ _fluid_density_qp

const MaterialProperty<std::vector<Real> >& PorousFlowAdvectiveFluxCalculatorBase::_fluid_density_qp
protected

Fluid density for each phase (at the qp)

Definition at line 97 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by computeVelocity(), and executeOnElement().

◆ _flux_limiter_type

enum AdvectiveFluxCalculatorBase::FluxLimiterTypeEnum AdvectiveFluxCalculatorBase::_flux_limiter_type
protectedinherited

◆ _flux_out

std::vector<Real> AdvectiveFluxCalculatorBase::_flux_out
protectedinherited

_flux_out[i] = flux of "heat" from sequential node i

Definition at line 164 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::finalize(), AdvectiveFluxCalculatorBase::getFluxOut(), and AdvectiveFluxCalculatorBase::timestepSetup().

◆ _grad_p

const MaterialProperty<std::vector<RealGradient> >& PorousFlowAdvectiveFluxCalculatorBase::_grad_p
protected

Gradient of the pore pressure in each phase.

Definition at line 103 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by computeVelocity(), and executeOnElement().

◆ _grad_phi

const VariablePhiGradient& PorousFlowAdvectiveFluxCalculatorBase::_grad_phi
protected

grad(Kuzmin-Turek shape function)

Definition at line 118 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by computeVelocity(), and executeOnElement().

◆ _gravity

const RealVectorValue PorousFlowAdvectiveFluxCalculatorBase::_gravity
protected

◆ _kij

std::vector<std::vector<Real> > AdvectiveFluxCalculatorBase::_kij
protectedinherited

Kuzmin-Turek K_ij matrix.

Along with R+ and R-, this is the key quantity computed by this UserObject. _kij[i][j] = k_ij corresponding to the i-j node pair. Here i is a sequential node numbers according to the _connections object, and j represents the j^th node connected to i according to the _connections object.

Definition at line 161 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::executeOnElement(), AdvectiveFluxCalculatorBase::finalize(), AdvectiveFluxCalculatorBase::initialize(), AdvectiveFluxCalculatorBase::PQPlusMinus(), AdvectiveFluxCalculatorBase::threadJoin(), and AdvectiveFluxCalculatorBase::timestepSetup().

◆ _num_vars

const unsigned PorousFlowAdvectiveFluxCalculatorBase::_num_vars
protected

Number of PorousFlow variables.

Definition at line 79 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by execute(), executeOnElement(), finalize(), initialize(), threadJoin(), and timestepSetup().

◆ _number_of_nodes

std::size_t AdvectiveFluxCalculatorBase::_number_of_nodes
protectedinherited

◆ _permeability

const MaterialProperty<RealTensorValue>& PorousFlowAdvectiveFluxCalculatorBase::_permeability
protected

Permeability of porous material.

Definition at line 88 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by computeVelocity(), and executeOnElement().

◆ _phase

const unsigned int PorousFlowAdvectiveFluxCalculatorBase::_phase
protected

◆ _phi

const VariablePhiValue& PorousFlowAdvectiveFluxCalculatorBase::_phi
protected

Kuzmin-Turek shape function.

Definition at line 115 of file PorousFlowAdvectiveFluxCalculatorBase.h.

Referenced by computeVelocity(), and executeOnElement().

◆ _resizing_needed

bool AdvectiveFluxCalculatorBase::_resizing_needed
protectedinherited

whether _kij, etc, need to be sized appropriately (and valence recomputed) at the start of the timestep

Definition at line 115 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::meshChanged(), AdvectiveFluxCalculatorBase::timestepSetup(), and timestepSetup().

◆ _u_nodal

std::vector<Real> AdvectiveFluxCalculatorBase::_u_nodal
protectedinherited

◆ _u_nodal_computed_by_thread

std::vector<bool> AdvectiveFluxCalculatorBase::_u_nodal_computed_by_thread
protectedinherited

_u_nodal_computed_by_thread(i) = true if _u_nodal[i] has been computed in execute() by the thread on this processor

Definition at line 186 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::execute(), AdvectiveFluxCalculatorBase::initialize(), AdvectiveFluxCalculatorBase::threadJoin(), and AdvectiveFluxCalculatorBase::timestepSetup().

◆ _valence

std::vector<unsigned> AdvectiveFluxCalculatorBase::_valence
protectedinherited

_valence[i] = number of times, in a loop over elements seen by this processor (viz, including ghost elements) and are part of the block-restricted blocks of this UserObject, that the sequential node i is encountered

Definition at line 180 of file AdvectiveFluxCalculatorBase.h.

Referenced by AdvectiveFluxCalculatorBase::getValence(), and AdvectiveFluxCalculatorBase::timestepSetup().


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