https://mooseframework.inl.gov
RhieChowMassFlux.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
13 #include "CellCenteredMapFunctor.h"
14 #include "FaceCenteredMapFunctor.h"
15 #include "VectorComponentFunctor.h"
17 #include <unordered_map>
18 #include <set>
19 #include <unordered_set>
20 
21 #include "libmesh/petsc_vector.h"
22 
23 class MooseMesh;
26 namespace libMesh
27 {
28 class Elem;
29 class MeshBase;
30 }
31 
37 {
38 public:
40  RhieChowMassFlux(const InputParameters & params);
41 
43  Real getMassFlux(const FaceInfo & fi) const;
44 
46  Real getVolumetricFaceFlux(const FaceInfo & fi) const;
47 
49  const FaceInfo & fi,
50  const Moose::StateArg & time,
51  const THREAD_ID tid,
52  bool subtract_mesh_velocity) const override;
53 
55  void initFaceMassFlux();
57  void initCouplingField();
59  void computeFaceMassFlux();
61  void computeCellVelocity();
62 
63  virtual void meshChanged() override;
64  virtual void initialize() override;
65  virtual void execute() override {}
66  virtual void finalize() override {}
67  virtual void initialSetup() override;
68 
76  void linkMomentumPressureSystems(const std::vector<LinearSystem *> & momentum_systems,
77  const LinearSystem & pressure_system,
78  const std::vector<unsigned int> & momentum_system_numbers);
79 
84  void computeHbyA(const bool with_updated_pressure, const bool verbose);
85 
86 protected:
88  std::vector<std::unique_ptr<NumericVector<Number>>> &
89  selectPressureGradient(const bool updated_pressure);
90 
92  void setupMeshInformation();
93 
95  void
96  populateCouplingFunctors(const std::vector<std::unique_ptr<NumericVector<Number>>> & raw_hbya,
97  const std::vector<std::unique_ptr<NumericVector<Number>>> & raw_Ainv);
98 
102  template <typename VarType>
103  void checkBlocks(const VarType & var) const;
104 
105  virtual bool supportMeshVelocity() const override { return false; }
106 
109 
112 
114  const unsigned int _dim;
115 
118 
120  std::vector<const MooseLinearVariableFVReal *> _vel;
121 
124 
132 
137  std::vector<std::unique_ptr<NumericVector<Number>>> _HbyA_raw;
138 
144 
149  std::vector<std::unique_ptr<NumericVector<Number>>> _Ainv_raw;
150 
155 
160  std::vector<std::unique_ptr<NumericVector<Number>>> _grad_p_current;
161 
166 
168  std::vector<LinearSystem *> _momentum_systems;
169 
171  std::vector<unsigned int> _momentum_system_numbers;
172 
174  std::vector<unsigned int> _global_momentum_system_numbers;
175 
177  std::vector<libMesh::LinearImplicitSystem *> _momentum_implicit_systems;
178 
181 
184 
186  std::unique_ptr<NumericVector<Number>> _cell_volumes;
187 
190 
191 private:
194  std::vector<const FaceInfo *> _flow_face_info;
195 };
196 
197 template <typename VarType>
198 void
199 RhieChowMassFlux::checkBlocks(const VarType & var) const
200 {
201  const auto & var_blocks = var.blockIDs();
202  const auto & uo_blocks = blockIDs();
203 
204  // Error if this UO has any blocks that the variable does not
205  std::set<SubdomainID> uo_blocks_minus_var_blocks;
206  std::set_difference(uo_blocks.begin(),
207  uo_blocks.end(),
208  var_blocks.begin(),
209  var_blocks.end(),
210  std::inserter(uo_blocks_minus_var_blocks, uo_blocks_minus_var_blocks.end()));
211  if (uo_blocks_minus_var_blocks.size() > 0)
212  mooseError("Block restriction of interpolator user object '",
213  this->name(),
214  "' (",
216  ") includes blocks not in the block restriction of variable '",
217  var.name(),
218  "' (",
219  Moose::stringify(var.blocks()),
220  ")");
221 
222  // Get the blocks in the variable but not this UO
223  std::set<SubdomainID> var_blocks_minus_uo_blocks;
224  std::set_difference(var_blocks.begin(),
225  var_blocks.end(),
226  uo_blocks.begin(),
227  uo_blocks.end(),
228  std::inserter(var_blocks_minus_uo_blocks, var_blocks_minus_uo_blocks.end()));
229 
230  // For each block in the variable but not this UO, error if there is connection
231  // to any blocks on the UO.
232  for (auto & block_id : var_blocks_minus_uo_blocks)
233  {
234  const auto connected_blocks = _moose_mesh.getBlockConnectedBlocks(block_id);
235  std::set<SubdomainID> connected_blocks_on_uo;
236  std::set_intersection(connected_blocks.begin(),
237  connected_blocks.end(),
238  uo_blocks.begin(),
239  uo_blocks.end(),
240  std::inserter(connected_blocks_on_uo, connected_blocks_on_uo.end()));
241  if (connected_blocks_on_uo.size() > 0)
242  mooseError("Block restriction of interpolator user object '",
243  this->name(),
244  "' (",
245  Moose::stringify(uo_blocks),
246  ") doesn't match the block restriction of variable '",
247  var.name(),
248  "' (",
249  Moose::stringify(var_blocks),
250  ")");
251  }
252 }
std::vector< LinearSystem * > _momentum_systems
Pointers to the linear system(s) in moose corresponding to the momentum equation(s) ...
const MooseEnum _pressure_projection_method
Enumerator for the method used for pressure projection.
virtual void initialize() override
void setupMeshInformation()
Compute the cell volumes on the mesh.
User object responsible for determining the face fluxes using the Rhie-Chow interpolation in a segreg...
A functor whose evaluation relies on querying a map where the keys are face info ids and the values c...
void checkBlocks(const VarType &var) const
Check the block consistency between the passed in var and us.
std::vector< libMesh::LinearImplicitSystem * > _momentum_implicit_systems
Pointers to the momentum equation implicit system(s) from libmesh.
Real getMassFlux(const FaceInfo &fi) const
Get the face velocity times density (used in advection terms)
void computeFaceMassFlux()
Update the values of the face velocities in the containers.
const std::set< SubdomainID > & getBlockConnectedBlocks(const SubdomainID subdomain_id) const
std::vector< std::unique_ptr< NumericVector< Number > > > _Ainv_raw
We hold on to the cell-based 1/A vectors so that we can easily reconstruct the cell velocities as wel...
const LinearSystem * _pressure_system
Pointer to the pressure system.
FaceCenteredMapFunctor< RealVectorValue, std::unordered_map< dof_id_type, RealVectorValue > > _Ainv
A map functor from faces to $(1/A)_f$.
virtual const std::set< SubdomainID > & blockIDs() const
The following methods are specializations for using the Parallel::packed_range_* routines for a vecto...
std::vector< unsigned int > _momentum_system_numbers
Numbers of the momentum system(s)
virtual void initialSetup() override
virtual void meshChanged() override
virtual const std::string & name() const
std::vector< unsigned int > _global_momentum_system_numbers
Global numbers of the momentum system(s)
void computeHbyA(const bool with_updated_pressure, const bool verbose)
Computes the inverse of the diagonal (1/A) of the system matrix plus the H/A components for the press...
FaceCenteredMapFunctor< Real, std::unordered_map< dof_id_type, Real > > _HbyA_flux
A map functor from faces to $HbyA_{ij} = (A_{offdiag}*{(predicted~velocity)} - {Source})_{ij}/A_{ij}$...
virtual bool supportMeshVelocity() const override
Returns whether the UO can support mesh velocity advection.
std::vector< const MooseLinearVariableFVReal * > _vel
The thread 0 copy of the x-velocity variable.
const MooseLinearVariableFVReal *const _p
The thread 0 copy of the pressure variable.
unsigned int _global_pressure_system_number
Global number of the pressure system.
RhieChowMassFlux(const InputParameters &params)
void linkMomentumPressureSystems(const std::vector< LinearSystem *> &momentum_systems, const LinearSystem &pressure_system, const std::vector< unsigned int > &momentum_system_numbers)
Update the momentum system-related information.
const Moose::Functor< Real > & _rho
Functor describing the density of the fluid.
static InputParameters validParams()
void computeCellVelocity()
Update the cell values of the velocity variables.
const MooseMesh & _moose_mesh
The MooseMesh that this user object operates on.
Real getVolumetricFaceFlux(const FaceInfo &fi) const
Get the volumetric face flux (used in advection terms)
std::string stringify(const T &t)
std::vector< std::unique_ptr< NumericVector< Number > > > & selectPressureGradient(const bool updated_pressure)
Select the right pressure gradient field and return a reference to the container. ...
virtual void execute() override
FaceCenteredMapFunctor< Real, std::unordered_map< dof_id_type, Real > > & _face_mass_flux
A map functor from faces to mass fluxes which are used in the advection terms.
std::vector< std::unique_ptr< NumericVector< Number > > > _grad_p_current
for a PISO iteration we need to hold on to the original pressure gradient field.
const std::vector< SubdomainName > & blocks() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void mooseError(Args &&... args) const
std::vector< const FaceInfo * > _flow_face_info
The subset of the FaceInfo objects that actually cover the subdomains which the flow field is defined...
void populateCouplingFunctors(const std::vector< std::unique_ptr< NumericVector< Number >>> &raw_hbya, const std::vector< std::unique_ptr< NumericVector< Number >>> &raw_Ainv)
Populate the face values of the H/A and 1/A fields.
LinearFVAnisotropicDiffusion * _p_diffusion_kernel
Pointer to the pressure diffusion term in the pressure Poisson equation.
std::vector< std::unique_ptr< NumericVector< Number > > > _HbyA_raw
We hold on to the cell-based HbyA vectors so that we can easily reconstruct the cell velocities as we...
void initCouplingField()
Initialize the coupling fields (HbyA and Ainv)
unsigned int THREAD_ID
const libMesh::MeshBase & _mesh
The libMesh mesh that this object acts on.
void initFaceMassFlux()
Initialize the container for face velocities.
std::unique_ptr< NumericVector< Number > > _cell_volumes
We will hold a vector of cell volumes to make sure we can do volume corrections rapidly.
virtual void finalize() override
const unsigned int _dim
The dimension of the mesh, e.g. 3 for hexes and tets, 2 for quads and tris.