Line data Source code
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 :
12 : #include "MooseFunctor.h"
13 : #include "GreenGaussGradient.h"
14 : #include "MathFVUtils.h"
15 : #include "libmesh/utility.h"
16 : #include "libmesh/type_tensor.h"
17 : #include "libmesh/compare_types.h"
18 : #include "libmesh/threads.h"
19 :
20 : /**
21 : * A functor whose evaluation relies on querying a map where the keys are element ids and the values
22 : * correspond to the element/cell values. This is a very useful data type for storing the result of
23 : * (possibly repeated) interpolation and reconstruction
24 : */
25 : template <typename T, typename Map>
26 : class CellCenteredMapFunctor : public Moose::FunctorBase<T>, public Map
27 : {
28 : public:
29 : using typename Moose::FunctorBase<T>::ValueType;
30 : using typename Moose::FunctorBase<T>::GradientType;
31 : using typename Moose::FunctorBase<T>::DotType;
32 : using ElemArg = Moose::ElemArg;
33 : using FaceArg = Moose::FaceArg;
34 : using ElemQpArg = Moose::ElemQpArg;
35 : using ElemSideQpArg = Moose::ElemSideQpArg;
36 : using ElemPointArg = Moose::ElemPointArg;
37 : using StateArg = Moose::StateArg;
38 : using NodeArg = Moose::NodeArg;
39 :
40 : /**
41 : * Use this constructor when you want the object to live everywhere on the mesh
42 : */
43 : CellCenteredMapFunctor(const MooseMesh & mesh,
44 : const std::string & name,
45 : const bool extrapolated_boundary);
46 :
47 : /**
48 : * Use this constructor if you want to potentially restrict this object to a specified set of
49 : * subdomains/blocks
50 : */
51 : CellCenteredMapFunctor(const MooseMesh & mesh,
52 : const std::set<SubdomainID> & sub_ids,
53 : const std::string & name,
54 : const bool extrapolated_boundary);
55 :
56 : bool isExtrapolatedBoundaryFace(const FaceInfo & fi,
57 : const Elem * elem,
58 : const StateArg & state) const override;
59 : bool hasBlocks(SubdomainID sub_id) const override;
60 :
61 : /**
62 : * Checks whether we are defined on the provided element
63 : */
64 : bool hasBlocks(const Elem * elem) const;
65 :
66 0 : bool supportsFaceArg() const override final { return true; }
67 0 : bool supportsElemSideQpArg() const override final { return false; }
68 :
69 : private:
70 : /// The mesh that this functor lives on
71 : const MooseMesh & _mesh;
72 :
73 : /// The subdomain IDs that this functor lives on. If empty, then we consider the functor to live
74 : /// on all subdomains
75 : const std::set<SubdomainID> _sub_ids;
76 :
77 : const bool _extrapolated_boundary;
78 :
79 : ValueType evaluate(const ElemArg & elem_arg, const StateArg &) const override;
80 : ValueType evaluate(const ElemPointArg & elem_point, const StateArg & state) const override;
81 : ValueType evaluate(const FaceArg & face, const StateArg &) const override;
82 : ValueType evaluate(const ElemQpArg &, const StateArg &) const override;
83 : ValueType evaluate(const ElemSideQpArg &, const StateArg &) const override;
84 : ValueType evaluate(const NodeArg & elem_arg, const StateArg &) const override;
85 :
86 : using Moose::FunctorBase<T>::evaluateGradient;
87 : GradientType evaluateGradient(const ElemArg & elem_arg, const StateArg & state) const override;
88 : GradientType evaluateGradient(const FaceArg & face, const StateArg & state) const override;
89 : };
90 :
91 : template <typename T, typename Map>
92 2976 : CellCenteredMapFunctor<T, Map>::CellCenteredMapFunctor(const MooseMesh & mesh,
93 : const std::string & name,
94 : const bool extrapolated_boundary)
95 11904 : : Moose::FunctorBase<T>(name), _mesh(mesh), _extrapolated_boundary(extrapolated_boundary)
96 : {
97 5952 : }
98 :
99 : template <typename T, typename Map>
100 7974 : CellCenteredMapFunctor<T, Map>::CellCenteredMapFunctor(const MooseMesh & mesh,
101 : const std::set<SubdomainID> & sub_ids,
102 : const std::string & name,
103 : const bool extrapolated_boundary)
104 : : Moose::FunctorBase<T>(name),
105 7974 : _mesh(mesh),
106 7974 : _sub_ids(sub_ids == mesh.meshSubdomains() ? std::set<SubdomainID>() : sub_ids),
107 31896 : _extrapolated_boundary(extrapolated_boundary)
108 : {
109 15948 : }
110 :
111 : template <typename T, typename Map>
112 : bool
113 41090736 : CellCenteredMapFunctor<T, Map>::isExtrapolatedBoundaryFace(const FaceInfo & fi,
114 : const Elem *,
115 : const StateArg &) const
116 : {
117 : const bool defined_on_elem = hasBlocks(&fi.elem());
118 : const bool defined_on_neighbor = hasBlocks(fi.neighborPtr());
119 41090736 : const bool extrapolated = (defined_on_elem + defined_on_neighbor) == 1;
120 :
121 : mooseAssert(defined_on_elem || defined_on_neighbor,
122 : "This shouldn't be called if we aren't defined on either side.");
123 41090736 : return extrapolated;
124 : }
125 :
126 : template <typename T, typename Map>
127 : bool
128 : CellCenteredMapFunctor<T, Map>::hasBlocks(const Elem * const elem) const
129 : {
130 198167830 : if (!elem)
131 : return false;
132 :
133 198167830 : return hasBlocks(elem->subdomain_id());
134 : }
135 :
136 : template <typename T, typename Map>
137 : bool
138 217308282 : CellCenteredMapFunctor<T, Map>::hasBlocks(const SubdomainID sub_id) const
139 : {
140 217308282 : return _sub_ids.empty() || _sub_ids.count(sub_id);
141 : }
142 :
143 : template <typename T, typename Map>
144 : typename CellCenteredMapFunctor<T, Map>::ValueType
145 509577844 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemArg & elem_arg, const StateArg &) const
146 : {
147 509577844 : const Elem * const elem = elem_arg.elem;
148 :
149 : try
150 : {
151 806987818 : return libmesh_map_find(*this, elem->id());
152 : }
153 12 : catch (libMesh::LogicError &)
154 : {
155 12 : if (!_sub_ids.empty() && !_sub_ids.count(elem->subdomain_id()))
156 6 : mooseError("Attempted to evaluate CellCenteredMapFunctor '",
157 : this->functorName(),
158 : "' with an element subdomain id of '",
159 12 : elem->subdomain_id(),
160 : "' but that subdomain id is not one of the subdomain ids the functor is "
161 : "restricted to.");
162 : else
163 6 : mooseError("Attempted access into CellCenteredMapFunctor '",
164 : this->functorName(),
165 : "' with a key that does not yet exist in the map. Make sure to fill your "
166 : "CellCenteredMapFunctor for all elements you will attempt to access later.");
167 : }
168 : }
169 :
170 : template <typename T, typename Map>
171 : typename CellCenteredMapFunctor<T, Map>::ValueType
172 516096 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemPointArg & elem_point,
173 : const StateArg & state) const
174 : {
175 516096 : return (*this)(elem_point.makeElem(), state) +
176 516096 : (elem_point.point - elem_point.elem->vertex_average()) *
177 516096 : this->gradient(elem_point.makeElem(), state);
178 : }
179 :
180 : template <typename T, typename Map>
181 : typename CellCenteredMapFunctor<T, Map>::ValueType
182 58895146 : CellCenteredMapFunctor<T, Map>::evaluate(const FaceArg & face, const StateArg & state) const
183 : {
184 58895146 : const auto & fi = *face.fi;
185 : mooseAssert(face.limiter_type == Moose::FV::LimiterType::CentralDifference,
186 : "this implementation currently only supports linear interpolations");
187 :
188 : const bool defined_on_elem = hasBlocks(&fi.elem());
189 : const bool defined_on_neighbor = hasBlocks(fi.neighborPtr());
190 58895146 : if (defined_on_elem && defined_on_neighbor)
191 57967280 : return Moose::FV::linearInterpolation(*this, face, state);
192 :
193 927866 : if (defined_on_elem)
194 : {
195 927866 : const auto elem_arg = face.makeElem();
196 927866 : const auto elem_value = (*this)(elem_arg, state);
197 927866 : if (!_extrapolated_boundary)
198 526010 : return elem_value;
199 : // Two term expansion
200 785792 : return elem_value + this->gradient(elem_arg, state) * (fi.faceCentroid() - fi.elemCentroid());
201 : }
202 : else
203 : {
204 : mooseAssert(defined_on_neighbor, "We should be defined on one of the sides");
205 0 : const auto neighbor_arg = face.makeNeighbor();
206 0 : const auto neighbor_value = (*this)(neighbor_arg, state);
207 0 : if (!_extrapolated_boundary)
208 0 : return neighbor_value;
209 :
210 : // Two term expansion
211 : return neighbor_value +
212 0 : this->gradient(neighbor_arg, state) * (fi.faceCentroid() - fi.neighborCentroid());
213 : }
214 : }
215 :
216 : template <typename T, typename Map>
217 : typename CellCenteredMapFunctor<T, Map>::ValueType
218 0 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemQpArg &, const StateArg &) const
219 : {
220 0 : mooseError("not implemented");
221 : }
222 :
223 : template <typename T, typename Map>
224 : typename CellCenteredMapFunctor<T, Map>::ValueType
225 0 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemSideQpArg &, const StateArg &) const
226 : {
227 0 : mooseError("not implemented");
228 : }
229 :
230 : template <typename T, typename Map>
231 : typename CellCenteredMapFunctor<T, Map>::ValueType
232 0 : CellCenteredMapFunctor<T, Map>::evaluate(const NodeArg &, const StateArg &) const
233 : {
234 0 : mooseError("not implemented");
235 : }
236 :
237 : template <typename T, typename Map>
238 : typename CellCenteredMapFunctor<T, Map>::GradientType
239 990012 : CellCenteredMapFunctor<T, Map>::evaluateGradient(const ElemArg & elem_arg,
240 : const StateArg & state) const
241 : {
242 990012 : return Moose::FV::greenGaussGradient(elem_arg, state, *this, _extrapolated_boundary, _mesh);
243 : }
244 :
245 : template <typename T, typename Map>
246 : typename CellCenteredMapFunctor<T, Map>::GradientType
247 4107320 : CellCenteredMapFunctor<T, Map>::evaluateGradient(const FaceArg & face, const StateArg & state) const
248 : {
249 4107320 : return Moose::FV::greenGaussGradient(face, state, *this, _extrapolated_boundary, _mesh);
250 : }
|
