- 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.)
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
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
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@
.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)
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;
// 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;
}
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) ||
(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");
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)) ||
(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)
// 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
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;
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;
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;
} 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;
}
#! /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
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) {
}
// 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;
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) {
}
// 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;
--- /dev/null
+/** @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;
+}
--- /dev/null
+/** @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;
+}
++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) {
// 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
(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)
git://www.ginac.de / ginac.git / commitdiff