LagrangePolynomial.hpp
Go to the documentation of this file.
2 // See LICENSE.txt for details.
3
4 /// \file
5 /// Defines functions lagrange_polynomial
6
7 #pragma once
8
9 #include <iterator>
10
11 #include "ErrorHandling/Assert.hpp"
12
13 /// Evaluate the jth Lagrange interpolating polynomial with the given
14 /// control points where j is the index of index_point.
15 template <typename Iterator>
16 double lagrange_polynomial(const Iterator& index_point,
17  double x,
18  const Iterator& control_points_begin,
19  const Iterator& control_points_end) {
20  ASSERT(control_points_begin != control_points_end, "No control points");
21
22  const double x_j = *index_point;
23  double result = 1.;
24  for (auto m_it = control_points_begin; m_it != control_points_end; ++m_it) {
25  if (m_it == index_point) { continue; }
26  const double x_m = *m_it;
27  result *= (x_m - x) / (x_m - x_j);
28  }
29  return result;
30 }
31
32 /// Evaluate the jth (zero-indexed) Lagrange interpolating polynomial
33 /// with the given control points.
34 template <typename Iterator>
35 double lagrange_polynomial(size_t j, double x,
36  const Iterator& control_points_begin,
37  const Iterator& control_points_end) {
38  const auto j_diff =
39  static_cast<typename std::iterator_traits<Iterator>::difference_type>(j);
40  ASSERT(j_diff < std::distance(control_points_begin, control_points_end),
41  "Polynomial number out of range " << j << " > "
42  << (std::distance(control_points_begin, control_points_end) - 1));
43  return lagrange_polynomial(std::next(control_points_begin, j_diff),
44  x, control_points_begin, control_points_end);
45 }
double lagrange_polynomial(const Iterator &index_point, double x, const Iterator &control_points_begin, const Iterator &control_points_end)
Evaluate the jth Lagrange interpolating polynomial with the given control points where j is the index...
Definition: LagrangePolynomial.hpp:16
#define ASSERT(a, m)
Assert that an expression should be true.
Definition: Assert.hpp:49
Defines macro ASSERT.