Partage

[Exercices] Venez vous entraîner !

Ce mois: Parseur de fonctions mathématiques

27 avril 2009 à 21:58:55

Exercice supplémentaire !



Pour ceux qui en veulent encore plus, je vous invite à aller faire un tour sur le forum "Autres langages" qui, grâce à bluestorm (et d'autres) propose depuis aujourd'hui un atelier très intéressant. Il s'agit d'écrire un interpréteur pour un mini langage de programmation.

Il y a différentes difficultés et si vous allez jusqu'au bout, vous aurez déjà un vrai langage bien conséquent et capable de beaucoup plus que ce qu'on pourrait croire au premier abord.

[Atelier] P'tit langage .

Si vous avez des questions d'ordre C++, n'hésitez pas à les poser sur ce forum. Vous recevrez certainement plus d'aide. Pour le reste, postez là-bas pour garder la cohérence du sujet.

Bon courage à tous !
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
2 mai 2009 à 14:26:56

Exercice du mois de mai 2009



Nom : Diff
Sujet : Fichiers


Un des outils souvent utilisés par les administateurs de machines sous Linux est le petit programme appelé "diff". Il permet d'afficher la différence entre deux fichiers en spécifiant quelles lignes sont différentes dans les 2 fichiers passés en argument. Ce qui donne par exemple (sans aucune option) avec les deux fichiers suivants,

#include <iostream>
using namespace std;

int main()
{
        cout << "Hello, World !" << endl;
        return 0;
}


#include <iostream>


int main()
{
        std::cout << "Hello, World !" << std::endl;
        return 0;
}

//Fin du programme


le résultat suivant:

diff Fichier1 Fichier2
2c2
< using namespace std;
---
> 
6c6
<         cout << "Hello, World !" << endl;
---
>         std::cout << "Hello, World !" << std::endl;
8a9,10
> 
> //Fin du programme


Il faut lire le résultat de la manière suivante:
2c2 veut dire que la ligne 2 de chaque fichier est différente et ce qui suit correspond aux lignes dans chacun des fichiers séparées par "---".
Idem pour 6c6 qui indique que la ligne 6 a été changée. La ligne 8a9,10 indique que les lignes 9 et 10 ont été ajoutées.

Si on appelle le programme avec les fichiers dans l'ordre opposé, on obtient:
diff fichier2 fichier1
2c2
< 
---
> using namespace std;
6c6
<         std::cout << "Hello, World !" << std::endl;
---
>         cout << "Hello, World !" << endl;
9,10d8
< 
< //Fin du programme


Le début est le même, mais cette fois on voit que les lignes 9 et 10 ont été détruites dans le fichier de droite.

Vous remarquerez que les lignes identiques ne sont jamais affichées.

L'exercice



Le but de cet exercice est de réaliser un programme de diff. Fonctionnant sur le même principe que celui de Linux.

Votre programme devra recevoir deux fichiers (les chemins) en argument et montrer les différences entre les deux fichiers en indiquant si les lignes ont été ajoutées (a), changées (c) ou détruites (d) avec bien sûr le numéro des lignes différentes et le contenu de la ligne dans chacun des 2 fichiers en utilisant la notation < et > pour indiquer de quel fichier provient la ligne affichée.

Vous avez jusqu'au 31 mai pour soumettre vos réponses à Réponse_Exercices.

Bonne chance à tous et amusez-vous bien !
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
2 mai 2009 à 19:01:04


Citation : Nanoc

Votre programme devra recevoir deux fichiers en argument



Les chemins des fichiers ?
2 mai 2009 à 19:34:27

Oui.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
2 mai 2009 à 22:45:17

Ah, même si j'ai mis du temps à comprendre, c'est plus facile que les polynômes ^^ . Par contre, j'ai pas vraiment compris le truc des lignes détruites... c'est quand le fichier entré en premier a moins de lignes que l'autre ? J'ai pas vraiment compris l'affichage... pareil que 8a9,10...

Si les fichiers n'ont aucune différence, il n'y a rien d'affiché, si ?

Merci :) !
2 mai 2009 à 22:51:22

Si les fichiers sont identiques, il n'y a rien d'affiché.

8a9,10

il faut lire ça comme ça:

A partir de la ligne 8 du fichier de gauche, il y a des lignes en plus (9 à 10) dans le fichier de droite.

Pour:

9,10d8

Les lignes 9 à 10 ont été détruites dans le fichier de droite.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
2 mai 2009 à 23:19:54

Ok, merci, je m'y met tout de suite :) !
6 mai 2009 à 14:04:57

Si je comprend bien, on a qu'à récupérer le contenu des deux fichiers ligne par ligne, et après on les compare :
if(ligneFichier1 != ligneFichier2)
{
    cout << "< " << ligneFichier1 << endl;
    cout << "---" << endl;
    cout << "> " << ligneFichier2 << endl;
}

Pour ce qui est de savoir les lignes qui ont été ajoutées/détruites, je pense qu'il faut compter le nombre de lignes dans chacun des deux fichiers : s'ils n'en ont pas le même nombre, c'est qu'il y a des lignes qui ont été ajoutées ou détruites, sinon pas de souci. :D
7 mai 2009 à 1:17:47

Même si l'exercice est plus simple que les autres, merci de ne pas dévoiler directement le code.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
10 mai 2009 à 11:20:36

Solution du mois d'avril 2009



Je n'ai pas reçu beaucoup de réponses (3) pour l'exercice sur les polynômes.


La base de la classe



Comment re présenter un polynôme ? Le plus simple est certainement de stocker les coefficients. Le plus facile pour cela est d'utiliser un vector<double>.

On peut également utiliser une map ou une list si on veut optimiser certains cas particuliers du type <math>\(x^{42} + 2\)</math> où il y a beaucoup de termes intermédiaires nuls.


Nous avons donc la base suivante:

class Polynome{

    //...

private:

    std::vector<double> m_termes;

};



Constructeurs




Les constructeurs prenant en argument un paramètre de type double et le constructeur par défaut peuvent être créés ensembles.

Polynome(double a = 0.)
            :m_termes(1,a)  //On remplit le vector avec un terme valant a
    {}


Un constructeur recevant un vector peut également s'avérer utile. Comme l'utilisateur peut avoir entré un tableau de type 3,4,0,5,0,0,0 avecdes 0s inutiles, il faut les simplifier. On crée donc une fonction simplifie() qui fait cela.

Polynome(const std::vector<double>& tab)
            :m_termes(tab)
    {
        simplifie();
    }


Modification des termes



Le plus simple ici est d'utiliser les fonctions membres du vector. On a ainsi les versions [] et at() qui vérifie l'indice. Comme toujours, il y a une version const et non-const.

inline const double& operator[](size_t n) const
    {
        return m_termes[n];
    }

    inline double& operator[](size_t n)
    {
        return m_termes[n];
    }

    inline const double& at(size_t n) const
    {
        return m_termes.at(n);
    }

    inline double& at(size_t n)
    {
        return m_termes.at(n);
    }


Degré du polynôme



A nouveau, on peut utiliser les fonctions de vector pour connaître la taille du polynôme et donc son degré.

inline size_t degre() const
    {
        return m_termes.size() - 1;
    }

    inline size_t length() const
    {
        return m_termes.size();
    }


Evaluation du polynôme



L'évaluation du polynôme consiste simplement à appliquer la définition.

double evalue(double x) const
    {
        double resultat = 0.;

        for (size_t i(0); i<m_termes.size(); ++i)
        {
            resultat += m_termes[i] * puissance(x,i); //puissance renvoit x^i
        }

        return resultat;
    }


Opérateurs de comparaison



Vérifier si deux polynomes sont égaux consiste à vérifier l'égalité des termes. Pour se simplifier la vie, on utilise l'opérateur de comparaison de vector.

bool operator==(const Polynome& p, const Polynome& q)
{
    return p.m_termes == q.m_termes;
}

inline bool operator!=(const Polynome& p, const Polynome& q)
{
    return !(p == q);
}


Opérateur d'injection dans un flux



Il faut penser aux différents cas, signe, polynome nul, polynome de degré 0,...

std::ostream& operator<<(std::ostream& flux, const Polynome& poly)
{
    if (poly.degre() == 0 && poly[0] == 0.)
    {
        flux << 0.;
        return flux;
    }

    if (poly[poly.length() - 1] < 0.)
        flux << "- ";

    for (size_t i(1);i<poly.length();++i)
    {
        if (poly[poly.length() - i] != 0.)
            flux << std::abs(poly[poly.length() - i]) << "x^" << poly.length() - i;

        if (poly[poly.length() - (i+1)] > 0.)
            flux << " + ";
        else if (poly[poly.length() - (i+1)] < 0.)
            flux << " - ";
    }

    if (poly[0] != 0.)
        flux << std::abs(poly[0]);

    return flux;
}


Opérateurs arithmétiques



Comme toujours, on déclare les opérateurs +=,-= et *= dans la classe puis les versions +,- et * à l'extérieur. Au niveau des calculs, on utilise la définition (par exemple de Wikipédia).

Polynome& operator+=(const Polynome& q)
    {
        size_t i = 0;
        for (; i<q.length() && i<m_termes.size(); ++i)
        {
            m_termes[i] += q[i];
        }

        for (; i<q.length() ; ++i)
            m_termes.push_back(q[i]);

        return *this;
    }

    Polynome& operator-=(const Polynome& q)
    {
        size_t i = 0;
        for (; i<q.length() && i<m_termes.size(); ++i)
        {
            m_termes[i] -= q[i];
        }

        for (; i<q.length() ; ++i)
            m_termes.push_back(-q[i]);

        return *this;
    }

    Polynome& operator*=(const Polynome& q)
    {
        if (q == 0. || *this == 0.)
        {
            m_termes.assign(1,0.);
            return *this;
        }

        const size_t len(q.degre() + degre() + 1);
        std::vector<double> vec(len,0.);

        for (size_t k(0); k<len; ++k)
            for (size_t i(0); i< this->length();++i)
            {
                const size_t j = k-i;
                if (j<q.length())
                    vec[k] += m_termes[i]*q[j];
            }

        *this = Polynome(vec);
        return *this;
    }


La difficulté principale est de penser à ajouter des "cases" dans le vector si un des deux opérandes est de degré plus grand.

inline Polynome operator+(const Polynome& p, const Polynome& q)
{
    return Polynome(p) += q;
}

inline Polynome operator-(const Polynome& p, const Polynome& q)
{
    return Polynome(p) -= q;
}

inline Polynome operator*(const Polynome& p, const Polynome& q)
{
    return Polynome(p) *= q;
}


La division



La division euclidienne donne accès à 2 informations, le quotient et le reste. On peut ainsi créer un opérateur % et / en créant une seule fonction.
A nouveau, on utilise l'algorithme usuel. Comme la fonction n'a pas besoin d'un objet pour fonctionner, on peut la mettre en statique.

typedef std::pair<Polynome, Polynome >  Paire;

    static Paire division(const Polynome& dividende, const Polynome& diviseur)
    {
        if (diviseur == 0.)
            throw std::runtime_error("Polynome: Division par 0");

        if (dividende.degre() < diviseur.degre())
            return Paire(0, dividende);

        if (dividende == diviseur)
            return Paire(1,0);

        Polynome reste(dividende);
        Polynome quotient;

        while (reste.degre() > diviseur.degre())
        {
            std::vector<double> vec(reste.degre() - diviseur.degre() + 1,0);
            vec.back() = reste.m_termes.back() / diviseur.m_termes.back();

            Polynome p(vec);

            reste -= p*diviseur;
            reste.simplifie();

            quotient += p;
        }

        return Paire(quotient,reste);
    }



Avec cela, on a terminé les niveaux 1 et 2. Voici le code complet avec quelques fioritures en plus.

#ifndef __POLYNOMIAL_HPP_INCLUDED__
#define __POLYNOMIAL_HPP_INCLUDED__

#include <vector>
#include <cmath>
#include <stdexcept>
#include <ostream>

class Polynome;

bool operator==(const Polynome& p, const Polynome& q);
inline Polynome operator*(const Polynome& p, const Polynome& q);
double puissance(double x, size_t n);

class Polynome
{

public:

// Constructeurs -------------------------------------------------------------------------------

    Polynome(double a = 0.)
            :m_termes(1,a)
    {}

    Polynome(const std::vector<double>& tab)
            :m_termes(tab)
    {
        simplifie();
    }

    template<size_t N>
    Polynome(const double (&tab)[N])
           :m_termes(&tab[0], &tab[N])
    {
        simplifie();
    }

// Accès aux termes ------------------------------------------------------------------------------

    inline const double& operator[](size_t n) const
    {
        return m_termes[n];
    }

    inline double& operator[](size_t n)
    {
        return m_termes[n];
    }

    inline const double& at(size_t n) const
    {
        return m_termes.at(n);
    }

    inline double& at(size_t n)
    {
        return m_termes.at(n);
    }

// Degre -----------------------------------------------------------------------------------

    inline size_t degre() const
    {
        return m_termes.size() - 1;
    }

    inline size_t length() const
    {
        return m_termes.size();
    }

// Evaluation -----------------------------------------------------------------------------------

    double evalue(double x) const
    {
        double resultat = 0.;

        for (size_t i(0); i<m_termes.size(); ++i)
        {
            resultat += m_termes[i] * puissance(x,i);
        }

        return resultat;
    }

    inline double operator()(double& x) const
    {
        return evalue(x);
    }

//Operateurs ----------------------------------------------------------------------------

    Polynome& operator+=(const Polynome& q)
    {
        size_t i = 0;
        for (; i<q.length() && i<m_termes.size(); ++i)
        {
            m_termes[i] += q[i];
        }

        for (; i<q.length() ; ++i)
            m_termes.push_back(q[i]);

        return *this;
    }

    Polynome& operator-=(const Polynome& q)
    {
        size_t i = 0;
        for (; i<q.length() && i<m_termes.size(); ++i)
        {
            m_termes[i] -= q[i];
        }

        for (; i<q.length() ; ++i)
            m_termes.push_back(-q[i]);

        return *this;
    }

    Polynome& operator*=(const Polynome& q)
    {
        if (q == 0. || *this == 0.)
        {
            m_termes.assign(1,0.);
            return *this;
        }

        const size_t len(q.degre() + degre() + 1);
        std::vector<double> vec(len,0.);

        for (size_t k(0); k<len; ++k)
            for (size_t i(0); i< this->length();++i)
            {
                const size_t j = k-i;
                if (j<q.length())
                    vec[k] += m_termes[i]*q[j];
            }

        *this = Polynome(vec);
        return *this;
    }

    inline Polynome& operator/=(const Polynome& q)
    {
        *this = division(*this,q).first;
        return *this;
    }

    inline Polynome& operator%=(const Polynome& q)
    {
        *this = division(*this,q).second;
        return *this;
    }

//Derivée ------------------------------------------------------------------------------

    Polynome derivee() const
    {
        if (degre() == 0)
            return 0;

        Polynome ret(*this);
        for (size_t i(0); i<ret.degre(); ++i)
            ret[i] = ret[i+1] * (i+1);

        ret.m_termes.pop_back();

        return ret;
    }

private:

//Simplification -----------------------------------------------------------------------

    void simplifie()
    {
        size_t counter = 0;
        typedef std::vector<double>::reverse_iterator Iterator;
        for (Iterator it = m_termes.rbegin(); it != m_termes.rend() && *it == 0. ; ++it)
        {
            ++counter;
        }

        m_termes.resize(m_termes.size() - counter);
    }

    typedef std::pair<Polynome, Polynome >  Paire;

    static Paire division(const Polynome& dividende, const Polynome& diviseur)
    {
        if (diviseur == 0.)
            throw std::runtime_error("Polynome: Division par 0");

        if (dividende.degre() < diviseur.degre())
            return Paire(0, dividende);

        if (dividende == diviseur)
            return Paire(1,0);

        Polynome reste(dividende);
        Polynome quotient;

        while (reste.degre() > diviseur.degre())
        {
            std::vector<double> vec(reste.degre() - diviseur.degre() + 1,0);
            vec.back() = reste.m_termes.back() / diviseur.m_termes.back();

            Polynome p(vec);

            reste -= p*diviseur;
            reste.simplifie();

            quotient += p;
        }

        return Paire(quotient,reste);
    }

//Attributs-----------------------------------------------------------------------------

    std::vector<double> m_termes;

//Amities ------------------------------------------------------------------------------

    friend bool operator==(const Polynome& p, const Polynome& q);

};

// Flux -----------------------------------------------------------------------------------

std::ostream& operator<<(std::ostream& flux, const Polynome& poly)
{
    if (poly.degre() == 0 && poly[0] == 0.)
    {
        flux << 0.;
        return flux;
    }

    if (poly[poly.length() - 1] < 0.)
        flux << "- ";

    for (size_t i(1);i<poly.length();++i)
    {
        if (poly[poly.length() - i] != 0.)
            flux << std::abs(poly[poly.length() - i]) << "x^" << poly.length() - i;

        if (poly[poly.length() - (i+1)] > 0.)
            flux << " + ";
        else if (poly[poly.length() - (i+1)] < 0.)
            flux << " - ";
    }

    if (poly[0] != 0.)
        flux << std::abs(poly[0]);

    return flux;
}

//Opérateurs de comparaison ------------------------------------------------------------------------

bool operator==(const Polynome& p, const Polynome& q)
{
    return p.m_termes == q.m_termes;
}

inline bool operator!=(const Polynome& p, const Polynome& q)
{
    return !(p == q);
}

// Opérateurs binaires ---------------------------------------------------------------------------

inline Polynome operator+(const Polynome& p, const Polynome& q)
{
    return Polynome(p) += q;
}

inline Polynome operator-(const Polynome& p, const Polynome& q)
{
    return Polynome(p) -= q;
}

inline Polynome operator*(const Polynome& p, const Polynome& q)
{
    return Polynome(p) *= q;
}

inline Polynome operator/(const Polynome& p, const Polynome& q)
{
    return Polynome(p) /= q;
}

inline Polynome operator%(const Polynome& p, const Polynome& q)
{
    return Polynome(p) %= q;
}

//Opérateurs unaires -------------------------------------------------------------------------------

Polynome operator-(const Polynome& p)
{
    Polynome result(p);
    for (size_t i(0); i<result.length(); ++i)
        result[i] = -result[i];
    return result;
}

Polynome operator+(const Polynome& p)
{
    Polynome result(p);
    for (size_t i(0); i<result.length(); ++i)
        result[i] = +result[i];
    return result;
}

//Opérateurs logiques ----------------------------------------------------------------------------

inline bool operator!(const Polynome& p)
{
    return p == 0.;
}

//Dérivée -----------------------------------------------------------------------------------

inline Polynome derivee(const Polynome& p)
{
    return p.derivee();
}

//Fonction assistante -----------------------------------------------------------------------------

double puissance(double x, size_t n)
{
    double resultat = 1.;

    for (size_t i(0); i<n ; ++i)
        resultat *= x;
    return resultat;
}

#endif //__POLYNOM_HPP_INCLUDED__


Pour le niveau 3, il s'agissait d'ajouter des templates. Cela ne fait pas grand chose de nouveau. Voici une version complète avec des commentaires Doxygen (en anglais).

#ifndef __POLYNOMIAL_HPP_INCLUDED__
#define __POLYNOMIAL_HPP_INCLUDED__

/**
 * \file Polynom.hpp
 * \author Nanoc
 * \brief A Polynomial class
 * \date 19 April 2009
 * \version 1.1
 *
 * The operators are implemented naively. Could be improved
 * \todo Add a solve() function
 */

#include <vector>
#include <cmath>
#include <stdexcept>
#include <ostream>

/**
 * \class Polynomial<>
 *
 *  \brief A class to represent polynomials of any degree
 *
 *  All the usual operators are overloaded. One can evaluate the polynomial at any point or compute its derivative
 *  Access to all terms is also granted
 *
 */
template<typename T>
class Polynomial
{

public:

// Constructors -------------------------------------------------------------------------------

    /**
     *  \brief Constructor. Buils a Polynomial from a single value
     */
    Polynomial(const T& a = static_cast<T>(0))
            :m_terms(1,a)
    {}

    /**
     *  \brief Constructor. Builds a Polynomial from a vector of terms
     */
    Polynomial(const std::vector<T>& tab)
            :m_terms(tab)
    {
        simplify();
    }

    /**
     *  \brief Constructor. Buils a Polynomial from a "C-like" array
     */
    template<size_t N>
    Polynomial(const T (&tab)[N])
         :m_terms(&tab[0], &tab[N])
    {
       simplify();
    }

// Term access ------------------------------------------------------------------------------

    /**
     *  \brief Access to a term (read only)
     */
    inline const T& operator[](size_t n) const
    {
        return m_terms[n];
    }

    /**
     *  \brief Access to a term
     */
    inline T& operator[](size_t n)
    {
        return m_terms[n];
    }

    /**
     *  \brief Access to a term with rangecheck (read only)
     */
    inline const T& at(size_t n) const
    {
        return m_terms.at(n);
    }

    /**
     *  \brief Access to a term with rangecheck
     */
    inline T& at(size_t n)
    {
        return m_terms.at(n);
    }

// Degree -----------------------------------------------------------------------------------

    /**
     *  \brief Returns the degree of the polynomial
     */
    inline size_t degree() const
    {
        return m_terms.size() - 1;
    }

    /**
     *  \brief Returns the length of the internal array (should be degree()+1)
     */
    inline size_t length() const
    {
        return m_terms.size();
    }

// Evaluation -----------------------------------------------------------------------------------

    /**
     *  \brief Evaluates the polynom at a given point x
     */
    template <typename S>
    S evaluate(const S& x) const
    {
        S result = static_cast<S>(0);

        for (size_t i(0); i<m_terms.size(); ++i)
        {
            result += m_terms[i] * power(x,i);
        }

        return result;
    }

    /**
     *  \brief Evaluates the polynom at a given point x
     */
    template <typename S>
    inline S operator()(const S& x) const
    {
        return evaluate(x);
    }

//Operators ----------------------------------------------------------------------------

    /**
     *  \brief Sum term by term of two polynomials
     */
    Polynomial<T>& operator+=(const Polynomial<T>& q)
    {
        size_t i = 0;
        for (; i<q.length() && i<m_terms.size(); ++i)
        {
            m_terms[i] += q[i];
        }

        for (; i<q.length() ; ++i)
            m_terms.push_back(q[i]);

        return *this;
    }

    /**
     *  \brief Difference term by term of two polynomials
     */
    Polynomial<T>& operator-=(const Polynomial<T>& q)
    {
        size_t i = 0;
        for (; i<q.length() && i<m_terms.size(); ++i)
        {
            m_terms[i] -= q[i];
        }

        for (; i<q.length() ; ++i)
            m_terms.push_back(-q[i]);

        return *this;
    }

    /**
     *  \brief Product of two Polynomials
     */
    Polynomial<T>& operator*=(const Polynomial<T>& q)
    {
        if (q == Polynomial<T>(0) || *this == Polynomial<T>(0))
        {
            m_terms.assign(1,static_cast<T>(0));
            return *this;
        }

        const size_t len(q.degree() + degree() + 1);
        std::vector<T> vec(len,static_cast<T>(0));

        for (size_t k(0); k<len; ++k)
            for (size_t i(0); i< this->length();++i)
            {
                const size_t j = k-i;
                if (j<q.length())
                    vec[k] += m_terms[i]*q[j];
            }

        *this = Polynomial<T>(vec);
        return *this;
    }

    /**
     *  \brief Quotient of the Euclidian division
     */
    inline Polynomial<T>& operator/=(const Polynomial<T>& q)
    {
        *this = division(*this,q).first;
        return *this;
    }

    /**
     *  \brief Remainder of the Euclidian division
     */
    inline Polynomial<T>& operator%=(const Polynomial<T>& q)
    {
        *this = division(*this,q).second;
        return *this;
    }

//Derivative -------------------------------------------------------------------------

    /**
     *  \brief Returns the derivative of the polynom
     */
    Polynomial<T> derivative() const
    {
        if (degree() == 0)
            return 0;

        Polynomial<T> ret(*this);
        for (size_t i(0); i<ret.degree(); ++i)
            ret[i] = ret[i+1] * (i+1);

        ret.m_terms.pop_back();

        return ret;
    }

private:

//Simplification -----------------------------------------------------------------------

    /**
     *  \brief Simplifies the polynom by removing all 0-terms after the highest non-zero term
     */
    void simplify()
    {
        size_t counter = 0;
        typedef typename std::vector<T>::reverse_iterator Iterator;
        for (Iterator it = m_terms.rbegin(); it != m_terms.rend() && *it == static_cast<T>(0) ; ++it)
        {
            ++counter;
        }

        m_terms.resize(m_terms.size() - counter);
    }

    /**
     *  \brief Notation simplification
     */
    typedef std::pair<Polynomial<T>, Polynomial<T> >  Pair;

    /**
     *  \brief Returns a Pair containing the quotient and the remainder of a division
     */
    static Pair division(const Polynomial& dividend, const Polynomial& divisor)
    {
        if (divisor == Polynomial<T>(static_cast<T>(0)))
            throw std::runtime_error("Polnomial: Division by 0");

        if (dividend.degree() < divisor.degree())
            return Pair(0, dividend);

        if (dividend == divisor)
            return Pair(1,0);

        Polynomial<T> remainder(dividend);
        Polynomial quotient;

        while (remainder.degree() > divisor.degree())
        {
            std::vector<T> vec(remainder.degree() - divisor.degree() + 1,0);
            vec.back() = remainder.m_terms.back() / divisor.m_terms.back();

            Polynomial<T> p(vec);

            remainder -= p*divisor;
            remainder.simplify();

            quotient += p;
        }

        return Pair(quotient,remainder);
    }

//Attributes-----------------------------------------------------------------------------

    std::vector<T> m_terms; /*!< The coefficients */

//Friends--------------------------------------------------------------------------------

    friend bool operator==(const Polynomial<T>& p, const Polynomial<T>& q);

};

// Stream -----------------------------------------------------------------------------------
/**
 *  \brief Writes the polynom in a stream
 */
template <typename T>
std::ostream& operator<<(std::ostream& stream, const Polynomial<T>& poly)
{
    if (poly.degree() == 0 && poly[0] == 0)
    {
        stream << static_cast<T>(0);
        return stream;
    }

    if (poly[poly.length() - 1] < static_cast<T>(0))
        stream << "- ";

    for (size_t i(1);i<poly.length();++i)
    {
        if (poly[poly.length() - i] != static_cast<T>(0))
            stream << std::abs(poly[poly.length() - i]) << "x^" << poly.length() - i;

        if (poly[poly.length() - (i+1)] > static_cast<T>(0))
            stream << " + ";
        else if (poly[poly.length() - (i+1)] < static_cast<T>(0))
            stream << " - ";
    }

    if (poly[0] != static_cast<T>(0))
        stream << std::abs(poly[0]);

    return stream;
}

//Comparison operators ------------------------------------------------------------------------
/**
 *  \brief Compares two polynoms to see if they are equal
 */
template <typename T>
bool operator==(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return p.m_terms == q.m_terms;
}

/**
 *  \brief Compares two polynoms to see if they are different
 */
template <typename T>
inline bool operator!=(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return !(p == q);
}

// Binary Operators ---------------------------------------------------------------------------

template<typename T>
inline Polynomial<T> operator+(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return Polynomial<T>(p) += q;
}

template<typename T>
inline Polynomial<T> operator-(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return Polynomial<T>(p) -= q;
}

template<typename T>
inline Polynomial<T> operator*(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return Polynomial<T>(p) *= q;
}

template<typename T>
inline Polynomial<T> operator/(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return Polynomial<T>(p) /= q;
}

template<typename T>
inline Polynomial<T> operator%(const Polynomial<T>& p, const Polynomial<T>& q)
{
    return Polynomial<T>(p) %= q;
}

//Unary operators -------------------------------------------------------------------------------

template<typename T>
Polynomial<T> operator-(const Polynomial<T>& p)
{
    Polynomial<T> result(p);
    for (size_t i(0); i<result.length(); ++i)
        result[i] = -result[i];
    return result;
}

template<typename T>
Polynomial<T> operator+(const Polynomial<T>& p)
{
    Polynomial<T> result(p);
    for (size_t i(0); i<result.length(); ++i)
        result[i] = +result[i];
    return result;
}

//Logical operators ----------------------------------------------------------------------------

template<typename T>
inline bool operator!(const Polynomial<T>& p)
{
    return p == Polynomial<T>(0);
}

//Derivative -----------------------------------------------------------------------------------

/**
 *  \brief Returns the derivative of the Polynomial passed
 */
template<typename T>
inline Polynomial<T> derivative(const Polynomial<T>& p)
{
    return p.derivative();
}

//Useful functions -----------------------------------------------------------------------------

/**
 *  \brief Returns the n-th power of a given x
 */
template<typename T>
T power(const T& x, size_t n)
{
    T result = static_cast<T>(1);

    for (size_t i(0); i<n ; ++i)
        result *= x;
    return result;
}

//Simplification typedef -----------------------------------------------------------------------

/**
 *  \typedef Poly
 *  A simpler name for standard polynomial using double as coefficients
 */
typedef Polynomial<double> Poly;

/**
 *  \typedef Polyd
 *  A simpler name for standard polynomial using double as coefficients
 */
typedef Polynomial<double> Polyd;

/**
 *  \typedef Polyf
 *  A simpler name for standard polynomial using float as coefficients
 */
typedef Polynomial<float> Polyf;

#endif //__POLYNOM_HPP_INCLUDED__


A bientôt et bonne chance avec le nouvel exercice.

EDIT: Simplification de ==. Merci lmghs.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
10 mai 2009 à 15:20:03

Je ne t'avais pas envoyé ce que j'avais fait, j'avais beaucoup de choses que je voulais rédiger (expliquer comment on écrit les opérateurs), plus une expérience à mener au sujet des expressions template. Mais j'ai eu d'autres choses à faire en court de route.
Je peux t'envoyer le code, ou le poster ici, comme tu préfères.

PS:
- polynomes et matrices sont deux "exceptions", où il est préférable d'implémenter *= en fonction de *.
- il suffit d'utiliser operator==(vector, vector) pour la comparaison.
C++: Blog|FAQ C++ dvpz|FAQ fclc++|FAQ Comeau|FAQ C++lite|FAQ BS| Bons livres sur le C++| PS: Je ne réponds pas aux questions techniques par MP.
Anonyme
10 mai 2009 à 18:07:02

Tiens, j'avais pas vu cet exo ! Je suis un peu dans le même cas que lmghs, j'ai réalisé avec un ami un programme de résolution d'(in)équation ( enfin, il manque l'IHM en php ... ) . Si ça intéresse : http://gmath.servhome.org/index.php .

3 choses :
a) Le logiciel en ligne n'est pas activé donc chargement infini. Une autre version est en développement mais aucune date n'est annoncée pour sa mise en utilisation libre.
b) Dans la doc, ce qui est utile ici se réduit à deux classes sauf erreur de ma part ( ça fait un bail que je ne l'ai pas retouché ) : BasePolynome et Polynome.
c) La division polynomiale est prise en charge.
10 mai 2009 à 18:16:32

@lmghs: Volontiers. Poste simplement ici.

@hiura: intéressant, mais ça ne marche pas chez moi ;-) Ou alors j'ai pas compris comment enrter l'équation. Quel algorithme utilisez-vous ? Newton + Horner ? Bézout ? Tschirnhausen ?

Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
10 mai 2009 à 18:24:45

Citation : lmghs

- polynomes et matrices sont deux "exceptions", où il est préférable d'implémenter *= en fonction de *.



Mmmh. J'y ai pas pensé pour les polynômes. Je vais regarder ça.

Citation : lmghs


- il suffit d'utiliser operator==(vector, vector) pour la comparaison.



En effet. Je simplifie. Merci.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
Anonyme
10 mai 2009 à 18:35:04

Citation : Nanoc

@hiura: intéressant, mais ça ne marche pas chez moi ;-)

Oui c'est normal, on n'a pas les droits pour exécuter l'application actuellement.

Citation : Nanoc

Quel algorithme utilisez-vous ? Newton + Horner ? Bézout ? Tschirnhausen ?

Je ne suis pas responsable de cette partie de l'application mais je dirais Newton. ( En trois ans de développement, et 6 mois d'inactivisme j'ai oublié certains détails :s . ) Je vais demander. Je peux simplement affirmé que pour les degrés 0 ( x € R ), et 1 ( x = -b / a ) et 2 ( Viète ) c'est d'autres méthodes plus simples qui sont employés.

[ Comme ça je dirais Newton car c'est le seul nom que je reconnais avec Bézout, mais comme quand on l'a codé on ne connaissait même pas les dérivées ... ]

PS : la conversion string -> polynôme n'est pas terrible, il aurait fallu qu'on connaisse les techniques de grammaire pour faire mieux.
10 mai 2009 à 18:41:52

Puis-je poster mon code aussi ?
Je viens de me rendre compte qu'il n'était pas si mal que ça, il permet quelques trucs en plus :
  • évaluation par la méthode de Horner (numériquement stable et moins lourd).
  • Calcul de la dérivée/primitive N-ième.
  • Calcul du pgcd.
  • Calcul du nombre de racines réelles comprises dans un segment (via le théorème de Sturm).
  • Et le plus gros : calcul des racines réelles par un mélange du théorème de Sturm, de la méthode de Newton, d'une approche dichotomique.
10 mai 2009 à 18:55:14

C'est en trois parties, et je vais juste vous donner la solution template, la non template est sans grand intérêt (et celle reposant sur les expressions template tout juste démarrée).

Pour pouvoir écrire: "const Polynomial A = (X^5) - (X^4) - (X^3) + (3*X^2) - (2*X);", je suis passé par une classe utilitaire pour les monômes. Et je continue avec deux très vilains hacks:
- le premier avec l'opérateur xor que je détourne pour les puissances -- opérateur dont la priorité est inadaptée à ce que l'on veut faire, d'où le parenthésage systématique car nécessaire dans l'expression ci-dessus.
- un second pour utiliser juste X sans avoir de .cpp, et en jouant avec de très vilains #define -- bref, les enfants ne faites pas ça au travail, ni même à la maison.
- l'évaluation repose sur la méthode de Hörner comme il se doit.

Ensuite, on peut remarquer que :
- pour des raisons de DRY et de symétrie, + et - sont définis en termes de += et -=, tandis que *= et /= sont définis en terme de * et /.
- les retours par valeur se font constants (pour refuser (2*X^2)*=42 -- à l'image de ce que fait le compilo pour les types entiers)
- je ne fournis pas de version non constante de l'opérateur[], ainsi nul besoin de sortir des proxys pour différencier l'accès RW (read-write) de l'accès RO (read-only) et s'enquiquiner à faire des trim à la volée.

Dans les dépendances, il faut rajouter array_helpers.hpp:
// array_helpers.hpp 
namespace luc_lib {
    // ...
    template <typename T, int N>
	inline std::reverse_iterator<T*> rbegin(T (&array)[N])
	{ return std::reverse_iterator<T*>(array + N); }

    template <typename T, int N>
	inline std::reverse_iterator<T*> rend  (T (&array)[N])
	{ return std::reverse_iterator<T*>(array + 0); }
}


Donc, sans plus attendre, voici le code principal.
#ifndef POLYNOMIAL_HPP
#define POLYNOMIAL_HPP

#include "../array_helpers.hpp"
#include <vector>
#include <iosfwd>
#include <cmath>

namespace luc_lib
{
    typedef size_t                 degree_type;

/*===========================================================================*/
/*================================[ is_null ]================================*/
/*===========================================================================*/
    template <typename N>
    inline 
    bool is_null(N n) {
        static const N k_err = 1e5;
        const bool null =
            (std::abs(n) <= k_err * std::numeric_limits<N>::epsilon());
        return null;
    }

    template <>
    inline 
    bool is_null<int>(int n) {
        return n == 0;
    }

/*===========================================================================*/
/*===============================[ Monomial ]================================*/
/*===========================================================================*/
   /** Monomial class.
    *
    * This class serves as an optimisation hack when defining new
    * polynomials.
    * @note default copy operations, and destructor are fine.
    */
    template <typename real_type>
    class Monomial
    {
    public:
	Monomial(real_type c, degree_type d) : m_coefficient(c), m_degree(d) {}
	explicit Monomial(degree_type d) : m_degree(d) {}

	real_type   coefficient() const { return m_coefficient; }
	degree_type degree     () const { return m_degree;      }

	Monomial& operator*=(real_type c) {
	    m_coefficient *= c;
	    return *this;
	}
	Monomial& operator^=(degree_type d) {
	    m_degree *= d;
	    return *this;
	}
    private:
        real_type      m_coefficient;
        degree_type    m_degree;
    };

    template <typename real_type>
    inline
    std::ostream& operator<<(std::ostream & os, Monomial<real_type> const& m)
    {
	return os << m.coefficient() << "X^" << m.degree() ;
    }

/*===========================================================================*/
/*==============================[ Polynomial ]===============================*/
/*===========================================================================*/
   /** Polynomial class.
    *
    * @note default copy operations, and destructor, provided thanks to
    * \c std::vector<> are fine.
    */
    template <typename real_type>
    class Polynomial 
    {
    public:
        typedef std::vector<real_type> coeffs_type;

	Polynomial(real_type c = real_type())
	    : m_coefficients(1)
	    { m_coefficients[0] = c; }

#if 1
        template <typename F>
            Polynomial(F c)
            : m_coefficients(1)
            { m_coefficients[0] = c; }
#endif

	/*explicit*/ Polynomial(Monomial<real_type> const& m)
	    : m_coefficients(m.degree()+1)
	    { m_coefficients.back() = m.coefficient();}

	explicit Polynomial(coeffs_type const& c)
	    : m_coefficients(c)
	    {}

	template <degree_type N> Polynomial(real_type const (&coeffs)[N])
	    : m_coefficients(luc_lib::rbegin(coeffs), luc_lib::rend(coeffs))
	    {
	    }

        void swap(Polynomial & other) {
            m_coefficients.swap(other.m_coefficients);
        }

	degree_type const degree() const { return m_coefficients.size()-1; }

	/** Evaluates the polynomial function.
	 * Implements Horner algorithm.
	 */
	real_type const operator()(real_type x) const
	{
	    real_type res = real_type(); // 0-init
	    degree_type i = degree()+1;
	    do {
		res = res * x + m_coefficients[--i];
	    } while (i) ;
	    return res;
	}

	/** Coefficients accessor.
	 */
	real_type operator[](degree_type d) const { return m_coefficients[d]; }

	Polynomial & operator+=(Polynomial const& rhs)
	{
	    const degree_type N = rhs.degree()+1;
	    if (m_coefficients.size() < N) { m_coefficients.resize(N); }
	    for (degree_type i=0; i!=N ; ++i)
		m_coefficients[i] += rhs[i];
	    do_trim();
	    return *this;
	}

	Polynomial & operator+=(Monomial<real_type> const& rhs)
	{
	    const degree_type N = rhs.degree()+1;
	    if (m_coefficients.size() < N) { m_coefficients.resize(N); }
	    // std::cout << *this << "(" << m_coefficients.size() << ") PLUS ";
	    m_coefficients[N-1] += rhs.coefficient();
	    // std::cout << rhs << " = " << *this<<"\n";
	    do_trim();
	    return *this;
	}

	Polynomial const operator-() const
	{
	    const degree_type N = degree()+1;
	    Polynomial res(N); // not very efficient
	    for (degree_type i=0; i!=N ; ++i)
		res.m_coefficients[i] = - m_coefficients[i];
	    return res;
	}

	Polynomial & operator-=(Polynomial const& rhs)
	{
	    const degree_type N = rhs.degree()+1;
	    if (m_coefficients.size() < N) { m_coefficients.resize(N); }
	    for (degree_type i=0; i!=N ; ++i)
		m_coefficients[i] -= rhs[i];
	    do_trim();
	    return *this;
	}

	Polynomial & operator-=(Monomial<real_type> const& rhs)
	{
	    const degree_type N = rhs.degree()+1;
	    if (m_coefficients.size() < N) { m_coefficients.resize(N); }
	    m_coefficients[N-1] -= rhs.coefficient();
	    do_trim();
	    return *this;
	}

	Polynomial & operator*=(Polynomial<real_type> const& rhs)
	{
	    Polynomial tmp = (*this) * rhs;
	    tmp.m_coefficients.swap(m_coefficients);
	    return *this;
	}

	friend Polynomial const operator*(Polynomial const& lhs, Polynomial const& rhs)
	{
	    const degree_type I = lhs.degree();
	    const degree_type J = rhs.degree();
	    /*Polynomial::*/coeffs_type res(I + J + 1);
	    for (degree_type i=0; i<=I ; ++i)
		for (degree_type j=0; j<=J ; ++j)
                    res[i+j] += lhs[i] * rhs[j];
	    return Polynomial(res);
	}

	Polynomial & operator*=(Monomial<real_type> const& rhs)
	{
	    (*this) *= rhs.coefficient();
	    m_coefficients.insert(m_coefficients.begin(), rhs.degree(), real_type());
	    return *this;
	}

	static void div(Polynomial lhs, Polynomial const& rhs, Polynomial &Q, Polynomial &R)
	    ;

	Polynomial & operator/=(Monomial<real_type> const& rhs) {
            Polynomial Q, R;
            div(*this, rhs, Q, R);
            this->swap(Q);
            return *this;
        }

	Polynomial & operator%=(Monomial<real_type> const& rhs) {
            Polynomial Q, R;
            div(*this, rhs, Q, R);
            this->swap(R);
            return *this;
        }

	friend bool operator==(Polynomial const& lhs, Polynomial const& rhs) {
	    return lhs.m_coefficients == rhs.m_coefficients;
	}

    private:
	void do_trim() {
	    while (m_coefficients.size()>1 && is_null(m_coefficients.back()))
	    {
                // std::cout << "trim:" << m_coefficients.back() << "\n";
		m_coefficients.pop_back();
	    }
	}
    private:
	coeffs_type m_coefficients;

	friend std::ostream& operator<<(std::ostream & os, Polynomial const& p)
	{
	    degree_type i = p.degree()+1;
	    if (i>1)
	    {
                bool needPlus = false;
		do
		{
		    if (! is_null(p[--i]))
		    {
			if (p[i] < 0) {
			    os << " - ";
			} else if (needPlus) {
			    os << " + ";
			}
			needPlus = true;
			if (p[i] != 1 || i==0) {
			    os << std::abs(p[i]) ; 
			}
			if (i>0) {
			    os << 'X';
			    if (i>1) {
				os << '^' << (i);
			    }
			}
		    }
		}
		while (i) ;
	    } else {
		os << (i==0 ? 0 : p[0]);
	    }
	    return os;
	}
    };


/*===========================================================================*/
/*===============================[ operators ]===============================*/
/*===========================================================================*/
    template <typename F, typename real_type>
    inline
    Monomial<real_type> const operator*(F c, Monomial<real_type> const& rhs)  {
        return Monomial<real_type>(rhs) *= c;
    }
    template <typename real_type>
    inline
    Monomial<real_type> const operator^(Monomial<real_type> const& lhs, degree_type d)  {
        return Monomial<real_type>(lhs) ^= d;
    }

/*===========================================================================*/

    template <typename real_type, typename convertible_real_type>
    inline
    Polynomial<real_type> const operator+(Polynomial<real_type> const& lhs, convertible_real_type const& rhs)  {
        return Polynomial<real_type>(lhs) += rhs;
    }
    template <typename real_type, typename convertible_real_type>
    inline
    Polynomial<real_type> const operator+(convertible_real_type const& lhs, Polynomial<real_type> const& rhs)  {
        return Polynomial<real_type>(rhs) += lhs;
    }

    template <typename real_type>
    inline
    Polynomial<real_type> const operator+(Polynomial<real_type> const& lhs, Monomial<real_type> const& rhs)  {
        return Polynomial<real_type>(lhs) += rhs;
    }
    template <typename real_type>
    inline
    Polynomial<real_type> const operator+(Monomial<real_type> const& lhs, Polynomial<real_type> const& rhs)  {
        return Polynomial<real_type>(rhs) += lhs;
    }

    template <typename real_type>
    inline
    Polynomial<real_type> const operator+(Polynomial<real_type> const& lhs, Polynomial<real_type> const& rhs)  {
        return Polynomial<real_type>(lhs) += rhs;
    }

    template <typename real_type>
    inline
    Polynomial<real_type> const operator+(Monomial<real_type> const& lhs, Monomial<real_type> const& rhs)  {
        return Polynomial<real_type>(lhs) += rhs;
    }
    template <typename real_type, typename convertible_real_type>
    inline
    Polynomial<real_type> const operator+(Monomial<real_type> const& lhs, convertible_real_type const& rhs)  {
        return Polynomial<real_type>(lhs) += rhs;
    }
    template <typename real_type, typename convertible_real_type>
    inline
    Polynomial<real_type> const operator+(convertible_real_type const& lhs, Monomial<real_type> const& rhs)  {
        return Polynomial<real_type>(rhs) += lhs;
    }

/*===========================================================================*/
    template <typename real_type>
    inline
    Polynomial<real_type> const operator-(Polynomial<real_type> const& lhs, Monomial<real_type> const& rhs)  {
        return Polynomial<real_type>(lhs) -= rhs;
    }
    template <typename real_type>
    inline
    Polynomial<real_type> const operator-(Monomial<real_type> const& lhs, Polynomial<real_type> const& rhs)  {
        return Polynomial<real_type>(-rhs) += lhs;
    }

    template <typename real_type>
    inline
    Polynomial<real_type> const operator-(Polynomial<real_type> const& lhs, Polynomial<real_type> const& rhs)  {
        return Polynomial<real_type>(lhs) -= rhs;
    }

    template <typename real_type, typename convertible_real_type>
    inline
    Polynomial<real_type> const operator-(Monomial<real_type> const& lhs, convertible_real_type const& rhs)  {
        return Polynomial<real_type>(lhs) -= rhs;
    }

    template <typename real_type>
    inline
    Polynomial<real_type> const operator-(Monomial<real_type> const& lhs, Monomial<real_type> const& rhs)  {
        return Polynomial<real_type>(lhs) -= rhs;
    }

    template <typename real_type, typename convertible_real_type>
    inline
    Polynomial<real_type> const operator-(convertible_real_type const& lhs, Monomial<real_type> const& rhs)  {
        return Polynomial<real_type>(rhs) -= lhs;
    }

/*===========================================================================*/
    template <typename real_type>
    inline
    Polynomial<real_type> const operator/(Polynomial<real_type> const& lhs, Polynomial<real_type> const& rhs)  {
        Polynomial<real_type> Q,R;
        Polynomial<real_type>::div(lhs, rhs, Q, R);
        return Q;
    }

    template <typename real_type>
    inline
    Polynomial<real_type> const operator%(Polynomial<real_type> const& lhs, Polynomial<real_type> const& rhs)  {
        Polynomial<real_type> Q,R;
        Polynomial<real_type>::div(lhs, rhs, Q, R);
        return R;
    }

#if  0
    inline
        std::ostream& operator<<(std::ostream & os, Polynomial const& p)
        {
        }
#endif

    template <typename real_type>
    inline
    void Polynomial<real_type>::div(
            Polynomial<real_type> lhs, Polynomial<real_type> const& rhs,
            Polynomial<real_type> &Q, Polynomial<real_type> &R)
    {
        Q.m_coefficients.resize(0); Q.m_coefficients.reserve(lhs.degree()-rhs.degree()+1);
        // R.m_coefficients.resize(0); R.m_coefficients.reserve(rhs.degree());

        while (lhs.degree() >= rhs.degree()) {
            const real_type   c = lhs[lhs.degree()] / rhs[rhs.degree()];
            const degree_type d = lhs.degree() - rhs.degree(); 
            const Monomial<real_type> m(c,d);
            lhs -= m * rhs;
            Q += m;
        }
        R.swap(lhs);
    }

    // const Monomial X(1,1);
#define X Monomial(1,1)    
} // namespace 


#endif // POLYNOMIAL_HPP
// Vim: let g:BTW_project='test_polynomial'

Et pour voir comment il s'utilise voici une baterie de test unitaires (reposant sur le framwork CxxTest -- tiens, le projet a trouvé des re-preneurs)
#ifndef ESSAIS_POLYNOMIAL_TEST_POLYNOMIAL_HPP__
#define ESSAIS_POLYNOMIAL_TEST_POLYNOMIAL_HPP__
#include <cxxtest/TestSuite.h>

#include <sstream>

// #define DOUBLE
#define TEMPLATE
// #define EXPR_TEMPLATE

#if defined(DOUBLE)
#  include "polynomial.hpp"
#elif defined(TEMPLATE)
#  include "polynomial_T.hpp"
#  define Polynomial Polynomial<double>
#  define Monomial   Monomial<double>
#elif defined(EXPR_TEMPLATE)
#  include "polynomial_ET.hpp"
#  define Polynomial Polynomial<double>
#  define Monomial   Monomial<double>
#endif
using namespace luc_lib;


namespace CxxTest 
{
    CXXTEST_TEMPLATE_INSTANTIATION
        class ValueTraits<Polynomial>
        {
            std::string value;
        public:
            ValueTraits(Polynomial const& p) {
                std::ostringstream oss;
                oss << p;
                value = oss.str();
            }
            const char *asString( void ) const { return value.c_str(); }
        };
}

/*===========================================================================*/
/*==============================[ ValueTraits ]==============================*/
/*===========================================================================*/

class test_polunomial : public CxxTest::TestSuite
{
public:
    void test_epsilon()
    {
	TS_TRACE(is_null(0));
	TS_TRACE(is_null(1e-9));
	TS_TRACE(is_null(1e-10));
	TS_TRACE(is_null(1e-11));

    }

/*===========================================================================*/
/*===============================[ Monomials ]===============================*/
/*===========================================================================*/
    void test_X()
    {
        Monomial m = X;
        TS_ASSERT_EQUALS(m.degree(), 1);
        TS_ASSERT_EQUALS(m.coefficient(), 1);
    }
    void test_X3()
    {
        Monomial m = X^3;
        TS_ASSERT_EQUALS(m.degree(), 3);
        TS_ASSERT_EQUALS(m.coefficient(), 1);
    }
    void test_2X()
    {
        Monomial m = 2 * X;
        TS_ASSERT_EQUALS(m.degree(), 1);
        TS_ASSERT_EQUALS(m.coefficient(), 2);
    }
    void test_2X3()
    {
        Monomial m = 2 * X^3;
        TS_ASSERT_EQUALS(m.degree(), 3);
        TS_ASSERT_EQUALS(m.coefficient(), 2);
    }

/*===========================================================================*/
/*==============================[ Polynomials ]==============================*/
/*===========================================================================*/
    void test_XpX()
    {
        Polynomial p = X + X;
        TS_ASSERT_EQUALS(p.degree(), 1);
        TS_ASSERT_EQUALS(p[1], 2);
        TS_ASSERT_EQUALS(p[0], 0);

        TS_ASSERT_EQUALS(p(0), 0);
        TS_ASSERT_EQUALS(p(1), 2);
        TS_ASSERT_EQUALS(p(2), 4);
        TS_ASSERT_EQUALS(p(3), 6);
    }
    void test_Xp5()
    {
        Polynomial p = X + 5;
        TS_ASSERT_EQUALS(p.degree(), 1);
        TS_ASSERT_EQUALS(p[1], 1);
        TS_ASSERT_EQUALS(p[0], 5);

        TS_ASSERT_EQUALS(p(0), 5);
        TS_ASSERT_EQUALS(p(1), 6);
        TS_ASSERT_EQUALS(p(2), 7);
        TS_ASSERT_EQUALS(p(3), 8);
    }
    void test_22X2p11Xp5()
    {
        Polynomial p = (22*X^2) + (11*X) + 5;
        p = (22*X^2) + (11*X) + 5;
        TS_ASSERT_EQUALS(p.degree(), 2);
        TS_ASSERT_EQUALS(p[2], 22);
        TS_ASSERT_EQUALS(p[1], 11);
        TS_ASSERT_EQUALS(p[0], 5);

        TS_ASSERT_EQUALS(p(0), 5);
        TS_ASSERT_EQUALS(p(1), 22*1*1 + 11*1 + 5);
        TS_ASSERT_EQUALS(p(2), 22*2*2 + 11*2 + 5);
        TS_ASSERT_EQUALS(p(3), 22*3*3 + 11*3 + 5);
    }
    void test_33X322X2p11Xp5()
    {
        Polynomial p = (33*X^3) + (22*X^2) + (11*X) + 5;
        p = (33*X^3) + (22*X^2) + (11*X) + 5;
        TS_ASSERT_EQUALS(p.degree(), 3);
        TS_ASSERT_EQUALS(p[3], 33);
        TS_ASSERT_EQUALS(p[2], 22);
        TS_ASSERT_EQUALS(p[1], 11);
        TS_ASSERT_EQUALS(p[0], 5);

        TS_ASSERT_EQUALS(p(0), 5);
        TS_ASSERT_EQUALS(p(1), 33*1*1*1 + 22*1*1 + 11*1 + 5);
        TS_ASSERT_EQUALS(p(2), 33*2*2*2 + 22*2*2 + 11*2 + 5);
        TS_ASSERT_EQUALS(p(3), 33*3*3*3 + 22*3*3 + 11*3 + 5);
    }

    void test_mult_Xp1_sq()
    {
        const Polynomial p1 = (X+1);
        const Polynomial m = p1 * p1;

        TS_ASSERT_EQUALS(m.degree(), 2);
        TS_ASSERT_EQUALS(m[2], 1);
        TS_ASSERT_EQUALS(m[1], 2);
        TS_ASSERT_EQUALS(m[0], 1);
    }
    void test_mult_Xp1_dot_Xm1()
    {
        const Polynomial p1 = (X+1);
        const Polynomial p2 = (X-1);
        const Polynomial m = p1 * p2;

        TS_ASSERT_EQUALS(m.degree(), 2);
        TS_ASSERT_EQUALS(m[2], 1);
        TS_ASSERT_EQUALS(m[1], 0);
        TS_ASSERT_EQUALS(m[0], -1);
    }

    void test_div()
    {
	// const Polynomial A = (-7*X^3) + (2*X^2) + (3*X) + 1;
	// const Polynomial B = (-2*X^2) + (1*X) + 1;
	// const Polynomial eQ = (-5*X^3) + (2*X^2) + (2*X) + 1;
	// const Polynomial eR = -10*X + 9;

	const Polynomial A = (X^5) - (X^4) - (X^3) + (3*X^2) - (2*X);
	const Polynomial B = (X^2) - (X) + 1;
	const Polynomial eQ = (X^3) - (2*X) + 1;
	const Polynomial eR = X - 1;

	// const Polynomial A = (X^5) + (2*X^3) - (3*X) - 2;
	// const Polynomial B = (X^3) + (X) + 1;
	// const Polynomial eQ = (X^2) + 1;
	// const Polynomial eR = (-1*X^2) - (4*X) - 3;

        Polynomial Q, R;
	Polynomial::div(A, B, Q, R);
	TS_ASSERT_EQUALS(Q, eQ);
	TS_ASSERT_EQUALS(R, eR);
	TS_ASSERT_EQUALS(A, B*Q+R);
    }
};


#endif // ESSAIS_POLYNOMIAL_TEST_POLYNOMIAL_HPP__

et le Makefile qui va bien:
# ---------------------------------------------------------------------------
# Development tools directories definition
# ---------------------------------------------------------------------------
# CC_OPT_XARCH    = -xarch=v9
# LINK_OPT_XARCH  = -xarch=v9
# CC_DIR          = /opt/SUNWspro/bin
# CXX             = ${CC_DIR}/CC
# CXXFLAGS        = -g ${CC_OPT_XARCH} -DUNIX -DSUN 
# INCLUDE_PATH    = -I${HOME}/include

INCLUDE_PATH    = -I/cygdrive/e/Dev/include
CXXFLAGS        = -g -pedantic -Wall ${INCLUDE_PATH}
# ---------------------------------------------------------------------------
# Définition de la construction de la cible principale
# ---------------------------------------------------------------------------
all: test_polynomial

# ---------------------------------------------------------------------------
test_polynomial.cpp: test_polynomial.hpp polynomial.hpp polynomial_T.hpp polynomial_ET.hpp
        cxxtestgen.pl -error-printer -o $@ $^

# test_polynomial: test_polynomial.cpp 
        # $(CPP_COMPILER) ${CC_FLAGS} ${CC_OPT_FLAGS} ${INCLUDE_PATH} -o $@ $^ ${LIBRARIES}  
        # ./$@


# ---------------------------------------------------------------------------
# cibles de nettoyage
# ---------------------------------------------------------------------------
clean: 
        rm test_polynomial test_polynomial.o test_polynomial.cpp
C++: Blog|FAQ C++ dvpz|FAQ fclc++|FAQ Comeau|FAQ C++lite|FAQ BS| Bons livres sur le C++| PS: Je ne réponds pas aux questions techniques par MP.
10 mai 2009 à 20:24:41

Bon bah je m'installe :D , mon code n'est pas aussi pro que celui de lmghs, ne gère pas les templates (bah oui j'en suis pas encore là dans le cours ^^ ) mais je le préfère pour les outils mathématiques élémentaires qu'il apporte.

Le header

#ifndef POLYNOME_H
#define POLYNOME_H

#include <vector>
#include <iostream>
#include <stdarg.h>
#include <stdexcept>

#define POLYNOME_NB_ITERATIONS_RECHERCHE 10
#define POLYNOME_NB_ITERATIONS_NEWTON 1000
#define POLYNOME_NB_ITERATIONS_DICHO 10
#define POLYNOME_BIG_EPSILON 0.00000001
#define POLYNOME_EPSILON 0.000000000000000001

class Polynome{
    public:
    // Constructeurs
        Polynome();
        Polynome(const Polynome &);
        Polynome(const std::vector<double> &);
        Polynome(int, const double, ...);
        Polynome(const int, const double[]);
    // Destructeur
        ~Polynome();

    // Opérateurs
        Polynome operator=(const Polynome &);

        Polynome operator+=(const Polynome &P);
        Polynome operator-=(const Polynome &P);
        Polynome operator*=(const Polynome &P);
        Polynome operator/=(const Polynome &P) throw(std::domain_error);
        Polynome operator%=(const Polynome &P) throw(std::domain_error);

        Polynome operator+(const Polynome &) const;
        Polynome operator-(const Polynome &) const;
        Polynome operator-() const;
        Polynome operator*(const Polynome &) const;
        Polynome operator/(const Polynome &) const throw(std::domain_error);
        Polynome operator%(const Polynome &) const throw(std::domain_error);
        bool operator==(const Polynome &) const;
        bool operator!=(const Polynome &) const;
        double operator()(const double) const;

    // Arithmétique
        static void div_eucli(const Polynome &, const Polynome &, Polynome &, Polynome &) throw(std::domain_error);
        static Polynome pgcd(const Polynome &, const Polynome &);
        Polynome derive(unsigned int d=1) const;
        Polynome integre(unsigned int d=1) const;
        void unitaire();

    // Racines
        int nb_roots();
        int nb_roots(double, double);
        std::vector<double> roots();
        void sturm_compute();
        // Accesseurs
            bool get_sturmcomputed();
        double newton(double) throw(std::domain_error, std::runtime_error);
        void dicho(double &, double &, int);

    // Autres
        double eval(double) const;
        int degre() const;
        void affiche(std::ostream &) const;

    protected:
        int nb_chg_signe(double);

    // Attributs
        std::vector<double> m_coeffs;   // vecteur des coefficients
        std::vector<Polynome> Sturm;    // Suite de Sturm du polynôme
        bool sturm_computed;            // Si la suite de Sturm a été calculée
        double M;                       // Majorant de la valeur absolue de toutes les racines

};

std::ostream& operator<<(std::ostream &, const Polynome &);

#endif



le fichier cpp

#include "polynome.h"
#include <cmath>

#include <iomanip>

using namespace std;


/********************************************************
                    Constructeurs
 *******************************************************/

// Constructeur de copie
Polynome::Polynome(const Polynome &P){
    m_coeffs = P.m_coeffs;
    sturm_computed = P.sturm_computed;
    Sturm = P.Sturm;
}

// Par défault : Créer le polynome nul
Polynome::Polynome(){
    m_coeffs.push_back(0.);
    sturm_computed = false;
}

// Via un vecteur : les coefficients du polynôme sont recopié du vecteur
Polynome::Polynome(const std::vector<double> &coeffs){
    m_coeffs = coeffs;
    sturm_computed = false;
}

// Via une liste variable : Polynome(nb de coeffs, a0, a1, a2, ...)
// ex : Polynome(3, 1., -2., 4.) crée le polynôme 4x^2-2x+1
// La liste variable est constituée de double !!! Polynome(2, 1, 2) ne renverra pas 2x+1 !!!
Polynome::Polynome(int n, const double x, ...){
    sturm_computed = false;
    if (n <= 0)
        return;

    va_list varg; // Variable stockant la liste
    va_start(varg, x); // On initialise la variable
    m_coeffs.push_back(x);
    n--;
    do {
        m_coeffs.push_back(va_arg(varg, double));
    }while(--n!=0);
    va_end(varg); // On détruit la liste variable
}

// Via un tableau de double, ceci permet de profiter de l'auto-cast int -> double
// et donc de pouvoir écrire Polynome(3, (double[]) {1, 0, 1.}) pour le polynôme x^2+1
// ou bien double coef[] = {1, 0, 1.}; Polynome(3, coef)
Polynome::Polynome(const int n, const double *liste){
    sturm_computed = false;
    if (n <= 0)
        return;
    for(int i = 0 ; i < n ; i++)
        m_coeffs.push_back(liste[i]);
}



/********************************************************
 *******************  Destructeur  **********************
 *******************************************************/
Polynome::~Polynome(){
    // Rien de bien palpitant
}




/********************************************************
 ******************  Opérateurs  ************************
 *******************************************************/

// Opérateur d'affectation
Polynome Polynome::operator=(const Polynome &P){
    m_coeffs = P.m_coeffs;
    sturm_computed = P.sturm_computed;
    Sturm = P.Sturm;
    return *this;
}

Polynome Polynome::operator+=(const Polynome &P)
{
    int n = this->degre();
    int m = P.degre();
    int i;

    if (m <= n)
    {
        // On rajoute un polynome de degré plus petit ou égal
        for(i = 0 ; i <= m ; i++)
            m_coeffs[i] += P.m_coeffs[i];
    }
    else
    {
        // On rajoute un polynôme de degré plus grand
        m_coeffs.reserve(m+1);
        for(i = 0 ; i <= n ; i++)
            m_coeffs[i] += P.m_coeffs[i];
        for(i = n + 1 ; i<= m ; i++)
            m_coeffs[i] = P.m_coeffs[i];
    }
    return *this;
}

Polynome Polynome::operator-=(const Polynome &P)
{
    return (*this)+=-P;
}

Polynome Polynome::operator*=(const Polynome &P)
{
    *this = *this * P;
    return *this;
}

Polynome Polynome::operator/=(const Polynome &P) throw(std::domain_error)
{
    if(P.degre() == -1)
        throw std::domain_error("Division par zéro");
    *this = *this / P;
    return *this;
}

Polynome Polynome::operator%=(const Polynome &P) throw(std::domain_error)
{
    if(P.degre() == -1)
        throw std::domain_error("Division par zéro");
    *this = *this % P;
    return *this;
}


// Ajoute 2 polynômes entre eux
Polynome Polynome::operator+(const Polynome &P) const{
    return Polynome(*this)+=P;
}

// Retranche 2 polynômes
Polynome Polynome::operator-(const Polynome &P) const{
    return Polynome(*this)+=-P;
}

// Renvoie l'opposée du polynôme
Polynome Polynome::operator-() const{
    int n = m_coeffs.size();

    Polynome res;
    res.m_coeffs.resize(n);

    for(int i = 0 ; i < n ; i++)
        res.m_coeffs[i] = - m_coeffs[i];
    return res;
}

// Multiplie 2 polynômes entre eux : X^i * X^j = X^(i+j)
Polynome Polynome::operator*(const Polynome &P) const{
    int n = this->degre();
    int m = P.degre();
    if (n < 0 || m < 0)
        return Polynome();
    vector<double> M(n+m+1,0);
    for(int i = 0 ; i<= n ; i++)
        for(int j = 0 ; j <= m ; j++)
            M[i+j] += m_coeffs[i]*P.m_coeffs[j];
    return Polynome(M);
}

// Renvoie le quotient de la division euclidienne
Polynome Polynome::operator/(const Polynome &P) const throw(std::domain_error){
    if(P.degre() == -1)
        throw std::domain_error("Division par zéro");
    Polynome Q,R;
    div_eucli(*this, P, Q, R);
    return Q;
}

// Renvoie le reste de la division euclidienne
Polynome Polynome::operator%(const Polynome &P) const  throw(std::domain_error){
    if(P.degre() == -1)
        throw std::domain_error("Division par zéro");
    Polynome Q;
    Polynome R;
    div_eucli(*this, P, Q, R);
    return R;
}

// Test l'égalité de 2 polynômes
bool Polynome::operator==(const Polynome &P) const{
    int n = this->degre();
    int m = P.degre();
    if (n != m)
        return false;
    for(int i = 0 ; i <=n ; i++)
        if(m_coeffs[i] != P.m_coeffs[i])
            return false;
    return true;
}

// Test la différence de 2 polynômes
bool Polynome::operator!=(const Polynome &P) const{
    return !(*this == P);
}

// Confie à Poly.affiche la tache cout << P
ostream &operator<<(ostream &out, const Polynome &Poly){
    Poly.affiche(out);
    return out;
}

// Evalue un polynôme en 1 point
double Polynome::operator()(const double x) const{
    return this->eval(x);
}






/********************************************************
 ******************  Arithmétique  **********************
 *******************************************************/

// Effectue la division euclidienne de A par B, le quotient va dans Q et le reste dans R
void Polynome::div_eucli(const Polynome &A, const Polynome &B, Polynome &Q, Polynome &R) throw(std::domain_error){
    if(B.degre() == -1)
        throw std::domain_error("Division par zéro");

    int n = A.degre();
    int m = B.degre();
    if (m < 0)
        return;
    if (n-m < 0){
        Q = Polynome();
        R = A;
        return;
    }

    R.m_coeffs = A.m_coeffs;
    Q.m_coeffs.resize(n-m+1);
    double pivot = 1/B.m_coeffs[m], r;

    for(int i = n-m ; i >= 0 ; i--){
        r = R.m_coeffs[i+m] * pivot;
        Q.m_coeffs[i] = r;
        if(r != 0)
            for(int j = i ; j<=i+m ; j++){
                R.m_coeffs[j] -= r * B.m_coeffs[j-i];
                if ( abs(R.m_coeffs[j]) < POLYNOME_BIG_EPSILON)
                    R.m_coeffs[j] = 0;
            }
        R.m_coeffs[i+m] = 0;
    }
}

// Renvoie le pgcd de 2 polynômes
Polynome Polynome::pgcd(const Polynome &A, const Polynome &B){
    Polynome A1 = A, B1 = B, Q, R;
    do{
        R = A1 % B1;
        A1 = B1;
        B1 = R;
    } while(R.degre() != -1);

    return A1;
}

// Renvoie la dérivée du polynôme
Polynome::Polynome Polynome::derive(unsigned int d) const{
    int n = this->degre();
    if (n+1-(int)d <= 0)
        return Polynome();
    vector<double> res(n+1-d);
    int p = 1;
    for(unsigned int i = 1 ; i <= d ; i++)
        p *= i;
    for(int i = 0; i <= n-(int)d ; i++){
        res[i] = p*m_coeffs[i+d];
        p = p * (d+i+1)/(i+1);
    }
    return Polynome(res);
}

// Renvoie la primitive du polynôme nulle en 0
Polynome::Polynome Polynome::integre(unsigned int d) const{
    int n = this->degre();
    if (n+1+d <= 0)
        return Polynome();
    vector<double> res(n+1+d);
    int p = 1;
    for(unsigned int i = 1 ; i <= d ; i++)
        p *= i;
    for(unsigned int i = 0 ; i < d ; i++)
        res[i] = 0;
    for(int i = 0 ; i <= n ; i++){
        res[i+d] = m_coeffs[i]/p;
        p = p * (d+i+1)/(i+1);
    }
    return Polynome(res);
}

// Rend le polynôme unitaire (son coefficient de plus haut degré vaudra 1)
void Polynome::unitaire(){
    int n = this->degre();
    if (n == -1)
        return;
    double div = 1 / m_coeffs[n];
    for(int i = 0 ; i <= n ; i++)
        m_coeffs[i] *= div;
    return;
}





/********************************************************
 ********************* Racines  *************************
 *******************************************************/

// Renvoie le nombre de racines réelles du polynômes (Théorème de Sturm)
int Polynome::nb_roots(){

    if( !sturm_computed )
        sturm_compute();

    return abs(nb_chg_signe(-M) - nb_chg_signe(M));
}

// Renvoie le nombre de racines réelles du polynômes dans l'intervalle [a,b] (Théorème de Sturm)
int Polynome::nb_roots(double a, double b){

    if( !sturm_computed )
        sturm_compute();

    return abs(nb_chg_signe(a) - nb_chg_signe(b));
}

// Calcul les racines réelles d'un polynôme
vector<double> Polynome::roots(){
    int nb_racines = nb_roots();

    if (nb_racines == 0)
        return vector<double>();

// Recherches les intervalles sur lesquels les racines sont isolés
// Cette recherche se fait via une dichotomie utilisée sur le théorème de Sturm
    double *bornes = new double[nb_racines + 1]; // Pour n racines, il y a n+1 intervalles
    bool *borne_fixe = new bool[nb_racines + 1];
    for(int i=0;i<nb_racines+1;i++)
        borne_fixe[i] = false;
    vector<double> racines;

    bornes[0] = -M; bornes[nb_racines] = M;        // On fixe les bornes limites [-M,M]
    borne_fixe[0] = true; borne_fixe[nb_racines] = true;

    for(int i = 1 ; i < nb_racines ; i++){

        // On recherche le premier intervalle qui contient plusieurs racines pour le scindé en 2
        // Derniere borne fixe
        int gauche = 0, droite, rac_gauche, rac_droite, rac_milieu, nb_rac;
        double x_g, x_d, milieu;
        bool continuer;

        while(borne_fixe[gauche])
            gauche++;
        gauche--;
        // Cette erreur ne devrait jamais se produire
        if (gauche > nb_racines){
            std::cerr << "Erreur lors du scindage de l'intervalle (recherche de la borne gauche\n";
            std::cerr << left << setw(20) << "Valeur des bornes" << "Fixé\n";
            for(int i = 0 ; i <= nb_racines ; i++)
                std::cerr << left << setw(20) << bornes[i] << borne_fixe[i] << "\n";
            std::cerr << flush;
            return vector<double>(0);
        }

        // Prochaine borne fixe
        droite = gauche + 1;
        while(!borne_fixe[droite])
            droite++;
        // Cette erreur ne devrait jamais se produire
        if (droite > nb_racines){
            std::cerr << "Erreur lors du scindage de l'intervalle (recherche de la borne droite\n";
            std::cerr << left << setw(20) << "Valeur des bornes" << "Fixé\n";
            for(int i = 0 ; i <= nb_racines ; i++)
                std::cerr << left << setw(20) << bornes[i] << borne_fixe[i] << "\n";
            std::cerr << flush;
            return vector<double>(0);
        }

        x_g = bornes[gauche];
        x_d = bornes[droite];
        rac_gauche = nb_chg_signe(bornes[gauche]);
        rac_droite = nb_chg_signe(bornes[droite]);
        // On cherche à couper l'intervalle [x_g, x_d] en 2 intervalles [x_g, milieu] et [x_d, milieu]
        // où chacun des 2 contient au moins une racine
        continuer = true;
        while(continuer){
            milieu = (x_g + x_d)/2;
            rac_milieu = nb_chg_signe(milieu);
            if( rac_milieu == rac_gauche)       // Il n'y a pas de racines entre x_g et milieu
                x_g = milieu;                       // donc on décale la borne gauche
            else if( rac_milieu == rac_droite)  // Il n'y a pas de racines entre milieu et x_d
                x_d = milieu;                       // donc on décale la borne droite
            else
                continuer = false;              // On a réussi à couper l'intervalle
        }

        // On scinde l'intervalle
        nb_rac = abs(rac_gauche - rac_milieu);
        bornes[gauche + nb_rac] = milieu;
        borne_fixe[gauche + nb_rac] = true;

        // Réduit la taille de l'intervalle total [-M, M] pour accélerer la convergence de la méthode de newton
        if (gauche == 0)
            bornes[gauche] = x_g;
        if (droite == nb_racines)
            bornes[droite] = x_d;
    }

    Polynome calc = *this / Polynome::pgcd(*this, this->derive()); // Permet de transformer toutes les racines de multiplicité > 1 en racines de multiplicité 1
    for(int i = 0 ; i < nb_racines ; i++){               // Ce qui permet une convergence quadratique de la méthode de newton
        int iter = 0;       // Nombre d'itérations pour le calcul d'une racine
        double rac;
        double borne_gauche = bornes[i], borne_droite = bornes[i+1];    // Bornes de recherche de la racine
        bool calcul = true; // On effectue le calcul
        while(calcul){
            try{
                iter++;
                rac = calc.newton((borne_gauche + borne_droite)/2);
                if (rac >= bornes[i] && rac <= bornes[i+1]){ // Si la racine se situe dans l'intervalle de recherche
                    racines.push_back(rac); // On l'empile sur le vecteur de solutions
                    calcul = false;         // Le calcul s'est effectué jusqu'au bout
                    calc = calc / Polynome(2, -racines[i], 1.); // Divise le polynôme de recherche par (X-rac) pour supprimer la racine (enlève un bassin d'attraction pour la méthode de newton)
                }
            }
            catch(const std::domain_error &e){
                // La dérivée s'est annulée, laisse calcul à true
            }
            catch(const std::runtime_error &e){
                // L'algorithme tourne en boucle, laisse calcul à true
            }
            if (calcul){ // Si l'on doit continuer à chercher la racine
                calc.dicho(borne_gauche, borne_droite, POLYNOME_NB_ITERATIONS_DICHO); // On divise l'intervalle de recherche par 2^POLYNOME_NB_ITERATIONS_DICHO
            }
            if (iter > POLYNOME_NB_ITERATIONS_RECHERCHE) // Si l'on a effectué trop d'itérations, on abandonne la recherche
                calcul = false;
        }
    }

    delete[] bornes;
    delete[] borne_fixe;

    return racines;
}

// Construit la suite de Sturm
void Polynome::sturm_compute(){
    Polynome A = *this;
    Polynome B = A.derive();
    Polynome R;
    A.unitaire();
    Sturm.clear();
    Sturm.push_back(A);
    int r = 1;

    do{
        Sturm.push_back(B);
        r++;
        R = - (A % B);
        A = B;
        B = R;
    } while(R.degre() != -1);

 // Si le polynôme a des racines multiples on réduit leur multiplicité à 1
    if (Sturm[r-1].degre() != 0)
        for(int i = 0 ; i < r ; i++)
            Sturm[i] = Sturm[i] / Sturm[r-1];

// Toutes les racines seront dans [-M, M]
    int degre = this->degre();
    M = 0;
    double coeff = m_coeffs[degre];
    for(int i = 0 ; i <= degre ; i++)
        M += abs(m_coeffs[i] / coeff);
    if (M < 1)
        M = 1;

    sturm_computed = true;
}

// Indique si la suite de Sturm est construite ou non
bool Polynome::get_sturmcomputed(){
    return sturm_computed;
}

// Calcul le nombre de changement de signe utilisé dans le théorème de Sturm
int Polynome::nb_chg_signe(double x){
    int nb = 0;
    int n = Sturm.size();

    // On recherche le premier terme non nul
    int i = 0;
    while(i < n && Sturm[i].eval(x) == 0)
        i++;
    if (i == n)
        return 0; // Tout les termes sont nuls

    double terme_precedent = Sturm[i].eval(x), calc;
    for(i = i+1 ; i < n ; i++){
        calc = Sturm[i].eval(x);
        if (terme_precedent * calc < 0){
            nb++;
            terme_precedent = calc;
        }
    }

    return nb;
}

// Approxime la racine via l'algorithme de Newton
double Polynome::newton(double x) throw(std::domain_error, std::runtime_error){
    Polynome DP = this->derive();
    int iter = 0;
    double f_x = eval(x), fp_x;
    while(abs(f_x) > POLYNOME_EPSILON){
        fp_x = DP.eval(x);

        // Exceptions
        if (abs(fp_x) < POLYNOME_EPSILON)    // La dérivée est nulle
            throw std::domain_error("Division par zéro");
        if (iter > POLYNOME_NB_ITERATIONS_NEWTON){                   // On a effectué trop d'itérations
            if (abs(f_x) < POLYNOME_BIG_EPSILON)            // Mais on est presque sur un zéro
                return x;
            throw std::runtime_error("Trop d'itérations");
        }

        x = x - f_x/fp_x;
        f_x = eval(x);
        iter++;
    }
    return x;
}

// Effectue un dichotomie sur un nombre d'itérations fixés
// Elle réduit l'intervalle [b_g, b_d] autour de la racine
void Polynome::dicho(double &b_g, double &b_d, int iter){
    double val_g = eval(b_g), val_d = eval(b_d);
    double milieu, val_m;

    if (val_g * val_d > 0)      // On ne peut pas lancer la dichotomie
        return;

    while (iter > 0){
        milieu = (b_g + b_d) / 2;
        val_m = eval(milieu);
        if (val_g * val_m <= 0){ // La racine est située dans [b_g, milieu]
            b_d = milieu; // On se place sur le segment [b_g, milieu]
        }
        else{                   // La racine est située dans [milieu, b_d]
            b_g = milieu; // On se place sur le segment [milieu, b_d]
            val_g = val_m;
        }
        iter--;
    }
}


/********************************************************
 ********************** Autres  *************************
 *******************************************************/

// Evalue le polynôme en un point, méthode de Horner
double Polynome::eval(double x) const{
    double res = 0;
    int n = this->degre();
    for(int i = n ; i >= 0 ; i--){
        res*= x;
        res+= m_coeffs[i];
    }
    return res;
}

// Renvoie le degré du polynôme. Par convention le degré du polynôme nul sera -1 et non pas - l'infini
int Polynome::degre() const{
    int n = m_coeffs.size() - 1;
    while(n >= 0 && m_coeffs[n] == 0)
        n--;
    return n;
}

// Envoie l'expression litérale du polynôme sur out (utiliser par cout << P)
void Polynome::affiche(ostream &out) const{
    int i;
    double c;
    bool poly_nul = true;
    for(i=m_coeffs.size()-1 ; i!=-1 ; i--){
        c = m_coeffs[i];
        if (i == 0){
            if(poly_nul)
                out << noshowpos << c;
            else
                if(c != 0)
                    out << showpos << c;
        }
        else{
            if(poly_nul && c != 0){
                if(c == -1)
                    out << '-';
                else if(c != 1)
                    out << noshowpos << c;
                poly_nul = false;
            }
            else if(c != 0){
                if( c == 1)
                    out << '+';
                else if(c == -1)
                    out << '-';
                else
                    out << showpos << c;
            }
        }
        if(c != 0){
            if(i==1)
                out << 'X';
            else if(i > 1)
                out << "X^(" << noshowpos << i <<')';
        }
    }
    out << noshowpos;
}



et un petit main pour tester ça

#include <iostream>
#include "polynome.h"
#include <windows.h>
using namespace std;

// Affiche un vector<double>
std::ostream &operator<<(std::ostream &out, const std::vector<double> &vec){
    int n = vec.size();
    for(int i = 0 ; i < n ; i++)
        out << "[" << i << "] = " << vec[i] << endl;
    return out;
}

int main(int argc, char *argv[])
{
//    SetConsoleOutputCP( 1252 );

    Polynome A(6, (double[]) {1, 0, 0, -1, 4., 1.});
    cout << "A(X) = " << A << endl;

    Polynome B(3, 1., 0., 1.); // Création via une liste variable de double
    cout << "B(X) = " << B << endl;

    cout << endl << "Division euclidienne : A = B * Q + R" << endl;
    Polynome Q, R;
    Polynome::div_eucli(A,B,Q,R);

    cout << "Q(X) = " << Q << endl;
    cout << "R(X) = " << R << endl;

    cout << "Vérification (on teste en même temps l'opérateur d'égalité) : ";
    if (A == B * Q + R)
        cout << "ok" << endl;
    else
        cout << "erreur" << endl;

    cout << endl << "A(X) / B(X) = " << A / B << endl;
    cout << "A(X) % B(X) = " << A % B << endl;

    cout << endl << "Calcul de PGCD" << endl;
    double coeffs_A[] = {7, -7, 0, 8, -7, 0, 1};
    A = Polynome(7, coeffs_A);
    B = Polynome(6, (double[]) {-7, 0, 3, -7, 0, 3});
    cout << "A(X) = " << A << endl << "B(X) = " << B << endl;

    R = Polynome::pgcd(B,A);
    cout << "PGCD(A,B) : " << R << endl;
    R.unitaire();
    cout << "PGCD(A,B) unitarisé : " << R << endl;

    A = Polynome(6, -12., 20., -1., -9., 1., 1.);
    //A = Polynome(3, -1., 0., 1.);
    A = Polynome(11, 30., 30., -11., 35., -13., 36., 36., 20., 42., 15., -34.);
    cout << endl << "A(X) = " << A << endl;
    cout << "Nombre de racines réelles : " << A.nb_roots() << endl;

    cout << "Racines : " << endl << A.roots();

    return 0;
}



ce qui donne au final :
A(X) = X^(5)+4X^(4)-X^(3)+1
B(X) = X^(2)+1

Division euclidienne : A = B * Q + R
Q(X) = X^(3)+4X^(2)-2X-4
R(X) = 2X+5
Vérification (on teste en même temps l'opérateur d'égalité) : ok

A(X) / B(X) = X^(3)+4X^(2)-2X-4
A(X) % B(X) = 2X+5

Calcul de PGCD
A(X) = X^(6)-7X^(4)+8X^(3)-7X+7
B(X) = 3X^(5)-7X^(3)+3X^(2)-7
PGCD(A,B) : -0.25X^(3)-0.25
PGCD(A,B) unitarisé : X^(3)+1

A(X) = -34X^(10)+15X^(9)+42X^(8)+20X^(7)+36X^(6)+36X^(5)-13X^(4)+35X^(3)-11X^(2)
+30X+30
Nombre de racines réelles : 2
Racines :
[0] = -0.574217
[1] = 1.71694

Process returned 0 (0x0)   execution time : 0.044 s
Press any key to continue.

10 mai 2009 à 20:59:14

Rapidement,
a- ton += n'est pas bon, tu as utilisé reserve au lieu de resize (c'est reserve+push_back, ou resize + op[])
b- utilise la symétrie, tu as tout à y gagner
c- pas de iostream dans un .h (iosfwd si tu passes des stream en paramètre, ostream/istream pour implémenter des opérateurs sur les flux, et iostream pour cout/cin/clog/cerr -- en tirant moins de choses, on y gagne sur les temps de compilation)
d- les spécifications d'exception sont des attrape-couillons qu'il veut mieux bannir de nos codes
e- pas besoin de stdarg dans le .h non plus.
f- il n'y a pas besoin de destructeur, ni de constructeur de copie, ni d'opérateur d'affectation.
g- les opérateurs x= devraient renvoyer des références et non des copies de *this.
h- pour recevoir le tableau de doubles, cf mon constructeur ligne 104 qui permet de recevoir directement un tableau _statique_ de doubles.
i- Pour les racines, j'aurais géré ça de manière complètement externe : les racines ne font pas parti du cœur du polynôme -- tu aurais pû profiter de l'info pour la multiplication, voire même la division, mais cela n'a pas été le cas.
h- donne des noms à tes paramètres dans les signatures des fonctions que tu déclares. En particulier pour la division (
C++: Blog|FAQ C++ dvpz|FAQ fclc++|FAQ Comeau|FAQ C++lite|FAQ BS| Bons livres sur le C++| PS: Je ne réponds pas aux questions techniques par MP.
10 mai 2009 à 21:05:19

Merci pour tout ses conseils :)
Qu'entends tu par utilise la symétrie ?
Idem pour ton petit d (spécification d'exception) :D
10 mai 2009 à 22:04:49

Les problèmes de symétrie ont à plusieurs reprises été évoqués ici (ce fil) et ailleurs (le forum, et encore plus ailleurs). Je te laisse fouiller.

Pour les spécifications d'exception, je te laisse lire cette entrée du GOTW: http://www.gotw.ca/gotw/082.htm ...
Hum j'espérai que cela soit plus verbeux que cela.
Rapidement une spec d'exception n'est pas vérifiée à la compilation, mais à l'exécution. L'effet est que si une exception du mauvais type cherche à sortir, alors le programme sera arrêté. Oh on peut intervenir et faire 2-3 trucs au moment de l'arrêt, mais arrêt il y aura.
Et ce n'est pas du tout un comportement acceptable -- et en plus cela ralentit l'exécution vu que des vérifications sont réalisées.
C++: Blog|FAQ C++ dvpz|FAQ fclc++|FAQ Comeau|FAQ C++lite|FAQ BS| Bons livres sur le C++| PS: Je ne réponds pas aux questions techniques par MP.
11 mai 2009 à 1:13:56

Petite question au passage. Y a-t'il un avantage à déclarer le moins unaire comme fonction membre de la classe ?
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
11 mai 2009 à 2:07:03

Pas que je sache.
Enfin, il y a toujours ces histoires de conversions implicites qui doivent pouvoir s'appliquer.
Genre s'il existe un type Toto dépourvu du moins unaire, et un constructeur implicite qui prend un Toto. Je m'avance peut-être. C'est à tester.
C++: Blog|FAQ C++ dvpz|FAQ fclc++|FAQ Comeau|FAQ C++lite|FAQ BS| Bons livres sur le C++| PS: Je ne réponds pas aux questions techniques par MP.
11 mai 2009 à 9:47:37

Bonjour,
quelqu'un peut m'expliquer
typedef std::pair<Polynome, Polynome >  Paire;

    static Paire division(const Polynome& dividende, const Polynome& diviseur)

et en particulier le typedef std::pair<Polynome, Polynome > Paire;

ou dois je ouvrir un nouveau sujet?
Merci
11 mai 2009 à 10:53:33

Réponse rapide:

Pair<A,B> est une structure contenant deux valeurs de types différents A et B. A et B peuvent aussi être le même type. Dans notre cas, c'est long à écrire, je crée donc un typedef pour simplifier la notation. Le typedef sert à créer un alias sur un nom de type.

La ligne qui te pose problème définit donc un alias nommé "Paire" qui remplacera donc les std::pair<Polynome,Polynome>, autrement dit les couples composés de deux polynomes. Ce que l'on veut pour la division puisqu'il y a le quotient et le reste.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
11 mai 2009 à 11:06:19

non le pb c'est vraiment Pair<A,B>, je suppose que par "structure" tu n'entend pas:
struct poly
{
polynome A;
polynome B;
};
11 mai 2009 à 11:15:56

Non.

template<typename A, typename B>
struct pair
{
   A first;
   B second;
};


Avec dans notre cas, A=B=Polynome. Regarde le premier chapitre sur les templates pour comprendre ou regarde les docs de la STL.
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.
11 mai 2009 à 11:31:29

Donc ecrire :
std::pair<PolynomeA, PolynomeB >

est équivalent a

template<typename A, typename B>
struct pair
{
A first;
B second;
};
Sauf que je peut avoir n'importe qu'elle type pour A et B?
11 mai 2009 à 11:34:08

Non. Lis le tuto sur les templates.

<minicode type="cpp">std::pair<Polynome, Polynome></code>

est équivalent à:

struct pair
{
   Polynome first;
   Polynome second;
};
Co-auteur du cours de C++. ||| Posez vos questions sur le forum ||| Me contacter.