- Two new timimgs that are interesting for optimizing.
authorRichard Kreckel <Richard.Kreckel@uni-mainz.de>
Mon, 13 Mar 2000 17:14:21 +0000 (17:14 +0000)
committerRichard Kreckel <Richard.Kreckel@uni-mainz.de>
Mon, 13 Mar 2000 17:14:21 +0000 (17:14 +0000)
- Readjusted some parameters.
- Cleaned up syntax in everything involving matrices to reflect the
  policy: "normal if in field, expand otherwise".  (It makes many things
  much clearer.)

13 files changed:
check/Makefile.am
check/Makefile.in
check/exam_matrices.cpp
check/exam_numeric.cpp
check/exams.cpp
check/run_checks
check/time_dennyfliegner.cpp
check/time_gammaseries.cpp
check/time_toeplitz.cpp [new file with mode: 0644]
check/time_vandermonde.cpp [new file with mode: 0644]
check/times.cpp
check/times.h
check/times.ref

index 96ec4089f6b64111f6060fec582850b47742223b..b27d37ef5dabf5f9b132f055f6039585e7d753ba 100644 (file)
@@ -9,8 +9,8 @@ exams_LDADD = ../ginac/libginac.la
 checks_SOURCES =  check_numeric.cpp check_inifcns.cpp check_matrices.cpp \
   check_lsolve.cpp genex.cpp checks.cpp checks.h
 checks_LDADD = ../ginac/libginac.la
-times_SOURCES = time_dennyfliegner.cpp time_gammaseries.cpp timer.cpp \
-  times.cpp times.h
+times_SOURCES = time_dennyfliegner.cpp time_gammaseries.cpp \
+  time_vandermonde.cpp time_toeplitz.cpp timer.cpp times.cpp times.h
 times_LDADD = ../ginac/libginac.la
 INCLUDES = -I$(srcdir)/../ginac
 CLEANFILES = exams.out checks.out times.out
index 0410e062fffccac694b365c09a37c976693082ef..318e21916c8b6f8ca35b60337165cc5e1ccee901 100644 (file)
@@ -111,7 +111,7 @@ exams_LDADD = ../ginac/libginac.la
 checks_SOURCES = check_numeric.cpp check_inifcns.cpp check_matrices.cpp   check_lsolve.cpp genex.cpp checks.cpp checks.h
 
 checks_LDADD = ../ginac/libginac.la
-times_SOURCES = time_dennyfliegner.cpp time_gammaseries.cpp timer.cpp   times.cpp times.h
+times_SOURCES = time_dennyfliegner.cpp time_gammaseries.cpp   time_vandermonde.cpp time_toeplitz.cpp timer.cpp times.cpp times.h
 
 times_LDADD = ../ginac/libginac.la
 INCLUDES = -I$(srcdir)/../ginac
@@ -135,7 +135,8 @@ checks_OBJECTS =  check_numeric.o check_inifcns.o check_matrices.o \
 check_lsolve.o genex.o checks.o
 checks_DEPENDENCIES =  ../ginac/libginac.la
 checks_LDFLAGS = 
-times_OBJECTS =  time_dennyfliegner.o time_gammaseries.o timer.o times.o
+times_OBJECTS =  time_dennyfliegner.o time_gammaseries.o \
+time_vandermonde.o time_toeplitz.o timer.o times.o
 times_DEPENDENCIES =  ../ginac/libginac.la
 times_LDFLAGS = 
 CXXFLAGS = @CXXFLAGS@
@@ -161,7 +162,8 @@ DEP_FILES =  .deps/check_inifcns.P .deps/check_lsolve.P \
 .deps/exam_misc.P .deps/exam_noncommut.P .deps/exam_normalization.P \
 .deps/exam_numeric.P .deps/exam_paranoia.P .deps/exam_polygcd.P \
 .deps/exam_powerlaws.P .deps/exam_pseries.P .deps/exams.P .deps/genex.P \
-.deps/time_dennyfliegner.P .deps/time_gammaseries.P .deps/timer.P \
+.deps/time_dennyfliegner.P .deps/time_gammaseries.P \
+.deps/time_toeplitz.P .deps/time_vandermonde.P .deps/timer.P \
 .deps/times.P
 SOURCES = $(exams_SOURCES) $(checks_SOURCES) $(times_SOURCES)
 OBJECTS = $(exams_OBJECTS) $(checks_OBJECTS) $(times_OBJECTS)
index 81fed94ed2265fd49d8314e10eaa0d586bc6c0c6..a27d586ff97e09d7721ad80d7af2830c115a8afe 100644 (file)
@@ -55,7 +55,7 @@ static unsigned matrix_determinants(void)
     m3.set(0,0,a).set(0,1,b).set(0,2,c);
     m3.set(1,0,d).set(1,1,e).set(1,2,f);
     m3.set(2,0,g).set(2,1,h).set(2,2,i);
-    det = m3.determinant().expand();
+    det = m3.determinant();
     if (det != (a*e*i - a*f*h - d*b*i + d*c*h + g*b*f - g*c*e)) {
         clog << "determinant of 3x3 matrix " << m3
              << " erroneously returned " << det << endl;
@@ -76,10 +76,14 @@ static unsigned matrix_determinants(void)
     // check dense symbolic 2x2 matrix determinant
     m2.set(0,0,a/(a-b)).set(0,1,numeric(1));
     m2.set(1,0,b/(a-b)).set(1,1,numeric(1));
-    det = m2.determinant(true);
+    det = m2.determinant();
     if (det != 1) {
-        clog << "determinant of 2x2 matrix " << m2
-             << " erroneously returned " << det << endl;
+        if (det.normal() == 1)  // only half wrong
+            clog << "determinant of 2x2 matrix " << m2
+                 << " was returned unnormalized as " << det << endl;
+        else  // totally wrong
+            clog << "determinant of 2x2 matrix " << m2
+                 << " erroneously returned " << det << endl;
         ++result;
     }
 
@@ -108,28 +112,31 @@ static unsigned matrix_determinants(void)
 
 static unsigned matrix_invert1(void)
 {
+    unsigned result = 0;
     matrix m(1,1);
     symbol a("a");
-
+    
     m.set(0,0,a);
     matrix m_i = m.inverse();
     
     if (m_i(0,0) != pow(a,-1)) {
         clog << "inversion of 1x1 matrix " << m
              << " erroneously returned " << m_i << endl;
-        return 1;
+        ++result;
     }
-    return 0;
+    
+    return result;
 }
 
 static unsigned matrix_invert2(void)
 {
+    unsigned result = 0;
     matrix m(2,2);
     symbol a("a"), b("b"), c("c"), d("d");
     m.set(0,0,a).set(0,1,b);
     m.set(1,0,c).set(1,1,d);
     matrix m_i = m.inverse();
-    ex det = m.determinant().expand();
+    ex det = m.determinant();
     
     if ((normal(m_i(0,0)*det) != d) ||
         (normal(m_i(0,1)*det) != -b) ||
@@ -137,13 +144,15 @@ static unsigned matrix_invert2(void)
         (normal(m_i(1,1)*det) != a)) {
         clog << "inversion of 2x2 matrix " << m
              << " erroneously returned " << m_i << endl;
-        return 1;
+        ++result;
     }
-    return 0;
+    
+    return result;
 }
 
 static unsigned matrix_invert3(void)
 {
+    unsigned result = 0;
     matrix m(3,3);
     symbol a("a"), b("b"), c("c");
     symbol d("d"), e("e"), f("f");
@@ -152,7 +161,7 @@ static unsigned matrix_invert3(void)
     m.set(1,0,d).set(1,1,e).set(1,2,f);
     m.set(2,0,g).set(2,1,h).set(2,2,i);
     matrix m_i = m.inverse();
-    ex det = m.determinant().normal().expand();
+    ex det = m.determinant();
     
     if ((normal(m_i(0,0)*det) != (e*i-f*h)) ||
         (normal(m_i(0,1)*det) != (c*h-b*i)) ||
@@ -165,9 +174,10 @@ static unsigned matrix_invert3(void)
         (normal(m_i(2,2)*det) != (a*e-b*d))) {
         clog << "inversion of 3x3 matrix " << m
              << " erroneously returned " << m_i << endl;
-        return 1;
+        ++result;
     }
-    return 0;
+    
+    return result;
 }
 
 static unsigned matrix_misc(void)
@@ -206,11 +216,11 @@ static unsigned matrix_misc(void)
     // produce a runtime-error by inverting a singular matrix and catch it
     matrix m4(2,2);
     matrix m5;
-    bool caught=false;
+    bool caught = false;
     try {
         m5 = inverse(m4);
     } catch (std::runtime_error err) {
-        caught=true;
+        caught = true;
     }
     if (!caught) {
         cerr << "singular 2x2 matrix " << m4
@@ -227,7 +237,7 @@ unsigned exam_matrices(void)
     
     cout << "examining symbolic matrix manipulations" << flush;
     clog << "----------symbolic matrix manipulations:" << endl;
-
+    
     result += matrix_determinants();  cout << '.' << flush;
     result += matrix_invert1();  cout << '.' << flush;
     result += matrix_invert2();  cout << '.' << flush;
index d3c1264a3aa1fdde148631f3f9f18b4fcfd19dfa..e37c565d019e0f4faa0a40b46138b1eb0e3be367 100644 (file)
@@ -301,6 +301,25 @@ static unsigned exam_numeric4(void)
     return result;
 }
 
+/* This test examines that simplifications of the form 5^(3/2) -> 5*5^(1/2)
+ * are carried out properly. */
+static unsigned exam_numeric5(void)
+{
+    unsigned result = 0;
+    
+    // A variation of one of Ramanujan's wonderful identities must be
+    // verifiable with very primitive means:
+    ex e1 = pow(1 + pow(3,numeric(1,5)) - pow(3,numeric(2,5)),3);
+    ex e2 = expand(e1 - 10 + 5*pow(3,numeric(3,5)));
+    if (!e2.is_zero()) {
+        clog << "expand((1+3^(1/5)-3^(2/5))^3-10+5*3^(3/5)) returned "
+             << e2 << " instead of 0." << endl;
+        ++result;
+    }
+    
+    return result;
+}
+
 unsigned exam_numeric(void)
 {
     unsigned result = 0;
@@ -312,6 +331,7 @@ unsigned exam_numeric(void)
     result += exam_numeric2();  cout << '.' << flush;
     result += exam_numeric3();  cout << '.' << flush;
     result += exam_numeric4();  cout << '.' << flush;
+    result += exam_numeric5();  cout << '.' << flush;
     
     if (!result) {
         cout << " passed " << endl;
index 31b5cf7f145d36b9dfc849c4ecf9cb6ed3ba5ab4..01453c401cce10dc2674db4dfac91b94e739537b 100644 (file)
@@ -113,7 +113,7 @@ int main()
         } else {
             cout << "(" << result << " individual failures)" << endl;
         }
-        cout << "please check exam.out against exam.ref for more details."
+        cout << "please check exams.out against exams.ref for more details."
              << endl << "happy debugging!" << endl;
     }
     
index 92ba1a2bd6b08881a4c7d0028a842af805e764ac..6304d18f5b66b09086d9a690a36071bdd9058e9e 100755 (executable)
@@ -1,4 +1,4 @@
 #! /bin/sh
-echo "GiNaC will now run through some rather costly consistency checks:"
+echo "GiNaC will now run through some rather costly random consistency checks:"
 ./checks 2>checks.out
 cmp ${srcdir}/checks.ref checks.out
index 23c7a95351088ae9fe746db847624db8b835ff8d..591e7a4f6f4719bf853e3e12c33c0cee272f7786 100644 (file)
@@ -64,15 +64,16 @@ unsigned time_dennyfliegner(void)
     vector<double> times;
     timer rolex;
     
-    sizes.push_back(40);
-    sizes.push_back(60);
+    sizes.push_back(25);
+    sizes.push_back(50);
     sizes.push_back(100);
-    sizes.push_back(150);
+    sizes.push_back(200);
     
     for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i) {
         rolex.start();
-        result += expand_subs(*i);  cout << '.' << flush;
+        result += expand_subs(*i);
         times.push_back(rolex.read());
+        cout << '.' << flush;
     }
     
     if (!result) {
@@ -83,13 +84,11 @@ unsigned time_dennyfliegner(void)
     }
     // print the report:
     cout << endl << "    size:  ";
-    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i) {
+    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i)
         cout << '\t' << (*i);
-    }
     cout << endl << "    time/s:";
-    for (vector<double>::iterator i=times.begin(); i!=times.end(); ++i) {
-        cout << '\t' << (*i);
-    }
+    for (vector<double>::iterator i=times.begin(); i!=times.end(); ++i)
+        cout << '\t' << int(1000*(*i))*0.001;
     cout << endl;
     
     return result;
index afc9f530d5fd51a5fd8078189d9f4581dbec3308..6da460f66994f7432db4aeb5a0579f883551d553 100644 (file)
@@ -61,8 +61,9 @@ unsigned time_gammaseries(void)
     
     for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i) {
         omega.start();
-        result += gammaseries(*i);  cout << '.' << flush;
+        result += gammaseries(*i);
         times.push_back(omega.read());
+        cout << '.' << flush;
     }
     
     if (!result) {
@@ -73,13 +74,11 @@ unsigned time_gammaseries(void)
     }
     // print the report:
     cout << endl << "    order: ";
-    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i) {
+    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i)
         cout << '\t' << (*i);
-    }
     cout << endl << "    time/s:";
-    for (vector<double>::iterator i=times.begin(); i!=times.end(); ++i) {
-        cout << '\t' << (*i);
-    }
+    for (vector<double>::iterator i=times.begin(); i!=times.end(); ++i)
+        cout << '\t' << int(1000*(*i))*0.001;
     cout << endl;
     
     return result;
diff --git a/check/time_toeplitz.cpp b/check/time_toeplitz.cpp
new file mode 100644 (file)
index 0000000..6f138d3
--- /dev/null
@@ -0,0 +1,109 @@
+/** @file time_toeplitz.cpp
+ *
+ *  Calculates determinants of dense symbolic Toeplitz materices.
+ *  For 4x4 our matrix would look like this:
+ *  [[a,b,a+b,a^2+a*b+b^2], [b,a,b,a+b], [a+b,b,a,b], [a^2+a*b+b^2,a+b,b,a]]
+ */
+
+/*
+ *  GiNaC Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany
+ *
+ *  This program is free software; you can redistribute it and/or modify
+ *  it under the terms of the GNU General Public License as published by
+ *  the Free Software Foundation; either version 2 of the License, or
+ *  (at your option) any later version.
+ *
+ *  This program is distributed in the hope that it will be useful,
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ *  GNU General Public License for more details.
+ *
+ *  You should have received a copy of the GNU General Public License
+ *  along with this program; if not, write to the Free Software
+ *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+#include "times.h"
+
+static unsigned toeplitz_det(unsigned size)
+{
+    unsigned result = 0;
+    symbol a("a"), b("b");
+    ex p[8] = {a,
+               b,
+               a+b,
+               pow(a,2) + a*b + pow(b,2),
+               pow(a,3) + pow(a,2)*b - a*pow(b,2) + pow(b,3),
+               pow(a,4) + pow(a,3)*b + pow(a*b,2) + a*pow(b,3) + pow(b,4),
+               pow(a,5) + pow(a,4)*b + pow(a,3)*pow(b,2) - pow(a,2)*pow(b,3) + a*pow(b,4) + pow(b,5),
+               pow(a,6) + pow(a,5)*b + pow(a,4)*pow(b,2) + pow(a*b,3) + pow(a,2)*pow(b,4) + a*pow(b,5) + pow(b,6)
+    };
+    
+    // construct Toeplitz matrix:
+    matrix M(size,size);
+    for (unsigned ro=0; ro<size; ++ro)
+        for (unsigned nd=ro; nd<size; ++nd) {
+            M.set(nd-ro,nd,p[ro]);
+            M.set(nd,nd-ro,p[ro]);
+        }
+    
+    // compute determinant:
+    ex tdet = M.determinant();
+    
+    // dirty consistency check of result:
+    if (!tdet.subs(a==0).subs(b==0).is_zero()) {
+        clog << "Determaint of Toeplitz matrix " << endl
+             << "M==" << M << endl
+             << "was miscalculated: det(M)==" << tdet << endl;
+        ++result;
+    }
+    
+    return result;
+}
+
+unsigned time_toeplitz(void)
+{
+    unsigned result = 0;
+    
+    cout << "timing determinant of polyvariate symbolic Toeplitz matrices" << flush;
+    clog << "-------determinant of polyvariate symbolic Toeplitz matrices:" << endl;
+    
+    vector<unsigned> sizes;
+    vector<double> times;
+    timer longines;
+    
+    sizes.push_back(5);
+    sizes.push_back(6);
+    sizes.push_back(7);
+    sizes.push_back(8);
+    
+    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i) {
+        int count = 1;
+        longines.start();
+        result += toeplitz_det(*i);
+        // correct for very small times:
+        while (longines.read()<0.1) {
+            toeplitz_det(*i);
+            ++count;
+        }
+        times.push_back(longines.read()/count);
+        cout << '.' << flush;
+    }
+    
+    if (!result) {
+        cout << " passed ";
+        clog << "(no output)" << endl;
+    } else {
+        cout << " failed ";
+    }
+    // print the report:
+    cout << endl << "    dim:   ";
+    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i)
+        cout << '\t' << *i << 'x' << *i;
+    cout << endl << "    time/s:";
+    for (vector<double>::iterator i=times.begin(); i!=times.end(); ++i)
+        cout << '\t' << int(1000*(*i))*0.001;
+    cout << endl;
+    
+    return result;
+}
diff --git a/check/time_vandermonde.cpp b/check/time_vandermonde.cpp
new file mode 100644 (file)
index 0000000..364ca5a
--- /dev/null
@@ -0,0 +1,102 @@
+/** @file time_vandermonde.cpp
+ *
+ *  Calculates determinants of dense symbolic Vandermonde materices with
+ *  monomials in one single variable as entries.
+ *  For 4x4 our matrix would look like this:
+ *  [[1,a,a^2,a^3], [1,-a,a^2,-a^3], [1,a^2,a^4,a^6], [1,-a^2,a^4,-a^6]]
+ */
+
+/*
+ *  GiNaC Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany
+ *
+ *  This program is free software; you can redistribute it and/or modify
+ *  it under the terms of the GNU General Public License as published by
+ *  the Free Software Foundation; either version 2 of the License, or
+ *  (at your option) any later version.
+ *
+ *  This program is distributed in the hope that it will be useful,
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ *  GNU General Public License for more details.
+ *
+ *  You should have received a copy of the GNU General Public License
+ *  along with this program; if not, write to the Free Software
+ *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+#include "times.h"
+
+static unsigned vandermonde_det(unsigned size)
+{
+    unsigned result = 0;
+    symbol a("a");
+    
+    // construct Vandermonde matrix:
+    matrix M(size,size);
+    for (unsigned ro=0; ro<size; ++ro)
+        for (unsigned co=0; co<size; ++co)
+            if (ro%2)
+                M.set(ro,co,pow(-pow(a,1+ro/2),co));
+            else
+                M.set(ro,co,pow(pow(a,1+ro/2),co));
+    
+    // compute determinant:
+    ex vdet = M.determinant();
+    
+    // dirty consistency check of result:
+    if (!vdet.subs(a==1).is_zero()) {
+        clog << "Determaint of Vandermonde matrix " << endl
+             << "M==" << M << endl
+             << "was miscalculated: det(M)==" << vdet << endl;
+        ++result;
+    }
+    
+    return result;
+}
+
+unsigned time_vandermonde(void)
+{
+    unsigned result = 0;
+    
+    cout << "timing determinant of univariate symbolic Vandermonde matrices" << flush;
+    clog << "-------determinant of univariate symbolic Vandermonde matrices:" << endl;
+    
+    vector<unsigned> sizes;
+    vector<double> times;
+    timer swatch;
+    
+    sizes.push_back(4);
+    sizes.push_back(6);
+    sizes.push_back(8);
+    sizes.push_back(10);
+    
+    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i) {
+        int count = 1;
+        swatch.start();
+        result += vandermonde_det(*i);
+        // correct for very small times:
+        while (swatch.read()<0.02) {
+            vandermonde_det(*i);
+            ++count;
+        }
+        times.push_back(swatch.read()/count);
+        cout << '.' << flush;
+    }
+    
+    if (!result) {
+        cout << " passed ";
+        clog << "(no output)" << endl;
+    } else {
+        cout << " failed ";
+    }
+    // print the report:
+    cout << endl << "    dim:   ";
+    for (vector<unsigned>::iterator i=sizes.begin(); i!=sizes.end(); ++i)
+        cout << '\t' << *i << 'x' << *i;
+    cout << endl << "    time/s:";
+    for (vector<double>::iterator i=times.begin(); i!=times.end(); ++i)
+        cout << '\t' << int(1000*(*i))*0.001;
+    cout << endl;
+    
+    return result;
+}
index d2b40e81e9fea4dfb2451b7f678ae3d42ba05304..5ba754a07866bb76dc2aa64496a0081f19761b75 100644 (file)
@@ -45,6 +45,22 @@ int main()
         ++result;
     }
     
+    try {
+        for (int i=0; i<1; ++i)
+            result += time_vandermonde();
+    } catch (const exception &e) {
+        cout << "Error: caught exception " << e.what() << endl;
+        ++result;
+    }
+    
+    try {
+        for (int i=0; i<1; ++i)
+            result += time_toeplitz();
+    } catch (const exception &e) {
+        cout << "Error: caught exception " << e.what() << endl;
+        ++result;
+    }
+        
     if (result) {
         cout << "Error: something went wrong. ";
         if (result == 1) {
index bb54399e61ef4f968c18f1880f684f5ed8ca9cd9..dec00ace8342e30f1b5dc2c97bcd93501cc35b55 100644 (file)
@@ -49,5 +49,7 @@ private:
 // prototypes for all individual timings should be unsigned fcn():
 unsigned time_dennyfliegner();
 unsigned time_gammaseries();
+unsigned time_vandermonde();
+unsigned time_toeplitz();
 
 #endif // ndef CHECKS_H
index 4aef1c991f58ec86f0c6f07e1b5f6f9ebe618c9a..abadcbe419d67d68f45c179c25bf76091c4ae10b 100644 (file)
@@ -2,3 +2,7 @@
 (no output)
 -------Laurent series expansion of gamma function:
 (no output)
+-------determinant of univariate symbolic Vandermonde matrices:
+(no output)
+-------determinant of polyvariate symbolic Toeplitz matrices:
+(no output)