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 : #include "FunctorSmoother.h"
11 : #include "MooseUtils.h"
12 : #include "RelationshipManager.h"
13 :
14 : #include "metaphysicl/raw_type.h"
15 :
16 : registerMooseObject("MooseApp", FunctorSmoother);
17 :
18 : template <typename T>
19 : InputParameters
20 15068 : FunctorSmootherTempl<T>::validParams()
21 : {
22 15068 : InputParameters params = FunctorMaterial::validParams();
23 30136 : params.addClassDescription("Creates smoother functor(s) using various averaging techniques");
24 60272 : params.addParam<std::vector<MooseFunctorName>>("functors_in",
25 : "The name(s) of the functors to smooth");
26 60272 : params.addParam<std::vector<MooseFunctorName>>("functors_out",
27 : "The name(s) of the smooth output functors");
28 60272 : MooseEnum smoothing_techniques("face_average layered_elem_average remove_checkerboard");
29 45204 : params.addParam<MooseEnum>(
30 : "smoothing_technique", smoothing_techniques, "How to smooth the functor");
31 :
32 30136 : return params;
33 15068 : }
34 :
35 : template <typename T>
36 42 : FunctorSmootherTempl<T>::FunctorSmootherTempl(const InputParameters & parameters)
37 : : FunctorMaterial(parameters),
38 42 : _functors_in(getParam<std::vector<MooseFunctorName>>("functors_in")),
39 84 : _functors_out(getParam<std::vector<MooseFunctorName>>("functors_out")),
40 126 : _smoothing_technique(getParam<MooseEnum>("smoothing_technique"))
41 : {
42 42 : if (_tid == 0)
43 : {
44 : // We need one layer of ghosting at least to get neighbor values
45 39 : auto & factory = _app.getFactory();
46 39 : auto rm_params = factory.getValidParams("ElementSideNeighborLayers");
47 78 : rm_params.template set<std::string>("for_whom") = name();
48 78 : rm_params.template set<MooseMesh *>("mesh") = &const_cast<MooseMesh &>(_mesh);
49 78 : rm_params.template set<Moose::RelationshipManagerType>("rm_type") =
50 39 : Moose::RelationshipManagerType::GEOMETRIC | Moose::RelationshipManagerType::ALGEBRAIC |
51 : Moose::RelationshipManagerType::COUPLING;
52 78 : rm_params.template set<unsigned short>("layers") = 1;
53 39 : rm_params.template set<bool>("use_point_neighbors") = false;
54 39 : rm_params.template set<bool>("use_displaced_mesh") =
55 39 : parameters.template get<bool>("use_displaced_mesh");
56 : mooseAssert(rm_params.areAllRequiredParamsValid(),
57 : "All relationship manager parameters should be valid.");
58 117 : auto rm_obj = factory.template create<RelationshipManager>(
59 39 : "ElementSideNeighborLayers", name() + "_functor_smoothing", rm_params);
60 :
61 : // Delete the resources created on behalf of the RM if it ends up not being added to the
62 : // App.
63 39 : if (!_app.addRelationshipManager(rm_obj))
64 39 : factory.releaseSharedObjects(*rm_obj);
65 39 : }
66 :
67 42 : const std::set<ExecFlagType> clearance_schedule(_execute_enum.begin(), _execute_enum.end());
68 :
69 84 : for (const auto i : index_range(_functors_in))
70 : {
71 : // This potentially upcasts functors from non-AD to AD
72 42 : const auto & functor_in = getFunctor<typename Moose::ADType<T>::type>(_functors_in[i]);
73 :
74 : // We always add the AD type
75 42 : addFunctorProperty<typename Moose::ADType<T>::type>(
76 : _functors_out[i],
77 2700 : [this, &functor_in](const auto & r, const auto & t) -> typename Moose::ADType<T>::type
78 : {
79 2700 : typename Moose::ADType<T>::type base_value = functor_in(r, t);
80 2700 : typename Moose::ADType<T>::type average = 0;
81 : using std::abs;
82 :
83 2700 : const Elem * r_elem = nullptr;
84 : if constexpr (std::is_same_v<const Moose::ElemArg &, decltype(r)>)
85 2700 : r_elem = r.elem;
86 :
87 : // Handle the single isolated element case
88 2700 : if (r_elem && r_elem->n_neighbors() == 0)
89 0 : return base_value;
90 :
91 : // If these techniques are of interest elsewhere, they should be moved to a
92 : // FunctorAveragingUtils file and namespaced
93 2700 : if (_smoothing_technique == FACE_AVERAGE)
94 : {
95 : if constexpr (std::is_same_v<const Moose::ElemArg &, decltype(r)>)
96 : {
97 900 : unsigned int n_faces = 0;
98 4500 : for (const auto side_index : r_elem->side_index_range())
99 : {
100 3600 : auto fi = _mesh.faceInfo(r_elem, side_index);
101 3600 : if (!fi)
102 : fi =
103 891 : _mesh.faceInfo(r_elem->neighbor_ptr(side_index),
104 : r_elem->neighbor_ptr(side_index)->which_neighbor_am_i(r_elem));
105 3600 : Moose::FaceArg face_arg{
106 : fi, Moose::FV::LimiterType::CentralDifference, true, false, nullptr, nullptr};
107 3600 : if (face_arg.fi)
108 : {
109 3600 : average += functor_in(face_arg, t);
110 3600 : n_faces++;
111 : }
112 : }
113 900 : average /= n_faces;
114 : }
115 : else
116 : // a conceptually simple option here would be to use the smoothed version of the
117 : // ElemArg to compute all the other functor arguments. Maybe with the same smoothing
118 : // on the gradient calls
119 0 : mooseError("Face averaging smoothing has only been defined for the ElemArg functor "
120 : "argument, not for ",
121 : MooseUtils::prettyCppType(&r),
122 : ". Please contact a MOOSE developer or implement it yourself.");
123 : }
124 1800 : else if (_smoothing_technique == LAYERED_AVERAGE)
125 : {
126 : if constexpr (std::is_same_v<const Moose::ElemArg &, decltype(r)>)
127 : {
128 900 : unsigned int n_neighbors = 0;
129 4500 : for (const auto neighbor : r_elem->neighbor_ptr_range())
130 3600 : if (neighbor)
131 : {
132 1782 : n_neighbors++;
133 1782 : average += functor_in(Moose::ElemArg{neighbor, false}, t);
134 : }
135 900 : average /= n_neighbors;
136 : }
137 : else
138 0 : mooseError("Element layered averaging smoothing has only been defined for the "
139 : "ElemArg functor argument, not for ",
140 : MooseUtils::prettyCppType(&r),
141 : ". Please contact a MOOSE developer or implement it yourself.");
142 : }
143 : else
144 : {
145 : if constexpr (std::is_same_v<const Moose::ElemArg &, decltype(r)>)
146 : {
147 : // We average the local value with the average of the two values furthest away
148 : // among the neighbors. For a checkerboarding but overall linear evolution, this will
149 : // compute the average of the "high" and "low" profiles
150 900 : average += base_value / 2;
151 :
152 : // Find the neighbor with the furthest value
153 900 : typename Moose::ADType<T>::type extreme_value = 0;
154 900 : typename Moose::ADType<T>::type delta = 0;
155 :
156 900 : unsigned int furthest_one = 0;
157 900 : unsigned int num_neighbors = 0;
158 4500 : for (const auto side_index : r_elem->side_index_range())
159 : {
160 3600 : auto neighbor = r_elem->neighbor_ptr(side_index);
161 3600 : if (neighbor)
162 : {
163 1782 : num_neighbors++;
164 1782 : auto neighbor_value = functor_in(Moose::ElemArg{neighbor, false}, t);
165 1782 : if (abs(neighbor_value - base_value) > delta)
166 : {
167 900 : furthest_one = side_index;
168 900 : extreme_value = neighbor_value;
169 900 : delta = abs(neighbor_value - base_value);
170 : }
171 1782 : }
172 : }
173 :
174 : // We're on a boundary in 1D, or maybe an odd shape corner in 2D
175 900 : if (num_neighbors == 1)
176 : {
177 18 : average += extreme_value / 2;
178 18 : return average;
179 : }
180 : else
181 882 : average += extreme_value / 4;
182 :
183 : // Get the value from the neighbor opposite the furthest-value neighbor
184 : try
185 : {
186 882 : auto opposite_side = r_elem->opposite_side(furthest_one);
187 882 : auto neighbor = r_elem->neighbor_ptr(opposite_side);
188 882 : if (neighbor)
189 : {
190 882 : average += functor_in(Moose::ElemArg{neighbor, false}, t) / 4;
191 : }
192 : else
193 : // We're probably at a boundary
194 0 : average += extreme_value / 4;
195 : }
196 : // if there is no opposite side (for example a triangle)
197 0 : catch (libMesh::LogicError & e)
198 : {
199 : // find the second furthest
200 0 : delta = 0;
201 0 : typename Moose::ADType<T>::type second_extreme_value = 0;
202 :
203 0 : for (const auto side_index : r_elem->side_index_range())
204 : {
205 0 : auto neighbor = r_elem->neighbor_ptr(side_index);
206 0 : auto neighbor_value =
207 0 : neighbor ? functor_in(Moose::ElemArg{neighbor, false}, t) : base_value;
208 0 : if (abs(neighbor_value - base_value) > delta && side_index != furthest_one)
209 : {
210 0 : second_extreme_value = neighbor_value;
211 0 : delta = abs(neighbor_value - base_value);
212 : }
213 : }
214 0 : average += second_extreme_value / 4;
215 0 : }
216 918 : }
217 : else
218 0 : mooseError("Checkerboard removal smoothing has only been defined for the "
219 : "ElemArg functor argument, not for ",
220 : MooseUtils::prettyCppType(&r),
221 : ". Please contact a MOOSE developer or implement it yourself.");
222 : }
223 2682 : return average;
224 2700 : },
225 : clearance_schedule);
226 : }
227 42 : }
228 :
229 : template class FunctorSmootherTempl<Real>;
|
