https://mooseframework.inl.gov
RhieChowInterpolatorBase.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 "TaggingInterface.h"
14 #include "INSFVPressureVariable.h"
15 #include "ADReal.h"
16 #include "ADFunctorInterface.h"
17 #include "INSFVPressureVariable.h"
18 #include "libmesh/vector_value.h"
19 #include "libmesh/id_types.h"
20 #include "libmesh/stored_range.h"
21 
22 class MooseMesh;
25 namespace libMesh
26 {
27 class Elem;
28 class MeshBase;
29 }
30 
32  public TaggingInterface,
33  public ADFunctorInterface
34 {
35 public:
38 
47  virtual void addToA(const libMesh::Elem * elem, unsigned int component, const ADReal & value) = 0;
48 
57  virtual VectorValue<ADReal> getVelocity(const Moose::FV::InterpMethod m,
58  const FaceInfo & fi,
59  const Moose::StateArg & time,
60  const THREAD_ID tid,
61  bool subtract_mesh_velocity) const = 0;
62 
64  const FaceInfo & fi,
65  const Moose::StateArg & time,
66  const THREAD_ID tid,
67  bool subtract_mesh_velocity) const override;
68 
71 
76  virtual void ghostADataOnBoundary(const BoundaryID /*boundary_id*/) {}
77 
81  const INSFVPressureVariable & pressure(THREAD_ID tid) const;
82 
83  const INSFVVelocityVariable * vel() const { return _u; }
84 
86  virtual bool segregated() const = 0;
87 
88 protected:
92  virtual const Moose::FunctorBase<ADReal> & epsilon(THREAD_ID tid) const;
93 
97  template <typename Container>
98  void fillContainer(const std::string & var_name, Container & container);
99 
103  template <typename VarType>
104  void checkBlocks(const VarType & var) const;
105 
108 
111 
113  const unsigned int _dim;
114 
117 
120 
123 
126 
128  std::vector<MooseVariableFVReal *> _ps;
129 
131  std::vector<MooseVariableFVReal *> _us;
132 
134  std::vector<MooseVariableFVReal *> _vs;
135 
137  std::vector<MooseVariableFVReal *> _ws;
138 
140  std::vector<unsigned int> _var_numbers;
141 
144 
146  const bool _displaced;
147 
148 private:
151 };
152 
153 inline const Moose::FunctorBase<ADReal> &
155 {
156  return _unity_functor;
157 }
158 
159 inline const INSFVPressureVariable &
161 {
162  mooseAssert(tid < _ps.size(), "Attempt to access out-of-bounds in pressure variable container");
163  return *static_cast<INSFVPressureVariable *>(_ps[tid]);
164 }
165 
166 template <typename Container>
167 void
168 RhieChowInterpolatorBase::fillContainer(const std::string & name, Container & container)
169 {
170  typedef typename Container::value_type ContainedType;
171  for (const auto tid : make_range(libMesh::n_threads()))
172  {
173  auto * const var = static_cast<ContainedType>(
174  &UserObject::_subproblem.getVariable(tid, getParam<VariableName>(name)));
175  container[tid] = var;
176  }
177 }
178 
179 template <typename VarType>
180 void
181 RhieChowInterpolatorBase::checkBlocks(const VarType & var) const
182 {
183  const auto & var_blocks = var.blockIDs();
184  const auto & uo_blocks = blockIDs();
185 
186  // Error if this UO has any blocks that the variable does not
187  std::set<SubdomainID> uo_blocks_minus_var_blocks;
188  std::set_difference(uo_blocks.begin(),
189  uo_blocks.end(),
190  var_blocks.begin(),
191  var_blocks.end(),
192  std::inserter(uo_blocks_minus_var_blocks, uo_blocks_minus_var_blocks.end()));
193  if (uo_blocks_minus_var_blocks.size() > 0)
194  mooseError("Block restriction of interpolator user object '",
195  this->name(),
196  "' (",
198  ") includes blocks not in the block restriction of variable '",
199  var.name(),
200  "' (",
201  Moose::stringify(var.blocks()),
202  ")");
203 
204  // Get the blocks in the variable but not this UO
205  std::set<SubdomainID> var_blocks_minus_uo_blocks;
206  std::set_difference(var_blocks.begin(),
207  var_blocks.end(),
208  uo_blocks.begin(),
209  uo_blocks.end(),
210  std::inserter(var_blocks_minus_uo_blocks, var_blocks_minus_uo_blocks.end()));
211 
212  // For each block in the variable but not this UO, error if there is connection
213  // to any blocks on the UO.
214  for (auto & block_id : var_blocks_minus_uo_blocks)
215  {
216  const auto connected_blocks = _moose_mesh.getBlockConnectedBlocks(block_id);
217  std::set<SubdomainID> connected_blocks_on_uo;
218  std::set_intersection(connected_blocks.begin(),
219  connected_blocks.end(),
220  uo_blocks.begin(),
221  uo_blocks.end(),
222  std::inserter(connected_blocks_on_uo, connected_blocks_on_uo.end()));
223  if (connected_blocks_on_uo.size() > 0)
224  mooseError("Block restriction of interpolator user object '",
225  this->name(),
226  "' (",
227  Moose::stringify(uo_blocks),
228  ") doesn't match the block restriction of variable '",
229  var.name(),
230  "' (",
231  Moose::stringify(var_blocks),
232  ")");
233  }
234 }
INSFVPressureVariable *const _p
The thread 0 copy of the pressure variable.
std::vector< MooseVariableFVReal * > _ps
All the thread copies of the pressure variable.
const INSFVPressureVariable & pressure(THREAD_ID tid) const
unsigned int n_threads()
MooseMesh & _moose_mesh
The MooseMesh that this user object operates on.
virtual const Moose::FunctorBase< ADReal > & epsilon(THREAD_ID tid) const
A virtual method that allows us to only implement getVelocity once for free and porous flows...
INSFVVelocityVariable *const _u
The thread 0 copy of the x-velocity variable.
INSFVVelocityVariable *const _v
The thread 0 copy of the y-velocity variable (null if the problem is 1D)
static const std::string component
Definition: NS.h:153
const std::set< SubdomainID > & getBlockConnectedBlocks(const SubdomainID subdomain_id) const
const bool _displaced
Whether this object is operating on the displaced mesh.
Moose::FV::InterpMethod velocityInterpolationMethod() const
Return the interpolation method used for velocity.
virtual const std::set< SubdomainID > & blockIDs() const
The following methods are specializations for using the Parallel::packed_range_* routines for a vecto...
virtual void ghostADataOnBoundary(const BoundaryID)
makes sure coefficient data gets communicated on both sides of a given boundary.
DualNumber< Real, DNDerivativeType, true > ADReal
SubProblem & _subproblem
virtual VectorValue< ADReal > getVelocity(const Moose::FV::InterpMethod m, const FaceInfo &fi, const Moose::StateArg &time, const THREAD_ID tid, bool subtract_mesh_velocity) const =0
Retrieve a face velocity.
std::vector< MooseVariableFVReal * > _us
All the thread copies of the x-velocity variable.
const std::string & name() const
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
RhieChowInterpolatorBase(const InputParameters &params)
virtual bool segregated() const =0
Bool of the Rhie Chow user object is used in monolithic/segregated approaches.
boundary_id_type BoundaryID
const std::string name
Definition: Setup.h:20
const unsigned int _dim
The dimension of the mesh, e.g. 3 for hexes and tets, 2 for quads and tris.
virtual const MooseVariableFieldBase & getVariable(const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const=0
std::string stringify(const T &t)
Moose::FV::InterpMethod _velocity_interp_method
The interpolation method to use for the velocity.
void fillContainer(const std::string &var_name, Container &container)
Fill the passed-in variable container with the thread copies of var_name.
const libMesh::MeshBase & _mesh
The libMesh mesh that this object acts on.
virtual Real getVolumetricFaceFlux(const Moose::FV::InterpMethod m, const FaceInfo &fi, const Moose::StateArg &time, const THREAD_ID tid, bool subtract_mesh_velocity) const override
Retrieve the volumetric face flux, will not include derivatives.
void checkBlocks(const VarType &var) const
Check the block consistency between the passed in var and us.
const INSFVVelocityVariable * vel() const
const std::vector< SubdomainName > & blocks() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void addToA(const libMesh::Elem *elem, unsigned int component, const ADReal &value)=0
API for momentum residual objects that have on-diagonals for velocity call.
INSFVVelocityVariable *const _w
The thread 0 copy of the z-velocity variable (null if the problem is not 3D)
std::vector< unsigned int > _var_numbers
The velocity variable numbers.
IntRange< T > make_range(T beg, T end)
void mooseError(Args &&... args) const
std::vector< MooseVariableFVReal * > _ws
All the thread copies of the z-velocity variable.
const Moose::ConstantFunctor< ADReal > _unity_functor
A unity functor used in the epsilon virtual method.
std::vector< MooseVariableFVReal * > _vs
All the thread copies of the y-velocity variable.
static InputParameters validParams()
unsigned int THREAD_ID