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 1828 : CellCenteredMapFunctor<T, Map>::CellCenteredMapFunctor(const MooseMesh & mesh,
93 : const std::string & name,
94 : const bool extrapolated_boundary)
95 7312 : : Moose::FunctorBase<T>(name), _mesh(mesh), _extrapolated_boundary(extrapolated_boundary)
96 : {
97 3656 : }
98 :
99 : template <typename T, typename Map>
100 4722 : 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 4722 : _mesh(mesh),
106 4722 : _sub_ids(sub_ids == mesh.meshSubdomains() ? std::set<SubdomainID>() : sub_ids),
107 18888 : _extrapolated_boundary(extrapolated_boundary)
108 : {
109 9444 : }
110 :
111 : template <typename T, typename Map>
112 : bool
113 27331760 : 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 27331760 : 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 27331760 : return extrapolated;
124 : }
125 :
126 : template <typename T, typename Map>
127 : bool
128 : CellCenteredMapFunctor<T, Map>::hasBlocks(const Elem * const elem) const
129 : {
130 131010620 : if (!elem)
131 : return false;
132 :
133 131010620 : return hasBlocks(elem->subdomain_id());
134 : }
135 :
136 : template <typename T, typename Map>
137 : bool
138 143612264 : CellCenteredMapFunctor<T, Map>::hasBlocks(const SubdomainID sub_id) const
139 : {
140 143612264 : return _sub_ids.empty() || _sub_ids.count(sub_id);
141 : }
142 :
143 : template <typename T, typename Map>
144 : typename CellCenteredMapFunctor<T, Map>::ValueType
145 357135312 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemArg & elem_arg, const StateArg &) const
146 : {
147 357135312 : const Elem * const elem = elem_arg.elem;
148 :
149 357135312 : auto it = this->find(elem->id());
150 357135312 : if (it == this->end())
151 : {
152 12 : if (!_sub_ids.empty() && !_sub_ids.count(elem->subdomain_id()))
153 6 : mooseError("Attempted to evaluate CellCenteredMapFunctor '",
154 : this->functorName(),
155 : "' with an element subdomain id of '",
156 6 : elem->subdomain_id(),
157 : "' but that subdomain id is not one of the subdomain ids the functor is "
158 : "restricted to.");
159 : else
160 6 : mooseError("Attempted access into CellCenteredMapFunctor '",
161 : this->functorName(),
162 : "' with a key that does not yet exist in the map. Make sure to fill your "
163 : "CellCenteredMapFunctor for all elements you will attempt to access later.");
164 : }
165 :
166 357135300 : return it->second;
167 : }
168 :
169 : template <typename T, typename Map>
170 : typename CellCenteredMapFunctor<T, Map>::ValueType
171 516096 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemPointArg & elem_point,
172 : const StateArg & state) const
173 : {
174 516096 : return (*this)(elem_point.makeElem(), state) +
175 516096 : (elem_point.point - elem_point.elem->vertex_average()) *
176 516096 : this->gradient(elem_point.makeElem(), state);
177 : }
178 :
179 : template <typename T, typename Map>
180 : typename CellCenteredMapFunctor<T, Map>::ValueType
181 38788216 : CellCenteredMapFunctor<T, Map>::evaluate(const FaceArg & face, const StateArg & state) const
182 : {
183 38788216 : const auto & fi = *face.fi;
184 : mooseAssert(face.limiter_type == Moose::FV::LimiterType::CentralDifference,
185 : "this implementation currently only supports linear interpolations");
186 :
187 : const bool defined_on_elem = hasBlocks(&fi.elem());
188 : const bool defined_on_neighbor = hasBlocks(fi.neighborPtr());
189 38788216 : if (defined_on_elem && defined_on_neighbor)
190 38120632 : return Moose::FV::linearInterpolation(*this, face, state);
191 :
192 667584 : if (defined_on_elem)
193 : {
194 667584 : const auto elem_arg = face.makeElem();
195 667584 : const auto elem_value = (*this)(elem_arg, state);
196 667584 : if (!_extrapolated_boundary)
197 418752 : return elem_value;
198 : // Two term expansion
199 479744 : return elem_value + this->gradient(elem_arg, state) * (fi.faceCentroid() - fi.elemCentroid());
200 : }
201 : else
202 : {
203 : mooseAssert(defined_on_neighbor, "We should be defined on one of the sides");
204 0 : const auto neighbor_arg = face.makeNeighbor();
205 0 : const auto neighbor_value = (*this)(neighbor_arg, state);
206 0 : if (!_extrapolated_boundary)
207 0 : return neighbor_value;
208 :
209 : // Two term expansion
210 : return neighbor_value +
211 0 : this->gradient(neighbor_arg, state) * (fi.faceCentroid() - fi.neighborCentroid());
212 : }
213 : }
214 :
215 : template <typename T, typename Map>
216 : typename CellCenteredMapFunctor<T, Map>::ValueType
217 0 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemQpArg &, const StateArg &) const
218 : {
219 0 : mooseError("not implemented");
220 : }
221 :
222 : template <typename T, typename Map>
223 : typename CellCenteredMapFunctor<T, Map>::ValueType
224 0 : CellCenteredMapFunctor<T, Map>::evaluate(const ElemSideQpArg &, const StateArg &) const
225 : {
226 0 : mooseError("not implemented");
227 : }
228 :
229 : template <typename T, typename Map>
230 : typename CellCenteredMapFunctor<T, Map>::ValueType
231 0 : CellCenteredMapFunctor<T, Map>::evaluate(const NodeArg &, const StateArg &) const
232 : {
233 0 : mooseError("not implemented");
234 : }
235 :
236 : template <typename T, typename Map>
237 : typename CellCenteredMapFunctor<T, Map>::GradientType
238 808284 : CellCenteredMapFunctor<T, Map>::evaluateGradient(const ElemArg & elem_arg,
239 : const StateArg & state) const
240 : {
241 808284 : return Moose::FV::greenGaussGradient(elem_arg, state, *this, _extrapolated_boundary, _mesh);
242 : }
243 :
244 : template <typename T, typename Map>
245 : typename CellCenteredMapFunctor<T, Map>::GradientType
246 2478312 : CellCenteredMapFunctor<T, Map>::evaluateGradient(const FaceArg & face, const StateArg & state) const
247 : {
248 2478312 : return Moose::FV::greenGaussGradient(face, state, *this, _extrapolated_boundary, _mesh);
249 : }
|