HyperHDG
one_dimensional.hxx
Go to the documentation of this file.
1 #pragma once // Ensure that file is included only once in a single compilation.
2 
4 #include <tpp/tpp_assert.hxx>
5 
6 namespace TPP
7 {
8 namespace ShapeType
9 {
10 /*!*************************************************************************************************
11  * \brief Struct that handles the evaluation of one-dimensional Legendre polynomials.
12  *
13  * We use Clenshaw's algorithm (cf. https://en.wikipedia.org/wiki/Clenshaw_algorithm; date: Jan. 09,
14  * 2021) to evaluate the three-term recusion formula defining Legendre polynomials as defined on
15  * https://en.wikipedia.org/wiki/Legendre_polynomials; date Jan. 09, 2021).
16  *
17  * The shape functions in this struct, however, are defined with respect to the unit interval [0,1]
18  * and normalized to be L^2 orthonormal (not just orthogonal).
19  *
20  * \tparam poly_deg The maximum degree of the polynomial.
21  *
22  * \authors Andreas Rupp, Heidelberg University, 2021.
23  **************************************************************************************************/
24 template <unsigned int poly_deg>
25 struct Legendre
26 {
27  /*!***********************************************************************************************
28  * \brief Number of shape functions that span the local polynomial space.
29  ************************************************************************************************/
30  static constexpr unsigned int n_fun() { return degree() + 1; }
31  /*!***********************************************************************************************
32  * \brief Dimension of abscissas of the polynomials.
33  ************************************************************************************************/
34  static constexpr unsigned int dim() { return 1U; }
35  /*!***********************************************************************************************
36  * \brief Maximum degree of all polynomials.
37  ************************************************************************************************/
38  static constexpr unsigned int degree() { return poly_deg; }
39 
40  /*!***********************************************************************************************
41  * \brief Evaluate variable alpha of the Clenshaw algorithm.
42  ************************************************************************************************/
43  template <typename return_t, typename input_t>
44  static constexpr return_t clenshaw_alpha(const unsigned int k, const input_t& x)
45  {
46  return (((return_t)(2 * k + 1)) / ((return_t)(k + 1))) * ((return_t)x);
47  }
48  /*!***********************************************************************************************
49  * \brief Evaluate variable beta of the Clenshaw algorithm.
50  ************************************************************************************************/
51  template <typename return_t, typename input_t>
52  static constexpr return_t clenshaw_beta(const unsigned int k, const input_t&)
53  {
54  return -((return_t)k) / ((return_t)(k + 1));
55  }
56  /*!***********************************************************************************************
57  * \brief Evaluate variable b of the Clenshaw algorithm.
58  ************************************************************************************************/
59  template <typename return_t, typename input_t>
60  static constexpr return_t clenshaw_b(const unsigned int k, const unsigned int n, const input_t& x)
61  {
62  if (k > n)
63  return (return_t)0.;
64  else if (k == n)
65  return (return_t)1.;
66  else
67  return clenshaw_alpha<return_t>(k, x) * clenshaw_b<return_t>(k + 1, n, x) +
68  clenshaw_beta<return_t>(k + 1, x) * clenshaw_b<return_t>(k + 2, n, x);
69  }
70 
71  /*!***********************************************************************************************
72  * \brief Evaluate value of one shape function.
73  *
74  * Evaluates value of the \c index one-dimensional shape function on the reference interval
75  * \f$[0,1]\f$ at abscissas \c x_val.
76  *
77  * \tparam return_t Floating type specification for return value.
78  * \tparam input_t Floating type specification for input value.
79  * \param index Index of evaluated shape function.
80  * \param x_val Abscissa of evaluated shape function.
81  * \param normalized Decide whether L^2 normalization is conducted. Defaults to true.
82  * \retval fct_value Evaluated value of shape function.
83  ************************************************************************************************/
84  template <typename return_t, typename input_t>
85  static constexpr return_t fct_val(const unsigned int index,
86  const input_t& x_val,
87  const bool normalized = true)
88  {
89  // Transform x value from reference interval [0,1] -> [-1,1].
90  const return_t x = 2. * (return_t)x_val - 1.;
91 
92  // Evaluate the value of the Legendre polynomial at x.
93  return_t legendre_val;
94  switch (index)
95  {
96  case 0:
97  legendre_val = 1.;
98  break;
99  case 1:
100  legendre_val = x;
101  break;
102  default:
103  legendre_val =
104  x * clenshaw_b<return_t>(1, index, x) - 0.5 * clenshaw_b<return_t>(2, index, x);
105  }
106 
107  // Return L^2 normalized value.
108  if (normalized)
109  legendre_val *= TPP::heron_root((return_t)(2. * index + 1.));
110 
111  return legendre_val;
112  }
113  /*!***********************************************************************************************
114  * \brief Evaluate value of the derivative of orthonormal shape function.
115  *
116  * Evaluates the value of the derivative of the \c index orthonormal, one-dimensional shape
117  * function on the reference interval \f$[0,1]\f$ at abscissa \c x_val.
118  *
119  * \tparam return_t Floating type specification for return value.
120  * \tparam input_t Floating type specification for input value.
121  * \param index Index of evaluated shape function.
122  * \param x_val Abscissa of evaluated shape function.
123  * \param normalized Decide whether L^2 normalization is conducted. Defaults to true.
124  * \retval fct_value Evaluated value of shape function's derivative.
125  ************************************************************************************************/
126  template <typename return_t, typename input_t>
127  static constexpr return_t der_val(const unsigned int index,
128  const input_t& x_val,
129  const bool normalized = true)
130  {
131  if (index == 0)
132  return (return_t)0.;
133 
134  // Transform x value from reference interval [0,1] -> [-1,1].
135  const return_t x = 2. * (return_t)x_val - 1.;
136 
137  // Non-recusrive version of:
138  // return_t legendre_val = ((return_t)index) * fct_val<return_t>(index - 1, x_val, false) +
139  // x * der_val<return_t>(index - 1, x_val, false, false);
140  return_t legendre_val = 0., x_pot = 1.;
141  for (unsigned int p = index; p > 0; --p)
142  {
143  legendre_val += ((return_t)p) * x_pot * fct_val<return_t>(p - 1, x_val, false);
144  x_pot *= x;
145  }
146 
147  // Return L^2 normalized value.
148  if (normalized)
149  legendre_val *= TPP::heron_root((return_t)(2. * index + 1.));
150 
151  return 2. * legendre_val;
152  }
153 }; // end of struct Legendre
154 
155 /*!*************************************************************************************************
156  * \brief Struct that handles the evaluation of one-dimensional Lobatto polynomials.
157  *
158  * The \f$n\f$-th Lobatto polynomial \f$L_n\f$ is defined to be the definite integral of the
159  * \f$n-1\f$-st Legendre polynomial \f$P_n$, i.e., \f$L_n(x) = \int_0^x P_{n-1}(t) \; dt\f$. Thus,
160  * the Lobatto polynomials turn out to be \f$H^1\f$ orthonormal. They can be calculated according to
161  * \f$L_n(x) = \frac{1}{n} (x P_{n-1}(x) - P_{n-2}(x))\f$, which can be deduced from the derivation
162  * rules of the Legendre polynomials. The zeroth polynomial is constant.
163  *
164  * \tparam poly_deg The maximum degree of the polynomial.
165  *
166  * \authors Andreas Rupp, Heidelberg University, 2021.
167  **************************************************************************************************/
168 template <unsigned int poly_deg>
169 struct Lobatto
170 {
171  /*!***********************************************************************************************
172  * \brief Number of shape functions that span the local polynomial space.
173  ************************************************************************************************/
174  static constexpr unsigned int n_fun() { return degree() + 1; }
175  /*!***********************************************************************************************
176  * \brief Dimension of abscissas of the polynomials.
177  ************************************************************************************************/
178  static constexpr unsigned int dim() { return 1U; }
179  /*!***********************************************************************************************
180  * \brief Maximum degree of all polynomials.
181  ************************************************************************************************/
182  static constexpr unsigned int degree() { return poly_deg; }
183 
184  /*!***********************************************************************************************
185  * \brief Evaluate value of one shape function.
186  *
187  * Evaluates value of the \c index one-dimensional shape function on the reference interval
188  * \f$[0,1]\f$ at abscissas \c x_val.
189  *
190  * \tparam return_t Floating type specification for return value.
191  * \tparam input_t Floating type specification for input value.
192  * \param index Index of evaluated shape function.
193  * \param x_val Abscissa of evaluated shape function.
194  * \param normalized Decide whether L^2 normalization is conducted. Defaults to true.
195  * \retval fct_value Evaluated value of shape function.
196  ************************************************************************************************/
197  template <typename return_t, typename input_t>
198  static constexpr return_t fct_val(const unsigned int index,
199  const input_t& x_val,
200  const bool normalized = true)
201  {
202  if (index == 0)
203  return 1.;
204  else if (index == 1)
205  return x_val;
206 
207  // Transform x value from reference interval [0,1] -> [-1,1].
208  const return_t x = 2. * (return_t)x_val - 1.;
209 
210  // Evaluate the value of the Lobatto polynomial at x.
211  return_t lobatto_val =
212  (x * Legendre<poly_deg>::template fct_val<return_t>(index - 1, x_val, false) -
213  Legendre<poly_deg>::template fct_val<return_t>(index - 2, x_val, false)) /
214  (return_t)(index);
215 
216  // Return L^2 normalized value.
217  if (normalized)
218  lobatto_val *= TPP::heron_root((return_t)(2. * index - 1.));
219 
220  return 0.5 * lobatto_val;
221  }
222  /*!***********************************************************************************************
223  * \brief Evaluate value of the derivative of orthonormal shape function.
224  *
225  * Evaluates the value of the derivative of the \c index orthonormal, one-dimensional shape
226  * function on the reference interval \f$[0,1]\f$ at abscissa \c x_val.
227  *
228  * \tparam return_t Floating type specification for return value.
229  * \tparam input_t Floating type specification for input value.
230  * \param index Index of evaluated shape function.
231  * \param x_val Abscissa of evaluated shape function.
232  * \param normalized Decide whether L^2 normalization is conducted. Defaults to true.
233  * \retval fct_value Evaluated value of shape function's derivative.
234  ************************************************************************************************/
235  template <typename return_t, typename input_t>
236  static constexpr return_t der_val(const unsigned int index,
237  const input_t& x_val,
238  const bool normalized = true)
239  {
240  if (index == 0)
241  return 0.;
242 
243  return Legendre<poly_deg>::template fct_val<return_t>(index - 1, x_val, normalized);
244  }
245 }; // end of struct Lobatto
246 
247 } // end of namespace ShapeType
248 
249 } // end of namespace TPP
TPP::ShapeType::Legendre::clenshaw_beta
static constexpr return_t clenshaw_beta(const unsigned int k, const input_t &)
Evaluate variable beta of the Clenshaw algorithm.
Definition: one_dimensional.hxx:52
check_push_test.index
index
Definition: check_push_test.py:10
TPP::ShapeType::Legendre
Struct that handles the evaluation of one-dimensional Legendre polynomials.
Definition: one_dimensional.hxx:25
TPP::ShapeType::Legendre::clenshaw_alpha
static constexpr return_t clenshaw_alpha(const unsigned int k, const input_t &x)
Evaluate variable alpha of the Clenshaw algorithm.
Definition: one_dimensional.hxx:44
TPP::heron_root
static constexpr float_t heron_root(const float_t square)
Calculate the non-negative square-root of a non-negative number at compile time.
Definition: compile_time_tricks.hxx:17
tpp_assert.hxx
This file provides the function tpp_assert.
TPP
Contains the functionalities to evaluate and integrate tensor product type polynomials.
Definition: compile_time_tricks.hxx:7
TPP::ShapeType::Legendre::degree
static constexpr unsigned int degree()
Maximum degree of all polynomials.
Definition: one_dimensional.hxx:38
TPP::ShapeType::Legendre::fct_val
static constexpr return_t fct_val(const unsigned int index, const input_t &x_val, const bool normalized=true)
Evaluate value of one shape function.
Definition: one_dimensional.hxx:85
TPP::ShapeType::Legendre::der_val
static constexpr return_t der_val(const unsigned int index, const input_t &x_val, const bool normalized=true)
Evaluate value of the derivative of orthonormal shape function.
Definition: one_dimensional.hxx:127
TPP::ShapeType::Legendre::n_fun
static constexpr unsigned int n_fun()
Number of shape functions that span the local polynomial space.
Definition: one_dimensional.hxx:30
TPP::ShapeType::Legendre::clenshaw_b
static constexpr return_t clenshaw_b(const unsigned int k, const unsigned int n, const input_t &x)
Evaluate variable b of the Clenshaw algorithm.
Definition: one_dimensional.hxx:60
TPP::ShapeType::Legendre::dim
static constexpr unsigned int dim()
Dimension of abscissas of the polynomials.
Definition: one_dimensional.hxx:34
compile_time_tricks.hxx