1 files changed, 569 insertions, 0 deletions
diff --git a/include/2geom/symbolic/polynomial.h b/include/2geom/symbolic/polynomial.h
new file mode 100644
index 0000000..fea7e6c
--- /dev/null
+++ b/include/2geom/symbolic/polynomial.h
@@ -0,0 +1,569 @@
+/*
+ * Polynomial<CoeffT> class template
+ *
+ * Authors:
+ *      Marco Cecchetti <mrcekets at gmail.com>
+ *
+ * Copyright 2008  authors
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it either under the terms of the GNU Lesser General Public
+ * License version 2.1 as published by the Free Software Foundation
+ * (the "LGPL") or, at your option, under the terms of the Mozilla
+ * Public License Version 1.1 (the "MPL"). If you do not alter this
+ * notice, a recipient may use your version of this file under either
+ * the MPL or the LGPL.
+ *
+ * You should have received a copy of the LGPL along with this library
+ * in the file COPYING-LGPL-2.1; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+ * You should have received a copy of the MPL along with this library
+ * in the file COPYING-MPL-1.1
+ *
+ * The contents of this file are subject to the Mozilla Public License
+ * Version 1.1 (the "License"); you may not use this file except in
+ * compliance with the License. You may obtain a copy of the License at
+ * http://www.mozilla.org/MPL/
+ *
+ * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY
+ * OF ANY KIND, either express or implied. See the LGPL or the MPL for
+ * the specific language governing rights and limitations.
+ */
+
+#ifndef _GEOM_SL_POLYNOMIAL_H_
+#define _GEOM_SL_POLYNOMIAL_H_
+
+
+#include <2geom/symbolic/unity-builder.h>
+
+#include <vector>
+#include <string>
+
+#include <2geom/exception.h>
+
+
+
+
+namespace Geom { namespace SL {
+
+/*
+ * Polynomial<CoeffT> class template
+ *
+ * It represents a generic univariate polynomial with coefficients
+ * of type CoeffT. One way to get a multi-variate polynomial is
+ * to utilize a Polynomial instantiation as coefficient type
+ * in a recursive style.
+ *
+ */
+
+template< typename CoeffT >
+class Polynomial
+{
+  public:
+    typedef CoeffT coeff_type;
+    typedef std::vector<coeff_type> coeff_container_t;
+    typedef typename coeff_container_t::iterator iterator;
+    typedef typename coeff_container_t::const_iterator const_iterator;
+
+    /*
+     * a Polynomial should be never empty
+     */
+    Polynomial()
+    {
+        m_coeff.push_back(zero_coeff);
+    }
+
+    explicit
+    Polynomial(CoeffT const& c, size_t i = 0)
+    {
+        m_coeff.resize(i, zero_coeff);
+        m_coeff.push_back(c);
+    }
+
+    /*
+     *  forwarding of some std::vector methods
+     */
+
+    size_t size() const
+    {
+        return m_coeff.size();
+    }
+
+    const_iterator begin() const
+    {
+        return m_coeff.begin();
+    }
+
+    const_iterator end() const
+    {
+        return m_coeff.end();
+    }
+
+    iterator begin()
+    {
+        return m_coeff.begin();
+    }
+
+    iterator end()
+    {
+        return m_coeff.end();
+    }
+
+    void reserve(size_t n)
+    {
+        m_coeff.reserve(n);
+    }
+
+    size_t capacity() const
+    {
+        return m_coeff.capacity();
+    }
+
+    /*
+     *  degree of the term with the highest degree
+     *  and an initialized coefficient (even if zero)
+     */
+    size_t max_degree() const
+    {
+        if (size() == 0)
+            THROW_INVARIANTSVIOLATION (0);
+
+        return (size() - 1);
+    }
+
+    void max_degree(size_t n)
+    {
+        m_coeff.resize(n+1, zero_coeff);
+    }
+
+    /*
+     *  degree of the term with the highest degree
+     *  and an initialized coefficient that is not null
+     */
+    size_t real_degree() const
+    {
+        if (size() == 0)
+            THROW_INVARIANTSVIOLATION (0);
+
+        const_iterator it = end() - 1;
+        for (; it != begin(); --it)
+        {
+            if (*it != zero_coeff) break;
+        }
+        size_t i = static_cast<size_t>(it - begin());
+        return i;
+    }
+
+    bool is_zero() const
+    {
+        if (size() == 0)
+            THROW_INVARIANTSVIOLATION (0);
+
+        if (real_degree() != 0) return false;
+        if (m_coeff[0] != zero_coeff) return false;
+        return true;
+    }
+
+    /*
+     * trim leading zero coefficients
+     * after calling normalize max_degree == real_degree
+     */
+    void normalize()
+    {
+        size_t rd = real_degree();
+        if (rd != max_degree())
+        {
+            m_coeff.erase(begin() + rd + 1, end());
+        }
+    }
+
+    coeff_type const& operator[] (size_t i) const
+    {
+        return m_coeff[i];
+    }
+
+    coeff_type & operator[] (size_t i)
+    {
+        return m_coeff[i];
+    }
+
+    // safe coefficient getter routine
+    coeff_type const& coefficient(size_t i) const
+    {
+        if (i > max_degree())
+        {
+            return zero_coeff;
+        }
+        else
+        {
+            return m_coeff[i];
+        }
+    }
+
+    // safe coefficient setter routine
+    void coefficient(size_t i, coeff_type const& c)
+    {
+        //std::cerr << "i: " << i << " c: " << c << std::endl;
+        if (i > max_degree())
+        {
+            if (c == zero_coeff) return;
+            reserve(i+1);
+            m_coeff.resize(i, zero_coeff);
+            m_coeff.push_back(c);
+        }
+        else
+        {
+            m_coeff[i] = c;
+        }
+    }
+
+    coeff_type const& leading_coefficient() const
+    {
+        return m_coeff[real_degree()];
+    }
+
+    coeff_type & leading_coefficient()
+    {
+        return m_coeff[real_degree()];
+    }
+
+    /*
+     * polynomail evaluation:
+     * T can be any type that is able to be + and * with the coefficient type
+     */
+    template <typename T>
+    T operator() (T const& x) const
+    {
+        T r = zero<T>()();
+        for(size_t i = max_degree(); i > 0; --i)
+        {
+            r += (*this)[i];
+            r *= x;
+        }
+        r += (*this)[0];
+        return r;
+    }
+
+    // opposite polynomial
+    Polynomial operator-() const
+    {
+        Polynomial r;
+        // we need r.m_coeff to be empty so we can utilize push_back
+        r.m_coeff.pop_back();
+        r.reserve(size());
+        for(size_t i = 0; i < size(); ++i)
+        {
+            r.m_coeff.push_back( -(*this)[i] );
+        }
+        return r;
+    }
+
+    /*
+     *  polynomial-polynomial mutating operators
+     */
+
+    Polynomial& operator+=(Polynomial const& p)
+    {
+        size_t sz = std::min(size(), p.size());
+        for (size_t i = 0; i < sz; ++i)
+        {
+            (*this)[i] += p[i];
+        }
+        if (size() < p.size())
+        {
+            m_coeff.insert(end(), p.begin() + size(), p.end());
+        }
+        return (*this);
+    }
+
+    Polynomial& operator-=(Polynomial const& p)
+    {
+        size_t sz = std::min(size(), p.size());
+        for (size_t i = 0; i < sz; ++i)
+        {
+            (*this)[i] -= p[i];
+        }
+        reserve(p.size());
+        for(size_t i = sz; i < p.size(); ++i)
+        {
+            m_coeff.push_back( -p[i] );
+        }
+        return (*this);
+    }
+
+    Polynomial& operator*=(Polynomial const& p)
+    {
+        Polynomial r;
+        r.m_coeff.resize(size() + p.size() - 1, zero_coeff);
+
+        for (size_t i = 0; i < size(); ++i)
+        {
+            for (size_t j = 0; j < p.size(); ++j)
+            {
+                r[i+j] += (*this)[i] * p[j];
+            }
+        }
+        (*this) = r;
+        return (*this);
+    }
+
+    /*
+     *  equivalent to multiply by x^n
+     */
+    Polynomial& operator<<=(size_t n)
+    {
+        m_coeff.insert(begin(), n, zero_coeff);
+        return (*this);
+    }
+
+    /*
+     *  polynomial-coefficient mutating operators
+     */
+
+    Polynomial& operator=(coeff_type const& c)
+    {
+        m_coeff[0] = c;
+        return (*this);
+    }
+
+    Polynomial& operator+=(coeff_type const& c)
+    {
+        (*this)[0] += c;
+        return (*this);
+    }
+
+    Polynomial& operator-=(coeff_type const& c)
+    {
+        (*this)[0] -= c;
+        return (*this);
+    }
+
+    Polynomial& operator*=(coeff_type const& c)
+    {
+        for (size_t i = 0; i < size(); ++i)
+        {
+            (*this)[i] *= c;
+        }
+        return (*this);
+    }
+
+    // return the poly in a string form
+    std::string str() const;
+
+  private:
+    // with zero_coeff defined as a static data member
+    // coefficient(size_t i) safe get method can always
+    // return a (const) reference
+    static const coeff_type zero_coeff;
+    coeff_container_t m_coeff;
+
+}; // end class Polynomial
+
+
+/*
+ *  zero and one element spezcialization for Polynomial
+ */
+
+template< typename CoeffT >
+struct zero<Polynomial<CoeffT>, false>
+{
+    Polynomial<CoeffT> operator() () const
+    {
+        CoeffT zc = zero<CoeffT>()();
+        Polynomial<CoeffT> z(zc);
+        return z;
+    }
+};
+
+template< typename CoeffT >
+struct one<Polynomial<CoeffT>, false>
+{
+    Polynomial<CoeffT> operator() ()
+    {
+        CoeffT _1c = one<CoeffT>()();
+        Polynomial<CoeffT> _1(_1c);
+        return _1;
+    }
+};
+
+
+/*
+ * initialization of Polynomial static data members
+ */
+
+template< typename CoeffT >
+const typename Polynomial<CoeffT>::coeff_type Polynomial<CoeffT>::zero_coeff
+    = zero<typename Polynomial<CoeffT>::coeff_type>()();
+
+/*
+ * Polynomial - Polynomial binary mathematical operators
+ */
+
+template< typename CoeffT >
+inline
+bool operator==(Polynomial<CoeffT> const& p, Polynomial<CoeffT> const& q)
+{
+    size_t d = p.real_degree();
+    if (d != q.real_degree()) return false;
+    for (size_t i = 0; i <= d; ++i)
+    {
+        if (p[i] != q[i]) return false;
+    }
+    return true;
+}
+
+template< typename CoeffT >
+inline
+bool operator!=(Polynomial<CoeffT> const& p, Polynomial<CoeffT> const& q)
+{
+    return !(p == q);
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator+( Polynomial<CoeffT> const& p, Polynomial<CoeffT> const& q )
+{
+    Polynomial<CoeffT> r(p);
+    r += q;
+    return r;
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator-( Polynomial<CoeffT> const& p, Polynomial<CoeffT> const& q )
+{
+    Polynomial<CoeffT> r(p);
+    r -= q;
+    return r;
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator*( Polynomial<CoeffT> const& p, Polynomial<CoeffT> const& q )
+{
+    Polynomial<CoeffT> r(p);
+    r *= q;
+    return r;
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT> operator<<(Polynomial<CoeffT> const& p, size_t n)
+{
+    Polynomial<CoeffT> r(p);
+    r <<= n;
+    return r;
+}
+
+
+/*
+ *  polynomial-coefficient and coefficient-polynomial mathematical operators
+ */
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator+( Polynomial<CoeffT> const& p, CoeffT const& c )
+{
+    Polynomial<CoeffT> r(p);
+    r += c;
+    return r;
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator+( CoeffT const& c, Polynomial<CoeffT> const& p)
+{
+    return (p + c);
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator-( Polynomial<CoeffT> const& p, CoeffT const& c )
+{
+    Polynomial<CoeffT> r(p);
+    r -= c;
+    return r;
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator-( CoeffT const& c, Polynomial<CoeffT> const& p)
+{
+    return (p - c);
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator*( Polynomial<CoeffT> const& p, CoeffT const& c )
+{
+    Polynomial<CoeffT> r(p);
+    r *= c;
+    return r;
+}
+
+template< typename CoeffT >
+inline
+Polynomial<CoeffT>
+operator*( CoeffT const& c, Polynomial<CoeffT> const& p)
+{
+    return (p * c);
+}
+
+
+/*
+ *  operator<< extension for printing Polynomial
+ *  and str() method for transforming a Polynomial into a string
+ */
+
+template< typename charT, typename CoeffT >
+inline
+std::basic_ostream<charT> &
+operator<< (std::basic_ostream<charT> & os, const Polynomial<CoeffT> & p)
+{
+    if (p.size() == 0) return os;
+    os << "{" << p[0];
+    for (size_t i = 1; i < p.size(); ++i)
+    {
+        os << ", " << p[i];
+    }
+    os << "}";
+    return os;
+}
+
+
+template< typename CoeffT >
+inline
+std::string Polynomial<CoeffT>::str() const
+{
+    std::ostringstream oss;
+    oss << (*this);
+    return oss.str();
+}
+
+
+} /*end namespace Geom*/  } /*end namespace SL*/
+
+
+
+
+#endif // _GEOM_SL_POLYNOMIAL_H_
+
+
+/*
+  Local Variables:
+  mode:c++
+  c-file-style:"stroustrup"
+  c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+  indent-tabs-mode:nil
+  fill-column:99
+  End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :