1 files changed, 156 insertions, 0 deletions
diff --git a/tests/root-find-test.cpp b/tests/root-find-test.cpp
new file mode 100644
index 0000000..b866f66
--- /dev/null
+++ b/tests/root-find-test.cpp
@@ -0,0 +1,156 @@
+#include <2geom/polynomial.h>
+#include <vector>
+#include <iterator>
+
+#include <2geom/sbasis.h>
+#include <2geom/sbasis-poly.h>
+#include <2geom/sbasis-to-bezier.h>
+#include <2geom/solver.h>
+#include <time.h>
+
+using namespace std;
+using namespace Geom;
+
+Poly lin_poly(double a, double b) { // ax + b
+    Poly p;
+    p.push_back(b);
+    p.push_back(a);
+    return p;
+}
+
+Linear linear(double ax, double b) {
+    return Linear(b, ax+b);
+}
+
+double uniform() {
+    return double(rand()) / RAND_MAX;
+}
+
+int main() {
+    Poly a, b, r;
+    double timer_precision = 0.01;
+    double units = 1e6; // us
+    
+    a = Poly::linear(1, -0.3)*Poly::linear(1, -0.25)*Poly::linear(1, -0.2);
+    
+    std::cout << a <<std::endl;
+    SBasis B = poly_to_sbasis(a);
+    std::cout << B << std::endl;
+    Bezier bez;
+    sbasis_to_bezier(bez, B);
+    cout << bez << endl;
+    //copy(bez.begin(), bez.end(), ostream_iterator<double>(cout, ", "));
+    cout << endl;
+    cout << endl;
+    cout << endl;
+    
+    std::vector<std::vector<double> > trials;
+    
+    // evenly spaced roots
+    for(int N = 2; N <= 5; N++)
+    {
+        std::vector<double> r;
+        for(int i = 0; i < N; i++)
+            r.push_back(double(i)/(N-1));
+        trials.push_back(r);
+    }
+    // sort of evenish
+    for(int N = 0; N <= 5; N++)
+    {
+        std::vector<double> r;
+        for(int i = 0; i < N; i++)
+            r.push_back(double(i+0.5)/(2*N));
+        trials.push_back(r);
+    }
+    // one at 0.1
+    for(int N = 0; N <= 5; N++)
+    {
+        std::vector<double> r;
+        for(int i = 0; i < N; i++)
+            r.push_back(i+0.1);
+        trials.push_back(r);
+    }
+    for(int N = 0; N <= 6; N++)
+    {
+        std::vector<double> r;
+        for(int i = 0; i < N; i++)
+            r.push_back(i*0.8+0.1);
+        trials.push_back(r);
+    }
+    for(int N = 0; N <= 20; N++)
+    {
+        std::vector<double> r;
+        for(int i = 0; i < N/2; i++) {
+            r.push_back(0.1);
+            r.push_back(0.9);
+        }
+        trials.push_back(r);
+    }
+    for(int i = 0; i <= 20; i++)
+    {
+        std::vector<double> r;
+        for(int i = 0; i < 4; i++) {
+            r.push_back(uniform()*5 - 2.5);
+            r.push_back(0.9);
+        }
+        trials.push_back(r);
+    }
+    double ave_left = 0;
+    cout << "err from exact\n";
+    for(auto & trial : trials) {
+        SBasis B = Linear(1.,1);
+        sort(trial.begin(), trial.end());
+        for(double j : trial) {
+            B = B*linear(1, -j);
+        }
+        double left_time;
+        clock_t end_t = clock()+clock_t(timer_precision*CLOCKS_PER_SEC);
+        unsigned iterations = 0;
+        while(end_t > clock()) {
+            roots(B);
+            iterations++;
+        }
+        left_time = timer_precision*units/iterations;
+        vector<double> rt = roots(B);
+        double err = 0;
+        for(double r : rt) {
+            double best = fabs(r - trial[0]);
+            for(unsigned j = 1; j < trial.size(); j++) {
+                if(fabs(r - trial[j]) < best)
+                    best = fabs(r - trial[j]);
+            }
+            err += best;
+        }
+        if(err > 1e-8){
+            for(double j : trial) {
+                cout << j << ", ";
+            }
+            cout << endl;
+        }
+        cout << " e: " << err << std::endl;
+        ave_left += left_time;
+    }
+    cout << "average time = " << ave_left/trials.size() << std::endl;
+    
+    for(int i = 10; i >= 0; i--) {
+        vector<double> rt = roots(Linear(i,-1));
+        for(double j : rt) {
+            cout << j << ", ";
+        }
+        cout << endl;
+    }
+
+    return 0;
+}
+
+
+/*
+  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 :