From: Richard Kreckel Date: Mon, 13 Dec 1999 18:43:11 +0000 (+0000) Subject: - introduced info_flags::cinteger, info_flags::crational, X-Git-Tag: release_0-5-0~95 X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=commitdiff_plain;h=26741891dadf23162799009b6fd57b4984bd4ce5 - introduced info_flags::cinteger, info_flags::crational, info_flags::cinteger_polynomial, info_flags::crational_polynomial with intuitive behaviour. - extended documentation - made things like evlaf(zeta(3)) work - killed several bugs --- diff --git a/Makefile.in b/Makefile.in index 561c21bd..a6e7630f 100644 --- a/Makefile.in +++ b/Makefile.in @@ -352,7 +352,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \ diff --git a/README b/README index 5a97502c..880be62b 100644 --- a/README +++ b/README @@ -1,8 +1,8 @@ General Information =================== -GiNaC (which stands for "GiNaC is Not a CAS (computer algebra system)) is a -C++ library for symbolic mathematical calculations. It is designed to allow +GiNaC (which stands for "GiNaC is Not a CAS" (computer algebra system)) is a +C++ library for symbolic mathematical calculations. It is designed to allow the creation of integrated systems that embed symbolic manipulations together with more established areas of computer science (like computation-intense numeric applications, graphical interfaces, etc.) under one roof. @@ -36,7 +36,7 @@ How to report bugs ================== If you have identified a bug in GiNaC you are welcome to send a detailed -bug report to . Please think about your bug! This +bug report to . Please think about your bug! This means that you should include * Information about your system @@ -58,5 +58,5 @@ means that you should include together with the output you get and the output you expect will help us to reproduce it quickly. -Patches are most welcome. If possible please make them with diff -c and +Patches are most welcome. If possible please make them with diff -c and include ChangeLog entries. diff --git a/check/Makefile.in b/check/Makefile.in index f9e5c2cd..9189758a 100644 --- a/check/Makefile.in +++ b/check/Makefile.in @@ -245,7 +245,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \ diff --git a/check/inifcns_consist.cpp b/check/inifcns_consist.cpp index 25e34197..f48db0ab 100644 --- a/check/inifcns_consist.cpp +++ b/check/inifcns_consist.cpp @@ -277,7 +277,7 @@ unsigned inifcns_consist(void) result += inifcns_consist_psi(); result += inifcns_consist_zeta(); - if ( !result ) { + if (!result) { cout << " passed "; clog << "(no output)" << endl; } else { diff --git a/check/matrix_checks.cpp b/check/matrix_checks.cpp index 7bce7872..4443654d 100644 --- a/check/matrix_checks.cpp +++ b/check/matrix_checks.cpp @@ -126,10 +126,10 @@ static unsigned matrix_invert2(void) matrix m_i = m.inverse(); ex det = m.determinant().expand(); - if ( (normal(m_i(0,0)*det) != d) || - (normal(m_i(0,1)*det) != -b) || - (normal(m_i(1,0)*det) != -c) || - (normal(m_i(1,1)*det) != a) ) { + if ((normal(m_i(0,0)*det) != d) || + (normal(m_i(0,1)*det) != -b) || + (normal(m_i(1,0)*det) != -c) || + (normal(m_i(1,1)*det) != a)) { clog << "inversion of 2x2 matrix " << m << " erroneously returned " << m_i << endl; return 1; @@ -149,15 +149,15 @@ static unsigned matrix_invert3(void) matrix m_i = m.inverse(); ex det = m.determinant().normal().expand(); - if ( (normal(m_i(0,0)*det) != (e*i-f*h)) || - (normal(m_i(0,1)*det) != (c*h-b*i)) || - (normal(m_i(0,2)*det) != (b*f-c*e)) || - (normal(m_i(1,0)*det) != (f*g-d*i)) || - (normal(m_i(1,1)*det) != (a*i-c*g)) || - (normal(m_i(1,2)*det) != (c*d-a*f)) || - (normal(m_i(2,0)*det) != (d*h-e*g)) || - (normal(m_i(2,1)*det) != (b*g-a*h)) || - (normal(m_i(2,2)*det) != (a*e-b*d)) ) { + if ((normal(m_i(0,0)*det) != (e*i-f*h)) || + (normal(m_i(0,1)*det) != (c*h-b*i)) || + (normal(m_i(0,2)*det) != (b*f-c*e)) || + (normal(m_i(1,0)*det) != (f*g-d*i)) || + (normal(m_i(1,1)*det) != (a*i-c*g)) || + (normal(m_i(1,2)*det) != (c*d-a*f)) || + (normal(m_i(2,0)*det) != (d*h-e*g)) || + (normal(m_i(2,1)*det) != (b*g-a*h)) || + (normal(m_i(2,2)*det) != (a*e-b*d))) { clog << "inversion of 3x3 matrix " << m << " erroneously returned " << m_i << endl; return 1; diff --git a/check/numeric_consist.cpp b/check/numeric_consist.cpp index 60088493..d7edb3b6 100644 --- a/check/numeric_consist.cpp +++ b/check/numeric_consist.cpp @@ -297,7 +297,7 @@ static unsigned numeric_consist4(void) // square roots of squares of integers: passed = true; for (int i=0; i<42; ++i) { - if ( !sqrt(numeric(i*i)).is_integer() ) { + if (!sqrt(numeric(i*i)).is_integer()) { passed = false; } } @@ -310,7 +310,7 @@ static unsigned numeric_consist4(void) passed = true; for (int num=0; num<41; ++num) { for (int den=1; den<42; ++den) { - if ( !sqrt(numeric(num*num)/numeric(den*den)).is_rational() ) { + if (!sqrt(numeric(num*num)/numeric(den*den)).is_rational()) { passed = false; } } diff --git a/check/numeric_output.cpp b/check/numeric_output.cpp index 50dd8c1f..b98d61f0 100644 --- a/check/numeric_output.cpp +++ b/check/numeric_output.cpp @@ -34,7 +34,7 @@ unsigned numeric_output(void) clog << "---------output of numeric types:" << endl; unsigned long Digits_before = Digits; - Digits = 200; + Digits = 222; clog << "Using " << Digits << " digits" << endl; clog << Pi << " evalfs to: " << Pi.evalf() << endl; clog << Catalan << " evalfs to: " << Catalan.evalf() << endl; diff --git a/check/paranoia_check.cpp b/check/paranoia_check.cpp index 3eb4185c..2e0c40df 100644 --- a/check/paranoia_check.cpp +++ b/check/paranoia_check.cpp @@ -1,9 +1,9 @@ /** @file paranoia_check.cpp * * This set of tests checks for some of GiNaC's oopses which showed up during - * development. Things were evaluated wrongly and so. It should not find such - * a sick behaviour again. But since we are paranoic and we want to exclude - * that behaviour for good... */ + * development. Things were evaluated wrongly and so. Such a sick behaviour + * shouldn't occur any more. But we are paranoic and we want to exclude these + * these oopses for good, so we run those stupid tests... */ /* * GiNaC Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany @@ -43,7 +43,7 @@ static unsigned paranoia_check1(void) g = e / f; // In the first one expand did not do any job at all: - if ( !g.expand().is_equal(x) ) { + if (!g.expand().is_equal(x)) { clog << "e = x*y*z; f = y*z; expand(e/f) erroneously returned " << g.expand() << endl; ++result; diff --git a/check/result.ref b/check/result.ref index a3280584..0e3a9296 100644 --- a/check/result.ref +++ b/check/result.ref @@ -1,10 +1,10 @@ ---------several ex-bugs just out of pure paranoia: (no output) ---------output of numeric types: -Using 200 digits -Pi evalfs to: 3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679821480865132823066470938446095505822317253594081284811174502841027019385211055596446229489549303819644L0 -Catalan evalfs to: 0.9159655941772190150546035149323841107741493742816721342664981196217630197762547694793565129261151062485744226191961995790358988033258590594315947374811584069953320287733194605190387274781640878659090247L0 -EulerGamma evalfs to: 0.577215664901532860606512090082402431042159335939923598805767234884867726777664670936947063291746749514631447249807082480960504014486542836224173997644923536253500333742937337737673942792595258247094916L0 +Using 222 digits +Pi evalfs to: 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648235L0 +Catalan evalfs to: 0.915965594177219015054603514932384110774149374281672134266498119621763019776254769479356512926115106248574422619196199579035898803325859059431594737481158406995332028773319460519038727478164087865909024706484152163000228727640942388L0 +EulerGamma evalfs to: 0.577215664901532860606512090082402431042159335939923598805767234884867726777664670936947063291746749514631447249807082480960504014486542836224173997644923536253500333742937337737673942792595258247094916008735203948165670853233151777L0 Complex integers: {(0,0)=0} {(1,0)=1} {(1,1)=1+I} {(0,1)=I} {(-1,1)=-1+I} {(-1,0)=-1} {(-1,-1)=-1-I} {(0,-1)=-I} {(1,-1)=1-I} ---------consistency of numeric types: (no output) diff --git a/check/series_expansion.cpp b/check/series_expansion.cpp index b648f9d7..b67705d7 100644 --- a/check/series_expansion.cpp +++ b/check/series_expansion.cpp @@ -28,125 +28,172 @@ using namespace GiNaC; static symbol x("x"); -static unsigned check_series(const ex &e, const ex &point, const ex &d) +static unsigned check_series(const ex &e, const ex &point, const ex &d, int order = 8) { - ex es = e.series(x, point, 8); - ex ep = static_cast(es.bp)->convert_to_poly(); - if ((ep - d).compare(exZERO()) != 0) { - clog << "series expansion of " << e << " at " << point + ex es = e.series(x, point, order); + ex ep = static_cast(es.bp)->convert_to_poly(); + if ((ep - d).compare(exZERO()) != 0) { + clog << "series expansion of " << e << " at " << point << " erroneously returned " << ep << " (instead of " << d << ")" << endl; - (ep-d).printtree(clog); - return 1; - } - return 0; + (ep-d).printtree(clog); + return 1; + } + return 0; } // Series expansion static unsigned series1(void) { - unsigned result = 0; - ex e, d; - - e = sin(x); - d = x - pow(x, 3) / 6 + pow(x, 5) / 120 - pow(x, 7) / 5040 + Order(pow(x, 8)); - result += check_series(e, exZERO(), d); - - e = cos(x); - d = 1 - pow(x, 2) / 2 + pow(x, 4) / 24 - pow(x, 6) / 720 + Order(pow(x, 8)); - result += check_series(e, exZERO(), d); - - e = exp(x); - d = 1 + x + pow(x, 2) / 2 + pow(x, 3) / 6 + pow(x, 4) / 24 + pow(x, 5) / 120 + pow(x, 6) / 720 + pow(x, 7) / 5040 + Order(pow(x, 8)); - result += check_series(e, exZERO(), d); - - e = pow(1 - x, -1); - d = 1 + x + pow(x, 2) + pow(x, 3) + pow(x, 4) + pow(x, 5) + pow(x, 6) + pow(x, 7) + Order(pow(x, 8)); - result += check_series(e, exZERO(), d); - - e = x + pow(x, -1); - d = x + pow(x, -1); - result += check_series(e, exZERO(), d); - - e = x + pow(x, -1); - d = 2 + pow(x-1, 2) - pow(x-1, 3) + pow(x-1, 4) - pow(x-1, 5) + pow(x-1, 6) - pow(x-1, 7) + Order(pow(x-1, 8)); - result += check_series(e, exONE(), d); - - e = pow(x + pow(x, 3), -1); - d = pow(x, -1) - x + pow(x, 3) - pow(x, 5) + Order(pow(x, 7)); - result += check_series(e, exZERO(), d); - - e = pow(pow(x, 2) + pow(x, 4), -1); - d = pow(x, -2) - 1 + pow(x, 2) - pow(x, 4) + Order(pow(x, 6)); - result += check_series(e, exZERO(), d); - - e = pow(sin(x), -2); - d = pow(x, -2) + numeric(1,3) + pow(x, 2) / 15 + pow(x, 4) * 2/189 + Order(pow(x, 5)); - result += check_series(e, exZERO(), d); - - e = sin(x) / cos(x); - d = x + pow(x, 3) / 3 + pow(x, 5) * 2/15 + pow(x, 7) * 17/315 + Order(pow(x, 8)); - result += check_series(e, exZERO(), d); - - e = cos(x) / sin(x); - d = pow(x, -1) - x / 3 - pow(x, 3) / 45 - pow(x, 5) * 2/945 + Order(pow(x, 6)); - result += check_series(e, exZERO(), d); - - e = pow(numeric(2), x); - ex t = log(ex(2)) * x; - d = 1 + t + pow(t, 2) / 2 + pow(t, 3) / 6 + pow(t, 4) / 24 + pow(t, 5) / 120 + pow(t, 6) / 720 + pow(t, 7) / 5040 + Order(pow(x, 8)); - result += check_series(e, exZERO(), d.expand()); - - e = pow(Pi, x); - t = log(Pi) * x; - d = 1 + t + pow(t, 2) / 2 + pow(t, 3) / 6 + pow(t, 4) / 24 + pow(t, 5) / 120 + pow(t, 6) / 720 + pow(t, 7) / 5040 + Order(pow(x, 8)); - result += check_series(e, exZERO(), d.expand()); - - return result; + unsigned result = 0; + ex e, d; + + e = sin(x); + d = x - pow(x, 3) / 6 + pow(x, 5) / 120 - pow(x, 7) / 5040 + Order(pow(x, 8)); + result += check_series(e, exZERO(), d); + + e = cos(x); + d = 1 - pow(x, 2) / 2 + pow(x, 4) / 24 - pow(x, 6) / 720 + Order(pow(x, 8)); + result += check_series(e, exZERO(), d); + + e = exp(x); + d = 1 + x + pow(x, 2) / 2 + pow(x, 3) / 6 + pow(x, 4) / 24 + pow(x, 5) / 120 + pow(x, 6) / 720 + pow(x, 7) / 5040 + Order(pow(x, 8)); + result += check_series(e, exZERO(), d); + + e = pow(1 - x, -1); + d = 1 + x + pow(x, 2) + pow(x, 3) + pow(x, 4) + pow(x, 5) + pow(x, 6) + pow(x, 7) + Order(pow(x, 8)); + result += check_series(e, exZERO(), d); + + e = x + pow(x, -1); + d = x + pow(x, -1); + result += check_series(e, exZERO(), d); + + e = x + pow(x, -1); + d = 2 + pow(x-1, 2) - pow(x-1, 3) + pow(x-1, 4) - pow(x-1, 5) + pow(x-1, 6) - pow(x-1, 7) + Order(pow(x-1, 8)); + result += check_series(e, exONE(), d); + + e = pow(x + pow(x, 3), -1); + d = pow(x, -1) - x + pow(x, 3) - pow(x, 5) + Order(pow(x, 7)); + result += check_series(e, exZERO(), d); + + e = pow(pow(x, 2) + pow(x, 4), -1); + d = pow(x, -2) - 1 + pow(x, 2) - pow(x, 4) + Order(pow(x, 6)); + result += check_series(e, exZERO(), d); + + e = pow(sin(x), -2); + d = pow(x, -2) + numeric(1,3) + pow(x, 2) / 15 + pow(x, 4) * 2/189 + Order(pow(x, 5)); + result += check_series(e, exZERO(), d); + + e = sin(x) / cos(x); + d = x + pow(x, 3) / 3 + pow(x, 5) * 2/15 + pow(x, 7) * 17/315 + Order(pow(x, 8)); + result += check_series(e, exZERO(), d); + + e = cos(x) / sin(x); + d = pow(x, -1) - x / 3 - pow(x, 3) / 45 - pow(x, 5) * 2/945 + Order(pow(x, 6)); + result += check_series(e, exZERO(), d); + + e = pow(numeric(2), x); + ex t = log(ex(2)) * x; + d = 1 + t + pow(t, 2) / 2 + pow(t, 3) / 6 + pow(t, 4) / 24 + pow(t, 5) / 120 + pow(t, 6) / 720 + pow(t, 7) / 5040 + Order(pow(x, 8)); + result += check_series(e, exZERO(), d.expand()); + + e = pow(Pi, x); + t = log(Pi) * x; + d = 1 + t + pow(t, 2) / 2 + pow(t, 3) / 6 + pow(t, 4) / 24 + pow(t, 5) / 120 + pow(t, 6) / 720 + pow(t, 7) / 5040 + Order(pow(x, 8)); + result += check_series(e, exZERO(), d.expand()); + + return result; } // Series addition static unsigned series2(void) { - unsigned result = 0; - ex e, d; - - e = pow(sin(x), -1).series(x, exZERO(), 8) + pow(sin(-x), -1).series(x, exZERO(), 12); - d = Order(pow(x, 6)); - result += check_series(e, exZERO(), d); - - return result; + unsigned result = 0; + ex e, d; + + e = pow(sin(x), -1).series(x, exZERO(), 8) + pow(sin(-x), -1).series(x, exZERO(), 12); + d = Order(pow(x, 6)); + result += check_series(e, exZERO(), d); + + return result; } // Series multiplication static unsigned series3(void) { - unsigned result = 0; - ex e, d; - - e = sin(x).series(x, exZERO(), 8) * pow(sin(x), -1).series(x, exZERO(), 12); - d = 1 + Order(pow(x, 7)); - result += check_series(e, exZERO(), d); + unsigned result = 0; + ex e, d; + + e = sin(x).series(x, exZERO(), 8) * pow(sin(x), -1).series(x, exZERO(), 12); + d = 1 + Order(pow(x, 7)); + result += check_series(e, exZERO(), d); + + return result; +} - return result; +// Series of special functions +static unsigned series4(void) +{ + unsigned result = 0; + ex e, d; + + e = gamma(2*x); + d = pow(x+1,-1)*numeric(1,4) + + pow(x+1,0)*(numeric(3,4) - + numeric(1,2)*EulerGamma) + + pow(x+1,1)*(numeric(7,4) - + numeric(3,2)*EulerGamma + + numeric(1,2)*pow(EulerGamma,2) + + numeric(1,12)*pow(Pi,2)) + + pow(x+1,2)*(numeric(15,4) - + numeric(7,2)*EulerGamma - + numeric(1,3)*pow(EulerGamma,3) + + numeric(1,4)*pow(Pi,2) + + numeric(3,2)*pow(EulerGamma,2) - + numeric(1,6)*pow(Pi,2)*EulerGamma - + numeric(2,3)*zeta(3)) + + pow(x+1,3)*(numeric(31,4) - pow(EulerGamma,3) - + numeric(15,2)*EulerGamma + + numeric(1,6)*pow(EulerGamma,4) + + numeric(7,2)*pow(EulerGamma,2) + + numeric(7,12)*pow(Pi,2) - + numeric(1,2)*pow(Pi,2)*EulerGamma - + numeric(2)*zeta(3) + + numeric(1,6)*pow(EulerGamma,2)*pow(Pi,2) + + numeric(1,40)*pow(Pi,4) + + numeric(4,3)*zeta(3)*EulerGamma) + + Order(pow(x+1,4)); + result += check_series(e, -1, d, 4); + + e = tan(x*Pi/2); + d = pow(x-1,-1)/Pi*(-2) + + pow(x-1,1)*Pi/6 + + pow(x-1,3)*pow(Pi,3)/360 + + pow(x-1,5)*pow(Pi,5)/15120 + + pow(x-1,7)*pow(Pi,7)/604800 + + Order(pow(x-1,8)); + result += check_series(e,1,d,8); + + return result; } unsigned series_expansion(void) { - unsigned result = 0; - - cout << "checking series expansion..." << flush; - clog << "---------series expansion:" << endl; - - result += series1(); - result += series2(); - result += series3(); - - if (!result) { - cout << " passed "; - clog << "(no output)" << endl; - } else { - cout << " failed "; - } - return result; + unsigned result = 0; + + cout << "checking series expansion..." << flush; + clog << "---------series expansion:" << endl; + + result += series1(); + result += series2(); + result += series3(); + result += series4(); + + if (!result) { + cout << " passed "; + clog << "(no output)" << endl; + } else { + cout << " failed "; + } + return result; } diff --git a/doc/Makefile.in b/doc/Makefile.in index 4a11c86a..2498724b 100644 --- a/doc/Makefile.in +++ b/doc/Makefile.in @@ -267,7 +267,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \ diff --git a/doc/reference/Makefile.in b/doc/reference/Makefile.in index 0931f3e0..152e777b 100644 --- a/doc/reference/Makefile.in +++ b/doc/reference/Makefile.in @@ -163,7 +163,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \ diff --git a/doc/tutorial/Makefile.in b/doc/tutorial/Makefile.in index 92f9ec15..c1715667 100644 --- a/doc/tutorial/Makefile.in +++ b/doc/tutorial/Makefile.in @@ -315,7 +315,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \ diff --git a/doc/tutorial/classhierarchy.fig b/doc/tutorial/classhierarchy.fig index 4df411c9..3f56ce5a 100644 --- a/doc/tutorial/classhierarchy.fig +++ b/doc/tutorial/classhierarchy.fig @@ -9,162 +9,170 @@ Single 1200 2 5 1 1 2 0 7 100 0 -1 4.000 0 0 1 0 1896.500 2807.128 791 579 1973 321 3002 579 1 1 1.00 68.57 137.14 -6 996 1144 2025 1401 -2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 - 2025 1144 996 1144 996 1401 2025 1401 2025 1144 -4 1 0 99 0 0 14 0.0000 4 195 840 1511 1350 expairseq\001 --6 6 3053 630 3876 887 2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 3053 630 3876 630 3876 887 3053 887 3053 630 4 1 0 99 0 0 14 0.0000 4 150 450 3465 836 basic\001 -6 -6 2340 1935 3195 2205 -2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 - 2340 1948 3195 1948 3195 2205 2340 2205 2340 1948 -4 1 0 99 0 0 14 0.0000 4 150 690 2764 2154 exprseq\001 +6 225 450 765 720 +2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 + 225 463 739 463 739 720 225 720 225 463 +4 1 0 99 0 0 14 0.0000 4 105 210 482 668 ex\001 -6 -6 2565 2925 3510 3195 -2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 - 2565 2925 3510 2925 3510 3195 2565 3195 2565 2925 -4 1 0 99 0 0 14 0.0000 4 150 675 3034 3131 indexed\001 +6 4770 1980 5760 2340 +1 2 0 1 0 30 100 0 20 0.000 1 0.0000 5287 2182 473 158 4814 2182 5760 2182 +4 1 0 99 0 0 14 0.0000 4 150 705 5310 2250 numeric\001 -6 -6 3060 3645 3420 3735 -1 4 0 1 0 0 100 0 20 0.000 1 0.0000 3127 3690 23 23 3104 3690 3150 3690 -1 4 0 1 0 0 100 0 20 0.000 1 0.0000 3263 3690 23 23 3240 3690 3286 3690 -1 4 0 1 0 0 100 0 20 0.000 1 0.0000 3397 3690 23 23 3374 3690 3420 3690 +6 4770 1530 5760 1890 +1 2 0 1 0 30 100 0 20 0.000 1 0.0000 5287 1732 473 158 4814 1732 5760 1732 +4 1 0 99 0 0 14 0.0000 4 135 735 5310 1800 constant\001 -6 -6 4410 3600 4770 3690 -1 4 0 1 0 0 100 0 20 0.000 1 0.0000 4477 3645 23 23 4454 3645 4500 3645 -1 4 0 1 0 0 100 0 20 0.000 1 0.0000 4613 3645 23 23 4590 3645 4636 3645 -1 4 0 1 0 0 100 0 20 0.000 1 0.0000 4747 3645 23 23 4724 3645 4770 3645 +6 4770 1080 5760 1440 +1 2 0 1 0 30 100 0 20 0.000 1 0.0000 5287 1282 473 158 4814 1282 5760 1282 +4 1 0 99 0 0 14 0.0000 4 195 615 5310 1350 symbol\001 -6 -6 585 1980 1125 2295 +6 996 1080 2025 1337 +2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 + 2025 1080 996 1080 996 1337 2025 1337 2025 1080 +4 1 0 99 0 0 14 0.0000 4 195 840 1511 1286 expairseq\001 +-6 +6 495 1800 1035 2115 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 585 2019 1099 2019 1099 2276 585 2276 585 2019 -4 1 0 99 0 0 14 0.0000 4 150 315 842 2224 add\001 + 495 1839 1009 1839 1009 2096 495 2096 495 1839 +4 1 0 99 0 0 14 0.0000 4 150 315 752 2044 add\001 -6 -6 225 450 765 720 +6 1260 1800 1800 2115 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 225 463 739 463 739 720 225 720 225 463 -4 1 0 99 0 0 14 0.0000 4 105 210 482 668 ex\001 + 1260 1839 1774 1839 1774 2096 1260 2096 1260 1839 +4 1 0 99 0 0 14 0.0000 4 150 315 1517 2044 mul\001 -6 -6 1350 1980 1890 2295 +6 2295 2295 3150 2565 +2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 + 2295 2308 3150 2308 3150 2565 2295 2565 2295 2308 +4 1 0 99 0 0 14 0.0000 4 150 690 2719 2514 exprseq\001 +-6 +6 2520 3285 3465 3555 +2 2 0 1 0 7 100 0 20 0.000 0 0 -1 0 0 5 + 2520 3285 3465 3285 3465 3555 2520 3555 2520 3285 +4 1 0 99 0 0 14 0.0000 4 150 675 2989 3491 indexed\001 +-6 +6 3015 4005 3375 4095 +1 4 0 1 0 0 100 0 20 0.000 1 0.0000 3082 4050 23 23 3059 4050 3105 4050 +1 4 0 1 0 0 100 0 20 0.000 1 0.0000 3218 4050 23 23 3195 4050 3241 4050 +1 4 0 1 0 0 100 0 20 0.000 1 0.0000 3352 4050 23 23 3329 4050 3375 4050 +-6 +6 1440 3240 2385 3510 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 1350 2019 1864 2019 1864 2276 1350 2276 1350 2019 -4 1 0 99 0 0 14 0.0000 4 150 315 1607 2224 mul\001 + 1440 3240 2385 3240 2385 3510 1440 3510 1440 3240 +4 1 0 99 0 0 14 0.0000 4 150 690 1909 3446 function\001 -6 -6 1305 2520 2025 2790 +6 1800 3870 2475 4140 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 1305 2533 2025 2533 2025 2790 1305 2790 1305 2533 -4 1 0 99 0 0 14 0.0000 4 150 525 1672 2739 ncmul\001 + 1800 3870 2475 3870 2475 4140 1800 4140 1800 3870 +4 1 0 99 0 0 14 0.0000 4 195 585 2134 4076 isospin\001 -6 -6 1485 2880 2430 3150 +6 2520 3870 3015 4140 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 1485 2880 2430 2880 2430 3150 1485 3150 1485 2880 -4 1 0 99 0 0 14 0.0000 4 150 690 1954 3086 function\001 + 2520 3870 3015 3870 3015 4140 2520 4140 2520 3870 +4 1 0 99 0 0 14 0.0000 4 150 435 2764 4076 color\001 -6 -6 1845 3510 2520 3780 +6 1350 2250 2070 2520 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 1845 3510 2520 3510 2520 3780 1845 3780 1845 3510 -4 1 0 99 0 0 14 0.0000 4 195 585 2179 3716 isospin\001 + 1350 2250 2070 2250 2070 2520 1350 2520 1350 2250 +4 1 0 99 0 0 14 0.0000 4 150 510 1684 2456 series\001 -6 -6 2565 3510 3060 3780 +6 1755 1440 2475 1710 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 2565 3510 3060 3510 3060 3780 2565 3780 2565 3510 -4 1 0 99 0 0 14 0.0000 4 150 435 2809 3716 color\001 + 1755 1453 2475 1453 2475 1710 1755 1710 1755 1453 +4 1 0 99 0 0 14 0.0000 4 150 555 2115 1659 power\001 -6 -6 3555 3465 4365 3780 -1 2 0 1 0 30 100 0 20 0.000 1 0.0000 3960 3622 405 157 3555 3622 4365 3622 -4 1 0 99 0 0 14 0.0000 4 150 690 3960 3690 coloridx\001 +6 1170 2835 1890 3105 +2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 + 1170 2848 1890 2848 1890 3105 1170 3105 1170 2848 +4 1 0 99 0 0 14 0.0000 4 150 525 1537 3054 ncmul\001 -6 -6 3600 2880 4275 3195 -1 2 0 1 0 30 100 0 20 0.000 1 0.0000 3937 3037 337 157 3600 3037 4274 3037 -4 1 0 99 0 0 14 0.0000 4 150 255 3960 3105 idx\001 +6 4455 3960 4815 4050 +1 4 0 1 0 0 100 0 20 0.000 1 0.0000 4522 4005 23 23 4499 4005 4545 4005 +1 4 0 1 0 0 100 0 20 0.000 1 0.0000 4658 4005 23 23 4635 4005 4681 4005 +1 4 0 1 0 0 100 0 20 0.000 1 0.0000 4792 4005 23 23 4769 4005 4815 4005 -6 -6 3780 2250 4320 2520 -2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 3806 2263 4320 2263 4320 2520 3806 2520 3806 2263 -4 1 0 99 0 0 14 0.0000 4 150 195 4063 2468 lst\001 +6 3600 3825 4410 4140 +1 2 0 1 0 30 100 0 20 0.000 1 0.0000 4005 3982 405 157 3600 3982 4410 3982 +4 1 0 99 0 0 14 0.0000 4 150 690 4005 4050 coloridx\001 -6 -6 2160 1485 2880 1755 -2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 2160 1498 2880 1498 2880 1755 2160 1755 2160 1498 -4 1 0 99 0 0 14 0.0000 4 150 555 2520 1704 power\001 +6 3645 3240 4320 3555 +1 2 0 1 0 30 100 0 20 0.000 1 0.0000 3982 3397 337 157 3645 3397 4319 3397 +4 1 0 99 0 0 14 0.0000 4 150 255 4005 3465 idx\001 -6 -6 4725 2970 5445 3240 +6 4770 2565 5760 2835 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 4725 2970 5445 2970 5445 3240 4725 3240 4725 2970 -4 1 0 99 0 0 14 0.0000 4 150 510 5059 3176 series\001 + 4770 2578 5760 2578 5760 2835 4770 2835 4770 2578 +4 1 0 99 0 0 14 0.0000 4 150 795 5271 2784 relational\001 -6 -6 4815 2520 5805 2790 +6 4455 2970 5265 3240 2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 - 4853 2533 5779 2533 5779 2790 4853 2790 4853 2533 -4 1 0 99 0 0 14 0.0000 4 150 555 5316 2739 matrix\001 --6 -6 4770 1980 5760 2340 -1 2 0 1 0 30 100 0 20 0.000 1 0.0000 5287 2182 473 158 4814 2182 5760 2182 -4 1 0 99 0 0 14 0.0000 4 150 705 5310 2250 numeric\001 --6 -6 4770 1530 5760 1890 -1 2 0 1 0 30 100 0 20 0.000 1 0.0000 5287 1732 473 158 4814 1732 5760 1732 -4 1 0 99 0 0 14 0.0000 4 135 735 5310 1800 constant\001 + 4455 2983 5265 2983 5265 3240 4455 3240 4455 2983 +4 1 0 99 0 0 14 0.0000 4 150 555 4866 3189 matrix\001 -6 -6 4770 1080 5760 1440 -1 2 0 1 0 30 100 0 20 0.000 1 0.0000 5287 1282 473 158 4814 1282 5760 1282 -4 1 0 99 0 0 14 0.0000 4 195 615 5310 1350 symbol\001 +6 3825 2385 4365 2655 +2 2 0 1 0 30 100 0 20 0.000 0 0 -1 0 0 5 + 3851 2398 4365 2398 4365 2655 3851 2655 3851 2398 +4 1 0 99 0 0 14 0.0000 4 150 195 4108 2603 lst\001 -6 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 3825 945 4802 2121 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3053 939 2076 1144 + 3053 939 2070 1080 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3156 939 2693 1453 + 3156 939 2430 1395 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3722 939 4802 2584 + 3722 939 4815 2565 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 3915 945 4815 1665 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3208 939 2925 1935 + 3285 945 2790 2250 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 1459 1453 894 1967 + 1305 1395 810 1800 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 1562 1453 1575 1980 + 1395 1395 1530 1800 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 3915 881 4815 1215 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3645 945 4770 2925 + 3645 945 4725 2925 +2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 + 1 1 1.00 60.00 120.00 + 3555 945 4050 2340 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 2430 2250 1935 2475 + 3420 945 3825 3240 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3599 937 4005 2205 + 2385 2610 1890 2835 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 2610 2295 2160 2835 + 2565 2610 2115 3195 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 2835 2295 2970 2880 + 2790 2610 2925 3240 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3330 945 3780 2880 + 2610 3600 2250 3825 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 2655 3240 2295 3465 + 2790 3600 2790 3825 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 2835 3240 2835 3465 + 3240 945 2115 2205 2 1 0 1 0 7 100 0 -1 0.000 0 0 -1 0 1 2 1 1 1.00 60.00 120.00 - 3870 3240 3870 3465 + 3915 3600 3915 3825 diff --git a/doc/tutorial/ginac.texi b/doc/tutorial/ginac.texi index 5a7ad9df..af43c973 100644 --- a/doc/tutorial/ginac.texi +++ b/doc/tutorial/ginac.texi @@ -173,7 +173,7 @@ leaves many open questions. @section How to use it from within C++ The GiNaC open framework for symbolic computation within the C++ programming -language does not try to define a language of it's own as conventional +language does not try to define a language of its own as conventional CAS do. Instead, it extends the capabilities of C++ by symbolic manipulations. Here is how to generate and print a simple (and rather pointless) bivariate polynomial with some large coefficients: @@ -256,7 +256,7 @@ convenient window into GiNaC's capabilities. @c node-name, next, previous, up @section What it can do for you -@cindex @code{ginsh} +@cindex @command{ginsh} After invoking @command{ginsh} one can test and experiment with GiNaC's features much like in other Computer Algebra Systems except that it does not provide programming constructs like loops or conditionals. For a @@ -364,6 +364,7 @@ You can differentiate functions and expand them as Taylor or Laurent series (the third argument of @code{series} is the evaluation point, the fourth defines the order): +@cindex Zeta function @example > diff(tan(x),x); tan(x)^2+1 @@ -371,15 +372,24 @@ tan(x)^2+1 x-1/6*x^3+Order(x^4) > series(1/tan(x),x,0,4); x^(-1)-1/3*x+Order(x^2) +> series(gamma(x),x,0,3); +x^(-1)-EulerGamma+(1/12*Pi^2+1/2*EulerGamma^2)*x ++(-1/3*zeta(3)-1/12*Pi^2*EulerGamma-1/6*EulerGamma^3)*x^2+Order(x^3) +> evalf("); +x^(-1.0)-0.5772156649015328606+(0.98905599532797255544)*x +-(0.90747907608088628905)*x^2+Order(x^(3.0)) > series(gamma(2*sin(x)-2),x,Pi/2,6); -(x-1/2*Pi)^(-2)+(-1/12*Pi^2-1/2*EulerGamma^2-1/240)*(x-1/2*Pi)^2 -EulerGamma-1/12+Order((x-1/2*Pi)^3) @end example -If you ever wanted to convert units in C or C++ and found this -is cumbersome, here is the solution. Symbolic types can always be -used as tags for different types of objects. Converting from wrong -units to the metric system is therefore easy: +Here we have made use of the @command{ginsh}-command @code{"} to pop the +previously evaluated element from @command{ginsh}'s internal stack. + +If you ever wanted to convert units in C or C++ and found this is +cumbersome, here is the solution. Symbolic types can always be used as +tags for different types of objects. Converting from wrong units to the +metric system is therefore easy: @example > in=.0254*m; @@ -455,7 +465,7 @@ when developing because it considerably speeds up compilation. @option{--prefix=@var{PREFIX}}: The directory where the compiled library and headers are installed. It defaults to @file{/usr/local} which means that the library is installed in the directory @file{/usr/local/lib}, -the header files in @file{/usr/local/include/GiNaC} and the documentation +the header files in @file{/usr/local/include/ginac} and the documentation (like this one) into @file{/usr/local/share/doc/GiNaC}. @item @@ -469,7 +479,7 @@ to have the header files installed in some other directory than @file{@var{PREFIX}/include/ginac/}. For instance, if you specify @option{--includedir=/usr/include} you will end up with the header files sitting in the directory @file{/usr/include/ginac/}. Note that the -subdirectory @file{GiNaC} is enforced by this process in order to +subdirectory @file{ginac} is enforced by this process in order to keep the header files separated from others. This avoids some clashes and allows for an easier deinstallation of GiNaC. This ought to be considered A Good Thing (tm). @@ -507,7 +517,7 @@ $ ./configure And here is a configuration for a private static GiNaC library with several components sitting in custom places (site-wide @acronym{GCC} and private @acronym{CLN}). The compiler is pursuaded to be picky and full -assertions and debugging are switched on: +assertions and debugging information are switched on: @example $ export CXX=/usr/local/gnu/bin/c++ @@ -623,8 +633,9 @@ meta-class for storing all mathematical objects. * Symbols:: Symbolic objects. * Numbers:: Numerical objects. * Constants:: Pre-defined constants. -* Fundamental operations:: The power, add and mul classes. +* Fundamental containers:: The power, add and mul classes. * Built-in functions:: Mathematical functions. +* Relations:: Equality, Inequality and all that. @end menu @@ -670,17 +681,19 @@ GiNaC's class hierarchy consists of several classes representing mathematical objects, all of which (except for @code{ex} and some helpers) are internally derived from one abstract base class called @code{basic}. You do not have to deal with objects of class -@code{basic}, instead you'll be dealing with symbols and functions of -symbols. You'll soon learn in this chapter how many of the functions on -symbols are really classes. This is because simple symbolic arithmetic -is not supported by languages like C++ so in a certain way GiNaC has to -implement its own arithmetic. - -To give an idea about what kinds of symbolic composits may be built we +@code{basic}, instead you'll be dealing with symbols, numbers, +containers of expressions and so on. You'll soon learn in this chapter +how many of the functions on symbols are really classes. This is +because simple symbolic arithmetic is not supported by languages like +C++ so in a certain way GiNaC has to implement its own arithmetic. + +@cindex container +@cindex atom +To get an idea about what kinds of symbolic composits may be built we have a look at the most important classes in the class hierarchy. The oval classes are atomic ones and the squared classes are containers. -The dashed line symbolizes a "points to" or "handles" relationship while -the solid lines stand for "inherits from" relationship in the class +The dashed line symbolizes a `points to' or `handles' relationship while +the solid lines stand for `inherits from' relationship in the class hierarchy: @image{classhierarchy} @@ -693,15 +706,17 @@ classes derived from them share certain features. An example would be consisting of one expression and a number (@code{numeric}). What @emph{is} visible to the user are the derived classes @code{add} and @code{mul}, representing sums of terms and products, respectively. -We'll come back later to some more details about these two classes and -motivate the use of pairs in sums and products here. +@xref{Internal Structures}, where these two classes are described in +more detail. @node Symbols, Numbers, The Class Hierarchy, Basic Concepts @c node-name, next, previous, up @section Symbols @cindex Symbols (class @code{symbol}) +@cindex hierarchy of classes +@cindex atom Symbols are for symbolic manipulation what atoms are for chemistry. You can declare objects of class @code{symbol} as any other object simply by saying @code{symbol x,y;}. There is, however, a catch in here having to @@ -727,6 +742,7 @@ function that declares a symbol with a name already existent in a symbol in the calling function. Again, comparing them (using @code{operator==} for instance) will always reveal their difference. Watch out, please. +@cindex @code{subs()} Although symbols can be assigned expressions for internal reasons, you should not do it (and we are not going to tell you how it is done). If you want to replace a symbol with something else in an expression, you @@ -738,7 +754,10 @@ can use the expression's @code{.subs()} method. @section Numbers @cindex numbers (class @code{numeric}) +@cindex GMP @cindex CLN +@cindex rational +@cindex fraction For storing numerical things, GiNaC uses Bruno Haible's library @acronym{CLN}. The classes therein serve as foundation classes for GiNaC. @acronym{CLN} stands for Class Library for Numbers or @@ -790,21 +809,22 @@ implicit constructor from C-integers directly to expressions handling numerics at work in most of our examples. This design really becomes convenient when one declares own functions having more than one parameter but it forbids using implicit constructors because that would -lead to ambiguities. +lead to compile-time ambiguities. It may be tempting to construct numbers writing @code{numeric r(3/2)}. This would, however, call C's built-in operator @code{/} for integers first and result in a numeric holding a plain integer 1. @strong{Never -use @code{/} on integers!} Use the constructor from two integers -instead, as shown in the example above. Writing @code{numeric(1)/2} may -look funny but works also. +use the operator @code{/} on integers} unless you know exactly what you +are doing! Use the constructor from two integers instead, as shown in +the example above. Writing @code{numeric(1)/2} may look funny but works +also. @cindex @code{Digits} @cindex accuracy We have seen now the distinction between exact numbers and floating point numbers. Clearly, the user should never have to worry about -dynamically created exact numbers, since their "exactness" always -determines how they ought to be handled, i.e. how "long" they are. The +dynamically created exact numbers, since their `exactness' always +determines how they ought to be handled, i.e. how `long' they are. The situation is different for floating point numbers. Their accuracy is controlled by one @emph{global} variable, called @code{Digits}. (For those readers who know about Maple: it behaves very much like Maple's @@ -897,50 +917,53 @@ zero. The full set of tests that can be applied is listed in the following table. @cartouche -@multitable @columnfractions .33 .67 -@item Method @tab Returns true if@dots{} +@multitable @columnfractions .30 .70 +@item @strong{Method} @tab @strong{Returns true if the object is@dots{}} @item @code{.is_zero()} -@tab object is equal to zero +@tab @dots{}equal to zero @item @code{.is_positive()} -@tab object is not complex and greater than 0 +@tab @dots{}not complex and greater than 0 @item @code{.is_integer()} -@tab object is a (non-complex) integer +@tab @dots{}a (non-complex) integer @item @code{.is_pos_integer()} -@tab object is an integer and greater than 0 +@tab @dots{}an integer and greater than 0 @item @code{.is_nonneg_integer()} -@tab object is an integer and greater equal 0 +@tab @dots{}an integer and greater equal 0 @item @code{.is_even()} -@tab object is an even integer +@tab @dots{}an even integer @item @code{.is_odd()} -@tab object is an odd integer +@tab @dots{}an odd integer @item @code{.is_prime()} -@tab object is a prime integer (probabilistic primality test) +@tab @dots{}a prime integer (probabilistic primality test) @item @code{.is_rational()} -@tab object is an exact rational number (integers are rational, too) +@tab @dots{}an exact rational number (integers are rational, too) @item @code{.is_real()} -@tab object is a real integer, rational or float (i.e. is not complex) +@tab @dots{}a real integer, rational or float (i.e. is not complex) @item @code{.is_cinteger()} -@tab object is a (complex) integer, such as @math{2-3*I} +@tab @dots{}a (complex) integer, such as @math{2-3*I} @item @code{.is_crational()} -@tab object is an exact (complex) rational number (such as @math{2/3+7/2*I}) +@tab @dots{}an exact (complex) rational number (such as @math{2/3+7/2*I}) @end multitable @end cartouche -@node Constants, Fundamental operations, Numbers, Basic Concepts +@node Constants, Fundamental containers, Numbers, Basic Concepts @c node-name, next, previous, up @section Constants @cindex constants (class @code{constant}) -@cindex @code{evalf()} +@cindex @code{Pi} +@cindex @code{Catalan} +@cindex @code{EulerGamma} +@cindex @code{evalf()} Constants behave pretty much like symbols except that they return some specific number when the method @code{.evalf()} is called. The predefined known constants are: @cartouche -@multitable @columnfractions .14 .29 .57 -@item Name @tab Common Name @tab Numerical Value (35 digits) +@multitable @columnfractions .14 .30 .56 +@item @strong{Name} @tab @strong{Common Name} @tab @strong{Numerical Value (to 35 digits)} @item @code{Pi} @tab Archimedes' constant @tab 3.14159265358979323846264338327950288 @@ -954,10 +977,10 @@ The predefined known constants are: @end cartouche -@node Fundamental operations, Built-in functions, Constants, Basic Concepts +@node Fundamental containers, Built-in functions, Constants, Basic Concepts @c node-name, next, previous, up -@section Fundamental operations: the @code{power}, @code{add} and @code{mul} classes -@cindex polynomials +@section Fundamental containers: the @code{power}, @code{add} and @code{mul} classes +@cindex polynomial @cindex @code{add} @cindex @code{mul} @cindex @code{power} @@ -983,6 +1006,7 @@ int main() @} @end example +@cindex @code{pow()} For exponentiation, you have already seen the somewhat clumsy (though C-ish) statement @code{pow(x,2);} to represent @code{x} squared. This direct construction is necessary since we cannot safely overload the constructor @@ -1038,16 +1062,22 @@ expression twice, either implicitly or explicitly, results in the same canonical form. -@node Built-in functions, Important Algorithms, Fundamental operations, Basic Concepts +@node Built-in functions, Relations, Fundamental containers, Basic Concepts @c node-name, next, previous, up @section Built-in functions +@cindex functions (class @code{function}) +@cindex trigonometric function +@cindex hyperbolic function -There are quite a number of useful functions built into GiNaC. They are -all objects of class @code{function}. They accept one or more +There are quite a number of useful functions hard-wired into GiNaC. For +instance, all trigonometric and hyperbolic functions are implemented. +They are all objects of class @code{function}. They accept one or more expressions as arguments and return one expression. If the arguments are not numerical, the evaluation of the function may be halted, as it does in the next example: +@cindex Gamma function +@cindex @code{subs()} @example #include using namespace GiNaC; @@ -1079,10 +1109,27 @@ Most of these functions can be differentiated, series expanded and so on. Read the next chapter in order to learn more about this. -@node Important Algorithms, Polynomial Expansion, Built-in functions, Top +@node Relations, Important Algorithms, Built-in functions, Basic Concepts +@c node-name, next, previous, up +@section Relations +@cindex relations (class @code{relational}) + +Sometimes, a relation holding between two expressions must be stored +somehow. The class @code{relational} is a convenient container for such +purposes. A relation is by definition a container for two @code{ex} and +a relation between them that signals equality, inequality and so on. +They are created by simply using the C++ operators @code{==}, @code{!=}, +@code{<}, @code{<=}, @code{>} and @code{>=} between two expressions. + +@xref{Built-in functions}, for examples where various applications of +the @code{.subs()} method show how objects of class relational are used +as arguments. There they provide an intuitive syntax for substitutions. + + +@node Important Algorithms, Polynomial Expansion, Relations, Top @c node-name, next, previous, up @chapter Important Algorithms -@cindex polynomials +@cindex polynomial In this chapter the most important algorithms provided by GiNaC will be described. Some of them are implemented as functions on expressions, @@ -1105,6 +1152,7 @@ int main() @} @end example +@cindex @code{subs()} The general rule is that wherever methods accept one or more parameters (@var{arg1}, @var{arg2}, @dots{}) the order of arguments the function wrapper accepts is the same but preceded by the object to act on @@ -1119,8 +1167,8 @@ here. Also, users of MuPAD will in most cases feel more comfortable with GiNaC's convention. All function wrappers are always implemented as simple inline functions which just call the corresponding method and are only provided for users uncomfortable with OO who are dead set to -avoid method invocations. Generally, a chain of function wrappers is -much harder to read than a chain of methods and should therefore be +avoid method invocations. Generally, nested function wrappers are much +harder to read than a sequence of methods and should therefore be avoided if possible. On the other hand, not everything in GiNaC is a method on class @code{ex} and sometimes calling a function cannot be avoided. @@ -1177,6 +1225,8 @@ Note that the original polynomial needs to be in expanded form in order to be able to find the coefficients properly. The range of occuring coefficients can be checked using the two methods +@cindex @code{degree()} +@cindex @code{ldegree()} @example int ex::degree(symbol const & s); int ex::ldegree(symbol const & s); @@ -1185,7 +1235,7 @@ int ex::ldegree(symbol const & s); where @code{degree()} returns the highest coefficient and @code{ldegree()} the lowest one. (These two methods work also reliably on non-expanded input polynomials). An application is illustrated in -the next example, where a multivariate polynomial is analysed: +the next example, where a multivariate polynomial is analyzed: @example #include @@ -1303,6 +1353,7 @@ normalized to @code{P_a/P_b} = @code{(4*y+z)/(y+3*z)}. @c node-name, next, previous, up @section Symbolic Differentiation @cindex differentiation +@cindex @code{diff()} @cindex chain rule @cindex product rule @@ -1500,8 +1551,8 @@ evaluate themselves numerically to a precision declared at runtime (using @code{Digits}). Some may be evaluated at certain points, but not generally. This ought to be fixed. However, doing numerical computations with GiNaC's quite abstract classes is doomed to be -inefficient. For this purpose, the underlying bignum-package -@acronym{CLN} is much better suited. +inefficient. For this purpose, the underlying foundation classes +provided by @acronym{CLN} are much better suited. @node Symbolic functions, A Comparison With Other CAS, What does not belong into GiNaC, Extending GiNaC @@ -1510,7 +1561,7 @@ inefficient. For this purpose, the underlying bignum-package The easiest and most instructive way to start with is probably to implement your own function. Objects of class @code{function} are -inserted into the system via a kind of "registry". They get a serial +inserted into the system via a kind of `registry'. They get a serial number that is used internally to identify them but you usually need not worry about this. What you have to care for are functions that are called when the user invokes certain methods. These are usual @@ -1523,14 +1574,16 @@ look something like this: @example static ex cos_eval_method(ex const & x) @{ - // if x%2*Pi return 1 - // if x%Pi return -1 - // if x%Pi/2 return 0 + // if (!x%(2*Pi)) return 1 + // if (!x%Pi) return -1 + // if (!x%Pi/2) return 0 // care for other cases... return cos(x).hold(); @} @end example +@cindex @code{hold()} +@cindex evaluation The last line returns @code{cos(x)} if we don't know what else to do and stops a potential recursive evaluation by saying @code{.hold()}. We should also implement a method for numerical evaluation and since we are @@ -1598,7 +1651,7 @@ have done our best to avoid them where we can.) That's it. May the source be with you! -@node A Comparison With Other CAS, Internal Structures, Symbolic functions, Top +@node A Comparison With Other CAS, Advantages, Symbolic functions, Top @c node-name, next, previous, up @chapter A Comparison With Other CAS @cindex advocacy @@ -1608,8 +1661,15 @@ other, traditional Computer Algebra Systems, like @emph{Maple}, @emph{Mathematica} or @emph{Reduce}, where it has advantages and disadvantages over these systems. +@menu +* Advantages:: Stengths of the GiNaC approach. +* Disadvantages:: Weaknesses of the GiNaC approach. +* Why C++?:: Attractiveness of C++. +@end menu -@heading Advantages +@node Advantages, Disadvantages, A Comparison With Other CAS, A Comparison With Other CAS +@c node-name, next, previous, up +@section Advantages GiNaC has several advantages over traditional Computer Algebra Systems, like @@ -1619,9 +1679,10 @@ Algebra Systems, like @item familiar language: all common CAS implement their own proprietary grammar which you have to learn first (and maybe learn again when your -vendor chooses to "enhance" it). With GiNaC you can write your program +vendor decides to `enhance' it). With GiNaC you can write your program in common C++, which is standardized. +@cindex STL @item structured data types: you can build up structured data types using @code{struct}s or @code{class}es together with STL features instead of @@ -1660,7 +1721,7 @@ programming language and vice versa. With GiNaC, your symbolic routines are part of your program. You can easily call third party libraries, e.g. for numerical evaluation or graphical interaction. All other approaches are much more cumbersome: they range from simply ignoring the -problem (i.e. @emph{Maple}) to providing a method for "embedding" the +problem (i.e. @emph{Maple}) to providing a method for `embedding' the system (i.e. @emph{Yacas}). @item @@ -1673,7 +1734,9 @@ CAS. @end itemize -@heading Disadvantages +@node Disadvantages, Why C++?, Advantages, A Comparison With Other CAS +@c node-name, next, previous, up +@section Disadvantages Of course it also has some disadvantages: @@ -1681,18 +1744,19 @@ Of course it also has some disadvantages: @item not interactive: GiNaC programs have to be written in an editor, -compiled and executed. You cannot play with expressions interactively. +compiled and executed. You cannot play with expressions interactively. However, such an extension is not inherently forbidden by design. In -fact, two interactive interfaces are possible: First, a simple shell -that exposes GiNaC's types to a command line can readily be written (and -has been written) and second, as a more consistent approach we plan an -integration with the @acronym{CINT} C++ interpreter. +fact, two interactive interfaces are possible: First, a shell that +exposes GiNaC's types to a command line can readily be written (the tiny +@command{ginsh} that is part of the distribution being an example) and +second, as a more consistent approach we are working on an integration +with the @acronym{CINT} C++ interpreter. @item advanced features: GiNaC cannot compete with a program like @emph{Reduce} which exists for more than 30 years now or @emph{Maple} which grows since 1981 by the work of dozens of programmers, with -respect to mathematical features. Integration, factorization, +respect to mathematical features. Integration, factorization, non-trivial simplifications, limits etc. are missing in GiNaC (and are not planned for the near future). @@ -1702,31 +1766,35 @@ platform dependent features (it should compile on any ANSI compliant C++ compiler), the currently used version of the CLN library (fast large integer and arbitrary precision arithmetics) can be compiled only on systems with a recently new C++ compiler from the GNU Compiler -Collection (@acronym{GCC}). GiNaC uses recent language features like -explicit constructors, mutable members, RTTI, @code{dynamic_cast}s and -STL, so ANSI compliance is meant literally. Recent @acronym{GCC} -versions starting at 2.95, although itself not yet ANSI compliant, -support all needed features. +Collection (@acronym{GCC}).@footnote{This is because CLN uses +PROVIDE/REQUIRE like macros to let the compiler gather all static +initializations, which works for GNU C++ only.} GiNaC uses recent +language features like explicit constructors, mutable members, RTTI, +@code{dynamic_cast}s and STL, so ANSI compliance is meant literally. +Recent @acronym{GCC} versions starting at 2.95, although itself not yet +ANSI compliant, support all needed features. @end itemize -@heading Why C++? +@node Why C++?, Internal Structures, Disadvantages, A Comparison With Other CAS +@c node-name, next, previous, up +@section Why C++? Why did we choose to implement GiNaC in C++ instead of Java or any other -language? C++ is not perfect: type checking is not strict (casting is +language? C++ is not perfect: type checking is not strict (casting is possible), separation between interface and implementation is not complete, object oriented design is not enforced. The main reason is -the often scolded feature of operator overloading in C++. While it may +the often scolded feature of operator overloading in C++. While it may be true that operating on classes with a @code{+} operator is rarely -meaningful, it is perfectly suited for algebraic expressions. Writing +meaningful, it is perfectly suited for algebraic expressions. Writing @math{3x+5y} as @code{3*x+5*y} instead of -@code{x.times(3).plus(y.times(5))} looks much more natural. Furthermore, -the main developers are more familiar with C++ than with any other -programming language. +@code{x.times(3).plus(y.times(5))} looks much more natural. +Furthermore, the main developers are more familiar with C++ than with +any other programming language. -@node Internal Structures, Expressions are reference counted, A Comparison With Other CAS, Top +@node Internal Structures, Expressions are reference counted, Why C++? , Top @c node-name, next, previous, up @appendix Internal Structures @@ -1849,10 +1917,10 @@ fashion: @cindex pair-wise representation However, doing so results in a rather deeply nested tree which will quickly become inefficient to manipulate. We can improve on this by -representing the sum instead as a sequence of terms, each one being a -pair of a purely numeric multiplicative coefficient and its rest. In -the same spirit we can store the multiplication as a sequence of terms, -each having a numeric exponent and a possibly complicated base, the tree +representing the sum as a sequence of terms, each one being a pair of a +purely numeric multiplicative coefficient and its rest. In the same +spirit we can store the multiplication as a sequence of terms, each +having a numeric exponent and a possibly complicated base, the tree becomes much more flat: @image{reppair} @@ -1871,15 +1939,15 @@ this is avoided by adding a field that carries an overall numeric coefficient. This results in the realistic picture of internal representation for @tex -$2d^3 \left( 4a + 5b - 3 \right)$ +$2d^3 \left( 4a + 5b - 3 \right)$: @end tex @ifnottex -@math{2*d^3*(4*a+5*b-3)} +@math{2*d^3*(4*a+5*b-3)}: @end ifnottex -: @image{repreal} +@cindex radical This also allows for a better handling of numeric radicals, since @code{sqrt(2)} can now be carried along calculations. Now it should be clear, why both classes @code{add} and @code{mul} are derived from the @@ -1923,16 +1991,17 @@ Prints '-I' flags pointing to the installed header files. Prints out the linker flags necessary to link a program against GiNaC. @item --prefix[=@var{PREFIX}] If @var{PREFIX} is specified, overrides the configured value of @env{$prefix}. -(And of exec-prefix, unless --exec-prefix is also specified) +(And of exec-prefix, unless @code{--exec-prefix} is also specified) Otherwise, prints out the configured value of @env{$prefix}. @item --exec-prefix[=@var{PREFIX}] If @var{PREFIX} is specified, overrides the configured value of @env{$exec_prefix}. Otherwise, prints out the configured value of @env{$exec_prefix}. @end table -Typically, @command{ginac-config} will be used within a configure script, -as described below. It, however, can also be used directly -from the command line to compile a simple program. For example: +Typically, @command{ginac-config} will be used within a configure +script, as described below. It, however, can also be used directly from +the command line using backquotes to compile a simple program. For +example: @example c++ -o simple `ginac-config --cppflags` simple.cpp `ginac-config --libs` @@ -1971,7 +2040,7 @@ either found in the user's path, or from the environment variable @env{GINACLIB_CONFIG}. @item -Tests the installed libraries to make sure that there version +Tests the installed libraries to make sure that their version is later than @var{MINIMUM-VERSION}. (A default version will be used if not specified) @@ -2118,7 +2187,7 @@ That command does the following: @display If a GiNaC version greater than 0.4.0 is found, adds @env{$GINACLIB_LIBS} to @env{$LIBS} and @env{$GINACLIB_CPPFLAGS} to @env{$CPPFLAGS}. Otherwise, dies -with the error message "need to have GiNaC installed" +with the error message `need to have GiNaC installed' @end display And the @file{Makefile.am}, which will be used to build the Makefile. diff --git a/doc/tutorial/version.texi b/doc/tutorial/version.texi index 027c1b4a..545a7822 100644 --- a/doc/tutorial/version.texi +++ b/doc/tutorial/version.texi @@ -1,3 +1,3 @@ -@set UPDATED 10 December 1999 +@set UPDATED 12 December 1999 @set EDITION 0.4.1 @set VERSION 0.4.1 diff --git a/ginac/Makefile.am b/ginac/Makefile.am index 0907823a..4f20370d 100644 --- a/ginac/Makefile.am +++ b/ginac/Makefile.am @@ -4,11 +4,10 @@ lib_LTLIBRARIES = libginac.la libginac_la_SOURCES = add.cpp basic.cpp constant.cpp diff.cpp ex.cpp \ expairseq.cpp exprseq.cpp fail.cpp function.cpp inifcns.cpp \ inifcns_trans.cpp inifcns_zeta.cpp inifcns_gamma.cpp matrix.cpp mul.cpp \ - normal.cpp numeric.cpp operators.cpp power.cpp print.cpp printraw.cpp \ - printtree.cpp printcsrc.cpp relational.cpp symbol.cpp utils.cpp series.cpp \ - ncmul.cpp clifford.cpp structure.cpp color.cpp indexed.cpp idx.cpp \ - isospin.cpp exprseq_suppl.cpp lst.cpp lst_suppl.cpp simp_lor.cpp \ - coloridx.cpp lorentzidx.cpp debugmsg.h utils.h + normal.cpp numeric.cpp operators.cpp power.cpp relational.cpp symbol.cpp \ + utils.cpp series.cpp ncmul.cpp clifford.cpp structure.cpp color.cpp \ + indexed.cpp idx.cpp isospin.cpp exprseq_suppl.cpp lst.cpp lst_suppl.cpp \ + simp_lor.cpp coloridx.cpp lorentzidx.cpp debugmsg.h utils.h libginac_la_LDFLAGS = -version-info $(LT_CURRENT):$(LT_REVISION):$(LT_AGE) \ -release $(LT_RELEASE) ginacincludedir = $(includedir)/ginac diff --git a/ginac/Makefile.in b/ginac/Makefile.in index 90ece770..4bd7009b 100644 --- a/ginac/Makefile.in +++ b/ginac/Makefile.in @@ -98,7 +98,7 @@ VERSION = @VERSION@ YACC = @YACC@ lib_LTLIBRARIES = libginac.la -libginac_la_SOURCES = add.cpp basic.cpp constant.cpp diff.cpp ex.cpp expairseq.cpp exprseq.cpp fail.cpp function.cpp inifcns.cpp inifcns_trans.cpp inifcns_zeta.cpp inifcns_gamma.cpp matrix.cpp mul.cpp normal.cpp numeric.cpp operators.cpp power.cpp print.cpp printraw.cpp printtree.cpp printcsrc.cpp relational.cpp symbol.cpp utils.cpp series.cpp ncmul.cpp clifford.cpp structure.cpp color.cpp indexed.cpp idx.cpp isospin.cpp exprseq_suppl.cpp lst.cpp lst_suppl.cpp simp_lor.cpp coloridx.cpp lorentzidx.cpp debugmsg.h utils.h +libginac_la_SOURCES = add.cpp basic.cpp constant.cpp diff.cpp ex.cpp expairseq.cpp exprseq.cpp fail.cpp function.cpp inifcns.cpp inifcns_trans.cpp inifcns_zeta.cpp inifcns_gamma.cpp matrix.cpp mul.cpp normal.cpp numeric.cpp operators.cpp power.cpp relational.cpp symbol.cpp utils.cpp series.cpp ncmul.cpp clifford.cpp structure.cpp color.cpp indexed.cpp idx.cpp isospin.cpp exprseq_suppl.cpp lst.cpp lst_suppl.cpp simp_lor.cpp coloridx.cpp lorentzidx.cpp debugmsg.h utils.h libginac_la_LDFLAGS = -version-info $(LT_CURRENT):$(LT_REVISION):$(LT_AGE) -release $(LT_RELEASE) @@ -120,10 +120,10 @@ libginac_la_LIBADD = libginac_la_OBJECTS = add.lo basic.lo constant.lo diff.lo ex.lo \ expairseq.lo exprseq.lo fail.lo function.lo inifcns.lo inifcns_trans.lo \ inifcns_zeta.lo inifcns_gamma.lo matrix.lo mul.lo normal.lo numeric.lo \ -operators.lo power.lo print.lo printraw.lo printtree.lo printcsrc.lo \ -relational.lo symbol.lo utils.lo series.lo ncmul.lo clifford.lo \ -structure.lo color.lo indexed.lo idx.lo isospin.lo exprseq_suppl.lo \ -lst.lo lst_suppl.lo simp_lor.lo coloridx.lo lorentzidx.lo +operators.lo power.lo relational.lo symbol.lo utils.lo series.lo \ +ncmul.lo clifford.lo structure.lo color.lo indexed.lo idx.lo isospin.lo \ +exprseq_suppl.lo lst.lo lst_suppl.lo simp_lor.lo coloridx.lo \ +lorentzidx.lo CXXFLAGS = @CXXFLAGS@ CXXCOMPILE = $(CXX) $(DEFS) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(AM_CXXFLAGS) $(CXXFLAGS) LTCXXCOMPILE = $(LIBTOOL) --mode=compile $(CXX) $(DEFS) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(AM_CXXFLAGS) $(CXXFLAGS) @@ -150,8 +150,7 @@ DEP_FILES = .deps/add.P .deps/basic.P .deps/clifford.P .deps/color.P \ .deps/inifcns_gamma.P .deps/inifcns_trans.P .deps/inifcns_zeta.P \ .deps/isospin.P .deps/lorentzidx.P .deps/lst.P .deps/lst_suppl.P \ .deps/matrix.P .deps/mul.P .deps/ncmul.P .deps/normal.P .deps/numeric.P \ -.deps/operators.P .deps/power.P .deps/print.P .deps/printcsrc.P \ -.deps/printraw.P .deps/printtree.P .deps/relational.P .deps/series.P \ +.deps/operators.P .deps/power.P .deps/relational.P .deps/series.P \ .deps/simp_lor.P .deps/structure.P .deps/symbol.P .deps/utils.P SOURCES = $(libginac_la_SOURCES) OBJECTS = $(libginac_la_OBJECTS) @@ -288,7 +287,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \ diff --git a/ginac/add.cpp b/ginac/add.cpp index c3f56191..804f07ae 100644 --- a/ginac/add.cpp +++ b/ginac/add.cpp @@ -159,15 +159,115 @@ basic * add::duplicate() const return new add(*this); } +void add::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("add print",LOGLEVEL_PRINT); + if (precedence<=upper_precedence) os << "("; + numeric coeff; + bool first=true; + for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) { + coeff = ex_to_numeric(cit->coeff); + if (!first) { + if (coeff.csgn()==-1) os << '-'; else os << '+'; + } else { + if (coeff.csgn()==-1) os << '-'; + first=false; + } + if (!coeff.is_equal(numONE()) && + !coeff.is_equal(numMINUSONE())) { + if (coeff.csgn()==-1) + (numMINUSONE()*coeff).print(os, precedence); + else + coeff.print(os, precedence); + os << '*'; + } + os << cit->rest; + } + // print the overall numeric coefficient, if present: + if (!overall_coeff.is_zero()) { + if (overall_coeff.info(info_flags::positive)) os << '+'; + os << overall_coeff; + } + if (precedence<=upper_precedence) os << ")"; +} + +void add::printraw(ostream & os) const +{ + debugmsg("add printraw",LOGLEVEL_PRINT); + + os << "+("; + for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { + os << "("; + (*it).rest.bp->printraw(os); + os << ","; + (*it).coeff.bp->printraw(os); + os << "),"; + } + os << ",hash=" << hashvalue << ",flags=" << flags; + os << ")"; +} + +void add::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const +{ + debugmsg("add print csrc", LOGLEVEL_PRINT); + if (precedence <= upper_precedence) + os << "("; + + // Print arguments, separated by "+" + epvector::const_iterator it = seq.begin(); + epvector::const_iterator itend = seq.end(); + while (it != itend) { + + // If the coefficient is -1, it is replaced by a single minus sign + if (it->coeff.compare(numONE()) == 0) { + it->rest.bp->printcsrc(os, type, precedence); + } else if (it->coeff.compare(numMINUSONE()) == 0) { + os << "-"; + it->rest.bp->printcsrc(os, type, precedence); + } else if (ex_to_numeric(it->coeff).numer().compare(numONE()) == 0) { + it->rest.bp->printcsrc(os, type, precedence); + os << "/"; + ex_to_numeric(it->coeff).denom().printcsrc(os, type, precedence); + } else if (ex_to_numeric(it->coeff).numer().compare(numMINUSONE()) == 0) { + os << "-"; + it->rest.bp->printcsrc(os, type, precedence); + os << "/"; + ex_to_numeric(it->coeff).denom().printcsrc(os, type, precedence); + } else { + it->coeff.bp->printcsrc(os, type, precedence); + os << "*"; + it->rest.bp->printcsrc(os, type, precedence); + } + + // Separator is "+", except if the following expression would have a leading minus sign + it++; + if (it != itend && !(it->coeff.compare(numZERO()) < 0 || (it->coeff.compare(numONE()) == 0 && is_ex_exactly_of_type(it->rest, numeric) && it->rest.compare(numZERO()) < 0))) + os << "+"; + } + + if (!overall_coeff.is_equal(exZERO())) { + if (overall_coeff.info(info_flags::positive)) os << '+'; + overall_coeff.bp->printcsrc(os,type,precedence); + } + + if (precedence <= upper_precedence) + os << ")"; +} + bool add::info(unsigned inf) const { // TODO: optimize - if (inf==info_flags::polynomial || inf==info_flags::integer_polynomial || inf==info_flags::rational_polynomial || inf==info_flags::rational_function) { + if (inf==info_flags::polynomial || + inf==info_flags::integer_polynomial || + inf==info_flags::cinteger_polynomial || + inf==info_flags::rational_polynomial || + inf==info_flags::crational_polynomial || + inf==info_flags::rational_function) { for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { if (!(recombine_pair_to_ex(*it).info(inf))) return false; } - return true; + return overall_coeff.info(inf); } else { return expairseq::info(inf); } diff --git a/ginac/add.h b/ginac/add.h index 1acb03d1..19438056 100644 --- a/ginac/add.h +++ b/ginac/add.h @@ -60,8 +60,8 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; - void printraw(ostream & os) const; void print(ostream & os, unsigned upper_precedence=0) const; + void printraw(ostream & os) const; void printcsrc(ostream & os, unsigned type, unsigned upper_precedence=0) const; bool info(unsigned inf) const; int degree(symbol const & s) const; diff --git a/ginac/basic.cpp b/ginac/basic.cpp index c2b4ce15..acacd234 100644 --- a/ginac/basic.cpp +++ b/ginac/basic.cpp @@ -110,6 +110,59 @@ basic::basic(unsigned ti) : flags(0), refcount(0), tinfo_key(ti) // public +/** Output to stream formatted to be useful as ginsh input. */ +void basic::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("basic print",LOGLEVEL_PRINT); + os << "[basic object]"; +} + +/** Output to stream in ugly raw format, so brave developers can have a look + * at the underlying structure. */ +void basic::printraw(ostream & os) const +{ + debugmsg("basic printraw",LOGLEVEL_PRINT); + os << "[basic object]"; +} + +/** Output to stream formatted in tree- (indented-) form, so developers can + * have a look at the underlying structure. */ +void basic::printtree(ostream & os, unsigned indent) const +{ + debugmsg("basic printtree",LOGLEVEL_PRINT); + os << string(indent,' ') << "type=" << typeid(*this).name() + << ", hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")" + << ", flags=" << flags + << ", nops=" << nops() << endl; + for (int i=0; idiff(s); - while ( nth>1 ) { + while (nth>1) { ndiff = ndiff.diff(s); --nth; } diff --git a/ginac/ex.cpp b/ginac/ex.cpp index a8b7b2f2..cdef6254 100644 --- a/ginac/ex.cpp +++ b/ginac/ex.cpp @@ -149,6 +149,7 @@ ex::ex(double const d) // public +/** Swap the contents of two expressions. */ void ex::swap(ex & other) { debugmsg("ex swap",LOGLEVEL_MEMBER_FUNCTION); @@ -163,6 +164,74 @@ void ex::swap(ex & other) other.bp=tmpbp; } +/** Output formatted to be useful as ginsh input. */ +void ex::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("ex print",LOGLEVEL_PRINT); + GINAC_ASSERT(bp!=0); + bp->print(os,upper_precedence); +} + +void ex::printraw(ostream & os) const +{ + debugmsg("ex printraw",LOGLEVEL_PRINT); + GINAC_ASSERT(bp!=0); + os << "ex("; + bp->printraw(os); + os << ")"; +} + +void ex::printtree(ostream & os, unsigned indent) const +{ + debugmsg("ex printtree",LOGLEVEL_PRINT); + GINAC_ASSERT(bp!=0); + // os << "refcount=" << bp->refcount << " "; + bp->printtree(os,indent); +} + +/** Print expression as a C++ statement. The output looks like + * " = ;". The "type" parameter has an effect + * on how number literals are printed. + * + * @param os output stream + * @param type variable type (one of the csrc_types) + * @param var_name variable name to be printed */ +void ex::printcsrc(ostream & os, unsigned type, const char *var_name) const +{ + debugmsg("ex print csrc", LOGLEVEL_PRINT); + GINAC_ASSERT(bp!=0); + switch (type) { + case csrc_types::ctype_float: + os << "float "; + break; + case csrc_types::ctype_double: + os << "double "; + break; + case csrc_types::ctype_cl_N: + os << "cl_N "; + break; + } + os << var_name << " = "; + bp->printcsrc(os, type, 0); + os << ";\n"; +} + +/** Little wrapper arount print to be called within a debugger. */ +void ex::dbgprint(void) const +{ + debugmsg("ex dbgprint",LOGLEVEL_PRINT); + GINAC_ASSERT(bp!=0); + bp->dbgprint(); +} + +/** Little wrapper arount printtree to be called within a debugger. */ +void ex::dbgprinttree(void) const +{ + debugmsg("ex dbgprinttree",LOGLEVEL_PRINT); + GINAC_ASSERT(bp!=0); + bp->dbgprinttree(); +} + bool ex::info(unsigned inf) const { if (inf == info_flags::normal_form) { @@ -433,7 +502,7 @@ void ex::makewriteable() void ex::construct_from_basic(basic const & other) { - if ( (other.flags & status_flags::evaluated)==0 ) { + if ((other.flags & status_flags::evaluated)==0) { // cf. copy constructor ex const & tmpex = other.eval(1); // evaluate only one (top) level bp = tmpex.bp; diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp index ec882394..209c7252 100644 --- a/ginac/expairseq.cpp +++ b/ginac/expairseq.cpp @@ -23,6 +23,7 @@ #include #include #include +#include #include "expairseq.h" #include "lst.h" @@ -164,6 +165,104 @@ basic * expairseq::duplicate() const return new expairseq(*this); } +void expairseq::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("expairseq print",LOGLEVEL_PRINT); + os << "[["; + printseq(os,',',precedence,upper_precedence); + os << "]]"; +} + +void expairseq::printraw(ostream & os) const +{ + debugmsg("expairseq printraw",LOGLEVEL_PRINT); + + os << "expairseq("; + for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) { + os << "("; + (*cit).rest.printraw(os); + os << ","; + (*cit).coeff.printraw(os); + os << "),"; + } + os << ")"; +} + +void expairseq::printtree(ostream & os, unsigned indent) const +{ + debugmsg("expairseq printtree",LOGLEVEL_PRINT); + + os << string(indent,' ') << "type=" << typeid(*this).name() + << ", hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")" + << ", flags=" << flags + << ", nops=" << nops() << endl; + for (unsigned i=0; i0) { + os << string(indent+delta_indent,' ') + << "bin " << i << " with entries "; + for (epplist::const_iterator it=hashtab[i].begin(); + it!=hashtab[i].end(); ++it) { + os << *it-seq.begin() << " "; + this_bin_fill++; + } + os << endl; + cum_fill += this_bin_fill; + cum_fill_sq += this_bin_fill*this_bin_fill; + } + if (this_bin_fill0) fact *= k; + double prob=pow(lambda,k)/fact*exp(-lambda); + cum_prob += prob; + os << string(indent+delta_indent,' ') << "bins with " << k << " entries: " + << int(1000.0*count[k]/hashtabsize)/10.0 << "% (expected: " + << int(prob*1000)/10.0 << ")" << endl; + } + os << string(indent+delta_indent,' ') << "bins with more entries: " + << int(1000.0*count[MAXCOUNT]/hashtabsize)/10.0 << "% (expected: " + << int((1-cum_prob)*1000)/10.0 << ")" << endl; + + os << string(indent+delta_indent,' ') << "variance: " + << 1.0/hashtabsize*cum_fill_sq-(1.0/hashtabsize*cum_fill)*(1.0/hashtabsize*cum_fill) + << endl; + os << string(indent+delta_indent,' ') << "average fill: " + << (1.0*cum_fill)/hashtabsize + << " (should be equal to " << (1.0*seq.size())/hashtabsize << ")" << endl; +#endif // def EXPAIRSEQ_USE_HASHTAB +} + bool expairseq::info(unsigned inf) const { return basic::info(inf); @@ -411,6 +510,33 @@ ex expairseq::thisexpairseq(epvector * vp, ex const & oc) const return expairseq(vp,oc); } +void expairseq::printpair(ostream & os, expair const & p, unsigned upper_precedence) const +{ + os << "[["; + p.rest.bp->print(os,precedence); + os << ","; + p.coeff.bp->print(os,precedence); + os << "]]"; +} + +void expairseq::printseq(ostream & os, char delim, unsigned this_precedence, + unsigned upper_precedence) const +{ + if (this_precedence<=upper_precedence) os << "("; + epvector::const_iterator it,it_last; + it_last=seq.end(); + --it_last; + for (it=seq.begin(); it!=it_last; ++it) { + printpair(os,*it,this_precedence); + os << delim; + } + printpair(os,*it,this_precedence); + if (!overall_coeff.is_equal(default_overall_coeff())) { + os << delim << overall_coeff; + } + if (this_precedence<=upper_precedence) os << ")"; +} + expair expairseq::split_ex_to_pair(ex const & e) const { return expair(e,exONE()); @@ -1421,7 +1547,7 @@ epvector * expairseq::evalchildren(int level) const ex const & evaled_ex=(*cit).rest.eval(level); if (!are_ex_trivially_equal((*cit).rest,evaled_ex)) { - // something changed, copy seq, eval and return it + // something changed, copy seq, eval and return it epvector *s=new epvector; s->reserve(seq.size()); diff --git a/ginac/expairseq.h b/ginac/expairseq.h index 72ac829c..50fba2ff 100644 --- a/ginac/expairseq.h +++ b/ginac/expairseq.h @@ -90,9 +90,9 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; + void print(ostream & os, unsigned upper_precedence=0) const; void printraw(ostream & os) const; void printtree(ostream & os, unsigned indent) const; - void print(ostream & os, unsigned upper_precedence=0) const; bool info(unsigned inf) const; int nops() const; ex op(int const i) const; diff --git a/ginac/fail.cpp b/ginac/fail.cpp index d4997196..bc75c47e 100644 --- a/ginac/fail.cpp +++ b/ginac/fail.cpp @@ -92,6 +92,18 @@ basic * fail::duplicate() const return new fail(*this); } +void fail::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("fail print",LOGLEVEL_PRINT); + os << "FAIL"; +} + +void fail::printraw(ostream & os) const +{ + debugmsg("fail printraw",LOGLEVEL_PRINT); + os << "FAIL"; +} + // protected int fail::compare_same_type(basic const & other) const diff --git a/ginac/fail.h b/ginac/fail.h index 292bf797..2a3617d2 100644 --- a/ginac/fail.h +++ b/ginac/fail.h @@ -50,8 +50,8 @@ protected: // functions overriding virtual functions from bases classes public: basic * duplicate() const; - void printraw(ostream & os) const; void print(ostream & os, unsigned upper_precedence=0) const; + void printraw(ostream & os) const; protected: int compare_same_type(basic const & other) const; unsigned return_type(void) const { return return_types::noncommutative_composite; }; diff --git a/ginac/flags.h b/ginac/flags.h index 6a1bfa4d..da2768d2 100644 --- a/ginac/flags.h +++ b/ginac/flags.h @@ -50,6 +50,8 @@ public: real, rational, integer, + crational, + cinteger, positive, negative, nonnegative, @@ -81,7 +83,9 @@ public: // answered by classes numeric, symbol, add, mul, power polynomial, integer_polynomial, + cinteger_polynomial, rational_polynomial, + crational_polynomial, rational_function, // answered by class ex diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp index 25eaae0d..8b7f327a 100644 --- a/ginac/inifcns.cpp +++ b/ginac/inifcns.cpp @@ -252,7 +252,7 @@ ex ncpower(ex const &basis, unsigned exponent) /** Force inclusion of functions from initcns_gamma and inifcns_zeta * for static lib (so ginsh will see them). */ unsigned force_include_gamma = function_index_gamma; -unsigned force_include_zeta = function_index_zeta; +unsigned force_include_zeta1 = function_index_zeta1; #ifndef NO_GINAC_NAMESPACE } // namespace GiNaC diff --git a/ginac/inifcns.h b/ginac/inifcns.h index a5a52700..d7642ba2 100644 --- a/ginac/inifcns.h +++ b/ginac/inifcns.h @@ -81,9 +81,17 @@ DECLARE_FUNCTION_1P(Li2) /** Trilogarithm. */ DECLARE_FUNCTION_1P(Li3) +// overloading at work: we cannot use the macros /** Riemann's Zeta-function. */ -DECLARE_FUNCTION_1P(zeta) -//DECLARE_FUNCTION_2P(zeta) +extern const unsigned function_index_zeta1; +inline function zeta(ex const & p1) { + return function(function_index_zeta1, p1); +} +/** Derivatives of Riemann's Zeta-function. */ +extern const unsigned function_index_zeta2; +inline function zeta(ex const & p1, ex const & p2) { + return function(function_index_zeta2, p1, p2); +} /** Gamma-function. */ DECLARE_FUNCTION_1P(gamma) @@ -91,18 +99,17 @@ DECLARE_FUNCTION_1P(gamma) /** Beta-function. */ DECLARE_FUNCTION_2P(beta) +// overloading at work: we cannot use the macros /** Psi-function (aka polygamma-function). */ -// overloading @ work: we cannot use the macros extern const unsigned function_index_psi1; inline function psi(ex const & p1) { return function(function_index_psi1, p1); } +/** Derivatives of Psi-function (aka polygamma-functions). */ extern const unsigned function_index_psi2; inline function psi(ex const & p1, ex const & p2) { return function(function_index_psi2, p1, p2); } -//DECLARE_FUNCTION_1P(psi) -//DECLARE_FUNCTION_2P(psi) /** Factorial function. */ DECLARE_FUNCTION_1P(factorial) @@ -119,7 +126,7 @@ ex ncpower(ex const &basis, unsigned exponent); inline bool is_order_function(ex const & e) { - return is_ex_the_function(e, Order); + return is_ex_the_function(e, Order); } #ifndef NO_GINAC_NAMESPACE diff --git a/ginac/inifcns_gamma.cpp b/ginac/inifcns_gamma.cpp index 7f3b9f5c..ecec37a6 100644 --- a/ginac/inifcns_gamma.cpp +++ b/ginac/inifcns_gamma.cpp @@ -42,6 +42,15 @@ namespace GiNaC { // Gamma-function ////////// +static ex gamma_evalf(ex const & x) +{ + BEGIN_TYPECHECK + TYPECHECK(x,numeric) + END_TYPECHECK(gamma(x)) + + return gamma(ex_to_numeric(x)); +} + /** Evaluation of gamma(x). Knows about integer arguments, half-integer * arguments and that's it. Somebody ought to provide some good numerical * evaluation some day... @@ -81,15 +90,6 @@ static ex gamma_eval(ex const & x) return gamma(x).hold(); } -static ex gamma_evalf(ex const & x) -{ - BEGIN_TYPECHECK - TYPECHECK(x,numeric) - END_TYPECHECK(gamma(x)) - - return gamma(ex_to_numeric(x)); -} - static ex gamma_diff(ex const & x, unsigned diff_param) { GINAC_ASSERT(diff_param==0); @@ -122,6 +122,17 @@ REGISTER_FUNCTION(gamma, gamma_eval, gamma_evalf, gamma_diff, gamma_series); // Beta-function ////////// +static ex beta_evalf(ex const & x, ex const & y) +{ + BEGIN_TYPECHECK + TYPECHECK(x,numeric) + TYPECHECK(y,numeric) + END_TYPECHECK(beta(x,y)) + + return gamma(ex_to_numeric(x))*gamma(ex_to_numeric(y)) + / gamma(ex_to_numeric(x+y)); +} + static ex beta_eval(ex const & x, ex const & y) { if (x.info(info_flags::numeric) && y.info(info_flags::numeric)) { @@ -155,17 +166,6 @@ static ex beta_eval(ex const & x, ex const & y) return beta(x,y).hold(); } -static ex beta_evalf(ex const & x, ex const & y) -{ - BEGIN_TYPECHECK - TYPECHECK(x,numeric) - TYPECHECK(y,numeric) - END_TYPECHECK(beta(x,y)) - - return gamma(ex_to_numeric(x))*gamma(ex_to_numeric(y)) - / gamma(ex_to_numeric(x+y)); -} - static ex beta_diff(ex const & x, ex const & y, unsigned diff_param) { GINAC_ASSERT(diff_param<2); @@ -184,6 +184,15 @@ REGISTER_FUNCTION(beta, beta_eval, beta_evalf, beta_diff, NULL); // Psi-function (aka polygamma-function) ////////// +static ex psi1_evalf(ex const & x) +{ + BEGIN_TYPECHECK + TYPECHECK(x,numeric) + END_TYPECHECK(psi(x)) + + return psi(ex_to_numeric(x)); +} + /** Evaluation of polygamma-function psi(x). * Somebody ought to provide some good numerical evaluation some day... */ static ex psi1_eval(ex const & x) @@ -214,15 +223,6 @@ static ex psi1_eval(ex const & x) return psi(x).hold(); } -static ex psi1_evalf(ex const & x) -{ - BEGIN_TYPECHECK - TYPECHECK(x,numeric) - END_TYPECHECK(psi(x)) - - return psi(ex_to_numeric(x)); -} - static ex psi1_diff(ex const & x, unsigned diff_param) { GINAC_ASSERT(diff_param==0); @@ -236,6 +236,16 @@ const unsigned function_index_psi1 = function::register_new("psi", psi1_eval, ps // Psi-functions (aka polygamma-functions) psi(0,x)==psi(x) ////////// +static ex psi2_evalf(ex const & n, ex const & x) +{ + BEGIN_TYPECHECK + TYPECHECK(n,numeric) + TYPECHECK(x,numeric) + END_TYPECHECK(psi(n,x)) + + return psi(ex_to_numeric(n), ex_to_numeric(x)); +} + /** Evaluation of polygamma-function psi(n,x). * Somebody ought to provide some good numerical evaluation some day... */ static ex psi2_eval(ex const & n, ex const & x) @@ -256,16 +266,6 @@ static ex psi2_eval(ex const & n, ex const & x) return psi(n, x).hold(); } -static ex psi2_evalf(ex const & n, ex const & x) -{ - BEGIN_TYPECHECK - TYPECHECK(n,numeric) - TYPECHECK(x,numeric) - END_TYPECHECK(psi(n,x)) - - return psi(ex_to_numeric(n), ex_to_numeric(x)); -} - static ex psi2_diff(ex const & n, ex const & x, unsigned diff_param) { GINAC_ASSERT(diff_param<2); diff --git a/ginac/inifcns_trans.cpp b/ginac/inifcns_trans.cpp index 01711a7f..8838360a 100644 --- a/ginac/inifcns_trans.cpp +++ b/ginac/inifcns_trans.cpp @@ -29,6 +29,9 @@ #include "constant.h" #include "numeric.h" #include "power.h" +#include "relational.h" +#include "symbol.h" +#include "utils.h" #ifndef NO_GINAC_NAMESPACE namespace GiNaC { @@ -71,11 +74,11 @@ static ex exp_eval(ex const & x) return x.op(0); // exp(float) - if (x.info(info_flags::numeric) && !x.info(info_flags::rational)) + if (x.info(info_flags::numeric) && !x.info(info_flags::crational)) return exp_evalf(x); return exp(x).hold(); -} +} static ex exp_diff(ex const & x, unsigned diff_param) { @@ -118,12 +121,18 @@ static ex log_eval(ex const & x) if (x.is_equal(exZERO())) throw(std::domain_error("log_eval(): log(0)")); // log(float) - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return log_evalf(x); } + // log(exp(t)) -> t (for real-valued t): + if (is_ex_the_function(x, exp)) { + ex t=x.op(0); + if (t.info(info_flags::real)) + return t; + } return log(x).hold(); -} +} static ex log_diff(ex const & x, unsigned diff_param) { @@ -149,17 +158,19 @@ static ex sin_evalf(ex const & x) static ex sin_eval(ex const & x) { - // sin(n*Pi) -> 0 ex xOverPi=x/Pi; - if (xOverPi.info(info_flags::integer)) - return exZERO(); - - // sin((2n+1)*Pi/2) -> {+|-}1 - ex xOverPiMinusHalf=xOverPi-exHALF(); - if (xOverPiMinusHalf.info(info_flags::even)) - return exONE(); - else if (xOverPiMinusHalf.info(info_flags::odd)) - return exMINUSONE(); + if (xOverPi.info(info_flags::numeric)) { + // sin(n*Pi) -> 0 + if (xOverPi.info(info_flags::integer)) + return exZERO(); + + // sin((2n+1)*Pi/2) -> {+|-}1 + ex xOverPiMinusHalf=xOverPi-exHALF(); + if (xOverPiMinusHalf.info(info_flags::even)) + return exONE(); + else if (xOverPiMinusHalf.info(info_flags::odd)) + return exMINUSONE(); + } if (is_ex_exactly_of_type(x, function)) { ex t=x.op(0); @@ -175,7 +186,7 @@ static ex sin_eval(ex const & x) } // sin(float) -> float - if (x.info(info_flags::numeric) && !x.info(info_flags::rational)) + if (x.info(info_flags::numeric) && !x.info(info_flags::crational)) return sin_evalf(x); return sin(x).hold(); @@ -205,17 +216,19 @@ static ex cos_evalf(ex const & x) static ex cos_eval(ex const & x) { - // cos(n*Pi) -> {+|-}1 ex xOverPi=x/Pi; - if (xOverPi.info(info_flags::even)) - return exONE(); - else if (xOverPi.info(info_flags::odd)) - return exMINUSONE(); - - // cos((2n+1)*Pi/2) -> 0 - ex xOverPiMinusHalf=xOverPi-exHALF(); - if (xOverPiMinusHalf.info(info_flags::integer)) - return exZERO(); + if (xOverPi.info(info_flags::numeric)) { + // cos(n*Pi) -> {+|-}1 + if (xOverPi.info(info_flags::even)) + return exONE(); + else if (xOverPi.info(info_flags::odd)) + return exMINUSONE(); + + // cos((2n+1)*Pi/2) -> 0 + ex xOverPiMinusHalf=xOverPi-exHALF(); + if (xOverPiMinusHalf.info(info_flags::integer)) + return exZERO(); + } if (is_ex_exactly_of_type(x, function)) { ex t=x.op(0); @@ -231,7 +244,7 @@ static ex cos_eval(ex const & x) } // cos(float) -> float - if (x.info(info_flags::numeric) && !x.info(info_flags::rational)) + if (x.info(info_flags::numeric) && !x.info(info_flags::crational)) return cos_evalf(x); return cos(x).hold(); @@ -292,7 +305,7 @@ static ex tan_eval(ex const & x) } // tan(float) -> float - if (x.info(info_flags::numeric) && !x.info(info_flags::rational)) { + if (x.info(info_flags::numeric) && !x.info(info_flags::crational)) { return tan_evalf(x); } @@ -306,7 +319,19 @@ static ex tan_diff(ex const & x, unsigned diff_param) return (1+power(tan(x),exTWO())); } -REGISTER_FUNCTION(tan, tan_eval, tan_evalf, tan_diff, NULL); +static ex tan_series(ex const & x, symbol const & s, ex const & point, int order) +{ + // method: + // Taylor series where there is no pole falls back to tan_diff. + // On a pole simply expand sin(x)/cos(x). + ex xpoint = x.subs(s==point); + if (!(2*xpoint/Pi).info(info_flags::odd)) + throw do_taylor(); + // if we got here we have to care for a simple pole + return (sin(x)/cos(x)).series(s, point, order+2); +} + +REGISTER_FUNCTION(tan, tan_eval, tan_evalf, tan_diff, tan_series); ////////// // inverse sine (arc sine) @@ -340,7 +365,7 @@ static ex asin_eval(ex const & x) if (x.is_equal(exMINUSONE())) return numeric(-1,2)*Pi; // asin(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return asin_evalf(x); } @@ -388,7 +413,7 @@ static ex acos_eval(ex const & x) if (x.is_equal(exMINUSONE())) return Pi; // acos(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return acos_evalf(x); } @@ -424,7 +449,7 @@ static ex atan_eval(ex const & x) if (x.is_equal(exZERO())) return exZERO(); // atan(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return atan_evalf(x); } @@ -456,8 +481,8 @@ static ex atan2_evalf(ex const & y, ex const & x) static ex atan2_eval(ex const & y, ex const & x) { - if (y.info(info_flags::numeric) && !y.info(info_flags::rational) && - x.info(info_flags::numeric) && !x.info(info_flags::rational)) { + if (y.info(info_flags::numeric) && !y.info(info_flags::crational) && + x.info(info_flags::numeric) && !x.info(info_flags::crational)) { return atan2_evalf(y,x); } @@ -498,10 +523,24 @@ static ex sinh_eval(ex const & x) if (x.is_zero()) return exZERO(); // sinh(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return sinh_evalf(x); } + ex xOverPiI=x/Pi/I; + if (xOverPiI.info(info_flags::numeric)) { + // sinh(n*Pi*I) -> 0 + if (xOverPiI.info(info_flags::integer)) + return exZERO(); + + // sinh((2n+1)*Pi*I/2) -> {+|-}I + ex xOverPiIMinusHalf=xOverPiI-exHALF(); + if (xOverPiIMinusHalf.info(info_flags::even)) + return I; + else if (xOverPiIMinusHalf.info(info_flags::odd)) + return -I; + } + if (is_ex_exactly_of_type(x, function)) { ex t=x.op(0); // sinh(asinh(x)) -> x @@ -547,10 +586,24 @@ static ex cosh_eval(ex const & x) if (x.is_zero()) return exONE(); // cosh(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return cosh_evalf(x); } + ex xOverPiI=x/Pi/I; + if (xOverPiI.info(info_flags::numeric)) { + // cosh(n*Pi*I) -> {+|-}1 + if (xOverPiI.info(info_flags::even)) + return exONE(); + else if (xOverPiI.info(info_flags::odd)) + return exMINUSONE(); + + // cosh((2n+1)*Pi*I/2) -> 0 + ex xOverPiIMinusHalf=xOverPiI-exHALF(); + if (xOverPiIMinusHalf.info(info_flags::integer)) + return exZERO(); + } + if (is_ex_exactly_of_type(x, function)) { ex t=x.op(0); // cosh(acosh(x)) -> x @@ -596,7 +649,7 @@ static ex tanh_eval(ex const & x) if (x.is_zero()) return exZERO(); // tanh(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return tanh_evalf(x); } @@ -623,7 +676,19 @@ static ex tanh_diff(ex const & x, unsigned diff_param) return exONE()-power(tanh(x),exTWO()); } -REGISTER_FUNCTION(tanh, tanh_eval, tanh_evalf, tanh_diff, NULL); +static ex tanh_series(ex const & x, symbol const & s, ex const & point, int order) +{ + // method: + // Taylor series where there is no pole falls back to tanh_diff. + // On a pole simply expand sinh(x)/cosh(x). + ex xpoint = x.subs(s==point); + if (!(2*I*xpoint/Pi).info(info_flags::odd)) + throw do_taylor(); + // if we got here we have to care for a simple pole + return (sinh(x)/cosh(x)).series(s, point, order+2); +} + +REGISTER_FUNCTION(tanh, tanh_eval, tanh_evalf, tanh_diff, tanh_series); ////////// // inverse hyperbolic sine (trigonometric function) @@ -645,7 +710,7 @@ static ex asinh_eval(ex const & x) if (x.is_zero()) return exZERO(); // asinh(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return asinh_evalf(x); } @@ -687,7 +752,7 @@ static ex acosh_eval(ex const & x) if (x.is_equal(exMINUSONE())) return Pi*I; // acosh(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return acosh_evalf(x); } @@ -726,7 +791,7 @@ static ex atanh_eval(ex const & x) if (x.is_equal(exONE()) || x.is_equal(exONE())) throw (std::domain_error("atanh_eval(): infinity")); // atanh(float) -> float - if (!x.info(info_flags::rational)) + if (!x.info(info_flags::crational)) return atanh_evalf(x); } diff --git a/ginac/inifcns_zeta.cpp b/ginac/inifcns_zeta.cpp index e343b6fb..6829d565 100644 --- a/ginac/inifcns_zeta.cpp +++ b/ginac/inifcns_zeta.cpp @@ -38,7 +38,16 @@ namespace GiNaC { // Riemann's Zeta-function ////////// -static ex zeta_eval(ex const & x) +static ex zeta1_evalf(ex const & x) +{ + BEGIN_TYPECHECK + TYPECHECK(x,numeric) + END_TYPECHECK(zeta(x)) + + return zeta(ex_to_numeric(x)); +} + +static ex zeta1_eval(ex const & x) { if (x.info(info_flags::numeric)) { numeric y = ex_to_numeric(x); @@ -64,16 +73,43 @@ static ex zeta_eval(ex const & x) return zeta(x).hold(); } -static ex zeta_evalf(ex const & x) +static ex zeta1_diff(ex const & x, unsigned diff_param) { - BEGIN_TYPECHECK - TYPECHECK(x,numeric) - END_TYPECHECK(zeta(x)) + GINAC_ASSERT(diff_param==0); - return zeta(ex_to_numeric(x)); + return zeta(exONE(), x); +} + +const unsigned function_index_zeta1 = function::register_new("zeta", zeta1_eval, zeta1_evalf, zeta1_diff, NULL); + +////////// +// Derivatives of Riemann's Zeta-function zeta(0,x)==zeta(x) +////////// + +static ex zeta2_eval(ex const & n, ex const & x) +{ + if (n.info(info_flags::numeric)) { + // zeta(0,x) -> zeta(x) + if (n.is_zero()) + return zeta(x); + } + + return zeta(n, x).hold(); +} + +static ex zeta2_diff(ex const & n, ex const & x, unsigned diff_param) +{ + GINAC_ASSERT(diff_param<2); + + if (diff_param==0) { + // d/dn zeta(n,x) + throw(std::logic_error("cannot diff zeta(n,x) with respect to n")); + } + // d/dx psi(n,x) + return zeta(n+1,x); } -REGISTER_FUNCTION(zeta, zeta_eval, zeta_evalf, NULL, NULL); +const unsigned function_index_zeta2 = function::register_new("zeta", zeta2_eval, NULL, zeta2_diff, NULL); #ifndef NO_GINAC_NAMESPACE } // namespace GiNaC diff --git a/ginac/matrix.cpp b/ginac/matrix.cpp index 43288df4..a4db03dd 100644 --- a/ginac/matrix.cpp +++ b/ginac/matrix.cpp @@ -119,6 +119,42 @@ basic * matrix::duplicate() const return new matrix(*this); } +void matrix::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("matrix print",LOGLEVEL_PRINT); + os << "[[ "; + for (int r=0; rprintraw(os); + os << ","; + (*it).coeff.bp->printraw(os); + os << "),"; + } + os << ",hash=" << hashvalue << ",flags=" << flags; + os << ")"; +} + +void mul::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const +{ + debugmsg("mul print csrc", LOGLEVEL_PRINT); + if (precedence <= upper_precedence) + os << "("; + + if (!overall_coeff.is_equal(exONE())) { + overall_coeff.bp->printcsrc(os,type,precedence); + os << "*"; + } + + // Print arguments, separated by "*" or "/" + epvector::const_iterator it = seq.begin(); + epvector::const_iterator itend = seq.end(); + while (it != itend) { + + // If the first argument is a negative integer power, it gets printed as "1.0/" + if (it == seq.begin() && ex_to_numeric(it->coeff).is_integer() && it->coeff.compare(numZERO()) < 0) { + if (type == csrc_types::ctype_cl_N) + os << "recip("; + else + os << "1.0/"; + } + + // If the exponent is 1 or -1, it is left out + if (it->coeff.compare(exONE()) == 0 || it->coeff.compare(numMINUSONE()) == 0) + it->rest.bp->printcsrc(os, type, precedence); + else + // outer parens around ex needed for broken gcc-2.95 parser: + (ex(power(it->rest, abs(ex_to_numeric(it->coeff))))).bp->printcsrc(os, type, upper_precedence); + + // Separator is "/" for negative integer powers, "*" otherwise + it++; + if (it != itend) { + if (ex_to_numeric(it->coeff).is_integer() && it->coeff.compare(numZERO()) < 0) + os << "/"; + else + os << "*"; + } + } + if (precedence <= upper_precedence) + os << ")"; +} + bool mul::info(unsigned inf) const { // TODO: optimize - if (inf==info_flags::polynomial || inf==info_flags::integer_polynomial || inf==info_flags::rational_polynomial || inf==info_flags::rational_function) { + if (inf==info_flags::polynomial || + inf==info_flags::integer_polynomial || + inf==info_flags::cinteger_polynomial || + inf==info_flags::rational_polynomial || + inf==info_flags::crational_polynomial || + inf==info_flags::rational_function) { for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { if (!(recombine_pair_to_ex(*it).info(inf))) return false; } - return true; + return overall_coeff.info(inf); } else { return expairseq::info(inf); } @@ -366,10 +458,10 @@ unsigned mul::return_type(void) const all_commutative=0; } if ((rt==return_types::noncommutative)&&(!all_commutative)) { - // another nc element found, compare type_infos + // another nc element found, compare type_infos if ((*cit_noncommutative_element).rest.return_type_tinfo()!=(*cit).rest.return_type_tinfo()) { - // diffent types -> mul is ncc - return return_types::noncommutative_composite; + // diffent types -> mul is ncc + return return_types::noncommutative_composite; } } } diff --git a/ginac/mul.h b/ginac/mul.h index 76546eb0..71b686bb 100644 --- a/ginac/mul.h +++ b/ginac/mul.h @@ -61,8 +61,8 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; - void printraw(ostream & os) const; void print(ostream & os, unsigned upper_precedence) const; + void printraw(ostream & os) const; void printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const; bool info(unsigned inf) const; int degree(symbol const & s) const; diff --git a/ginac/ncmul.cpp b/ginac/ncmul.cpp index 9e14edc0..d94e61c4 100644 --- a/ginac/ncmul.cpp +++ b/ginac/ncmul.cpp @@ -146,6 +146,39 @@ basic * ncmul::duplicate() const return new ncmul(*this); } +void ncmul::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("ncmul print",LOGLEVEL_PRINT); + printseq(os,'(','%',')',precedence,upper_precedence); +} + +void ncmul::printraw(ostream & os) const +{ + debugmsg("ncmul printraw",LOGLEVEL_PRINT); + + os << "%("; + for (exvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { + (*it).bp->printraw(os); + os << ","; + } + os << ",hash=" << hashvalue << ",flags=" << flags; + os << ")"; +} + +void ncmul::printcsrc(ostream & os, unsigned upper_precedence) const +{ + debugmsg("ncmul print csrc",LOGLEVEL_PRINT); + exvector::const_iterator it; + exvector::const_iterator itend = seq.end()-1; + os << "ncmul("; + for (it=seq.begin(); it!=itend; ++it) { + (*it).bp->printcsrc(os,precedence); + os << ","; + } + (*it).bp->printcsrc(os,precedence); + os << ")"; +} + bool ncmul::info(unsigned inf) const { throw(std::logic_error("which flags have to be implemented in ncmul::info()?")); diff --git a/ginac/ncmul.h b/ginac/ncmul.h index de48b59e..c4e77655 100644 --- a/ginac/ncmul.h +++ b/ginac/ncmul.h @@ -64,8 +64,8 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; - void printraw(ostream & os) const; void print(ostream & os, unsigned upper_precedence) const; + void printraw(ostream & os) const; void printcsrc(ostream & os, unsigned upper_precedence) const; bool info(unsigned inf) const; int degree(symbol const & s) const; diff --git a/ginac/normal.cpp b/ginac/normal.cpp index 3549dfee..9a3cd159 100644 --- a/ginac/normal.cpp +++ b/ginac/normal.cpp @@ -1308,21 +1308,20 @@ ex numeric::normal(lst &sym_lst, lst &repl_lst, int level) const if (is_real()) if (is_rational()) return *this; - else - return replace_with_symbol(*this, sym_lst, repl_lst); + else + return replace_with_symbol(*this, sym_lst, repl_lst); else { // complex numeric re = real(), im = imag(); - ex re_ex = re.is_rational() ? re : replace_with_symbol(re, sym_lst, repl_lst); - ex im_ex = im.is_rational() ? im : replace_with_symbol(im, sym_lst, repl_lst); - return re_ex + im_ex * replace_with_symbol(I, sym_lst, repl_lst); - } + ex re_ex = re.is_rational() ? re : replace_with_symbol(re, sym_lst, repl_lst); + ex im_ex = im.is_rational() ? im : replace_with_symbol(im, sym_lst, repl_lst); + return re_ex + im_ex * replace_with_symbol(I, sym_lst, repl_lst); + } } /* * Helper function for fraction cancellation (returns cancelled fraction n/d) */ - static ex frac_cancel(const ex &n, const ex &d) { ex num = n; diff --git a/ginac/numeric.cpp b/ginac/numeric.cpp index 9f7cb6ec..b5cc40b3 100644 --- a/ginac/numeric.cpp +++ b/ginac/numeric.cpp @@ -275,6 +275,49 @@ void numeric::print(ostream & os, unsigned upper_precedence) const } } +void numeric::printtree(ostream & os, unsigned indent) const +{ + debugmsg("numeric printtree", LOGLEVEL_PRINT); + os << string(indent,' ') << *value + << " (numeric): " + << "hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")" + << ", flags=" << flags << endl; +} + +void numeric::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const +{ + debugmsg("numeric print csrc", LOGLEVEL_PRINT); + ios::fmtflags oldflags = os.flags(); + os.setf(ios::scientific); + if (is_rational() && !is_integer()) { + if (compare(numZERO()) > 0) { + os << "("; + if (type == csrc_types::ctype_cl_N) + os << "cl_F(\"" << numer().evalf() << "\")"; + else + os << numer().to_double(); + } else { + os << "-("; + if (type == csrc_types::ctype_cl_N) + os << "cl_F(\"" << -numer().evalf() << "\")"; + else + os << -numer().to_double(); + } + os << "/"; + if (type == csrc_types::ctype_cl_N) + os << "cl_F(\"" << denom().evalf() << "\")"; + else + os << denom().to_double(); + os << ")"; + } else { + if (type == csrc_types::ctype_cl_N) + os << "cl_F(\"" << evalf() << "\")"; + else + os << to_double(); + } + os.flags(oldflags); +} + bool numeric::info(unsigned inf) const { switch (inf) { @@ -287,9 +330,15 @@ bool numeric::info(unsigned inf) const case info_flags::rational: case info_flags::rational_polynomial: return is_rational(); + case info_flags::crational: + case info_flags::crational_polynomial: + return is_crational(); case info_flags::integer: case info_flags::integer_polynomial: return is_integer(); + case info_flags::cinteger: + case info_flags::cinteger_polynomial: + return is_cinteger(); case info_flags::positive: return is_positive(); case info_flags::negative: @@ -317,7 +366,7 @@ bool numeric::info(unsigned inf) const * currently set. * * @param level ignored, but needed for overriding basic::evalf. - * @return an ex-handle to a numeric. */ + * @return an ex-handle to a numeric. */ ex numeric::evalf(int level) const { // level can safely be discarded for numeric objects. @@ -405,7 +454,7 @@ numeric numeric::mul(numeric const & other) const * @exception overflow_error (division by zero) */ numeric numeric::div(numeric const & other) const { - if (zerop(*other.value)) + if (::zerop(*other.value)) throw (std::overflow_error("division by zero")); return numeric((*value)/(*other.value)); } @@ -413,18 +462,17 @@ numeric numeric::div(numeric const & other) const numeric numeric::power(numeric const & other) const { static const numeric * numONEp=&numONE(); - if (&other==numONEp) { + if (&other==numONEp) return *this; - } - if (zerop(*value) && other.is_real() && minusp(realpart(*other.value))) + if (::zerop(*value) && other.is_real() && ::minusp(realpart(*other.value))) throw (std::overflow_error("division by zero")); - return numeric(expt(*value,*other.value)); + return numeric(::expt(*value,*other.value)); } /** Inverse of a number. */ numeric numeric::inverse(void) const { - return numeric(recip(*value)); // -> CLN + return numeric(::recip(*value)); // -> CLN } numeric const & numeric::add_dyn(numeric const & other) const @@ -453,7 +501,7 @@ numeric const & numeric::mul_dyn(numeric const & other) const numeric const & numeric::div_dyn(numeric const & other) const { - if (zerop(*other.value)) + if (::zerop(*other.value)) throw (std::overflow_error("division by zero")); return static_cast((new numeric((*value)/(*other.value)))-> setflag(status_flags::dynallocated)); @@ -462,29 +510,12 @@ numeric const & numeric::div_dyn(numeric const & other) const numeric const & numeric::power_dyn(numeric const & other) const { static const numeric * numONEp=&numONE(); - if (&other==numONEp) { + if (&other==numONEp) return *this; - } - // The ifs are only a workaround for a bug in CLN. It gets stuck otherwise: - if ( !other.is_integer() && - other.is_rational() && - (*this).is_nonneg_integer() ) { - if ( !zerop(*value) ) { - return static_cast((new numeric(exp(*other.value * log(*value))))-> - setflag(status_flags::dynallocated)); - } else { - if ( !zerop(*other.value) ) { // 0^(n/m) - return static_cast((new numeric(0))-> - setflag(status_flags::dynallocated)); - } else { // raise FPE (0^0 requested) - return static_cast((new numeric(1/(*other.value)))-> - setflag(status_flags::dynallocated)); - } - } - } else { // default -> CLN - return static_cast((new numeric(expt(*value,*other.value)))-> - setflag(status_flags::dynallocated)); - } + if (::zerop(*value) && other.is_real() && ::minusp(realpart(*other.value))) + throw (std::overflow_error("division by zero")); + return static_cast((new numeric(::expt(*value,*other.value)))-> + setflag(status_flags::dynallocated)); } numeric const & numeric::operator=(int i) @@ -526,13 +557,13 @@ int numeric::csgn(void) const { if (is_zero()) return 0; - if (!zerop(realpart(*value))) { - if (plusp(realpart(*value))) + if (!::zerop(realpart(*value))) { + if (::plusp(realpart(*value))) return 1; else return -1; } else { - if (plusp(imagpart(*value))) + if (::plusp(imagpart(*value))) return 1; else return -1; @@ -551,14 +582,14 @@ int numeric::compare(numeric const & other) const // Comparing two real numbers? if (is_real() && other.is_real()) // Yes, just compare them - return cl_compare(The(cl_R)(*value), The(cl_R)(*other.value)); + return ::cl_compare(The(cl_R)(*value), The(cl_R)(*other.value)); else { // No, first compare real parts - cl_signean real_cmp = cl_compare(realpart(*value), realpart(*other.value)); + cl_signean real_cmp = ::cl_compare(realpart(*value), realpart(*other.value)); if (real_cmp) return real_cmp; - return cl_compare(imagpart(*value), imagpart(*other.value)); + return ::cl_compare(imagpart(*value), imagpart(*other.value)); } } @@ -570,59 +601,53 @@ bool numeric::is_equal(numeric const & other) const /** True if object is zero. */ bool numeric::is_zero(void) const { - return zerop(*value); // -> CLN + return ::zerop(*value); // -> CLN } /** True if object is not complex and greater than zero. */ bool numeric::is_positive(void) const { - if (is_real()) { - return plusp(The(cl_R)(*value)); // -> CLN - } + if (is_real()) + return ::plusp(The(cl_R)(*value)); // -> CLN return false; } /** True if object is not complex and less than zero. */ bool numeric::is_negative(void) const { - if (is_real()) { - return minusp(The(cl_R)(*value)); // -> CLN - } + if (is_real()) + return ::minusp(The(cl_R)(*value)); // -> CLN return false; } /** True if object is a non-complex integer. */ bool numeric::is_integer(void) const { - return instanceof(*value, cl_I_ring); // -> CLN + return ::instanceof(*value, cl_I_ring); // -> CLN } /** True if object is an exact integer greater than zero. */ bool numeric::is_pos_integer(void) const { - return (is_integer() && - plusp(The(cl_I)(*value))); // -> CLN + return (is_integer() && ::plusp(The(cl_I)(*value))); // -> CLN } /** True if object is an exact integer greater or equal zero. */ bool numeric::is_nonneg_integer(void) const { - return (is_integer() && - !minusp(The(cl_I)(*value))); // -> CLN + return (is_integer() && !::minusp(The(cl_I)(*value))); // -> CLN } /** True if object is an exact even integer. */ bool numeric::is_even(void) const { - return (is_integer() && - evenp(The(cl_I)(*value))); // -> CLN + return (is_integer() && ::evenp(The(cl_I)(*value))); // -> CLN } /** True if object is an exact odd integer. */ bool numeric::is_odd(void) const { - return (is_integer() && - oddp(The(cl_I)(*value))); // -> CLN + return (is_integer() && ::oddp(The(cl_I)(*value))); // -> CLN } /** Probabilistic primality test. @@ -630,21 +655,20 @@ bool numeric::is_odd(void) const * @return true if object is exact integer and prime. */ bool numeric::is_prime(void) const { - return (is_integer() && - isprobprime(The(cl_I)(*value))); // -> CLN + return (is_integer() && ::isprobprime(The(cl_I)(*value))); // -> CLN } /** True if object is an exact rational number, may even be complex * (denominator may be unity). */ bool numeric::is_rational(void) const { - return instanceof(*value, cl_RA_ring); + return ::instanceof(*value, cl_RA_ring); // -> CLN } /** True if object is a real integer, rational or float (but not complex). */ bool numeric::is_real(void) const { - return instanceof(*value, cl_R_ring); // -> CLN + return ::instanceof(*value, cl_R_ring); // -> CLN } bool numeric::operator==(numeric const & other) const @@ -661,11 +685,11 @@ bool numeric::operator!=(numeric const & other) const * of the form a+b*I, where a and b are integers. */ bool numeric::is_cinteger(void) const { - if (instanceof(*value, cl_I_ring)) + if (::instanceof(*value, cl_I_ring)) return true; else if (!is_real()) { // complex case, handle n+m*I - if (instanceof(realpart(*value), cl_I_ring) && - instanceof(imagpart(*value), cl_I_ring)) + if (::instanceof(realpart(*value), cl_I_ring) && + ::instanceof(imagpart(*value), cl_I_ring)) return true; } return false; @@ -675,11 +699,11 @@ bool numeric::is_cinteger(void) const * (denominator may be unity). */ bool numeric::is_crational(void) const { - if (instanceof(*value, cl_RA_ring)) + if (::instanceof(*value, cl_RA_ring)) return true; else if (!is_real()) { // complex case, handle Q(i): - if (instanceof(realpart(*value), cl_RA_ring) && - instanceof(imagpart(*value), cl_RA_ring)) + if (::instanceof(realpart(*value), cl_RA_ring) && + ::instanceof(imagpart(*value), cl_RA_ring)) return true; } return false; @@ -690,9 +714,8 @@ bool numeric::is_crational(void) const * @exception invalid_argument (complex inequality) */ bool numeric::operator<(numeric const & other) const { - if ( is_real() && other.is_real() ) { + if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) < The(cl_R)(*other.value)); // -> CLN - } throw (std::invalid_argument("numeric::operator<(): complex inequality")); return false; // make compiler shut up } @@ -702,9 +725,8 @@ bool numeric::operator<(numeric const & other) const * @exception invalid_argument (complex inequality) */ bool numeric::operator<=(numeric const & other) const { - if ( is_real() && other.is_real() ) { + if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) <= The(cl_R)(*other.value)); // -> CLN - } throw (std::invalid_argument("numeric::operator<=(): complex inequality")); return false; // make compiler shut up } @@ -714,9 +736,8 @@ bool numeric::operator<=(numeric const & other) const * @exception invalid_argument (complex inequality) */ bool numeric::operator>(numeric const & other) const { - if ( is_real() && other.is_real() ) { + if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) > The(cl_R)(*other.value)); // -> CLN - } throw (std::invalid_argument("numeric::operator>(): complex inequality")); return false; // make compiler shut up } @@ -726,9 +747,8 @@ bool numeric::operator>(numeric const & other) const * @exception invalid_argument (complex inequality) */ bool numeric::operator>=(numeric const & other) const { - if ( is_real() && other.is_real() ) { + if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) >= The(cl_R)(*other.value)); // -> CLN - } throw (std::invalid_argument("numeric::operator>=(): complex inequality")); return false; // make compiler shut up } @@ -738,7 +758,7 @@ bool numeric::operator>=(numeric const & other) const int numeric::to_int(void) const { GINAC_ASSERT(is_integer()); - return cl_I_to_int(The(cl_I)(*value)); + return ::cl_I_to_int(The(cl_I)(*value)); // -> CLN } /** Converts numeric types to machine's double. You should check with is_real() @@ -746,19 +766,19 @@ int numeric::to_int(void) const double numeric::to_double(void) const { GINAC_ASSERT(is_real()); - return cl_double_approx(realpart(*value)); + return ::cl_double_approx(realpart(*value)); // -> CLN } /** Real part of a number. */ numeric numeric::real(void) const { - return numeric(realpart(*value)); // -> CLN + return numeric(::realpart(*value)); // -> CLN } /** Imaginary part of a number. */ numeric numeric::imag(void) const { - return numeric(imagpart(*value)); // -> CLN + return numeric(::imagpart(*value)); // -> CLN } #ifndef SANE_LINKER @@ -784,22 +804,22 @@ numeric numeric::numer(void) const return numeric(*this); } #ifdef SANE_LINKER - else if (instanceof(*value, cl_RA_ring)) { - return numeric(numerator(The(cl_RA)(*value))); + else if (::instanceof(*value, cl_RA_ring)) { + return numeric(::numerator(The(cl_RA)(*value))); } else if (!is_real()) { // complex case, handle Q(i): - cl_R r = realpart(*value); - cl_R i = imagpart(*value); - if (instanceof(r, cl_I_ring) && instanceof(i, cl_I_ring)) + cl_R r = ::realpart(*value); + cl_R i = ::imagpart(*value); + if (::instanceof(r, cl_I_ring) && ::instanceof(i, cl_I_ring)) return numeric(*this); - if (instanceof(r, cl_I_ring) && instanceof(i, cl_RA_ring)) - return numeric(complex(r*denominator(The(cl_RA)(i)), numerator(The(cl_RA)(i)))); - if (instanceof(r, cl_RA_ring) && instanceof(i, cl_I_ring)) - return numeric(complex(numerator(The(cl_RA)(r)), i*denominator(The(cl_RA)(r)))); - if (instanceof(r, cl_RA_ring) && instanceof(i, cl_RA_ring)) { - cl_I s = lcm(denominator(The(cl_RA)(r)), denominator(The(cl_RA)(i))); - return numeric(complex(numerator(The(cl_RA)(r))*(exquo(s,denominator(The(cl_RA)(r)))), - numerator(The(cl_RA)(i))*(exquo(s,denominator(The(cl_RA)(i)))))); + if (::instanceof(r, cl_I_ring) && ::instanceof(i, cl_RA_ring)) + return numeric(complex(r*::denominator(The(cl_RA)(i)), ::numerator(The(cl_RA)(i)))); + if (::instanceof(r, cl_RA_ring) && ::instanceof(i, cl_I_ring)) + return numeric(complex(::numerator(The(cl_RA)(r)), i*::denominator(The(cl_RA)(r)))); + if (::instanceof(r, cl_RA_ring) && ::instanceof(i, cl_RA_ring)) { + cl_I s = lcm(::denominator(The(cl_RA)(r)), ::denominator(The(cl_RA)(i))); + return numeric(complex(::numerator(The(cl_RA)(r))*(exquo(s,::denominator(The(cl_RA)(r)))), + ::numerator(The(cl_RA)(i))*(exquo(s,::denominator(The(cl_RA)(i)))))); } } #else @@ -836,19 +856,19 @@ numeric numeric::denom(void) const } #ifdef SANE_LINKER if (instanceof(*value, cl_RA_ring)) { - return numeric(denominator(The(cl_RA)(*value))); + return numeric(::denominator(The(cl_RA)(*value))); } if (!is_real()) { // complex case, handle Q(i): cl_R r = realpart(*value); cl_R i = imagpart(*value); - if (instanceof(r, cl_I_ring) && instanceof(i, cl_I_ring)) + if (::instanceof(r, cl_I_ring) && ::instanceof(i, cl_I_ring)) return numONE(); - if (instanceof(r, cl_I_ring) && instanceof(i, cl_RA_ring)) - return numeric(denominator(The(cl_RA)(i))); - if (instanceof(r, cl_RA_ring) && instanceof(i, cl_I_ring)) - return numeric(denominator(The(cl_RA)(r))); - if (instanceof(r, cl_RA_ring) && instanceof(i, cl_RA_ring)) - return numeric(lcm(denominator(The(cl_RA)(r)), denominator(The(cl_RA)(i)))); + if (::instanceof(r, cl_I_ring) && ::instanceof(i, cl_RA_ring)) + return numeric(::denominator(The(cl_RA)(i))); + if (::instanceof(r, cl_RA_ring) && ::instanceof(i, cl_I_ring)) + return numeric(::denominator(The(cl_RA)(r))); + if (::instanceof(r, cl_RA_ring) && ::instanceof(i, cl_RA_ring)) + return numeric(lcm(::denominator(The(cl_RA)(r)), ::denominator(The(cl_RA)(i)))); } #else if (instanceof(*value, cl_RA_ring)) { @@ -879,11 +899,10 @@ numeric numeric::denom(void) const * in two's complement if it is an integer, 0 otherwise. */ int numeric::int_length(void) const { - if (is_integer()) { - return integer_length(The(cl_I)(*value)); // -> CLN - } else { + if (is_integer()) + return ::integer_length(The(cl_I)(*value)); // -> CLN + else return 0; - } } @@ -1090,10 +1109,19 @@ numeric atanh(numeric const & x) * integer arguments. */ numeric zeta(numeric const & x) { - if (x.is_integer()) - return ::cl_zeta(x.to_int()); // -> CLN - else - clog << "zeta(): Does anybody know good way to calculate this numerically?" << endl; + // A dirty hack to allow for things like zeta(3.0), since CLN currently + // only knows about integer arguments and zeta(3).evalf() automatically + // cascades down to zeta(3.0).evalf(). The trick is to rely on 3.0-3 + // being an exact zero for CLN, which can be tested and then we can just + // pass the number casted to an int: + if (x.is_real()) { + int aux = (int)(::cl_double_approx(realpart(*x.value))); + if (zerop(*x.value-aux)) + return ::cl_zeta(aux); // -> CLN + } + clog << "zeta(" << x + << "): Does anybody know good way to calculate this numerically?" + << endl; return numeric(0); } @@ -1101,7 +1129,9 @@ numeric zeta(numeric const & x) * This is only a stub! */ numeric gamma(numeric const & x) { - clog << "gamma(): Does anybody know good way to calculate this numerically?" << endl; + clog << "gamma(" << x + << "): Does anybody know good way to calculate this numerically?" + << endl; return numeric(0); } @@ -1109,7 +1139,9 @@ numeric gamma(numeric const & x) * This is only a stub! */ numeric psi(numeric const & x) { - clog << "psi(): Does anybody know good way to calculate this numerically?" << endl; + clog << "psi(" << x + << "): Does anybody know good way to calculate this numerically?" + << endl; return numeric(0); } @@ -1117,7 +1149,9 @@ numeric psi(numeric const & x) * This is only a stub! */ numeric psi(numeric const & n, numeric const & x) { - clog << "psi(): Does anybody know good way to calculate this numerically?" << endl; + clog << "psi(" << n << "," << x + << "): Does anybody know good way to calculate this numerically?" + << endl; return numeric(0); } @@ -1126,10 +1160,8 @@ numeric psi(numeric const & n, numeric const & x) * @exception range_error (argument must be integer >= 0) */ numeric factorial(numeric const & nn) { - if ( !nn.is_nonneg_integer() ) { + if (!nn.is_nonneg_integer()) throw (std::range_error("numeric::factorial(): argument must be integer >= 0")); - } - return numeric(::factorial(nn.to_int())); // -> CLN } @@ -1213,7 +1245,7 @@ numeric binomial(numeric const & n, numeric const & k) return numZERO(); } else { return numMINUSONE().power(k)*binomial(k-n-numONE(),k); - } + } } // should really be gamma(n+1)/(gamma(r+1)/gamma(n-r+1) or a suitable limit @@ -1277,12 +1309,10 @@ numeric abs(numeric const & x) * integer, 0 otherwise. */ numeric mod(numeric const & a, numeric const & b) { - if (a.is_integer() && b.is_integer()) { + if (a.is_integer() && b.is_integer()) return ::mod(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN - } - else { + else return numZERO(); // Throw? - } } /** Modulus (in symmetric representation). @@ -1294,9 +1324,8 @@ numeric smod(numeric const & a, numeric const & b) if (a.is_integer() && b.is_integer()) { cl_I b2 = The(cl_I)(ceiling1(The(cl_I)(*b.value) / 2)) - 1; return ::mod(The(cl_I)(*a.value) + b2, The(cl_I)(*b.value)) - b2; - } else { + } else return numZERO(); // Throw? - } } /** Numeric integer remainder. @@ -1307,12 +1336,10 @@ numeric smod(numeric const & a, numeric const & b) * @return remainder of a/b if both are integer, 0 otherwise. */ numeric irem(numeric const & a, numeric const & b) { - if (a.is_integer() && b.is_integer()) { + if (a.is_integer() && b.is_integer()) return ::rem(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN - } - else { + else return numZERO(); // Throw? - } } /** Numeric integer remainder. @@ -1341,11 +1368,10 @@ numeric irem(numeric const & a, numeric const & b, numeric & q) * @return truncated quotient of a/b if both are integer, 0 otherwise. */ numeric iquo(numeric const & a, numeric const & b) { - if (a.is_integer() && b.is_integer()) { + if (a.is_integer() && b.is_integer()) return truncate1(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN - } else { + else return numZERO(); // Throw? - } } /** Numeric integer quotient. @@ -1382,12 +1408,12 @@ numeric sqrt(numeric const & z) /** Integer numeric square root. */ numeric isqrt(numeric const & x) { - if (x.is_integer()) { - cl_I root; - ::isqrt(The(cl_I)(*x.value), &root); // -> CLN - return root; - } else - return numZERO(); // Throw? + if (x.is_integer()) { + cl_I root; + ::isqrt(The(cl_I)(*x.value), &root); // -> CLN + return root; + } else + return numZERO(); // Throw? } /** Greatest Common Divisor. @@ -1397,7 +1423,7 @@ numeric isqrt(numeric const & x) numeric gcd(numeric const & a, numeric const & b) { if (a.is_integer() && b.is_integer()) - return ::gcd(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN + return ::gcd(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN else return numONE(); } @@ -1409,7 +1435,7 @@ numeric gcd(numeric const & a, numeric const & b) numeric lcm(numeric const & a, numeric const & b) { if (a.is_integer() && b.is_integer()) - return ::lcm(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN + return ::lcm(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN else return *a.value * *b.value; } diff --git a/ginac/numeric.h b/ginac/numeric.h index 168651d4..68c227d1 100644 --- a/ginac/numeric.h +++ b/ginac/numeric.h @@ -79,6 +79,7 @@ class numeric : public basic friend numeric asinh(numeric const & x); friend numeric acosh(numeric const & x); friend numeric atanh(numeric const & x); + friend numeric zeta(numeric const & x); friend numeric bernoulli(numeric const & n); friend numeric abs(numeric const & x); friend numeric mod(numeric const & a, numeric const & b); @@ -125,9 +126,9 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; + void print(ostream & os, unsigned precedence=0) const; void printraw(ostream & os) const; void printtree(ostream & os, unsigned indent) const; - void print(ostream & os, unsigned precedence=0) const; void printcsrc(ostream & os, unsigned type, unsigned precedence=0) const; bool info(unsigned inf) const; ex evalf(int level=0) const; @@ -309,6 +310,12 @@ inline bool is_rational(numeric const & x) inline bool is_real(numeric const & x) { return x.is_real(); } +inline bool is_cinteger(numeric const & x) +{ return x.is_cinteger(); } + +inline bool is_crational(numeric const & x) +{ return x.is_crational(); } + inline numeric real(numeric const & x) { return x.real(); } @@ -330,7 +337,7 @@ ex CatalanEvalf(void); // utility functions inline const numeric &ex_to_numeric(const ex &e) { - return static_cast(*e.bp); + return static_cast(*e.bp); } #ifndef NO_GINAC_NAMESPACE diff --git a/ginac/operators.h b/ginac/operators.h index 01210fc1..a60da5d6 100644 --- a/ginac/operators.h +++ b/ginac/operators.h @@ -40,24 +40,6 @@ ex operator*(ex const & lh, ex const & rh); ex operator/(ex const & lh, ex const & rh); ex operator%(ex const & lh, ex const & rh); // non-commutative multiplication -/* - -// binary arithmetic operators ex with numeric -ex operator+(ex const & lh, numeric const & rh); -ex operator-(ex const & lh, numeric const & rh); -ex operator*(ex const & lh, numeric const & rh); -ex operator/(ex const & lh, numeric const & rh); -ex operator%(ex const & lh, numeric const & rh); // non-commutative multiplication - -// binary arithmetic operators numeric with ex -ex operator+(numeric const & lh, ex const & rh); -ex operator-(numeric const & lh, ex const & rh); -ex operator*(numeric const & lh, ex const & rh); -ex operator/(numeric const & lh, ex const & rh); -ex operator%(numeric const & lh, ex const & rh); // non-commutative multiplication - -*/ - // binary arithmetic operators numeric with numeric numeric operator+(numeric const & lh, numeric const & rh); numeric operator-(numeric const & lh, numeric const & rh); @@ -71,17 +53,6 @@ ex const & operator*=(ex & lh, ex const & rh); ex const & operator/=(ex & lh, ex const & rh); ex const & operator%=(ex & lh, ex const & rh); // non-commutative multiplication -/* - -// binary arithmetic assignment operators with numeric -ex const & operator+=(ex & lh, numeric const & rh); -ex const & operator-=(ex & lh, numeric const & rh); -ex const & operator*=(ex & lh, numeric const & rh); -ex const & operator/=(ex & lh, numeric const & rh); -ex const & operator%=(ex & lh, numeric const & rh); // non-commutative multiplication - -*/ - // binary arithmetic assignment operators with numeric numeric const & operator+=(numeric & lh, numeric const & rh); numeric const & operator-=(numeric & lh, numeric const & rh); @@ -107,26 +78,6 @@ relational operator<=(ex const & lh, ex const & rh); relational operator>(ex const & lh, ex const & rh); relational operator>=(ex const & lh, ex const & rh); -/* - -// binary relational operators ex with numeric -relational operator==(ex const & lh, numeric const & rh); -relational operator!=(ex const & lh, numeric const & rh); -relational operator<(ex const & lh, numeric const & rh); -relational operator<=(ex const & lh, numeric const & rh); -relational operator>(ex const & lh, numeric const & rh); -relational operator>=(ex const & lh, numeric const & rh); - -// binary relational operators numeric with ex -relational operator==(numeric const & lh, ex const & rh); -relational operator!=(numeric const & lh, ex const & rh); -relational operator<(numeric const & lh, ex const & rh); -relational operator<=(numeric const & lh, ex const & rh); -relational operator>(numeric const & lh, ex const & rh); -relational operator>=(numeric const & lh, ex const & rh); - -*/ - // input/output stream operators ostream & operator<<(ostream & os, ex const & e); istream & operator>>(istream & is, ex & e); diff --git a/ginac/power.cpp b/ginac/power.cpp index 64dd1fdc..1c9a342d 100644 --- a/ginac/power.cpp +++ b/ginac/power.cpp @@ -116,9 +116,110 @@ basic * power::duplicate() const return new power(*this); } +void power::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("power print",LOGLEVEL_PRINT); + if (precedence<=upper_precedence) os << "("; + basis.print(os,precedence); + os << "^"; + exponent.print(os,precedence); + if (precedence<=upper_precedence) os << ")"; +} + +void power::printraw(ostream & os) const +{ + debugmsg("power printraw",LOGLEVEL_PRINT); + + os << "power("; + basis.printraw(os); + os << ","; + exponent.printraw(os); + os << ",hash=" << hashvalue << ",flags=" << flags << ")"; +} + +void power::printtree(ostream & os, unsigned indent) const +{ + debugmsg("power printtree",LOGLEVEL_PRINT); + + os << string(indent,' ') << "power: " + << "hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")" + << ", flags=" << flags << endl; + basis.printtree(os,indent+delta_indent); + exponent.printtree(os,indent+delta_indent); +} + +static void print_sym_pow(ostream & os, unsigned type, const symbol &x, int exp) +{ + // Optimal output of integer powers of symbols to aid compiler CSE + if (exp == 1) { + x.printcsrc(os, type, 0); + } else if (exp == 2) { + x.printcsrc(os, type, 0); + os << "*"; + x.printcsrc(os, type, 0); + } else if (exp & 1) { + x.printcsrc(os, 0); + os << "*"; + print_sym_pow(os, type, x, exp-1); + } else { + os << "("; + print_sym_pow(os, type, x, exp >> 1); + os << ")*("; + print_sym_pow(os, type, x, exp >> 1); + os << ")"; + } +} + +void power::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const +{ + debugmsg("power print csrc", LOGLEVEL_PRINT); + + // Integer powers of symbols are printed in a special, optimized way + if (exponent.info(info_flags::integer) && + (is_ex_exactly_of_type(basis, symbol) || + is_ex_exactly_of_type(basis, constant))) { + int exp = ex_to_numeric(exponent).to_int(); + if (exp > 0) + os << "("; + else { + exp = -exp; + if (type == csrc_types::ctype_cl_N) + os << "recip("; + else + os << "1.0/("; + } + print_sym_pow(os, type, static_cast(*basis.bp), exp); + os << ")"; + + // ^-1 is printed as "1.0/" or with the recip() function of CLN + } else if (exponent.compare(numMINUSONE()) == 0) { + if (type == csrc_types::ctype_cl_N) + os << "recip("; + else + os << "1.0/("; + basis.bp->printcsrc(os, type, 0); + os << ")"; + + // Otherwise, use the pow() or expt() (CLN) functions + } else { + if (type == csrc_types::ctype_cl_N) + os << "expt("; + else + os << "pow("; + basis.bp->printcsrc(os, type, 0); + os << ","; + exponent.bp->printcsrc(os, type, 0); + os << ")"; + } +} + bool power::info(unsigned inf) const { - if (inf==info_flags::polynomial || inf==info_flags::integer_polynomial || inf==info_flags::rational_polynomial) { + if (inf==info_flags::polynomial || + inf==info_flags::integer_polynomial || + inf==info_flags::cinteger_polynomial || + inf==info_flags::rational_polynomial || + inf==info_flags::crational_polynomial) { return exponent.info(info_flags::nonnegint); } else if (inf==info_flags::rational_function) { return exponent.info(info_flags::integer); @@ -239,17 +340,17 @@ ex power::eval(int level) const if (basis_is_numerical && exponent_is_numerical) { // ^(c1,c2) -> c1^c2 (c1, c2 numeric(), // except if c1,c2 are rational, but c1^c2 is not) - bool basis_is_rational = num_basis->is_rational(); - bool exponent_is_rational = num_exponent->is_rational(); + bool basis_is_crational = num_basis->is_crational(); + bool exponent_is_crational = num_exponent->is_crational(); numeric res = (*num_basis).power(*num_exponent); - if ((!basis_is_rational || !exponent_is_rational) - || res.is_rational()) { + if ((!basis_is_crational || !exponent_is_crational) + || res.is_crational()) { return res; } GINAC_ASSERT(!num_exponent->is_integer()); // has been handled by now // ^(c1,n/m) -> *(c1^q,c1^(n/m-q)), 0<(n/m-h)<1, q integer - if (basis_is_rational && exponent_is_rational + if (basis_is_crational && exponent_is_crational && num_exponent->is_real() && !num_exponent->is_integer()) { numeric r, q, n, m; diff --git a/ginac/power.h b/ginac/power.h index 25f4f8d7..b1136c35 100644 --- a/ginac/power.h +++ b/ginac/power.h @@ -59,9 +59,9 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; + void print(ostream & os, unsigned upper_precedence=0) const; void printraw(ostream & os) const; void printtree(ostream & os, unsigned indent) const; - void print(ostream & os, unsigned upper_precedence=0) const; void printcsrc(ostream & os, unsigned type, unsigned upper_precedence=0) const; bool info(unsigned inf) const; int nops() const; diff --git a/ginac/print.cpp b/ginac/print.cpp deleted file mode 100644 index c1e9841c..00000000 --- a/ginac/print.cpp +++ /dev/null @@ -1,266 +0,0 @@ -/** @file print.cpp - * - * The methods .print() are responsible for the nice default-output of - * objects. All related helper-functions go in here as well. */ - -/* - * GiNaC Copyright (C) 1999 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 - -#include "basic.h" -#include "ex.h" -#include "add.h" -#include "constant.h" -#include "expairseq.h" -#include "fail.h" -#include "indexed.h" -#include "inifcns.h" -#include "matrix.h" -#include "mul.h" -#include "ncmul.h" -#include "numeric.h" -#include "power.h" -#include "relational.h" -#include "series.h" -#include "symbol.h" -#include "debugmsg.h" - -#ifndef NO_GINAC_NAMESPACE -namespace GiNaC { -#endif // ndef NO_GINAC_NAMESPACE - -void ex::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("ex print",LOGLEVEL_PRINT); - GINAC_ASSERT(bp!=0); - bp->print(os,upper_precedence); -} - -void ex::dbgprint(void) const -{ - debugmsg("ex dbgprint",LOGLEVEL_PRINT); - GINAC_ASSERT(bp!=0); - bp->dbgprint(); -} - -void basic::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("basic print",LOGLEVEL_PRINT); - os << "[basic object]"; -} - -void basic::dbgprint(void) const -{ - print(cerr); - cerr << endl; -} - -void symbol::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("symbol print",LOGLEVEL_PRINT); - os << name; -} - -void constant::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("constant print",LOGLEVEL_PRINT); - os << name; -} - -void power::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("power print",LOGLEVEL_PRINT); - if (precedence<=upper_precedence) os << "("; - basis.print(os,precedence); - os << "^"; - exponent.print(os,precedence); - if (precedence<=upper_precedence) os << ")"; -} - -void fail::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("fail print",LOGLEVEL_PRINT); - os << "FAIL"; -} - -void expairseq::printpair(ostream & os, expair const & p, unsigned upper_precedence) const -{ - os << "[["; - p.rest.bp->print(os,precedence); - os << ","; - p.coeff.bp->print(os,precedence); - os << "]]"; -} - -void expairseq::printseq(ostream & os, char delim, unsigned this_precedence, - unsigned upper_precedence) const -{ - if (this_precedence<=upper_precedence) os << "("; - epvector::const_iterator it,it_last; - it_last=seq.end(); - --it_last; - for (it=seq.begin(); it!=it_last; ++it) { - printpair(os,*it,this_precedence); - os << delim; - } - printpair(os,*it,this_precedence); - if (!overall_coeff.is_equal(default_overall_coeff())) { - os << delim << overall_coeff; - } - if (this_precedence<=upper_precedence) os << ")"; -} - -void expairseq::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("expairseq print",LOGLEVEL_PRINT); - os << "[["; - printseq(os,',',precedence,upper_precedence); - os << "]]"; -} - -void add::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("add print",LOGLEVEL_PRINT); - if (precedence<=upper_precedence) os << "("; - numeric coeff; - bool first=true; - for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) { - coeff = ex_to_numeric(cit->coeff); - if (!first) { - if (coeff.csgn()==-1) os << '-'; else os << '+'; - } else { - if (coeff.csgn()==-1) os << '-'; - first=false; - } - if (!coeff.is_equal(numONE()) && - !coeff.is_equal(numMINUSONE())) { - if (coeff.csgn()==-1) - (numMINUSONE()*coeff).print(os, precedence); - else - coeff.print(os, precedence); - os << '*'; - } - os << cit->rest; - } - // print the overall numeric coefficient, if present: - if (!overall_coeff.is_zero()) { - if (overall_coeff > 0) os << '+'; - os << overall_coeff; - } - if (precedence<=upper_precedence) os << ")"; -} - -void mul::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("mul print",LOGLEVEL_PRINT); - if (precedence<=upper_precedence) os << "("; - bool first=true; - // first print the overall numeric coefficient: - if (ex_to_numeric(overall_coeff).csgn()==-1) os << '-'; - if (!overall_coeff.is_equal(exONE()) && - !overall_coeff.is_equal(exMINUSONE())) { - if (ex_to_numeric(overall_coeff).csgn()==-1) - (numMINUSONE()*overall_coeff).print(os, precedence); - else - overall_coeff.print(os, precedence); - os << '*'; - } - // then proceed with the remaining factors: - for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) { - if (!first) { - os << '*'; - } else { - first=false; - } - recombine_pair_to_ex(*cit).print(os,precedence); - } - if (precedence<=upper_precedence) os << ")"; -} - -void ncmul::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("ncmul print",LOGLEVEL_PRINT); - printseq(os,'(','%',')',precedence,upper_precedence); -} - -/*void function::print(ostream & os, unsigned upper_precedence) const - *{ - * debugmsg("function print",LOGLEVEL_PRINT); - * os << name; - * printseq(os,'(',',',')',exprseq::precedence,function::precedence); - *}*/ - -void series::print(ostream &os, unsigned upper_precedence) const -{ - debugmsg("symbol print", LOGLEVEL_PRINT); - convert_to_poly().print(os, upper_precedence); -} - -void relational::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("relational print",LOGLEVEL_PRINT); - if (precedence<=upper_precedence) os << "("; - lh.print(os,precedence); - switch (o) { - case equal: - os << "=="; - break; - case not_equal: - os << "!="; - break; - case less: - os << "<"; - break; - case less_or_equal: - os << "<="; - break; - case greater: - os << ">"; - break; - case greater_or_equal: - os << ">="; - break; - default: - os << "(INVALID RELATIONAL OPERATOR)"; - } - rh.print(os,precedence); - if (precedence<=upper_precedence) os << ")"; -} - -void matrix::print(ostream & os, unsigned upper_precedence) const -{ - debugmsg("matrix print",LOGLEVEL_PRINT); - os << "[[ "; - for (int r=0; r - -#include "basic.h" -#include "ex.h" -#include "add.h" -#include "constant.h" -#include "expairseq.h" -#include "indexed.h" -#include "inifcns.h" -#include "mul.h" -#include "ncmul.h" -#include "numeric.h" -#include "power.h" -#include "relational.h" -#include "series.h" -#include "symbol.h" -#include "debugmsg.h" - -#ifndef NO_GINAC_NAMESPACE -namespace GiNaC { -#endif // ndef NO_GINAC_NAMESPACE - -/** Print expression as a C++ statement. The output looks like - * " = ;". The "type" parameter has an effect - * on how number literals are printed. - * - * @param os output stream - * @param type variable type (one of the csrc_types) - * @param var_name variable name to be printed */ -void ex::printcsrc(ostream & os, unsigned type, const char *var_name) const -{ - debugmsg("ex print csrc", LOGLEVEL_PRINT); - GINAC_ASSERT(bp!=0); - switch (type) { - case csrc_types::ctype_float: - os << "float "; - break; - case csrc_types::ctype_double: - os << "double "; - break; - case csrc_types::ctype_cl_N: - os << "cl_N "; - break; - } - os << var_name << " = "; - bp->printcsrc(os, type, 0); - os << ";\n"; -} - -void basic::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("basic print csrc", LOGLEVEL_PRINT); -} - -void numeric::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("numeric print csrc", LOGLEVEL_PRINT); - ios::fmtflags oldflags = os.flags(); - os.setf(ios::scientific); - if (is_rational() && !is_integer()) { - if (compare(numZERO()) > 0) { - os << "("; - if (type == csrc_types::ctype_cl_N) - os << "cl_F(\"" << numer().evalf() << "\")"; - else - os << numer().to_double(); - } else { - os << "-("; - if (type == csrc_types::ctype_cl_N) - os << "cl_F(\"" << -numer().evalf() << "\")"; - else - os << -numer().to_double(); - } - os << "/"; - if (type == csrc_types::ctype_cl_N) - os << "cl_F(\"" << denom().evalf() << "\")"; - else - os << denom().to_double(); - os << ")"; - } else { - if (type == csrc_types::ctype_cl_N) - os << "cl_F(\"" << evalf() << "\")"; - else - os << to_double(); - } - os.flags(oldflags); -} - -void symbol::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("symbol print csrc", LOGLEVEL_PRINT); - os << name; -} - -void constant::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("constant print csrc",LOGLEVEL_PRINT); - os << name; -} - -static void print_sym_pow(ostream & os, unsigned type, const symbol &x, int exp) -{ - // Optimal output of integer powers of symbols to aid compiler CSE - if (exp == 1) { - x.printcsrc(os, type, 0); - } else if (exp == 2) { - x.printcsrc(os, type, 0); - os << "*"; - x.printcsrc(os, type, 0); - } else if (exp & 1) { - x.printcsrc(os, 0); - os << "*"; - print_sym_pow(os, type, x, exp-1); - } else { - os << "("; - print_sym_pow(os, type, x, exp >> 1); - os << ")*("; - print_sym_pow(os, type, x, exp >> 1); - os << ")"; - } -} - -void power::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("power print csrc", LOGLEVEL_PRINT); - - // Integer powers of symbols are printed in a special, optimized way - if (exponent.info(info_flags::integer) && - (is_ex_exactly_of_type(basis, symbol) || - is_ex_exactly_of_type(basis, constant))) { - int exp = ex_to_numeric(exponent).to_int(); - if (exp > 0) - os << "("; - else { - exp = -exp; - if (type == csrc_types::ctype_cl_N) - os << "recip("; - else - os << "1.0/("; - } - print_sym_pow(os, type, static_cast(*basis.bp), exp); - os << ")"; - - // ^-1 is printed as "1.0/" or with the recip() function of CLN - } else if (exponent.compare(numMINUSONE()) == 0) { - if (type == csrc_types::ctype_cl_N) - os << "recip("; - else - os << "1.0/("; - basis.bp->printcsrc(os, type, 0); - os << ")"; - - // Otherwise, use the pow() or expt() (CLN) functions - } else { - if (type == csrc_types::ctype_cl_N) - os << "expt("; - else - os << "pow("; - basis.bp->printcsrc(os, type, 0); - os << ","; - exponent.bp->printcsrc(os, type, 0); - os << ")"; - } -} - -void add::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("add print csrc", LOGLEVEL_PRINT); - if (precedence <= upper_precedence) - os << "("; - - // Print arguments, separated by "+" - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); - while (it != itend) { - - // If the coefficient is -1, it is replaced by a single minus sign - if (it->coeff.compare(numONE()) == 0) { - it->rest.bp->printcsrc(os, type, precedence); - } else if (it->coeff.compare(numMINUSONE()) == 0) { - os << "-"; - it->rest.bp->printcsrc(os, type, precedence); - } else if (ex_to_numeric(it->coeff).numer().compare(numONE()) == 0) { - it->rest.bp->printcsrc(os, type, precedence); - os << "/"; - ex_to_numeric(it->coeff).denom().printcsrc(os, type, precedence); - } else if (ex_to_numeric(it->coeff).numer().compare(numMINUSONE()) == 0) { - os << "-"; - it->rest.bp->printcsrc(os, type, precedence); - os << "/"; - ex_to_numeric(it->coeff).denom().printcsrc(os, type, precedence); - } else { - it->coeff.bp->printcsrc(os, type, precedence); - os << "*"; - it->rest.bp->printcsrc(os, type, precedence); - } - - // Separator is "+", except it the following expression would have a leading minus sign - it++; - if (it != itend && !(it->coeff.compare(numZERO()) < 0 || (it->coeff.compare(numONE()) == 0 && is_ex_exactly_of_type(it->rest, numeric) && it->rest.compare(numZERO()) < 0))) - os << "+"; - } - - if (!overall_coeff.is_equal(exZERO())) { - if (overall_coeff > 0) os << '+'; - overall_coeff.bp->printcsrc(os,type,precedence); - } - - if (precedence <= upper_precedence) - os << ")"; -} - -void mul::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("mul print csrc", LOGLEVEL_PRINT); - if (precedence <= upper_precedence) - os << "("; - - if (!overall_coeff.is_equal(exONE())) { - overall_coeff.bp->printcsrc(os,type,precedence); - os << "*"; - } - - // Print arguments, separated by "*" or "/" - epvector::const_iterator it = seq.begin(); - epvector::const_iterator itend = seq.end(); - while (it != itend) { - - // If the first argument is a negative integer power, it gets printed as "1.0/" - if (it == seq.begin() && ex_to_numeric(it->coeff).is_integer() && it->coeff.compare(numZERO()) < 0) { - if (type == csrc_types::ctype_cl_N) - os << "recip("; - else - os << "1.0/"; - } - - // If the exponent is 1 or -1, it is left out - if (it->coeff.compare(exONE()) == 0 || it->coeff.compare(numMINUSONE()) == 0) - it->rest.bp->printcsrc(os, type, precedence); - else - // outer parens around ex needed for broken gcc-2.95 parser: - (ex(power(it->rest, abs(ex_to_numeric(it->coeff))))).bp->printcsrc(os, type, upper_precedence); - - // Separator is "/" for negative integer powers, "*" otherwise - it++; - if (it != itend) { - if (ex_to_numeric(it->coeff).is_integer() && it->coeff.compare(numZERO()) < 0) - os << "/"; - else - os << "*"; - } - } - if (precedence <= upper_precedence) - os << ")"; -} - -void ncmul::printcsrc(ostream & os, unsigned upper_precedence) const -{ - debugmsg("ncmul print csrc",LOGLEVEL_PRINT); - exvector::const_iterator it; - exvector::const_iterator itend = seq.end()-1; - os << "ncmul("; - for (it=seq.begin(); it!=itend; ++it) { - (*it).bp->printcsrc(os,precedence); - os << ","; - } - (*it).bp->printcsrc(os,precedence); - os << ")"; -} - -void relational::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const -{ - debugmsg("relational print csrc", LOGLEVEL_PRINT); - if (precedence<=upper_precedence) - os << "("; - - // Print left-hand expression - lh.bp->printcsrc(os, type, precedence); - - // Print relational operator - switch (o) { - case equal: - os << "=="; - break; - case not_equal: - os << "!="; - break; - case less: - os << "<"; - break; - case less_or_equal: - os << "<="; - break; - case greater: - os << ">"; - break; - case greater_or_equal: - os << ">="; - break; - default: - os << "(INVALID RELATIONAL OPERATOR)"; - break; - } - - // Print right-hand operator - rh.bp->printcsrc(os, type, precedence); - - if (precedence <= upper_precedence) - os << ")"; -} - -#ifndef NO_GINAC_NAMESPACE -} // namespace GiNaC -#endif // ndef NO_GINAC_NAMESPACE diff --git a/ginac/printraw.cpp b/ginac/printraw.cpp deleted file mode 100644 index 331622a7..00000000 --- a/ginac/printraw.cpp +++ /dev/null @@ -1,232 +0,0 @@ -/** @file printraw.cpp - * - * print in ugly raw format, so brave developers can have a look at the - * underlying structure. */ - -/* - * GiNaC Copyright (C) 1999 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 - */ - -/* We are cheating here, because we don't want to include the underlying - * bignum package's headers again, so in this file we omit the definition of - * void numeric::printraw(ostream & os) const; */ - -#include - -#include "basic.h" -#include "ex.h" -#include "add.h" -#include "constant.h" -#include "expairseq.h" -#include "fail.h" -#include "indexed.h" -#include "inifcns.h" -#include "matrix.h" -#include "mul.h" -#include "ncmul.h" -#include "numeric.h" -#include "power.h" -#include "relational.h" -#include "series.h" -#include "symbol.h" -#include "debugmsg.h" - -#ifndef NO_GINAC_NAMESPACE -namespace GiNaC { -#endif // ndef NO_GINAC_NAMESPACE - -void ex::printraw(ostream & os) const -{ - debugmsg("ex printraw",LOGLEVEL_PRINT); - GINAC_ASSERT(bp!=0); - os << "ex("; - bp->printraw(os); - os << ")"; -} - -void basic::printraw(ostream & os) const -{ - debugmsg("basic printraw",LOGLEVEL_PRINT); - os << "[basic object]"; -} - -void symbol::printraw(ostream & os) const -{ - debugmsg("symbol printraw",LOGLEVEL_PRINT); - os << "symbol(" << "name=" << name << ",serial=" << serial - << ",hash=" << hashvalue << ",flags=" << flags << ")"; -} - -void constant::printraw(ostream & os) const -{ - debugmsg("constant printraw",LOGLEVEL_PRINT); - os << "constant(" << name << ")"; -} - -void power::printraw(ostream & os) const -{ - debugmsg("power printraw",LOGLEVEL_PRINT); - - os << "power("; - basis.printraw(os); - os << ","; - exponent.printraw(os); - os << ",hash=" << hashvalue << ",flags=" << flags << ")"; -} - -void fail::printraw(ostream & os) const -{ - debugmsg("fail printraw",LOGLEVEL_PRINT); - os << "FAIL"; -} - -void expairseq::printraw(ostream & os) const -{ - debugmsg("expairseq printraw",LOGLEVEL_PRINT); - - os << "expairseq("; - for (epvector::const_iterator cit=seq.begin(); cit!=seq.end(); ++cit) { - os << "("; - (*cit).rest.printraw(os); - os << ","; - (*cit).coeff.printraw(os); - os << "),"; - } - os << ")"; -} - -void add::printraw(ostream & os) const -{ - debugmsg("add printraw",LOGLEVEL_PRINT); - - os << "+("; - for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { - os << "("; - (*it).rest.bp->printraw(os); - os << ","; - (*it).coeff.bp->printraw(os); - os << "),"; - } - os << ",hash=" << hashvalue << ",flags=" << flags; - os << ")"; -} - -void mul::printraw(ostream & os) const -{ - debugmsg("mul printraw",LOGLEVEL_PRINT); - - os << "*("; - for (epvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { - os << "("; - (*it).rest.bp->printraw(os); - os << ","; - (*it).coeff.bp->printraw(os); - os << "),"; - } - os << ",hash=" << hashvalue << ",flags=" << flags; - os << ")"; -} - -void ncmul::printraw(ostream & os) const -{ - debugmsg("ncmul printraw",LOGLEVEL_PRINT); - - os << "%("; - for (exvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { - (*it).bp->printraw(os); - os << ","; - } - os << ",hash=" << hashvalue << ",flags=" << flags; - os << ")"; -} - -/*void function::printraw(ostream & os) const - *{ - * debugmsg("function printraw",LOGLEVEL_PRINT); - * - * os << "function." << name << "("; - * for (exvector::const_iterator it=seq.begin(); it!=seq.end(); ++it) { - * (*it).bp->print(os); - * os << ","; - * } - * os << ")"; - *}*/ - -void series::printraw(ostream &os) const -{ - debugmsg("symbol printraw", LOGLEVEL_PRINT); - os << "series(" << var << ";" << point << ";"; - for (epvector::const_iterator i=seq.begin(); i!=seq.end(); i++) { - os << "(" << (*i).rest << "," << (*i).coeff << "),"; - } - os << ")"; -} - -void relational::printraw(ostream & os) const -{ - debugmsg("relational printraw",LOGLEVEL_PRINT); - os << "RELATIONAL("; - lh.printraw(os); - os << ","; - rh.printraw(os); - os << ","; - switch (o) { - case equal: - os << "=="; - break; - case not_equal: - os << "!="; - break; - case less: - os << "<"; - break; - case less_or_equal: - os << "<="; - break; - case greater: - os << ">"; - break; - case greater_or_equal: - os << ">="; - break; - default: - os << "(INVALID RELATIONAL OPERATOR)"; - } - os << ")"; -} - -void matrix::printraw(ostream & os) const -{ - debugmsg("matrix printraw",LOGLEVEL_PRINT); - os << "matrix(" << row << "," << col <<","; - for (int r=0; r -#include - -#include "basic.h" -#include "ex.h" -#include "add.h" -#include "constant.h" -#include "expairseq.h" -#include "indexed.h" -#include "inifcns.h" -#include "mul.h" -#include "ncmul.h" -#include "numeric.h" -#include "power.h" -#include "relational.h" -#include "series.h" -#include "symbol.h" -#include "debugmsg.h" - -#ifndef NO_GINAC_NAMESPACE -namespace GiNaC { -#endif // ndef NO_GINAC_NAMESPACE - -void ex::printtree(ostream & os, unsigned indent) const -{ - debugmsg("ex printtree",LOGLEVEL_PRINT); - GINAC_ASSERT(bp!=0); - // os << "refcount=" << bp->refcount << " "; - bp->printtree(os,indent); -} - -void ex::dbgprinttree(void) const -{ - debugmsg("ex dbgprinttree",LOGLEVEL_PRINT); - GINAC_ASSERT(bp!=0); - bp->dbgprinttree(); -} - -void basic::printtree(ostream & os, unsigned indent) const -{ - debugmsg("basic printtree",LOGLEVEL_PRINT); - os << string(indent,' ') << "type=" << typeid(*this).name() - << ", hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")" - << ", flags=" << flags - << ", nops=" << nops() << endl; - for (int i=0; i0) { - os << string(indent+delta_indent,' ') - << "bin " << i << " with entries "; - for (epplist::const_iterator it=hashtab[i].begin(); - it!=hashtab[i].end(); ++it) { - os << *it-seq.begin() << " "; - this_bin_fill++; - } - os << endl; - cum_fill += this_bin_fill; - cum_fill_sq += this_bin_fill*this_bin_fill; - } - if (this_bin_fill0) fact *= k; - double prob=pow(lambda,k)/fact*exp(-lambda); - cum_prob += prob; - os << string(indent+delta_indent,' ') << "bins with " << k << " entries: " - << int(1000.0*count[k]/hashtabsize)/10.0 << "% (expected: " - << int(prob*1000)/10.0 << ")" << endl; - } - os << string(indent+delta_indent,' ') << "bins with more entries: " - << int(1000.0*count[MAXCOUNT]/hashtabsize)/10.0 << "% (expected: " - << int((1-cum_prob)*1000)/10.0 << ")" << endl; - - os << string(indent+delta_indent,' ') << "variance: " - << 1.0/hashtabsize*cum_fill_sq-(1.0/hashtabsize*cum_fill)*(1.0/hashtabsize*cum_fill) - << endl; - os << string(indent+delta_indent,' ') << "average fill: " - << (1.0*cum_fill)/hashtabsize - << " (should be equal to " << (1.0*seq.size())/hashtabsize << ")" << endl; -#endif // def EXPAIRSEQ_USE_HASHTAB -} - -#ifndef NO_GINAC_NAMESPACE -} // namespace GiNaC -#endif // ndef NO_GINAC_NAMESPACE diff --git a/ginac/relational.cpp b/ginac/relational.cpp index 76bd7da2..19a2e3a2 100644 --- a/ginac/relational.cpp +++ b/ginac/relational.cpp @@ -104,6 +104,111 @@ basic * relational::duplicate() const return new relational(*this); } +void relational::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("relational print",LOGLEVEL_PRINT); + if (precedence<=upper_precedence) os << "("; + lh.print(os,precedence); + switch (o) { + case equal: + os << "=="; + break; + case not_equal: + os << "!="; + break; + case less: + os << "<"; + break; + case less_or_equal: + os << "<="; + break; + case greater: + os << ">"; + break; + case greater_or_equal: + os << ">="; + break; + default: + os << "(INVALID RELATIONAL OPERATOR)"; + } + rh.print(os,precedence); + if (precedence<=upper_precedence) os << ")"; +} + +void relational::printraw(ostream & os) const +{ + debugmsg("relational printraw",LOGLEVEL_PRINT); + os << "RELATIONAL("; + lh.printraw(os); + os << ","; + rh.printraw(os); + os << ","; + switch (o) { + case equal: + os << "=="; + break; + case not_equal: + os << "!="; + break; + case less: + os << "<"; + break; + case less_or_equal: + os << "<="; + break; + case greater: + os << ">"; + break; + case greater_or_equal: + os << ">="; + break; + default: + os << "(INVALID RELATIONAL OPERATOR)"; + } + os << ")"; +} + +void relational::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const +{ + debugmsg("relational print csrc", LOGLEVEL_PRINT); + if (precedence<=upper_precedence) + os << "("; + + // Print left-hand expression + lh.bp->printcsrc(os, type, precedence); + + // Print relational operator + switch (o) { + case equal: + os << "=="; + break; + case not_equal: + os << "!="; + break; + case less: + os << "<"; + break; + case less_or_equal: + os << "<="; + break; + case greater: + os << ">"; + break; + case greater_or_equal: + os << ">="; + break; + default: + os << "(INVALID RELATIONAL OPERATOR)"; + break; + } + + // Print right-hand operator + rh.bp->printcsrc(os, type, precedence); + + if (precedence <= upper_precedence) + os << ")"; +} + bool relational::info(unsigned inf) const { switch (inf) { diff --git a/ginac/relational.h b/ginac/relational.h index 2767c89b..750b36f0 100644 --- a/ginac/relational.h +++ b/ginac/relational.h @@ -64,8 +64,8 @@ public: // functions overriding virtual functions from bases classes public: basic * duplicate() const; - void printraw(ostream & os) const; void print(ostream & os, unsigned upper_precedence=0) const; + void printraw(ostream & os) const; void printcsrc(ostream & os, unsigned type, unsigned upper_precedence=0) const; bool info(unsigned inf) const; int nops() const; diff --git a/ginac/series.cpp b/ginac/series.cpp index e21b2ba4..0cc2be34 100644 --- a/ginac/series.cpp +++ b/ginac/series.cpp @@ -112,6 +112,22 @@ basic *series::duplicate() const return new series(*this); } +void series::print(ostream &os, unsigned upper_precedence) const +{ + debugmsg("symbol print", LOGLEVEL_PRINT); + convert_to_poly().print(os, upper_precedence); +} + +void series::printraw(ostream &os) const +{ + debugmsg("symbol printraw", LOGLEVEL_PRINT); + os << "series(" << var << ";" << point << ";"; + for (epvector::const_iterator i=seq.begin(); i!=seq.end(); i++) { + os << "(" << (*i).rest << "," << (*i).coeff << "),"; + } + os << ")"; +} + // Highest degree of variable int series::degree(symbol const &s) const { diff --git a/ginac/series.h b/ginac/series.h index 08f859ff..a8a7ebf8 100644 --- a/ginac/series.h +++ b/ginac/series.h @@ -55,8 +55,8 @@ public: // functions overriding virtual functions from base classes public: basic *duplicate() const; - void printraw(ostream &os) const; void print(ostream &os, unsigned upper_precedence=0) const; + void printraw(ostream &os) const; int degree(symbol const &s) const; int ldegree(symbol const &s) const; ex coeff(symbol const &s, int const n=1) const; diff --git a/ginac/symbol.cpp b/ginac/symbol.cpp index 1bea2cec..2ea8ed51 100644 --- a/ginac/symbol.cpp +++ b/ginac/symbol.cpp @@ -110,10 +110,43 @@ basic * symbol::duplicate() const return new symbol(*this); } +void symbol::print(ostream & os, unsigned upper_precedence) const +{ + debugmsg("symbol print",LOGLEVEL_PRINT); + os << name; +} + +void symbol::printraw(ostream & os) const +{ + debugmsg("symbol printraw",LOGLEVEL_PRINT); + os << "symbol(" << "name=" << name << ",serial=" << serial + << ",hash=" << hashvalue << ",flags=" << flags << ")"; +} + +void symbol::printtree(ostream & os, unsigned indent) const +{ + debugmsg("symbol printtree",LOGLEVEL_PRINT); + os << string(indent,' ') << name << " (symbol): " + << "serial=" << serial + << ", hash=" << hashvalue << " (0x" << hex << hashvalue << dec << ")" + << ", flags=" << flags << endl; +} + +void symbol::printcsrc(ostream & os, unsigned type, unsigned upper_precedence) const +{ + debugmsg("symbol print csrc", LOGLEVEL_PRINT); + os << name; +} + bool symbol::info(unsigned inf) const { if (inf==info_flags::symbol) return true; - if (inf==info_flags::polynomial || inf==info_flags::integer_polynomial || inf==info_flags::rational_polynomial || inf==info_flags::rational_function) { + if (inf==info_flags::polynomial || + inf==info_flags::integer_polynomial || + inf==info_flags::cinteger_polynomial || + inf==info_flags::rational_polynomial || + inf==info_flags::crational_polynomial || + inf==info_flags::rational_function) { return true; } else { return basic::info(inf); diff --git a/ginac/symbol.h b/ginac/symbol.h index d315efc5..21c004c9 100644 --- a/ginac/symbol.h +++ b/ginac/symbol.h @@ -67,9 +67,9 @@ public: // functions overriding virtual functions from base classes public: basic * duplicate() const; + void print(ostream & os, unsigned upper_precedence=0) const; void printraw(ostream & os) const; void printtree(ostream & os, unsigned indent) const; - void print(ostream & os, unsigned upper_precedence=0) const; void printcsrc(ostream & os, unsigned type, unsigned upper_precedence=0) const; bool info(unsigned inf) const; ex expand(unsigned options=0) const; diff --git a/ginsh/Makefile.in b/ginsh/Makefile.in index 342ed462..8bb6e4f5 100644 --- a/ginsh/Makefile.in +++ b/ginsh/Makefile.in @@ -309,7 +309,7 @@ distdir: $(DISTFILES) @for file in $(DISTFILES); do \ d=$(srcdir); \ if test -d $$d/$$file; then \ - cp -pr $$d/$$file $(distdir)/$$file; \ + cp -pr $$/$$file $(distdir)/$$file; \ else \ test -f $(distdir)/$$file \ || ln $$d/$$file $(distdir)/$$file 2> /dev/null \