Branch 1.3.5 tree to 1.4.x tree, goodbye 1.3.x

git-svn-id: svn://scribus.net/branches/Version14x/Scribus@17163 11d20701-8431-0410-a711-e3c959e3b870

1 files changed, 197 insertions, 0 deletions

diff --git a/scribus/plugins/tools/2geomtools/lib2geom/poly.cpp b/scribus/plugins/tools/2geomtools/lib2geom/poly.cpp
new file mode 100644
index 0000000..17f2864
--- /dev/null
+++ b/scribus/plugins/tools/2geomtools/lib2geom/poly.cpp
@@ -0,0 +1,197 @@
+#include "poly.h"
+
+Poly Poly::operator*(const Poly& p) const {
+    Poly result; 
+    result.resize(degree() +  p.degree()+1);
+    
+    for(unsigned i = 0; i < size(); i++) {
+        for(unsigned j = 0; j < p.size(); j++) {
+            result[i+j] += (*this)[i] * p[j];
+        }
+    }
+    return result;
+}
+#ifdef HAVE_GSL
+#include <gsl/gsl_poly.h>
+#endif
+
+/*double Poly::eval(double x) const {
+    return gsl_poly_eval(&coeff[0], size(), x);
+    }*/
+
+void Poly::normalize() {
+    while(back() == 0)
+        pop_back();
+}
+
+void Poly::monicify() {
+    normalize();
+    
+    double scale = 1./back(); // unitize
+    
+    for(unsigned i = 0; i < size(); i++) {
+        (*this)[i] *= scale;
+    }
+}
+
+
+#ifdef HAVE_GSL
+std::vector<std::complex<double> > solve(Poly const & pp) {
+    Poly p(pp);
+    p.normalize();
+    gsl_poly_complex_workspace * w 
+        = gsl_poly_complex_workspace_alloc (p.size());
+       
+    gsl_complex_packed_ptr z = new double[p.degree()*2];
+    double* a = new double[p.size()];
+    for(int i = 0; i < p.size(); i++)
+        a[i] = p[i];
+    std::vector<std::complex<double> > roots;
+    //roots.resize(p.degree());
+    
+    gsl_poly_complex_solve (a, p.size(), w, z);
+    delete[]a;
+     
+    gsl_poly_complex_workspace_free (w);
+     
+    for (int i = 0; i < p.degree(); i++) {
+        roots.push_back(std::complex<double> (z[2*i] ,z[2*i+1]));
+        //printf ("z%d = %+.18f %+.18f\n", i, z[2*i], z[2*i+1]);
+    }    
+    delete[] z;
+    return roots;
+}
+
+std::vector<double > solve_reals(Poly const & p) {
+    std::vector<std::complex<double> > roots = solve(p);
+    std::vector<double> real_roots;
+    
+    for(int i = 0; i < roots.size(); i++) {
+        if(roots[i].imag() == 0) // should be more lenient perhaps
+            real_roots.push_back(roots[i].real());
+    }
+    return real_roots;
+}
+#endif
+
+double polish_root(Poly const & p, double guess, double tol) {
+    Poly dp = derivative(p);
+    
+    double fn = p(guess);
+    while(fabs(fn) > tol) {
+        guess -= fn/dp(guess);
+        fn = p(guess);
+    }
+    return guess;
+}
+
+Poly integral(Poly const & p) {
+    Poly result;
+    
+    result.reserve(p.size()+1);
+    result.push_back(0); // arbitrary const
+    for(unsigned i = 0; i < p.size(); i++) {
+        result.push_back(p[i]/(i+1));
+    }
+    return result;
+
+}
+
+Poly derivative(Poly const & p) {
+    Poly result;
+    
+    if(p.size() <= 1)
+        return Poly(0);
+    result.reserve(p.size()-1);
+    for(unsigned i = 1; i < p.size(); i++) {
+        result.push_back(i*p[i]);
+    }
+    return result;
+}
+
+Poly compose(Poly const & a, Poly const & b) {
+    Poly result;
+    
+    for(unsigned i = a.size(); i > 0; i--) {
+        result = Poly(a[i-1]) + result * b;
+    }
+    return result;
+    
+}
+
+/* This version is backwards - dividing taylor terms
+Poly divide(Poly const &a, Poly const &b, Poly &r) {
+    Poly c;
+    r = a; // remainder
+    
+    const unsigned k = a.size();
+    r.resize(k, 0);
+    c.resize(k, 0);
+
+    for(unsigned i = 0; i < k; i++) {
+        double ci = r[i]/b[0];
+        c[i] += ci;
+        Poly bb = ci*b;
+        std::cout << ci <<"*" << b << ", r= " << r << std::endl;
+        r -= bb.shifted(i);
+    }
+    
+    return c;
+}
+*/
+
+Poly divide(Poly const &a, Poly const &b, Poly &r) {
+    Poly c;
+    r = a; // remainder
+    assert(b.size() > 0);
+    
+    const unsigned k = a.degree();
+    const unsigned l = b.degree();
+    c.resize(k, 0.);
+    
+    for(unsigned i = k; i >= l; i--) {
+        assert(i >= 0);
+        double ci = r.back()/b.back();
+        c[i-l] += ci;
+        Poly bb = ci*b;
+        //std::cout << ci <<"*(" << b.shifted(i-l) << ") = " 
+        //          << bb.shifted(i-l) << "     r= " << r << std::endl;
+        r -= bb.shifted(i-l);
+        r.pop_back();
+    }
+    //std::cout << "r= " << r << std::endl;
+    r.normalize();
+    c.normalize();
+    
+    return c;
+}
+
+Poly gcd(Poly const &a, Poly const &b, const double tol) {
+    if(a.size() < b.size())
+        return gcd(b, a);
+    if(b.size() <= 0)
+        return a;
+    if(b.size() == 1)
+        return a;
+    Poly r;
+    divide(a, b, r);
+    return gcd(b, r);
+}
+
+
+
+/*Poly divide_out_root(Poly const & p, double x) {
+    assert(1);
+    }*/
+
+
+/*
+  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:encoding=utf-8:textwidth=99 :