debut du codage de la classe polynome

src/analysis/Polynomial.hpp

+#ifndef POLYNOMIAL_HPP
+#define POLYNOMIAL_HPP
+
+#include <algebra/TinyVector.hpp>
+
+template <size_t N>
+class Polynomial
+{
+ private:
+  using Coefficients = TinyVector<N + 1, double>;
+  Coefficients m_coefficients;
+
+ public:
+  PUGS_INLINE
+  constexpr const TinyVector<N + 1>&
+  coefficients() const
+  {
+    return m_coefficients;
+  }
+
+  PUGS_INLINE
+  constexpr TinyVector<N + 1>&
+  coefficients()
+  {
+    return m_coefficients;
+  }
+
+  PUGS_INLINE
+  constexpr bool
+  operator==(const Polynomial<N>& q) const
+  {
+    if (m_coefficients != q.coefficients()) {
+      return false;
+    }
+    return true;
+  }
+
+  PUGS_INLINE
+  constexpr bool
+  operator!=(const Polynomial<N>& q) const
+  {
+    return not this->operator==(q);
+  }
+
+  PUGS_INLINE
+  constexpr Polynomial<N>
+  operator+(const Polynomial<N>& q) const
+  {
+    TinyVector<N + 1> sum_coefs = m_coefficients + q.coefficients();
+    Polynomial<N> S(sum_coefs);
+    return S;
+  }
+
+  PUGS_INLINE
+  constexpr Polynomial<N>
+  operator-(const Polynomial<N>& q) const
+  {
+    TinyVector<N + 1> diff_coefs = m_coefficients - q.coefficients();
+    Polynomial<N> D(diff_coefs);
+    return D;
+  }
+
+  PUGS_INLINE
+  constexpr Polynomial<N> operator*(const double& lambda) const
+  {
+    TinyVector<N + 1> mult_coefs = lambda * m_coefficients;
+    Polynomial<N> M(mult_coefs);
+    return M;
+  }
+
+  PUGS_INLINE
+  constexpr friend Polynomial<N> operator*(const double& lambda, const Polynomial<N> P)
+  {
+    return P * lambda;
+  }
+  // evaluation using Horner's method https://en.wikipedia.org/wiki/Horner's_method
+  PUGS_INLINE
+  constexpr double
+  evaluate(const double& x) const
+  {
+    TinyVector<N + 1> coefs = this->coefficients();
+    double bcoef            = coefs[N];
+    for (size_t i = N; i > 0; --i) {
+      bcoef *= x;
+      bcoef += coefs[i - 1];
+    }
+    return bcoef;
+  }
+
+  PUGS_INLINE
+  constexpr double
+  operator()(const double x) const
+  {
+    TinyVector<N + 1> coefs = this->coefficients();
+    double bcoef            = coefs[N];
+    for (size_t i = N; i > 0; --i) {
+      bcoef *= x;
+      bcoef += coefs[i - 1];
+    }
+    return bcoef;
+  }
+
+  PUGS_INLINE
+  constexpr friend Polynomial<N + 1>
+  primitive(const Polynomial<N>& P)
+  {
+    TinyVector<N + 2> coefs;
+    for (size_t i = 0; i < N + 1; ++i) {
+      coefs[i + 1] = P.coefficients()[i] / double(i + 1);
+    }
+    coefs[0] = 0;
+    return Polynomial<N + 1>{coefs};
+  }
+
+  PUGS_INLINE
+  constexpr friend std::ostream&
+  operator<<(std::ostream& os, const Polynomial<N>& P)
+  {
+    //    os << "P(x) = ";
+    bool all_coef_zero = true;
+    for (size_t i = N; i > 1; --i) {
+      if (P.coefficients()[i] != 0.) {
+        if (P.coefficients()[i] != 1 && P.coefficients()[i] != -1) {
+          os << std::abs(P.coefficients()[i]);
+        }
+        os << "x^" << i;
+        if (P.coefficients()[i - 1] > 0.)
+          os << " + ";
+        else
+          os << " - ";
+        all_coef_zero = false;
+      } else {
+        os << " ";
+      }
+    }
+    if (P.coefficients()[1] != 0.) {
+      if (P.coefficients()[1] != 1 && P.coefficients()[1] != -1) {
+        os << std::abs(P.coefficients()[1]);
+      }
+      os << "x";
+      if (P.coefficients()[0] > 0.)
+        os << " + ";
+      else if (P.coefficients()[0] < 0.)
+        os << " - ";
+      all_coef_zero = false;
+    } else {
+      os << " ";
+    }
+    if (P.coefficients()[0] != 0. || all_coef_zero)
+      os << std::abs(P.coefficients()[0]);
+    return os;
+  }
+
+  PUGS_INLINE
+  constexpr Polynomial(const TinyVector<N + 1>& coefficients) noexcept : m_coefficients{coefficients} {}
+
+  PUGS_INLINE
+  constexpr Polynomial(TinyVector<N + 1>&& coefficients) noexcept : m_coefficients{coefficients} {}
+
+  PUGS_INLINE
+  constexpr Polynomial() noexcept = default;
+  ~Polynomial()                   = default;
+};
+
+#endif   // POLYNOMIAL_HPP

tests/test_Polynomial.cpp

+#include <catch2/catch.hpp>
+
+#include <Kokkos_Core.hpp>
+
+#include <utils/PugsAssert.hpp>
+#include <utils/Types.hpp>
+
+#include <algebra/TinyMatrix.hpp>
+#include <analysis/Polynomial.hpp>
+
+// Instantiate to ensure full coverage is performed
+template class Polynomial<1>;
+template class Polynomial<2>;
+template class Polynomial<3>;
+template class Polynomial<4>;
+
+// clazy:excludeall=non-pod-global-static
+
+TEST_CASE("Polynomial", "[analysis]")
+{
+  SECTION("construction")
+  {
+    REQUIRE_NOTHROW(Polynomial<2>{{2, 3, 4}});
+  }
+  SECTION("equality")
+  {
+    Polynomial<2> P({2, 3, 4});
+    Polynomial<2> Q({2, 3, 4});
+
+    REQUIRE(P == Q);
+  }
+  SECTION("addition")
+  {
+    Polynomial<2> P({2, 3, 4});
+    Polynomial<2> Q({-1, -3, 2});
+    Polynomial<2> S({1, 0, 6});
+    REQUIRE(S == (P + Q));
+  }
+  SECTION("difference")
+  {
+    Polynomial<2> P({2, 3, 4});
+    Polynomial<2> Q({3, 4, 5});
+    Polynomial<2> D({-1, -1, -1});
+    REQUIRE(D == (P - Q));
+  }
+  SECTION("product_by_scalar")
+  {
+    Polynomial<2> P({2, 3, 4});
+    Polynomial<2> M({6, 9, 12});
+    REQUIRE(M == (P * 3));
+    REQUIRE(M == (3 * P));
+  }
+  SECTION("evaluation")
+  {
+    Polynomial<2> P({2, -3, 4});
+    double result = P.evaluate(3);
+    REQUIRE(result == 29);
+    result = P(3);
+    REQUIRE(result == 29);
+  }
+  SECTION("primitive")
+  {
+    Polynomial<2> P({2, -3, 4});
+    TinyVector<4> coefs = zero;
+    Polynomial<3> Q(coefs);
+    Q = primitive(P);
+    Polynomial<3> R({0, 2, -3. / 2, 4. / 3});
+    REQUIRE(Q == R);
+  }
+}