Libreria polinomiale che utilizza modelli di espressione

Oct 18 2020

In un recente progetto di giocattoli, ho fatto un uso massiccio di calcoli che coinvolgono polinomi in Z [x]. Come esperimento, ho deciso di provare l'implementazione di modelli di espressione e sembrava fornire drastiche accelerazioni nel mio caso d'uso principale, con modifiche minime al codice consumer. (Nel caso in cui possa essere utile vedere il contesto in cui sto usando questa libreria, vedi:https://github.com/dschepler/groebner-zx .)

Qui, la mia preoccupazione principale è con la parte di implementazione del modello di espressione della libreria. Alcuni dei miei pensieri attuali (anche se queste domande particolari sono troppo aperte per fornire buone risposte, è giusto; d'altra parte, se ci sono risposte standard, in particolare per la domanda di riferimento o spostamento di oggetti di sottoespressione, allora io sono tutto orecchie):

  • Ho scelto qui di catturare riferimenti a sottoespressioni ovunque; così, ad esempio, auto myexpr = 2 * p + q;risulterebbe immediatamente in riferimenti penzolanti alla 2 * psottoespressione. Potrebbe anche essere possibile catturare in movimento sottoespressioni intermedie; sebbene ciò aumenterebbe la complessità dell'implementazione (in particolare, vorrei comunque evitare di fare copie di polynomialespressioni foglia che sono riferimenti a variabili esistenti).
  • Una possibilità a cui stavo pensando sarebbe: fornire una classe speciale per le sottoespressioni in cui è facile determinare in fase di compilazione che il risultato della sottoespressione è un monomiale; e quindi ciò consentirebbe di fornire implementazioni ottimizzate della moltiplicazione per un monomio. Ciò renderebbe il codice consumer un po 'più leggibile: ad esempio, invece di p -= 2 * times_x_to(q, 3);leggerlo p -= 2 * x_to(3) * q;. Ma ancora una volta, al prezzo di una maggiore complessità di implementazione.

Le cose di cui sono a conoscenza ma che non rappresentano una grande preoccupazione ai fini di questa recensione includono:

  • Una possibilità di templatizzare il tipo polinomiale in modo che possa essere utilizzato anche per polinomi su Q, su double, su reali a precisione arbitraria mpfr, ecc.
  • Scelte di rappresentazione interna - in particolare, se sia meglio memorizzare i coefficienti in ordine crescente di grado.

polynomial.h:

#pragma once

#include <gmpxx.h>
#include <initializer_list>
#include <iostream>
#include <string>
#include <vector>

using Z = mpz_class;
inline Z operator""_Z(const char* s)
{
    return operator""_mpz(s);
}

// This is a library for manipulating polynomials with integer
// coefficients.  Note that because the library uses expression
// templates behind the scenes, the following general restrictions
// apply:
//
// * Use of "auto" keyword to initialize a variable with the result
//   of a calculation will generally not work as expected.  Instead,
//   explicitly set the result type to "polynomial".  e.g. instead of
//     auto p = 2 * q + r;
//   write
//     polynomial p = 2 * q + r;
//
// * Assignment expressions aliasing the destination variable will not
//   work as expected.  To work around this, use a polynomial constructor
//   to materialize the intermediate result before assignment.  e.g.
//   instead of
//     p += 2 * times_x_to(p + q, 3);
//   write
//     p += polynomial{2 * times_x_to(p + q, 3)};



// polynomial expression templates: each type should declare a member type
// is_polynomial_expr, and implement methods compatible with:
//    int degree_bound() const;
//      upper bound on degree of the result (does not need to be exact in
//      cases such as sum or difference of two polynomials)
//    Z coefficient(int d) const;
//      return the coefficient of x^d - where d can still be greater than
//      degree_bound()

template <typename CoeffCallable>
class polynomial_expr {
public:
    using is_polynomial_expr = std::true_type;

    polynomial_expr(int degree_bound, CoeffCallable&& coeff_callable)
        : m_degree_bound(degree_bound)
        , m_coeff_callable(std::forward<CoeffCallable>(coeff_callable))
    {
    }
    int degree_bound() const { return m_degree_bound; }
    Z coefficient(int d) const
    {
        return m_coeff_callable(d);
    }

private:
    int m_degree_bound;
    typename std::decay_t<CoeffCallable> m_coeff_callable;
};

template <typename PolyExpr1, typename PolyExpr2,
    typename = typename std::decay_t<PolyExpr1>::is_polynomial_expr,
    typename = typename std::decay_t<PolyExpr2>::is_polynomial_expr>
auto operator+(PolyExpr1&& p, PolyExpr2&& q)
{
    return polynomial_expr {
        std::max(p.degree_bound(), q.degree_bound()),
        [&p, &q](int d) -> Z { return p.coefficient(d) + q.coefficient(d); }
    };
}

template <typename PolyExpr1, typename PolyExpr2,
    typename = typename std::decay_t<PolyExpr1>::is_polynomial_expr,
    typename = typename std::decay_t<PolyExpr2>::is_polynomial_expr>
auto operator-(PolyExpr1&& p, PolyExpr2&& q)
{
    return polynomial_expr {
        std::max(p.degree_bound(), q.degree_bound()),
        [&p, &q](int d) -> Z { return p.coefficient(d) - q.coefficient(d); }
    };
}

template <typename PolyExpr,
    typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
auto operator-(PolyExpr&& p)
{
    return polynomial_expr {
        p.degree_bound(),
        [&p](int d) -> Z { return -(p.coefficient(d)); }
    };
}

template <typename PolyExpr,
    typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
auto operator*(const Z& n, PolyExpr&& p)
{
    return polynomial_expr {
        n == 0 ? -1 : p.degree_bound(),
        [&n, &p](int d) -> Z { return n * p.coefficient(d); }
    };
}

template <typename PolyExpr,
    typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
auto operator*(PolyExpr&& p, const Z& n)
{
    return polynomial_expr {
        n == 0 ? -1 : p.degree_bound(),
        [&n, &p](int d) -> Z { return p.coefficient(d) * n; }
    };
}

class polynomial;
polynomial operator*(const polynomial& p, const polynomial& q);

template <typename PolyExpr,
    typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
auto times_x_to(PolyExpr&& p, int d)
{
    return polynomial_expr {
        p.degree_bound() < 0 ? -1 : p.degree_bound() + d,
        [&p, d](int e) -> Z { return e >= d ? p.coefficient(e - d) : 0; }
    };
}

template <typename PolyExpr1, typename PolyExpr2,
    typename = typename std::decay_t<PolyExpr1>::is_polynomial_expr,
    typename = typename std::decay_t<PolyExpr2>::is_polynomial_expr>
bool operator==(PolyExpr1&& p, PolyExpr2&& q)
{
    auto d = std::max(p.degree_bound(), q.degree_bound());
    for (int i = 0; i <= d; ++i)
        if (p.coefficient(i) != q.coefficient(i))
            return false;
    return true;
}

template <typename PolyExpr1, typename PolyExpr2,
    typename = typename std::decay_t<PolyExpr1>::is_polynomial_expr,
    typename = typename std::decay_t<PolyExpr2>::is_polynomial_expr>
bool operator!=(PolyExpr1&& p, PolyExpr2&& q)
{
    auto d = std::max(p.degree_bound(), q.degree_bound());
    for (int i = 0; i <= d; ++i)
        if (p.coefficient(i) != q.coefficient(i))
            return true;
    return false;
}

class polynomial {
public:
    using is_polynomial_expr = std::true_type;

    polynomial() = default;
    polynomial(std::initializer_list<Z> coeffs);
    explicit polynomial(std::vector<Z> coeffs);

    polynomial(const polynomial&) = default;
    polynomial(polynomial&&) = default;
    template <typename PolyExpr,
        typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
    polynomial(PolyExpr&& p)
    {
        int d = p.degree_bound();
        if (d >= 0) {
            m_coeffs.reserve(d + 1);
            for (; d >= 0; --d)
                m_coeffs.push_back(p.coefficient(d));
            normalize();
        }
    }

    polynomial& operator=(const polynomial&) = default;
    polynomial& operator=(polynomial&&) = default;
    template <typename PolyExpr,
        typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
    polynomial& operator=(PolyExpr&& p)
    {
        auto deg = p.degree_bound();
        m_coeffs.resize(deg + 1);
        for (int d = 0; d <= deg; ++d)
            m_coeffs[deg - d] = p.coefficient(d);
        normalize();
        return *this;
    }

    int degree() const { return m_coeffs.size() - 1; }
    int degree_bound() const { return degree(); }
    const Z& coefficient(int d) const
    {
        static Z static_zero = 0;
        return d > degree() ? static_zero : m_coeffs[degree() - d];
    }

    // leading_coefficient has as a precondition that the polynomial must not be 0
    const Z& leading_coefficient() const
    {
        return m_coeffs.front();
    }

    void negate();

    template <typename PolyExpr,
        typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
    polynomial& operator+=(PolyExpr&& p)
    {
        if (p.degree_bound() > degree())
            m_coeffs.insert(m_coeffs.begin(), p.degree_bound() - degree(), 0);
        for (int d = p.degree_bound(); d >= 0; --d)
            m_coeffs[m_coeffs.size() - d - 1] += p.coefficient(d);
        normalize();
        return *this;
    }
    template <typename PolyExpr,
        typename = typename std::decay_t<PolyExpr>::is_polynomial_expr>
    polynomial& operator-=(PolyExpr&& p)
    {
        if (p.degree_bound() > degree())
            m_coeffs.insert(m_coeffs.begin(), p.degree_bound() - degree(), 0);
        for (int d = p.degree_bound(); d >= 0; --d)
            m_coeffs[m_coeffs.size() - d - 1] -= p.coefficient(d);
        normalize();
        return *this;
    }
    polynomial& operator*=(Z n);
    polynomial& operator*=(const polynomial& p)
    {
        return *this = (*this * p);
    }

    std::string to_string() const;

    friend std::ostream& operator<<(std::ostream& os, const polynomial& p)
    {
        return os << p.to_string();
    }

private:
    std::vector<Z> m_coeffs;

    void normalize();
    static std::string monomial_to_string(const Z& coeff, int d);
};

polynomial.cpp:

#include "polynomial.h"
#include <algorithm>

polynomial::polynomial(std::initializer_list<Z> coeffs)
    : m_coeffs(coeffs)
{
    normalize();
}

polynomial::polynomial(std::vector<Z> coeffs)
    : m_coeffs(std::move(coeffs))
{
    normalize();
}

void polynomial::normalize()
{
    auto first_nonzero = find_if(m_coeffs.begin(), m_coeffs.end(),
        [](const Z& coeff) { return coeff != 0; });
    m_coeffs.erase(m_coeffs.begin(), first_nonzero);
}

void polynomial::negate()
{
    for (auto& coeff : m_coeffs)
        coeff = -coeff;
}

polynomial& polynomial::operator*=(Z n)
{
    if (n == 0)
        m_coeffs.clear();
    else {
        for (Z& coeff : m_coeffs)
            coeff *= n;
    }
    return *this;
}

namespace polynomial_mult_details {

auto evenpart(const polynomial& p)
{
    return polynomial_expr {
        p.degree() / 2,
        [&p](int d) -> const Z& { return p.coefficient(d * 2); }
    };
}

auto oddpart(const polynomial& p)
{
    return polynomial_expr {
        (p.degree() - 1) / 2,
        [&p](int d) -> const Z& { return p.coefficient(d * 2 + 1); }
    };
}

template <typename PolyExpr1, typename PolyExpr2,
    typename = typename PolyExpr1::is_polynomial_expr,
    typename = typename PolyExpr2::is_polynomial_expr>
auto interleave(PolyExpr1&& p, PolyExpr2&& q)
{
    return polynomial_expr {
        std::max(2 * p.degree_bound(), 2 * q.degree_bound() + 1),
        [&p, &q](int d) -> Z {
            if (d % 2 == 0)
                return p.coefficient(d / 2);
            else
                return q.coefficient(d / 2);
        }
    };
}

} // namespace polynomial_mult_details

polynomial operator*(const polynomial& p, const polynomial& q)
{
    // Following the classic recursive algorithm with O(d^lg(3)) multiplications of Z values

    if (p == polynomial {} || q == polynomial {})
        return polynomial {};
    if (p.degree() == 0)
        return p.coefficient(0) * q;
    if (q.degree() == 0)
        return p * q.coefficient(0);

    // In the following, we'll be using each coefficient of p and q
    // multiple times, which is why we have designed the interface to let
    // the caller materialize p and q for us.

    auto pe = polynomial_mult_details::evenpart(p);
    auto po = polynomial_mult_details::oddpart(p);
    auto qe = polynomial_mult_details::evenpart(q);
    auto qo = polynomial_mult_details::oddpart(q);

    polynomial pe_qe = pe * qe;
    polynomial po_qo = po * qo;
    polynomial pepo_qeqo = (pe + po) * (qe + qo);

    return polynomial_mult_details::interleave(
        pe_qe + times_x_to(po_qo, 1),
        pepo_qeqo - pe_qe - po_qo);
}

std::string polynomial::monomial_to_string(const Z& coeff, int d)
{
    std::string result;
    if (coeff == 1 && d > 0)
        result = "";
    else if (coeff == -1 && d > 0)
        result = "-";
    else {
        result = coeff.get_str();
        if (d > 0)
            result += " ";
    }
    if (d == 1)
        result += "x";
    else if (d > 1) {
        result += "x^";
        result += std::to_string(d);
    }
    return result;
}

std::string polynomial::to_string() const
{
    if (m_coeffs.empty())
        return "0";
    std::string result = monomial_to_string(m_coeffs.front(), degree());
    for (int d = degree() - 1; d >= 0; --d) {
        auto coeff = coefficient(d);
        if (coeff > 0) {
            result += " + ";
            result += monomial_to_string(coeff, d);
        } else if (coeff < 0) {
            result += " - ";
            result += monomial_to_string(-coeff, d);
        }
    }

    return result;
}

E, nel caso in cui possa essere utile vedere la mia batteria di prova corrente (anche se sono consapevole che attualmente mancano i test per operator+=,-=,*=): polynomial_test.cpp:

#include "polynomial.h"
#include <gtest/gtest.h>

TEST(Polynomial, Equality)
{
    EXPECT_EQ((polynomial {}), (polynomial {}));
    EXPECT_EQ((polynomial { 0, 0 }), (polynomial {}));
    EXPECT_EQ((polynomial {}), (polynomial { 0, 0 }));
    EXPECT_EQ((polynomial { 1, 2 }), (polynomial { 1, 2 }));
    EXPECT_EQ((polynomial { 0, 1, 2 }), (polynomial { 1, 2 }));
    EXPECT_EQ((polynomial { 1, 2 }), (polynomial { 0, 1, 2 }));
    EXPECT_NE((polynomial { 1, 2 }), (polynomial { 1, 2, 0 }));
}

TEST(Polynomial, Degree)
{
    EXPECT_EQ((polynomial {}.degree()), -1);
    EXPECT_EQ((polynomial { 5 }.degree()), 0);
    EXPECT_EQ((polynomial { 1, 2, 3 }.degree()), 2);
}

TEST(Polynomial, Coefficient)
{
    EXPECT_EQ((polynomial {}.coefficient(0)), 0);
    EXPECT_EQ((polynomial {}.coefficient(1)), 0);
    EXPECT_EQ((polynomial { 5 }.coefficient(0)), 5);
    EXPECT_EQ((polynomial { 5 }.coefficient(1)), 0);
    EXPECT_EQ((polynomial { 5 }.coefficient(2)), 0);
    EXPECT_EQ((polynomial { 1, 2, 3 }.coefficient(0)), 3);
    EXPECT_EQ((polynomial { 1, 2, 3 }.coefficient(1)), 2);
    EXPECT_EQ((polynomial { 1, 2, 3 }.coefficient(2)), 1);
    EXPECT_EQ((polynomial { 1, 2, 3 }.coefficient(3)), 0);
    EXPECT_EQ((polynomial { 1, 2, 3 }.coefficient(4)), 0);
}

TEST(Polynomial, Negate)
{
    EXPECT_EQ((-polynomial {}), (polynomial {}));
    EXPECT_EQ((-polynomial { 1, 2 }), (polynomial { -1, -2 }));
    EXPECT_EQ((-polynomial { 1, -3, -2 }), (polynomial { -1, 3, 2 }));
}

TEST(Polynomial, TimesXTo)
{
    EXPECT_EQ(times_x_to(polynomial {}, 5), (polynomial {}));
    EXPECT_EQ(times_x_to(polynomial { 1, 2, 3 }, 0), (polynomial { 1, 2, 3 }));
    EXPECT_EQ(times_x_to(polynomial { 1, 2, 3 }, 3),
        (polynomial { 1, 2, 3, 0, 0, 0 }));
}

TEST(Polynomial, ScalarMult)
{
    EXPECT_EQ((2 * polynomial {}), (polynomial {}));
    EXPECT_EQ((0 * polynomial { 1, 2, 3 }), (polynomial {}));
    EXPECT_EQ((3 * polynomial { 1, -2, 3 }), (polynomial { 3, -6, 9 }));
    EXPECT_EQ((-3 * polynomial { 1, -2, 3 }), (polynomial { -3, 6, -9 }));
}

TEST(Polynomial, Add)
{
    EXPECT_EQ((polynomial {} + polynomial {}), (polynomial {}));
    EXPECT_EQ((polynomial { 1, 2, 3 } + polynomial {}), (polynomial { 1, 2, 3 }));
    EXPECT_EQ((polynomial {} + polynomial { 1, 2, 3 }), (polynomial { 1, 2, 3 }));
    EXPECT_EQ((polynomial { 1, 2, 3 } + polynomial { 1, 2, 3 }), (polynomial { 2, 4, 6 }));
    EXPECT_EQ((polynomial { 1, 2, 3 } + polynomial { -1, -2, -3 }), (polynomial {}));
    EXPECT_EQ((polynomial { 1, 2, 3, 4 } + polynomial { -1, -2, 5, 6 }),
        (polynomial { 8, 10 }));
    EXPECT_EQ((polynomial { 1, 2, 3, 4 } + polynomial { 5, 6 }),
        (polynomial { 1, 2, 8, 10 }));
    EXPECT_EQ((polynomial { 5, 6 } + polynomial { 1, 2, 3, 4 }),
        (polynomial { 1, 2, 8, 10 }));
}

TEST(Polynomial, Subtract)
{
    EXPECT_EQ((polynomial {} - polynomial {}), (polynomial {}));
    EXPECT_EQ((polynomial { 1, 2, 3 } - polynomial {}), (polynomial { 1, 2, 3 }));
    EXPECT_EQ((polynomial {} - polynomial { 1, 2, 3 }), (polynomial { -1, -2, -3 }));
    EXPECT_EQ((polynomial { 1, 2, 3 } - polynomial { 1, 2, 3 }), (polynomial {}));
    EXPECT_EQ((polynomial { 1, 2, 3, 4 } - polynomial { 1, 2, 5, 6 }),
        (polynomial { -2, -2 }));
    EXPECT_EQ((polynomial { 1, 2, 3, 4 } - polynomial { 5, 6 }),
        (polynomial { 1, 2, -2, -2 }));
    EXPECT_EQ((polynomial { 5, 6 } - polynomial { 1, 2, 3, 4 }),
        (polynomial { -1, -2, 2, 2 }));
}

TEST(Polynomial, Multiply)
{
    EXPECT_EQ((polynomial {} * polynomial {}), (polynomial {}));
    EXPECT_EQ((polynomial {} * polynomial { 1, 2, 3 }), (polynomial {}));
    EXPECT_EQ((polynomial { 1, 2, 3 } * polynomial {}), (polynomial {}));
    EXPECT_EQ((polynomial { 2 } * polynomial { 1, 2, 3 }), (polynomial { 2, 4, 6 }));
    EXPECT_EQ((polynomial { 1, 2, 3 } * polynomial { 2 }), (polynomial { 2, 4, 6 }));
    EXPECT_EQ((polynomial { 1, 5 } * polynomial { 1, 6 }), (polynomial { 1, 11, 30 }));
    EXPECT_EQ((polynomial { 1, 2 } * polynomial { 1, -2 }), (polynomial { 1, 0, -4 }));
    EXPECT_EQ((polynomial { 1, 1, 1, 1, 1 } * polynomial { 1, 1, 1 }), (polynomial { 1, 2, 3, 3, 3, 2, 1 }));
}

TEST(Polynomial, ExprTemplates)
{
    polynomial p { 2, 3, 4 };
    polynomial q { 5, 6 };
    EXPECT_EQ(p * q - p, p * (q - polynomial { 1 }));
    EXPECT_EQ(p + q + p, 2 * p + q);
    EXPECT_EQ(3 * p + 4 * times_x_to(q, 2), (polynomial { 20, 30, 9, 12 }));
    EXPECT_EQ(5 * p - 2 * times_x_to(q, 1), (polynomial { 3, 20 }));
}

TEST(Polynomial, ToString)
{
    EXPECT_EQ((polynomial {}.to_string()), "0");
    EXPECT_EQ((polynomial { 13 }.to_string()), "13");
    EXPECT_EQ((polynomial { -3 }.to_string()), "-3");
    EXPECT_EQ((polynomial { 1, 0 }.to_string()), "x");
    EXPECT_EQ((polynomial { 1, 4 }.to_string()), "x + 4");
    EXPECT_EQ((polynomial { 1, -4 }.to_string()), "x - 4");
    EXPECT_EQ((polynomial { 2, 3 }.to_string()), "2 x + 3");
    EXPECT_EQ((polynomial { -2, -3 }.to_string()), "-2 x - 3");
    EXPECT_EQ((polynomial { -1, 5 }.to_string()), "-x + 5");
    EXPECT_EQ((polynomial { 1, 3, 2 }.to_string()), "x^2 + 3 x + 2");
    EXPECT_EQ((polynomial { 1, 0, -3, 0, 0, 1, 5 }.to_string()), "x^6 - 3 x^4 + x + 5");
    EXPECT_EQ((polynomial { -1, 0, 0, 0 }.to_string()), "-x^3");
    EXPECT_EQ((polynomial { 1, 0, 1 }.to_string()), "x^2 + 1");
    EXPECT_EQ((polynomial { 1, 0, -1 }.to_string()), "x^2 - 1");
    EXPECT_EQ((polynomial { -1 }.to_string()), "-1");
}

Risposte

2 G.Sliepen Oct 20 2020 at 02:10

Considera il principio del minimo stupore

A meno che tu non conosca i dettagli di implementazione delle tue classi, è difficile capire perché funziona quanto segue:

std::cout << polynomial{1} << "\n";

Ma questo non:

std::cout << polynomial{1} + polynomial{2} << "\n";

Hai già spiegato perché, ma doverlo ricordare e racchiudere l'espressione in un'altra polynomial{}è frustrante e noioso. Alcuni utenti, di fronte a qualcosa di simile, aggiungeranno polynomial{}tutto, anche se non ne ha bisogno, e quindi qualsiasi vantaggio in termini di prestazioni andrà perso.

Dovresti fare in modo che l'interfaccia pubblica segua il principio del minimo stupore e che funzioni esattamente come ti aspetti che le espressioni regolari funzionino, per quanto possibile. Le prestazioni dovrebbero essere di secondaria importanza qui. Tuttavia, ciò non significa che non dovresti averlo polynomial_expraffatto. Significa tuttavia che dovresti provare a:

  • Rendersi polynomialpiù efficiente, o
  • Effettua polynomial_expril cast implicito di s polynomialdove necessario, o
  • Fornisci sovraccarichi per polynomial_exprfarlo comportare come polynomial.

Ad esempio, posso compilare ed eseguire la seconda istruzione che ho scritto sopra aggiungendo:

template <typename CoeffCallable>
class polynomial_expr {
    ...
    friend std::ostream &operator<<(const std::ostream &out, const polynomial_expr &expr);
};

... // after definition of class polynomial:

template <typename CoeffCallable>
std::ostream &operator<<(std::ostream &out, const polynomial_expr<CoeffCallable> &expr)
{
    return out << polynomial{expr};
}

Perché polynomial_exprs sono più veloci di polynomials?

Penso che il motivo per cui restituire polynomial_exprs sia più veloce è perché eviti di costruire un nuovo polynomialcon il risultato, il che implica la costruzione di a std::vectorcon i coefficienti. Poiché ciò richiede l'allocazione della memoria, che ha effetti collaterali e potrebbe cambiare il comportamento se, ad esempio, il globale operator newè sovraccarico, i compilatori potrebbero non essere in grado di ottimizzarlo, anche se possono vedere che il risultato polynomialstesso è solo temporaneo.

Ma gli utenti C ++ sono abituati a questo tipo di problemi di prestazioni; è il motivo per cui, ad esempio, modifichiamo spesso i contenitori sul posto. Quindi, invece di scrivere:

polynomial p1, p2, p3 = ...;
polynomial result = p1 + p2 * p3;

Se sai che non è più necessario mantenere i valori originali in giro, puoi scrivere:

polynomial p1, p2, p3 = ...;
p1 += p2 *= p3;

Anche se questo non vince alcun concorso di bellezza, è almeno meno sorprendente, ed è un modo in cui un utente attento alle prestazioni può ottenere codice ad alte prestazioni solo con la tua polynomialclasse.

Facendo polynomialandare più veloce

Come ho detto prima, il fatto che polynomialutilizzi un std::vectorsignifica che deve fare allocazioni di memoria nell'heap. Considerare l'utilizzo di un contenitore diverso per memorizzare i coefficienti che implementa una piccola ottimizzazione del vettore .

Inoltre, hai già menzionato la possibilità di rendere i monomi un caso speciale. Penso che sia una buona idea, ma puoi forse generalizzare questa ottimizzazione e fare in modo che memorizzi solo i coefficienti tra il coefficiente più alto e quello più basso diverso da zero. Ad esempio, fai in modo che il polinomio \$x^{102} + x^{101} + x^{100}\$memorizza solo tre coefficienti. L'implementazione può essere semplice come l'aggiunta di una variabile membro class polynomialche memorizza l'offset al coefficiente più basso diverso da zero. Naturalmente, questo potrebbe significare un po 'più di lavoro qua e là, ma alcune operazioni ne trarranno grandi benefici. Ad esempio, moltiplicare o dividere per un'unità monomiale sarà banale.

Rendendolo ancora più facile da usare

Hai già detto che x_to()sarebbe un bel modo breve per scrivere un monomiale. Anche più breve sarebbe X(), ma poiché si tratta di un nome molto breve, probabilmente vorrai inserirlo in uno spazio dei nomi, per evitare di inquinare lo spazio dei nomi globale e in modo che gli utenti possano partecipare a quello:

namespace polynomial_utilities {
polynomial_expr X(int d = 1) {
    return {d, [&p, d](int e) -> Z { return int(e == d); }};
}
}

...

using polynomial_utilies::X;
auto expr = 3 * X(2) - 5 * X() + 4;

In alternativa puoi creare dun parametro del modello in modo da utilizzare le parentesi angolari.

Aggiungi un modo per valutare un polynomialin un dato punto

È molto carino che tu possa costruire polinomi e stamparli, ma di solito vorrai valutare il polinomio in un dato punto. Quindi aggiungerei un operator()()per fare questo:

class polynomial {
    ...
    template<typename T>
    T operator()(T x) {
        T result{0};
        T multiplier{1};

        for (auto &coeff: m_coeffs) {
            result += coeff * multiplier;
            multiplier *= x;
        }

        return result;
    }
}

Quindi puoi scrivere:

polynomial func{3, -5, 4};
std::cout << func(1.5) << "\n"; // should print: 3.25