www.mooseframework.org
PorousFlowAdvectiveFluxCalculatorBase.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
13 #include "PorousFlowDictator.h"
14 
16 
17 template <>
19 
33 {
34 public:
35  PorousFlowAdvectiveFluxCalculatorBase(const InputParameters & parameters);
36 
42  const std::map<dof_id_type, std::vector<Real>> & getdFluxOut_dvars(unsigned node_id) const;
43 
44 protected:
45  virtual void timestepSetup() override;
46  virtual void initialize() override;
47  virtual void execute() override;
48  virtual void finalize() override;
49  virtual void threadJoin(const UserObject & uo) override;
50 
51  virtual Real computeVelocity(unsigned i, unsigned j, unsigned qp) const override;
52 
53  virtual void executeOnElement(dof_id_type global_i,
54  dof_id_type global_j,
55  unsigned local_i,
56  unsigned local_j,
57  unsigned qp) override;
58 
64  virtual Real computedU_dvar(unsigned i, unsigned pvar) const = 0;
65 
71  const std::map<dof_id_type, std::vector<Real>> & getdK_dvar(dof_id_type node_i,
72  dof_id_type node_j) const;
73 
74  virtual void buildCommLists() override;
75  virtual void exchangeGhostedInfo() override;
76 
79 
81  const unsigned _num_vars;
82 
84  const RealVectorValue _gravity;
85 
87  const unsigned int _phase;
88 
90  const MaterialProperty<RealTensorValue> & _permeability;
91 
93  const MaterialProperty<std::vector<RealTensorValue>> & _dpermeability_dvar;
94 
96  const MaterialProperty<std::vector<std::vector<RealTensorValue>>> & _dpermeability_dgradvar;
97 
99  const MaterialProperty<std::vector<Real>> & _fluid_density_qp;
100 
102  const MaterialProperty<std::vector<std::vector<Real>>> & _dfluid_density_qp_dvar;
103 
105  const MaterialProperty<std::vector<RealGradient>> & _grad_p;
106 
108  const MaterialProperty<std::vector<std::vector<Real>>> & _dgrad_p_dgrad_var;
109 
111  const MaterialProperty<std::vector<std::vector<RealGradient>>> & _dgrad_p_dvar;
112 
114  const FEType _fe_type;
115 
117  const VariablePhiValue & _phi;
118 
120  const VariablePhiGradient & _grad_phi;
121 
123  std::vector<std::vector<Real>> _du_dvar;
124 
126  std::vector<bool> _du_dvar_computed_by_thread;
127 
133  std::vector<std::vector<std::map<dof_id_type, std::vector<Real>>>> _dkij_dvar;
134 
136  std::vector<std::map<dof_id_type, std::vector<Real>>> _dflux_out_dvars;
137 
147  std::map<processor_id_type, std::vector<dof_id_type>> _triples_to_receive;
148 
158  std::map<processor_id_type, std::vector<dof_id_type>> _triples_to_send;
159 };
virtual Real computedU_dvar(unsigned i, unsigned pvar) const =0
Compute d(u)/d(porous_flow_variable)
virtual void buildCommLists() override
When using multiple processors, other processors will compute:
const MaterialProperty< RealTensorValue > & _permeability
Permeability of porous material.
const VariablePhiGradient & _grad_phi
grad(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) override
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< std::vector< RealTensorValue > > > & _dpermeability_dgradvar
d(permeabiity)/d(grad(PorousFlow variable))
const std::map< dof_id_type, std::vector< Real > > & getdFluxOut_dvars(unsigned node_id) const
Returns d(flux_out)/d(porous_flow_variables.
const VariablePhiValue & _phi
Kuzmin-Turek shape function.
const MaterialProperty< std::vector< RealTensorValue > > & _dpermeability_dvar
d(permeabiity)/d(PorousFlow variable)
Base class to compute the advective flux of fluid in PorousFlow situations using the Kuzmin-Turek FEM...
const MaterialProperty< std::vector< RealGradient > > & _grad_p
Gradient of the pore pressure in each phase.
std::map< processor_id_type, std::vector< dof_id_type > > _triples_to_receive
_triples_to_receive[proc_id] indicates the dk(i, j)/du_nodal information that we will receive from pr...
PorousFlowAdvectiveFluxCalculatorBase(const InputParameters &parameters)
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 cur...
virtual void threadJoin(const UserObject &uo) override
virtual void exchangeGhostedInfo() override
Sends and receives multi-processor information regarding u_nodal and k_ij.
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...
Base class to compute Advective fluxes.
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]) ...
std::map< processor_id_type, std::vector< dof_id_type > > _triples_to_send
_triples_to_send[proc_id] indicates the dk(i, j)/du_nodal information that we will send to proc_id...
InputParameters validParams< PorousFlowAdvectiveFluxCalculatorBase >()
const MaterialProperty< std::vector< std::vector< Real > > > & _dgrad_p_dgrad_var
Derivative of Grad porepressure in each phase wrt grad(PorousFlow variables)
This holds maps between the nonlinear variables used in a PorousFlow simulation and the variable numb...
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...
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...