www.mooseframework.org
Legendre.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 "Legendre.h"
11 #include <functional>
12 
17 #define MAX_DIRECT_CALCULATION_LEGENDRE 12
18 
19 Legendre::Legendre() : SingleSeriesBasisInterface() {}
20 
21 Legendre::Legendre(const std::vector<MooseEnum> & domain,
22  const std::vector<std::size_t> & order,
23  MooseEnum expansion_type,
24  MooseEnum generation_type)
25  : SingleSeriesBasisInterface(domain, order, calculatedNumberOfTermsBasedOnOrder(order))
26 {
27  if (expansion_type == "orthonormal")
28  _evaluateExpansionWrapper = [this]() { this->evaluateOrthonormal(); };
29  else if (expansion_type == "sqrt_mu")
30  _evaluateExpansionWrapper = [this]() { this->evaluateSqrtMu(); };
31  else if (expansion_type == "standard")
32  _evaluateExpansionWrapper = [this]() { this->evaluateStandard(); };
33  else
34  mooseError("The specified type of normalization for expansion does not exist");
35 
36  if (generation_type == "orthonormal")
37  _evaluateGenerationWrapper = [this]() { this->evaluateOrthonormal(); };
38  else if (generation_type == "sqrt_mu")
39  _evaluateGenerationWrapper = [this]() { this->evaluateSqrtMu(); };
40  else if (generation_type == "standard")
41  _evaluateGenerationWrapper = [this]() { this->evaluateStandard(); };
42  else
43  mooseError("The specified type of normalization for generation does not exist");
44 }
45 
46 std::size_t
47 Legendre::calculatedNumberOfTermsBasedOnOrder(const std::vector<std::size_t> & order) const
48 {
49  return order[0] + 1;
50 }
51 
52 void
53 Legendre::checkPhysicalBounds(const std::vector<Real> & bounds) const
54 {
55  // Each Legendre series should have a min and max bound
56  if (bounds.size() != 2)
57  mooseError("Legend: Invalid number of bounds specified for single series!");
58 }
59 
60 void
61 Legendre::evaluateOrthonormal()
62 {
63  std::size_t k;
64  const Real & x = _standardized_location[0];
65  const Real x2 = x * x;
66 
67  /*
68  * Use direct formula to efficiently evaluate the polynomials for n <= 12
69  *
70  * The performance benefit diminishes for higher n. It is expected that the cost of the direct
71  * calculation nears that of the recurrence relation in the neighborhood of n == 15, although this
72  * theory is untested due to only implementing the direct calculations up to n == 12.
73  *
74  * If you want to calculate the higher-order Legendre Coefficients and code them in then be my
75  * guest.
76  */
77  // clang-format off
78  switch (_orders[0])
79  {
80  default:
81  case MAX_DIRECT_CALCULATION_LEGENDRE: /* 12 */
82  save(12, ((((((676039 * x2 - 1939938) * x2 + 2078505) * x2 - 1021020) * x2 + 225225) * x2 - 18018) * x2 + 231) / 1024
83  * 12.5);
84  libmesh_fallthrough();
85 
86  case 11:
87  save(11, (((((88179 * x2 - 230945) * x2 + 218790) * x2 - 90090) * x2 + 15015) * x2 - 693) * x / 256
88  * 11.5);
89  libmesh_fallthrough();
90 
91  case 10:
92  save(10, (((((46189 * x2 - 109395) * x2 + 90090) * x2 - 30030) * x2 + 3465) * x2 - 63) / 256
93  * 10.5);
94  libmesh_fallthrough();
95 
96  case 9:
97  save(9, ((((12155 * x2 - 25740) * x2 + 18018) * x2 - 4620) * x2 + 315) * x / 128
98  * 9.5);
99  libmesh_fallthrough();
100 
101  case 8:
102  save(8, ((((6435 * x2 - 12012) * x2 + 6930) * x2 - 1260) * x2 + 35) / 128
103  * 8.5);
104  libmesh_fallthrough();
105 
106  case 7:
107  save(7, (((429 * x2 - 693) * x2 + 315) * x2 - 35) * x / 16
108  * 7.5);
109  libmesh_fallthrough();
110 
111  case 6:
112  save(6, (((231 * x2 - 315) * x2 + 105) * x2 - 5) / 16
113  * 6.5);
114  libmesh_fallthrough();
115 
116  case 5:
117  save(5, ((63 * x2 - 70) * x2 + 15) * x / 8
118  * 5.5);
119  libmesh_fallthrough();
120 
121  case 4:
122  save(4, ((35 * x2 - 30) * x2 + 3) / 8
123  * 4.5);
124  libmesh_fallthrough();
125 
126  case 3:
127  save(3, (5 * x2 - 3) * x / 2
128  * 3.5);
129  libmesh_fallthrough();
130 
131  case 2:
132  save(2, (3 * x2 - 1) / 2
133  * 2.5);
134  libmesh_fallthrough();
135 
136  case 1:
137  save(1, x
138  * 1.5);
139  libmesh_fallthrough();
140 
141  case 0:
142  save(0, 1
143  * 0.5);
144  }
145  // clang-format on
146 
147  /*
148  * Evaluate any remaining polynomials.
149  *
150  * The original recurrence relation is:
151  * (2 * k - 1) * x * L_(k-1) - (k - 1) * L_(k-2)
152  * L_k = ---------------------------------------------
153  * k
154  *
155  * However, for FXs we are using a the orthonormalized version of the polynomials, so each
156  * polynomial L_k is multiplied by:
157  * (2 * k + 1)
158  * ----------- essentially: k + 0.5
159  * 2
160  * Reversing this in the previous polynomials and implementing for the current polynomial results
161  * in the orthonormalized recurrence:
162  * (2 * k + 1) / (k - 1) \
163  * L_k = ----------- * | x * L_(k-1) - ----------- * L_(k-2) |
164  * k \ (2 * k - 3) /
165  *
166  * The options are 1) to use this form, or 2) to not apply the orthonormalization at first, and
167  * then loop through all the values in a second loop and then apply the orthonormalization.
168  */
169  for (k = MAX_DIRECT_CALCULATION_LEGENDRE + 1; k <= _orders[0]; ++k)
170  save(k, ((k + k + 1) / Real(k)) * (x * load(k - 1) - ((k - 1) / (k + k - 3.0)) * load(k - 2)));
171 }
172 
173 void
174 Legendre::evaluateStandard()
175 {
176  std::size_t k;
177  const Real x = _standardized_location[0];
178  const Real x2 = x * x;
179 
180  /*
181  * Use direct formula to efficiently evaluate the polynomials for n <= 12
182  *
183  * The performance benefit diminishes for higher n. It is expected that the cost of the direct
184  * calculation nears that of the recurrence relation in the neighborhood of n == 15, although this
185  * theory is untested due to only implementing the direct calculations up to n == 12.
186  *
187  * If you want to calculate the higher-order Legendre Coefficients and
188  * code them in then be my guest.
189  */
190  // clang-format off
191  switch (_orders[0])
192  {
193  default:
194  case MAX_DIRECT_CALCULATION_LEGENDRE: /* 12 */
195  save(12, ((((((676039 * x2 - 1939938) * x2 + 2078505) * x2 - 1021020) * x2 + 225225) * x2 - 18018) * x2 + 231) / 1024);
196  libmesh_fallthrough();
197 
198  case 11:
199  save(11, (((((88179 * x2 - 230945) * x2 + 218790) * x2 - 90090) * x2 + 15015) * x2 - 693) * x / 256);
200  libmesh_fallthrough();
201 
202  case 10:
203  save(10, (((((46189 * x2 - 109395) * x2 + 90090) * x2 - 30030) * x2 + 3465) * x2 - 63) / 256);
204  libmesh_fallthrough();
205 
206  case 9:
207  save(9, ((((12155 * x2 - 25740) * x2 + 18018) * x2 - 4620) * x2 + 315) * x / 128);
208  libmesh_fallthrough();
209 
210  case 8:
211  save(8, ((((6435 * x2 - 12012) * x2 + 6930) * x2 - 1260) * x2 + 35) / 128);
212  libmesh_fallthrough();
213 
214  case 7:
215  save(7, (((429 * x2 - 693) * x2 + 315) * x2 - 35) * x / 16);
216  libmesh_fallthrough();
217 
218  case 6:
219  save(6, (((231 * x2 - 315) * x2 + 105) * x2 - 5) / 16);
220  libmesh_fallthrough();
221 
222  case 5:
223  save(5, ((63 * x2 - 70) * x2 + 15) * x / 8);
224  libmesh_fallthrough();
225 
226  case 4:
227  save(4, ((35 * x2 - 30) * x2 + 3) / 8);
228  libmesh_fallthrough();
229 
230  case 3:
231  save(3, (5 * x2 - 3) * x / 2);
232  libmesh_fallthrough();
233 
234  case 2:
235  save(2, (3 * x2 - 1) / 2);
236  libmesh_fallthrough();
237 
238  case 1:
239  save(1, x);
240  libmesh_fallthrough();
241 
242  case 0:
243  save(0, 1);
244  }
245  // clang-format on
246 
247  /*
248  * Evaluate any remaining polynomials.
249  * The recurrence relation is:
250  * (2 * k - 1) * x * L_(k-1) - (k - 1) * L_(k-2)
251  * L_k = ---------------------------------------------
252  * k
253  */
254  for (k = MAX_DIRECT_CALCULATION_LEGENDRE + 1; k <= _orders[0]; ++k)
255  save(k, (((2 * k - 1) * x * load(k - 1)) - ((k - 1) * load(k - 2))) / Real(k));
256 }
257 
258 void
259 Legendre::evaluateSqrtMu()
260 {
261  evaluateStandard();
262  for (size_t i = 0; i < getNumberOfTerms(); ++i)
263  save(i, load(i) * std::sqrt(i + 0.5));
264 }
265 
266 const std::vector<Real> &
267 Legendre::getStandardizedFunctionLimits() const
268 {
269  // Lazily instantiate the function limits array
270  static const std::vector<Real> standardizedFunctionLimits = {-1, 1};
271 
272  return standardizedFunctionLimits;
273 }
274 
275 Real
276 Legendre::getStandardizedFunctionVolume() const
277 {
278  return 2.0; // Span of [-1, 1]
279 }
280 
281 std::vector<Real>
282 Legendre::getStandardizedLocation(const std::vector<Real> & location) const
283 {
284  const Real difference = location[0] - _physical_bounds[0];
285  const Real span = _physical_bounds[1] - _physical_bounds[0];
286 
287  // Convert to [0, 1] (assuming that location[0] is within _physical_bounds)
288  const Real ratio = difference / span;
289 
290  // Legendre space is [-1, 1]
291  return {ratio * 2 - 1};
292 }
293 
294 bool
295 Legendre::isInPhysicalBounds(const Point & point) const
296 {
297  std::vector<Real> location = extractLocationFromPoint(point);
298 
299  if (location[0] < _physical_bounds[0] || _physical_bounds[1] < location[0])
300  return false;
301  else
302  return true;
303 }
SingleSeriesBasisInterface::_evaluateGenerationWrapper
std::function< void()> _evaluateGenerationWrapper
The generation evaluation wrapper.
Definition: SingleSeriesBasisInterface.h:91
SingleSeriesBasisInterface::_evaluateExpansionWrapper
std::function< void()> _evaluateExpansionWrapper
The expansion evaluation wrapper.
Definition: SingleSeriesBasisInterface.h:88
Legendre::isInPhysicalBounds
virtual bool isInPhysicalBounds(const Point &point) const override
Determines if the point provided is in within the physical bounds.
Definition: Legendre.C:295
Legendre::Legendre
Legendre()
Definition: Legendre.C:19
Legendre::getStandardizedFunctionLimits
virtual const std::vector< Real > & getStandardizedFunctionLimits() const override
Returns a vector of the lower and upper bounds of the standard functional space.
Definition: Legendre.C:267
FunctionalBasisInterface::load
Real load(std::size_t index) const
Helper function to load a value from #_series.
Definition: FunctionalBasisInterface.C:132
expansion_type
MooseEnum expansion_type
Legendre::evaluateOrthonormal
virtual void evaluateOrthonormal()
Evaluates the orthonormal form of the basis functions.
Definition: Legendre.C:61
Legendre.h
SingleSeriesBasisInterface::_standardized_location
std::vector< Real > _standardized_location
The standardized location of evaluation.
Definition: SingleSeriesBasisInterface.h:85
Legendre::evaluateSqrtMu
virtual void evaluateSqrtMu()
Evaluates the 1/sqrt(mu) normalized form of the basis functions.
Definition: Legendre.C:259
Legendre::calculatedNumberOfTermsBasedOnOrder
virtual std::size_t calculatedNumberOfTermsBasedOnOrder(const std::vector< std::size_t > &order) const override
Returns the number of terms in the single series given a value for the order.
Definition: Legendre.C:47
Legendre::getStandardizedLocation
virtual std::vector< Real > getStandardizedLocation(const std::vector< Real > &location) const override
Standardize the location according to the requirements of the underlying basis, which may actually co...
Definition: Legendre.C:282
SingleSeriesBasisInterface
This class is a simple wrapper around FunctionalBasisInterface, and intended for use by any single fu...
Definition: SingleSeriesBasisInterface.h:19
generation_type
MooseEnum generation_type
SingleSeriesBasisInterface::_physical_bounds
std::vector< Real > _physical_bounds
The physical bounds of the series.
Definition: SingleSeriesBasisInterface.h:82
Legendre::evaluateStandard
virtual void evaluateStandard()
Evaluates the standard form of the basis functions.
Definition: Legendre.C:174
FunctionalBasisInterface::save
void save(std::size_t index, Real value)
Helper function to store a value in #_series.
Definition: FunctionalBasisInterface.C:138
FunctionalBasisInterface::getNumberOfTerms
std::size_t getNumberOfTerms() const
Returns the number of terms in the series.
Definition: FunctionalBasisInterface.C:88
Legendre::checkPhysicalBounds
virtual void checkPhysicalBounds(const std::vector< Real > &bounds) const override
Checks the physical bounds according to the actual implementation.
Definition: Legendre.C:53
SingleSeriesBasisInterface::_orders
std::vector< std::size_t > _orders
The order of the series.
Definition: SingleSeriesBasisInterface.h:79
SingleSeriesBasisInterface::extractLocationFromPoint
std::vector< Real > extractLocationFromPoint(const Point &point) const
Convert a spatial point to a location that the series will use to determine the value at which to eva...
Definition: SingleSeriesBasisInterface.C:122
Legendre::getStandardizedFunctionVolume
virtual Real getStandardizedFunctionVolume() const override
Returns the volume within the standardized function local_limits.
Definition: Legendre.C:276
Generated on Sat Jan 25 2020 12:04:37 for www.mooseframework.org by   doxygen 1.8.16