From 4afb5bbe2c0b0a60928120a042997ba7d89e8f5c Mon Sep 17 00:00:00 2001 From: Richard Kreckel Date: Mon, 20 Dec 1999 20:12:36 +0000 Subject: [PATCH] - more logic on the trigonometric function stuff. - changed several occurences of numeric const & to const numeric &, which is the same, but doxygen kept being confused because declarations differed from implementations. --- ginac/add.cpp | 57 ++++++++- ginac/add.h | 2 +- ginac/basic.h | 2 +- ginac/ex.cpp | 17 ++- ginac/ex.h | 2 +- ginac/inifcns.cpp | 11 +- ginac/inifcns_gamma.cpp | 111 ++++++++++++------ ginac/inifcns_trans.cpp | 191 ++++++++++++++++-------------- ginac/inifcns_zeta.cpp | 2 +- ginac/mul.cpp | 22 ++-- ginac/mul.h | 2 +- ginac/numeric.cpp | 176 ++++++++++++++-------------- ginac/numeric.h | 230 ++++++++++++++++++------------------ ginac/operators.h | 20 ++-- ginac/power.cpp | 14 ++- ginac/utils.cpp | 252 ++++++++++++++++++++++++++++++++++++++++ ginac/utils.h | 96 ++++++++++----- 17 files changed, 817 insertions(+), 390 deletions(-) diff --git a/ginac/add.cpp b/ginac/add.cpp index ad0ac1fc..d376b5d1 100644 --- a/ginac/add.cpp +++ b/ginac/add.cpp @@ -160,7 +160,7 @@ basic * add::duplicate() const return new add(*this); } -void add::print(ostream & os, unsigned upper_precedence) const +/*void add::print(ostream & os, unsigned upper_precedence) const { debugmsg("add print",LOGLEVEL_PRINT); if (precedence<=upper_precedence) os << "("; @@ -176,10 +176,17 @@ void add::print(ostream & os, unsigned upper_precedence) const } if (!coeff.is_equal(_num1()) && !coeff.is_equal(_num_1())) { - if (coeff.csgn()==-1) - (_num_1()*coeff).print(os, precedence); - else - coeff.print(os, precedence); + if (coeff.is_rational()) { + if (coeff.is_negative()) + os << -coeff; + else + os << coeff; + } else { + if (coeff.csgn()==-1) + (-coeff).print(os, precedence); + else + coeff.print(os, precedence); + } os << '*'; } os << cit->rest; @@ -190,6 +197,46 @@ void add::print(ostream & os, unsigned upper_precedence) const os << overall_coeff; } if (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; + // first print the overall numeric coefficient, if present: + if (!overall_coeff.is_zero()) { + os << overall_coeff; + first = false; + } + // then proceed with the remaining factors: + 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(_num1()) && + !coeff.is_equal(_num_1())) { + if (coeff.is_rational()) { + if (coeff.is_negative()) + os << -coeff; + else + os << coeff; + } else { + if (coeff.csgn()==-1) + (-coeff).print(os, precedence); + else + coeff.print(os, precedence); + } + os << '*'; + } + os << cit->rest; + } + if (precedence<=upper_precedence) os << ")"; } void add::printraw(ostream & os) const diff --git a/ginac/add.h b/ginac/add.h index 19438056..40cc25a5 100644 --- a/ginac/add.h +++ b/ginac/add.h @@ -72,7 +72,7 @@ public: ex series(symbol const & s, ex const & point, int order) const; ex normal(lst &sym_lst, lst &repl_lst, int level=0) const; numeric integer_content(void) const; - ex smod(numeric const &xi) const; + ex smod(const numeric &xi) const; numeric max_coefficient(void) const; exvector get_indices(void) const; ex simplify_ncmul(exvector const & v) const; diff --git a/ginac/basic.h b/ginac/basic.h index 87dc1dee..d087186d 100644 --- a/ginac/basic.h +++ b/ginac/basic.h @@ -138,7 +138,7 @@ public: // only const functions please (may break reference counting) virtual ex subs(lst const & ls, lst const & lr) const; virtual ex normal(lst &sym_lst, lst &repl_lst, int level=0) const; virtual numeric integer_content(void) const; - virtual ex smod(numeric const &xi) const; + virtual ex smod(const numeric &xi) const; virtual numeric max_coefficient(void) const; virtual exvector get_indices(void) const; virtual ex simplify_ncmul(exvector const & v) const; diff --git a/ginac/ex.cpp b/ginac/ex.cpp index 6a7a0aa9..b853ced5 100644 --- a/ginac/ex.cpp +++ b/ginac/ex.cpp @@ -105,7 +105,22 @@ ex::ex(basic const & other) ex::ex(int const i) { debugmsg("ex constructor from int",LOGLEVEL_CONSTRUCT); - construct_from_basic(numeric(i)); + switch (i) { // some tiny efficiency-hack (FIXME: is this ok?) + case -1: + bp = _ex_1().bp; + ++bp->refcount; + break; + case 0: + bp = _ex0().bp; + ++bp->refcount; + break; + case 1: + bp = _ex1().bp; + ++bp->refcount; + break; + default: + construct_from_basic(numeric(i)); + } } ex::ex(unsigned int const i) diff --git a/ginac/ex.h b/ginac/ex.h index f8d2970c..a7191035 100644 --- a/ginac/ex.h +++ b/ginac/ex.h @@ -44,7 +44,7 @@ extern ex const & _ex0(void); /* FIXME: should this pollute headers? */ #define INLINE_EX_CONSTRUCTORS -/** Lightweight interface to GiNaC's symbolic objects. Basically all it does is +/** Lightweight wrapper for GiNaC's symbolic objects. Basically all it does is * to hold a pointer to the other objects, manage the reference counting and * provide methods for manipulation of these objects. */ class ex diff --git a/ginac/inifcns.cpp b/ginac/inifcns.cpp index 0a38c4b0..28cadfb1 100644 --- a/ginac/inifcns.cpp +++ b/ginac/inifcns.cpp @@ -50,9 +50,9 @@ static ex Li2_eval(ex const & x) if (x.is_zero()) return x; if (x.is_equal(_ex1())) - return power(Pi, 2) / 6; + return power(Pi, _ex2()) / _ex6(); if (x.is_equal(_ex_1())) - return -power(Pi, 2) / 12; + return -power(Pi, _ex2()) / _ex12(); return Li2(x).hold(); } @@ -142,7 +142,10 @@ static ex Order_series(ex const & x, symbol const & s, ex const & point, int ord REGISTER_FUNCTION(Order, Order_eval, NULL, NULL, Order_series); -/** linear solve. */ +////////// +// Solve linear system +////////// + ex lsolve(ex const &eqns, ex const &symbols) { // solve a system of linear equations @@ -180,7 +183,7 @@ ex lsolve(ex const &eqns, ex const &symbols) matrix sys(eqns.nops(),symbols.nops()); matrix rhs(eqns.nops(),1); matrix vars(symbols.nops(),1); - + for (int r=0; r Pi^(1/2)*(1*3*..*(2*n-1))/(2^n) - if ((x*2).info(info_flags::posint)) { + if ((x*_ex2()).info(info_flags::posint)) { numeric n = ex_to_numeric(x).sub(_num1_2()); numeric coefficient = doublefactorial(n.mul(_num2()).sub(_num1())); - coefficient = coefficient.div(_num2().power(n)); - return coefficient * pow(Pi,_num1_2()); + coefficient = coefficient.div(pow(_num2(),n)); + return coefficient * pow(Pi,_ex1_2()); } else { // trap negative x==(-n+1/2) // gamma(-n+1/2) -> Pi^(1/2)*(-2)^n/(1*3*..*(2*n-1)) numeric n = abs(ex_to_numeric(x).sub(_num1_2())); - numeric coefficient = numeric(-2).power(n); + numeric coefficient = pow(_num_2(), n); coefficient = coefficient.div(doublefactorial(n.mul(_num2()).sub(_num1())));; - return coefficient*sqrt(Pi); + return coefficient*power(Pi,_ex1_2()); } } } return gamma(x).hold(); } -static ex gamma_diff(ex const & x, unsigned diff_param) +static ex gamma_diff(const ex & x, unsigned diff_param) { GINAC_ASSERT(diff_param==0); @@ -99,10 +99,11 @@ static ex gamma_diff(ex const & x, unsigned diff_param) return psi(x)*gamma(x); } -static ex gamma_series(ex const & x, symbol const & s, ex const & point, int order) +static ex gamma_series(const ex & x, const symbol & s, const ex & point, int order) { // method: - // Taylor series where there is no pole falls back to psi function evaluation. + // Taylor series where there is no pole falls back to psi function + // evaluation. // On a pole at -m use the recurrence relation // gamma(x) == gamma(x+1) / x // from which follows @@ -126,7 +127,7 @@ REGISTER_FUNCTION(gamma, gamma_eval, gamma_evalf, gamma_diff, gamma_series); // Beta-function ////////// -static ex beta_evalf(ex const & x, ex const & y) +static ex beta_evalf(const ex & x, const ex & y) { BEGIN_TYPECHECK TYPECHECK(x,numeric) @@ -137,24 +138,25 @@ static ex beta_evalf(ex const & x, ex const & y) / gamma(ex_to_numeric(x+y)); } -static ex beta_eval(ex const & x, ex const & y) +static ex beta_eval(const ex & x, const ex & y) { if (x.info(info_flags::numeric) && y.info(info_flags::numeric)) { numeric nx(ex_to_numeric(x)); numeric ny(ex_to_numeric(y)); // treat all problematic x and y that may not be passed into gamma, - // because they would throw there although beta(x,y) is well-defined: + // because they would throw there although beta(x,y) is well-defined + // using the formula beta(x,y) == (-1)^y * beta(1-x-y, y) if (nx.is_real() && nx.is_integer() && ny.is_real() && ny.is_integer()) { if (nx.is_negative()) { if (nx<=-ny) - return _num_1().power(ny)*beta(1-x-y, y); + return pow(_num_1(), ny)*beta(1-x-y, y); else throw (std::domain_error("beta_eval(): simple pole")); } if (ny.is_negative()) { if (ny<=-nx) - return _num_1().power(nx)*beta(1-y-x, x); + return pow(_num_1(), nx)*beta(1-y-x, x); else throw (std::domain_error("beta_eval(): simple pole")); } @@ -164,13 +166,15 @@ static ex beta_eval(ex const & x, ex const & y) if ((nx+ny).is_real() && (nx+ny).is_integer() && !(nx+ny).is_positive()) - return _ex0(); + return _ex0(); + // everything is ok: return gamma(x)*gamma(y)/gamma(x+y); } + return beta(x,y).hold(); } -static ex beta_diff(ex const & x, ex const & y, unsigned diff_param) +static ex beta_diff(const ex & x, const ex & y, unsigned diff_param) { GINAC_ASSERT(diff_param<2); ex retval; @@ -179,18 +183,43 @@ static ex beta_diff(ex const & x, ex const & y, unsigned diff_param) if (diff_param==0) retval = (psi(x)-psi(x+y))*beta(x,y); // d/dy beta(x,y) -> (psi(y)-psi(x+y)) * beta(x,y) - if (diff_param==1) + if (diff_param==1) retval = (psi(y)-psi(x+y))*beta(x,y); return retval; } -REGISTER_FUNCTION(beta, beta_eval, beta_evalf, beta_diff, NULL); +static ex beta_series(const ex & x, const ex & y, const symbol & s, const ex & point, int order) +{ + // method: + // Taylor series where there is no pole falls back to beta function + // evaluation. + // On a pole at -m use the recurrence relation + // gamma(x) == gamma(x+1) / x + // from which follows + // series(gamma(x),x,-m,order) == + // series(gamma(x+m+1)/(x*(x+1)...*(x+m)),x,-m,order+1); + ex xpoint = x.subs(s==point); + ex ypoint = y.subs(s==point); + if ((!xpoint.info(info_flags::integer) || xpoint.info(info_flags::positive)) && + (!ypoint.info(info_flags::integer) || ypoint.info(info_flags::positive))) + throw do_taylor(); // caught by function::series() + // if we got here we have to care for a simple pole at -m: + throw (std::domain_error("beta_series(): please code me")); + /*numeric m = -ex_to_numeric(xpoint); + *ex ser_numer = gamma(x+m+_ex1()); + *ex ser_denom = _ex1(); + *for (numeric p; p<=m; ++p) + * ser_denom *= x+p; + *return (ser_numer/ser_denom).series(s, point, order+1);*/ +} + +REGISTER_FUNCTION(beta, beta_eval, beta_evalf, beta_diff, beta_series); ////////// // Psi-function (aka digamma-function) ////////// -static ex psi1_evalf(ex const & x) +static ex psi1_evalf(const ex & x) { BEGIN_TYPECHECK TYPECHECK(x,numeric) @@ -201,7 +230,7 @@ static ex psi1_evalf(ex const & x) /** Evaluation of digamma-function psi(x). * Somebody ought to provide some good numerical evaluation some day... */ -static ex psi1_eval(ex const & x) +static ex psi1_eval(const ex & x) { if (x.info(info_flags::numeric)) { if (x.info(info_flags::integer) && !x.info(info_flags::positive)) @@ -229,7 +258,7 @@ static ex psi1_eval(ex const & x) return psi(x).hold(); } -static ex psi1_diff(ex const & x, unsigned diff_param) +static ex psi1_diff(const ex & x, unsigned diff_param) { GINAC_ASSERT(diff_param==0); @@ -237,7 +266,7 @@ static ex psi1_diff(ex const & x, unsigned diff_param) return psi(_ex1(), x); } -static ex psi1_series(ex const & x, symbol const & s, ex const & point, int order) +static ex psi1_series(const ex & x, const symbol & s, const ex & point, int order) { // method: // Taylor series where there is no pole falls back to polygamma function @@ -264,7 +293,7 @@ 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) +static ex psi2_evalf(const ex & n, const ex & x) { BEGIN_TYPECHECK TYPECHECK(n,numeric) @@ -276,7 +305,7 @@ static ex psi2_evalf(ex const & n, ex const & 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) +static ex psi2_eval(const ex & n, const ex & x) { // psi(0,x) -> psi(x) if (n.is_zero()) @@ -288,13 +317,29 @@ static ex psi2_eval(ex const & n, ex const & x) x.info(info_flags::numeric)) { numeric nn = ex_to_numeric(n); numeric nx = ex_to_numeric(x); - if (x.is_equal(_ex1())) - return _num_1().power(nn+_num1())*factorial(nn)*zeta(ex(nn+_num1())); + if (nx.is_integer()) { + if (nx.is_equal(_num1())) + return pow(_num_1(), nn+_num1())*factorial(nn)*zeta(ex(nn+_num1())); + if (nx.is_positive()) { + // use the recurrence relation + // psi(n,m) == psi(n,m+1) - (-)^n * n! / m^(n+1) + // to relate psi(n,m) to psi(n,1): + // psi(n,m) == psi(n,1) + r + // where r == (-)^n * n! * (1^(-n-1) + ... + (m-1)^(-n-1)) + numeric recur; + for (numeric p(1); p psi(n+1,x) - return psi(n+1, x); + return psi(n+_ex1(), x); } -static ex psi2_series(ex const & n, ex const & x, symbol const & s, ex const & point, int order) +static ex psi2_series(const ex & n, const ex & x, const symbol & s, const ex & point, int order) { // method: // Taylor series where there is no pole falls back to polygamma function // evaluation. // On a pole at -m use the recurrence relation - // psi(n,x) == psi(n,x+1) - (-)^n * n! / z^(n+1) + // psi(n,x) == psi(n,x+1) - (-)^n * n! / x^(n+1) // from which follows // series(psi(x),x,-m,order) == - // series(psi(x+m+1) - (-1)^n * n! - // * ((x)^(-n-1) + (x+1)^(-n-1) + (x+m)^(-n-1))),x,-m,order); + // series(psi(x+m+1) - (-1)^n * n! * ((x)^(-n-1) + (x+1)^(-n-1) + ... + // ... + (x+m)^(-n-1))),x,-m,order); ex xpoint = x.subs(s==point); if (!xpoint.info(info_flags::integer) || xpoint.info(info_flags::positive)) throw do_taylor(); // caught by function::series() diff --git a/ginac/inifcns_trans.cpp b/ginac/inifcns_trans.cpp index 4473fe02..9937e305 100644 --- a/ginac/inifcns_trans.cpp +++ b/ginac/inifcns_trans.cpp @@ -122,7 +122,7 @@ static ex log_eval(ex const & x) } // log(exp(t)) -> t (for real-valued t): if (is_ex_the_function(x, exp)) { - ex t=x.op(0); + ex t = x.op(0); if (t.info(info_flags::real)) return t; } @@ -135,7 +135,7 @@ static ex log_diff(ex const & x, unsigned diff_param) GINAC_ASSERT(diff_param==0); // d/dx log(x) -> 1/x - return power(x, -1); + return power(x, _ex_1()); } REGISTER_FUNCTION(log, log_eval, log_evalf, log_diff, NULL); @@ -155,38 +155,42 @@ static ex sin_evalf(ex const & x) static ex sin_eval(ex const & x) { - // sin(n*Pi/6) -> { 0 | +/-1/2 | +/-sqrt(3)/2 | +/-1 } - ex SixExOverPi = _ex6()*x/Pi; - if (SixExOverPi.info(info_flags::integer)) { - numeric z = smod(ex_to_numeric(SixExOverPi),_num12()); - if (z.is_equal(_num_5())) // sin(7*Pi/6) -> -1/2 - return _ex_1_2(); - if (z.is_equal(_num_4())) // sin(8*Pi/6) -> -sqrt(3)/2 - return _ex_1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num_3())) // sin(9*Pi/6) -> -1 - return _ex_1(); - if (z.is_equal(_num_2())) // sin(10*Pi/6) -> -sqrt(3)/2 - return _ex_1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num_1())) // sin(11*Pi/6) -> -1/2 - return _ex_1_2(); - if (z.is_equal(_num0())) // sin(0) -> 0 - return _ex0(); - if (z.is_equal(_num1())) // sin(1*Pi/6) -> 1/2 - return _ex1_2(); - if (z.is_equal(_num2())) // sin(2*Pi/6) -> sqrt(3)/2 - return _ex1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num3())) // sin(3*Pi/6) -> 1 - return _ex1(); - if (z.is_equal(_num4())) // sin(4*Pi/6) -> sqrt(3)/2 - return _ex1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num5())) // sin(5*Pi/6) -> 1/2 - return _ex1_2(); - if (z.is_equal(_num6())) // sin(6*Pi/6) -> 0 + // sin(n/d*Pi) -> { all known non-nested radicals } + ex SixtyExOverPi = _ex60()*x/Pi; + ex sign = _ex1(); + if (SixtyExOverPi.info(info_flags::integer)) { + numeric z = mod(ex_to_numeric(SixtyExOverPi),_num120()); + if (z>=_num60()) { + // wrap to interval [0, Pi) + z -= _num60(); + sign = _ex_1(); + } + if (z>_num30()) { + // wrap to interval [0, Pi/2) + z = _num60()-z; + } + if (z.is_equal(_num0())) // sin(0) -> 0 return _ex0(); + if (z.is_equal(_num5())) // sin(Pi/12) -> sqrt(6)/4*(1-sqrt(3)/3) + return sign*_ex1_4()*power(_ex6(),_ex1_2())*(_ex1()+_ex_1_3()*power(_ex3(),_ex1_2())); + if (z.is_equal(_num6())) // sin(Pi/10) -> sqrt(5)/4-1/4 + return sign*(_ex1_4()*power(_ex5(),_ex1_2())+_ex_1_4()); + if (z.is_equal(_num10())) // sin(Pi/6) -> 1/2 + return sign*_ex1_2(); + if (z.is_equal(_num15())) // sin(Pi/4) -> sqrt(2)/2 + return sign*_ex1_2()*power(_ex2(),_ex1_2()); + if (z.is_equal(_num18())) // sin(3/10*Pi) -> sqrt(5)/4+1/4 + return sign*(_ex1_4()*power(_ex5(),_ex1_2())+_ex1_4()); + if (z.is_equal(_num20())) // sin(Pi/3) -> sqrt(3)/2 + return sign*_ex1_2()*power(_ex3(),_ex1_2()); + if (z.is_equal(_num25())) // sin(5/12*Pi) -> sqrt(6)/4*(1+sqrt(3)/3) + return sign*_ex1_4()*power(_ex6(),_ex1_2())*(_ex1()+_ex1_3()*power(_ex3(),_ex1_2())); + if (z.is_equal(_num30())) // sin(Pi/2) -> 1 + return sign*_ex1(); } if (is_ex_exactly_of_type(x, function)) { - ex t=x.op(0); + ex t = x.op(0); // sin(asin(x)) -> x if (is_ex_the_function(x, asin)) return t; @@ -230,38 +234,42 @@ static ex cos_evalf(ex const & x) static ex cos_eval(ex const & x) { - // cos(n*Pi/6) -> { 0 | +/-1/2 | +/-sqrt(3)/2 | +/-1 } - ex SixExOverPi = _ex6()*x/Pi; - if (SixExOverPi.info(info_flags::integer)) { - numeric z = smod(ex_to_numeric(SixExOverPi),_num12()); - if (z.is_equal(_num_5())) // cos(7*Pi/6) -> -sqrt(3)/2 - return _ex_1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num_4())) // cos(8*Pi/6) -> -1/2 - return _ex_1_2(); - if (z.is_equal(_num_3())) // cos(9*Pi/6) -> 0 - return _ex0(); - if (z.is_equal(_num_2())) // cos(10*Pi/6) -> 1/2 - return _ex1_2(); - if (z.is_equal(_num_1())) // cos(11*Pi/6) -> sqrt(3)/2 - return _ex1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num0())) // cos(0) -> 1 - return _ex1(); - if (z.is_equal(_num1())) // cos(1*Pi/6) -> sqrt(3)/2 - return _ex1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num2())) // cos(2*Pi/6) -> 1/2 - return _ex1_2(); - if (z.is_equal(_num3())) // cos(3*Pi/6) -> 0 - return _ex0(); - if (z.is_equal(_num4())) // cos(4*Pi/6) -> -1/2 - return _ex_1_2(); - if (z.is_equal(_num5())) // cos(5*Pi/6) -> -sqrt(3)/2 - return _ex_1_2()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num6())) // cos(6*Pi/6) -> -1 - return _ex_1(); + // cos(n/d*Pi) -> { all known non-nested radicals } + ex SixtyExOverPi = _ex60()*x/Pi; + ex sign = _ex1(); + if (SixtyExOverPi.info(info_flags::integer)) { + numeric z = mod(ex_to_numeric(SixtyExOverPi),_num120()); + if (z>=_num60()) { + // wrap to interval [0, Pi) + z = _num120()-z; + } + if (z>=_num30()) { + // wrap to interval [0, Pi/2) + z = _num60()-z; + sign = _ex_1(); + } + if (z.is_equal(_num0())) // cos(0) -> 1 + return sign*_ex1(); + if (z.is_equal(_num5())) // cos(Pi/12) -> sqrt(6)/4*(1+sqrt(3)/3) + return sign*_ex1_4()*power(_ex6(),_ex1_2())*(_ex1()+_ex1_3()*power(_ex3(),_ex1_2())); + if (z.is_equal(_num10())) // cos(Pi/6) -> sqrt(3)/2 + return sign*_ex1_2()*power(_ex3(),_ex1_2()); + if (z.is_equal(_num12())) // cos(Pi/5) -> sqrt(5)/4+1/4 + return sign*(_ex1_4()*power(_ex5(),_ex1_2())+_ex1_4()); + if (z.is_equal(_num15())) // cos(Pi/4) -> sqrt(2)/2 + return sign*_ex1_2()*power(_ex2(),_ex1_2()); + if (z.is_equal(_num20())) // cos(Pi/3) -> 1/2 + return sign*_ex1_2(); + if (z.is_equal(_num24())) // cos(2/5*Pi) -> sqrt(5)/4-1/4x + return sign*(_ex1_4()*power(_ex5(),_ex1_2())+_ex_1_4()); + if (z.is_equal(_num25())) // cos(5/12*Pi) -> sqrt(6)/4*(1-sqrt(3)/3) + return sign*_ex1_4()*power(_ex6(),_ex1_2())*(_ex1()+_ex_1_3()*power(_ex3(),_ex1_2())); + if (z.is_equal(_num30())) // cos(Pi/2) -> 0 + return sign*_ex0(); } if (is_ex_exactly_of_type(x, function)) { - ex t=x.op(0); + ex t = x.op(0); // cos(acos(x)) -> x if (is_ex_the_function(x, acos)) return t; @@ -305,26 +313,38 @@ static ex tan_evalf(ex const & x) static ex tan_eval(ex const & x) { - // tan(n*Pi/6) -> { 0 | +/-sqrt(3) | +/-sqrt(3)/2 } - ex SixExOverPi = _ex6()*x/Pi; - if (SixExOverPi.info(info_flags::integer)) { - numeric z = smod(ex_to_numeric(SixExOverPi),_num6()); - if (z.is_equal(_num_2())) // tan(4*Pi/6) -> -sqrt(3) - return _ex_1()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num_1())) // tan(5*Pi/6) -> -sqrt(3)/3 - return _ex_1_3()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num0())) // tan(0) -> 0 + // tan(n/d*Pi) -> { all known non-nested radicals } + ex SixtyExOverPi = _ex60()*x/Pi; + ex sign = _ex1(); + if (SixtyExOverPi.info(info_flags::integer)) { + numeric z = mod(ex_to_numeric(SixtyExOverPi),_num60()); + if (z>=_num60()) { + // wrap to interval [0, Pi) + z -= _num60(); + } + if (z>=_num30()) { + // wrap to interval [0, Pi/2) + z = _num60()-z; + sign = _ex_1(); + } + if (z.is_equal(_num0())) // tan(0) -> 0 return _ex0(); - if (z.is_equal(_num1())) // tan(1*Pi/6) -> sqrt(3)/3 - return _ex1_3()*power(_ex3(),_ex1_2()); - if (z.is_equal(_num2())) // tan(2*Pi/6) -> sqrt(3) - return power(_ex3(),_ex1_2()); - if (z.is_equal(_num3())) // tan(3*Pi/6) -> infinity + if (z.is_equal(_num5())) // tan(Pi/12) -> 2-sqrt(3) + return sign*(_ex2()-power(_ex3(),_ex1_2())); + if (z.is_equal(_num10())) // tan(Pi/6) -> sqrt(3)/3 + return sign*_ex1_3()*power(_ex3(),_ex1_2()); + if (z.is_equal(_num15())) // tan(Pi/4) -> 1 + return sign*_ex1(); + if (z.is_equal(_num20())) // tan(Pi/3) -> sqrt(3) + return sign*power(_ex3(),_ex1_2()); + if (z.is_equal(_num25())) // tan(5/12*Pi) -> 2+sqrt(3) + return sign*(power(_ex3(),_ex1_2())+_ex2()); + if (z.is_equal(_num30())) // tan(Pi/2) -> infinity throw (std::domain_error("tan_eval(): infinity")); } - + if (is_ex_exactly_of_type(x, function)) { - ex t=x.op(0); + ex t = x.op(0); // tan(atan(x)) -> x if (is_ex_the_function(x, atan)) return t; @@ -349,7 +369,7 @@ static ex tan_diff(ex const & x, unsigned diff_param) GINAC_ASSERT(diff_param==0); // d/dx tan(x) -> 1+tan(x)^2; - return (1+power(tan(x),_ex2())); + return (_ex1()+power(tan(x),_ex2())); } static ex tan_series(ex const & x, symbol const & s, ex const & point, int order) @@ -436,10 +456,10 @@ static ex acos_eval(ex const & x) return _ex0(); // acos(1/2) -> Pi/3 if (x.is_equal(_ex1_2())) - return numeric(1,3)*Pi; + return _ex1_3()*Pi; // acos(0) -> Pi/2 if (x.is_zero()) - return numeric(1,2)*Pi; + return _ex1_2()*Pi; // acos(-1/2) -> 2/3*Pi if (x.is_equal(_ex_1_2())) return numeric(2,3)*Pi; @@ -495,7 +515,8 @@ static ex atan_diff(ex const & x, unsigned diff_param) { GINAC_ASSERT(diff_param==0); - return power(1+x*x, -1); + // d/dx atan(x) -> 1/(1+x^2) + return power(_ex1()+power(x,_ex2()), _ex_1()); } REGISTER_FUNCTION(atan, atan_eval, atan_evalf, atan_diff, NULL); @@ -530,10 +551,10 @@ static ex atan2_diff(ex const & y, ex const & x, unsigned diff_param) if (diff_param==0) { // d/dy atan(y,x) - return x*pow(pow(x,2)+pow(y,2),-1); + return x*power(power(x,_ex2())+power(y,_ex2()),_ex_1()); } // d/dx atan(y,x) - return -y*pow(pow(x,2)+pow(y,2),-1); + return -y*power(power(x,_ex2())+power(y,_ex2()),_ex_1()); } REGISTER_FUNCTION(atan2, atan2_eval, atan2_evalf, atan2_diff, NULL); @@ -565,7 +586,7 @@ static ex sinh_eval(ex const & x) return I*sin(x/I); if (is_ex_exactly_of_type(x, function)) { - ex t=x.op(0); + ex t = x.op(0); // sinh(asinh(x)) -> x if (is_ex_the_function(x, asinh)) return t; @@ -617,7 +638,7 @@ static ex cosh_eval(ex const & x) return cos(x/I); if (is_ex_exactly_of_type(x, function)) { - ex t=x.op(0); + ex t = x.op(0); // cosh(acosh(x)) -> x if (is_ex_the_function(x, acosh)) return t; @@ -669,7 +690,7 @@ static ex tanh_eval(ex const & x) return I*tan(x/I); if (is_ex_exactly_of_type(x, function)) { - ex t=x.op(0); + ex t = x.op(0); // tanh(atanh(x)) -> x if (is_ex_the_function(x, atanh)) return t; @@ -738,7 +759,7 @@ static ex asinh_diff(ex const & x, unsigned diff_param) GINAC_ASSERT(diff_param==0); // d/dx asinh(x) -> 1/sqrt(1+x^2) - return power(1+power(x,_ex2()),_ex_1_2()); + return power(_ex1()+power(x,_ex2()),_ex_1_2()); } REGISTER_FUNCTION(asinh, asinh_eval, asinh_evalf, asinh_diff, NULL); diff --git a/ginac/inifcns_zeta.cpp b/ginac/inifcns_zeta.cpp index 9fdce1e8..dfd4638e 100644 --- a/ginac/inifcns_zeta.cpp +++ b/ginac/inifcns_zeta.cpp @@ -62,7 +62,7 @@ static ex zeta1_eval(ex const & x) if (x.info(info_flags::odd)) return zeta(x).hold(); else - return abs(bernoulli(y))*pow(Pi,x)*_num2().power(y-_num1())/factorial(y); + return abs(bernoulli(y))*pow(Pi,x)*pow(_num2(),y-_num1())/factorial(y); } else { if (x.info(info_flags::odd)) return -bernoulli(_num1()-y)/(_num1()-y); diff --git a/ginac/mul.cpp b/ginac/mul.cpp index 9bee46f2..b444b1e5 100644 --- a/ginac/mul.cpp +++ b/ginac/mul.cpp @@ -181,13 +181,21 @@ void mul::print(ostream & os, unsigned upper_precedence) const 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(_ex1()) && - !overall_coeff.is_equal(_ex_1())) { - if (ex_to_numeric(overall_coeff).csgn()==-1) - (_num_1()*overall_coeff).print(os, precedence); - else - overall_coeff.print(os, precedence); + numeric coeff = ex_to_numeric(overall_coeff); + if (coeff.csgn()==-1) os << '-'; + if (!coeff.is_equal(_num1()) && + !coeff.is_equal(_num_1())) { + if (coeff.is_rational()) { + if (coeff.is_negative()) + os << -coeff; + else + os << coeff; + } else { + if (coeff.csgn()==-1) + (-coeff).print(os, precedence); + else + coeff.print(os, precedence); + } os << '*'; } // then proceed with the remaining factors: diff --git a/ginac/mul.h b/ginac/mul.h index 71b686bb..1d950711 100644 --- a/ginac/mul.h +++ b/ginac/mul.h @@ -73,7 +73,7 @@ public: ex series(symbol const & s, ex const & point, int order) const; ex normal(lst &sym_lst, lst &repl_lst, int level=0) const; numeric integer_content(void) const; - ex smod(numeric const &xi) const; + ex smod(const numeric &xi) const; numeric max_coefficient(void) const; exvector get_indices(void) const; ex simplify_ncmul(exvector const & v) const; diff --git a/ginac/numeric.cpp b/ginac/numeric.cpp index 472d9c31..49374485 100644 --- a/ginac/numeric.cpp +++ b/ginac/numeric.cpp @@ -72,13 +72,13 @@ numeric::~numeric() destroy(0); } -numeric::numeric(numeric const & other) +numeric::numeric(const numeric & other) { debugmsg("numeric copy constructor", LOGLEVEL_CONSTRUCT); copy(other); } -numeric const & numeric::operator=(numeric const & other) +const numeric & numeric::operator=(const numeric & other) { debugmsg("numeric operator=", LOGLEVEL_ASSIGNMENT); if (this != &other) { @@ -90,7 +90,7 @@ numeric const & numeric::operator=(numeric const & other) // protected -void numeric::copy(numeric const & other) +void numeric::copy(const numeric & other) { basic::copy(other); value = new cl_N(*other.value); @@ -110,7 +110,7 @@ void numeric::destroy(bool call_parent) numeric::numeric(int i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from int",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from int",LOGLEVEL_CONSTRUCT); // Not the whole int-range is available if we don't cast to long // first. This is due to the behaviour of the cl_I-ctor, which // emphasizes efficiency: @@ -122,7 +122,7 @@ numeric::numeric(int i) : basic(TINFO_numeric) numeric::numeric(unsigned int i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from uint",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from uint",LOGLEVEL_CONSTRUCT); // Not the whole uint-range is available if we don't cast to ulong // first. This is due to the behaviour of the cl_I-ctor, which // emphasizes efficiency: @@ -134,7 +134,7 @@ numeric::numeric(unsigned int i) : basic(TINFO_numeric) numeric::numeric(long i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from long",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from long",LOGLEVEL_CONSTRUCT); value = new cl_I(i); calchash(); setflag(status_flags::evaluated| @@ -143,7 +143,7 @@ numeric::numeric(long i) : basic(TINFO_numeric) numeric::numeric(unsigned long i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from ulong",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from ulong",LOGLEVEL_CONSTRUCT); value = new cl_I(i); calchash(); setflag(status_flags::evaluated| @@ -155,7 +155,7 @@ numeric::numeric(unsigned long i) : basic(TINFO_numeric) * @exception overflow_error (division by zero) */ numeric::numeric(long numer, long denom) : basic(TINFO_numeric) { - debugmsg("numeric constructor from long/long",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from long/long",LOGLEVEL_CONSTRUCT); if (!denom) throw (std::overflow_error("division by zero")); value = new cl_I(numer); @@ -167,7 +167,7 @@ numeric::numeric(long numer, long denom) : basic(TINFO_numeric) numeric::numeric(double d) : basic(TINFO_numeric) { - debugmsg("numeric constructor from double",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from double",LOGLEVEL_CONSTRUCT); // We really want to explicitly use the type cl_LF instead of the // more general cl_F, since that would give us a cl_DF only which // will not be promoted to cl_LF if overflow occurs: @@ -180,7 +180,7 @@ numeric::numeric(double d) : basic(TINFO_numeric) numeric::numeric(char const *s) : basic(TINFO_numeric) { // MISSING: treatment of complex and ints and rationals. - debugmsg("numeric constructor from string",LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from string",LOGLEVEL_CONSTRUCT); if (strchr(s, '.')) value = new cl_LF(s); else @@ -194,7 +194,7 @@ numeric::numeric(char const *s) : basic(TINFO_numeric) * only. */ numeric::numeric(cl_N const & z) : basic(TINFO_numeric) { - debugmsg("numeric constructor from cl_N", LOGLEVEL_CONSTRUCT); + debugmsg("const numericructor from cl_N", LOGLEVEL_CONSTRUCT); value = new cl_N(z); calchash(); setflag(status_flags::evaluated| @@ -213,19 +213,11 @@ basic * numeric::duplicate() const return new numeric(*this); } -// The method printraw doesn't do much, it simply uses CLN's operator<<() for -// output, which is ugly but reliable. Examples: -// 2+2i -void numeric::printraw(ostream & os) const -{ - debugmsg("numeric printraw", LOGLEVEL_PRINT); - os << "numeric(" << *value << ")"; -} - -// The method print adds to the output so it blends more consistently together -// with the other routines and produces something compatible to Maple input. void numeric::print(ostream & os, unsigned upper_precedence) const { + // The method print adds to the output so it blends more consistently + // together with the other routines and produces something compatible to + // ginsh input. debugmsg("numeric print", LOGLEVEL_PRINT); if (is_real()) { // case 1, real: x or -x @@ -276,6 +268,14 @@ void numeric::print(ostream & os, unsigned upper_precedence) const } } + +void numeric::printraw(ostream & os) const +{ + // The method printraw doesn't do much, it simply uses CLN's operator<<() + // for output, which is ugly but reliable. e.g: 2+2i + debugmsg("numeric printraw", LOGLEVEL_PRINT); + os << "numeric(" << *value << ")"; +} void numeric::printtree(ostream & os, unsigned indent) const { debugmsg("numeric printtree", LOGLEVEL_PRINT); @@ -379,7 +379,7 @@ ex numeric::evalf(int level) const int numeric::compare_same_type(basic const & other) const { GINAC_ASSERT(is_exactly_of_type(other, numeric)); - numeric const & o = static_cast(const_cast(other)); + const numeric & o = static_cast(const_cast(other)); if (*value == *o.value) { return 0; @@ -391,7 +391,7 @@ int numeric::compare_same_type(basic const & other) const bool numeric::is_equal_same_type(basic const & other) const { GINAC_ASSERT(is_exactly_of_type(other,numeric)); - numeric const *o = static_cast(&other); + const numeric *o = static_cast(&other); return is_equal(*o); } @@ -424,21 +424,21 @@ unsigned numeric::calchash(void) const /** Numerical addition method. Adds argument to *this and returns result as * a new numeric object. */ -numeric numeric::add(numeric const & other) const +numeric numeric::add(const numeric & other) const { return numeric((*value)+(*other.value)); } /** Numerical subtraction method. Subtracts argument from *this and returns * result as a new numeric object. */ -numeric numeric::sub(numeric const & other) const +numeric numeric::sub(const numeric & other) const { return numeric((*value)-(*other.value)); } /** Numerical multiplication method. Multiplies *this and argument and returns * result as a new numeric object. */ -numeric numeric::mul(numeric const & other) const +numeric numeric::mul(const numeric & other) const { static const numeric * _num1p=&_num1(); if (this==_num1p) { @@ -453,14 +453,14 @@ numeric numeric::mul(numeric const & other) const * a new numeric object. * * @exception overflow_error (division by zero) */ -numeric numeric::div(numeric const & other) const +numeric numeric::div(const numeric & other) const { if (::zerop(*other.value)) throw (std::overflow_error("division by zero")); return numeric((*value)/(*other.value)); } -numeric numeric::power(numeric const & other) const +numeric numeric::power(const numeric & other) const { static const numeric * _num1p=&_num1(); if (&other==_num1p) @@ -476,19 +476,19 @@ numeric numeric::inverse(void) const return numeric(::recip(*value)); // -> CLN } -numeric const & numeric::add_dyn(numeric const & other) const +const numeric & numeric::add_dyn(const numeric & other) const { - return static_cast((new numeric((*value)+(*other.value)))-> + return static_cast((new numeric((*value)+(*other.value)))-> setflag(status_flags::dynallocated)); } -numeric const & numeric::sub_dyn(numeric const & other) const +const numeric & numeric::sub_dyn(const numeric & other) const { - return static_cast((new numeric((*value)-(*other.value)))-> + return static_cast((new numeric((*value)-(*other.value)))-> setflag(status_flags::dynallocated)); } -numeric const & numeric::mul_dyn(numeric const & other) const +const numeric & numeric::mul_dyn(const numeric & other) const { static const numeric * _num1p=&_num1(); if (this==_num1p) { @@ -496,55 +496,55 @@ numeric const & numeric::mul_dyn(numeric const & other) const } else if (&other==_num1p) { return *this; } - return static_cast((new numeric((*value)*(*other.value)))-> + return static_cast((new numeric((*value)*(*other.value)))-> setflag(status_flags::dynallocated)); } -numeric const & numeric::div_dyn(numeric const & other) const +const numeric & numeric::div_dyn(const numeric & other) const { if (::zerop(*other.value)) throw (std::overflow_error("division by zero")); - return static_cast((new numeric((*value)/(*other.value)))-> + return static_cast((new numeric((*value)/(*other.value)))-> setflag(status_flags::dynallocated)); } -numeric const & numeric::power_dyn(numeric const & other) const +const numeric & numeric::power_dyn(const numeric & other) const { static const numeric * _num1p=&_num1(); if (&other==_num1p) return *this; 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)))-> + return static_cast((new numeric(::expt(*value,*other.value)))-> setflag(status_flags::dynallocated)); } -numeric const & numeric::operator=(int i) +const numeric & numeric::operator=(int i) { return operator=(numeric(i)); } -numeric const & numeric::operator=(unsigned int i) +const numeric & numeric::operator=(unsigned int i) { return operator=(numeric(i)); } -numeric const & numeric::operator=(long i) +const numeric & numeric::operator=(long i) { return operator=(numeric(i)); } -numeric const & numeric::operator=(unsigned long i) +const numeric & numeric::operator=(unsigned long i) { return operator=(numeric(i)); } -numeric const & numeric::operator=(double d) +const numeric & numeric::operator=(double d) { return operator=(numeric(d)); } -numeric const & numeric::operator=(char const * s) +const numeric & numeric::operator=(char const * s) { return operator=(numeric(s)); } @@ -553,7 +553,7 @@ numeric const & numeric::operator=(char const * s) * csgn(x)==0 for x==0, csgn(x)==1 for Re(x)>0 or Re(x)=0 and Im(x)>0, * csgn(x)==-1 for Re(x)<0 or Re(x)=0 and Im(x)<0. * - * @see numeric::compare(numeric const & other) */ + * @see numeric::compare(const numeric & other) */ int numeric::csgn(void) const { if (is_zero()) @@ -578,7 +578,7 @@ int numeric::csgn(void) const * * @return csgn(*this-other) * @see numeric::csgn(void) */ -int numeric::compare(numeric const & other) const +int numeric::compare(const numeric & other) const { // Comparing two real numbers? if (is_real() && other.is_real()) @@ -594,7 +594,7 @@ int numeric::compare(numeric const & other) const } } -bool numeric::is_equal(numeric const & other) const +bool numeric::is_equal(const numeric & other) const { return (*value == *other.value); } @@ -672,12 +672,12 @@ bool numeric::is_real(void) const return ::instanceof(*value, cl_R_ring); // -> CLN } -bool numeric::operator==(numeric const & other) const +bool numeric::operator==(const numeric & other) const { return (*value == *other.value); // -> CLN } -bool numeric::operator!=(numeric const & other) const +bool numeric::operator!=(const numeric & other) const { return (*value != *other.value); // -> CLN } @@ -713,7 +713,7 @@ bool numeric::is_crational(void) const /** Numerical comparison: less. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator<(numeric const & other) const +bool numeric::operator<(const numeric & other) const { if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) < The(cl_R)(*other.value)); // -> CLN @@ -724,7 +724,7 @@ bool numeric::operator<(numeric const & other) const /** Numerical comparison: less or equal. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator<=(numeric const & other) const +bool numeric::operator<=(const numeric & other) const { if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) <= The(cl_R)(*other.value)); // -> CLN @@ -735,7 +735,7 @@ bool numeric::operator<=(numeric const & other) const /** Numerical comparison: greater. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator>(numeric const & other) const +bool numeric::operator>(const numeric & other) const { if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) > The(cl_R)(*other.value)); // -> CLN @@ -746,7 +746,7 @@ bool numeric::operator>(numeric const & other) const /** Numerical comparison: greater or equal. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator>=(numeric const & other) const +bool numeric::operator>=(const numeric & other) const { if (is_real() && other.is_real()) return (bool)(The(cl_R)(*value) >= The(cl_R)(*other.value)); // -> CLN @@ -928,7 +928,7 @@ const numeric I = numeric(complex(cl_I(0),cl_I(1))); /** Exponential function. * * @return arbitrary precision numerical exp(x). */ -numeric exp(numeric const & x) +numeric exp(const numeric & x) { return ::exp(*x.value); // -> CLN } @@ -938,7 +938,7 @@ numeric exp(numeric const & x) * @param z complex number * @return arbitrary precision numerical log(x). * @exception overflow_error (logarithmic singularity) */ -numeric log(numeric const & z) +numeric log(const numeric & z) { if (z.is_zero()) throw (std::overflow_error("log(): logarithmic singularity")); @@ -948,7 +948,7 @@ numeric log(numeric const & z) /** Numeric sine (trigonometric function). * * @return arbitrary precision numerical sin(x). */ -numeric sin(numeric const & x) +numeric sin(const numeric & x) { return ::sin(*x.value); // -> CLN } @@ -956,7 +956,7 @@ numeric sin(numeric const & x) /** Numeric cosine (trigonometric function). * * @return arbitrary precision numerical cos(x). */ -numeric cos(numeric const & x) +numeric cos(const numeric & x) { return ::cos(*x.value); // -> CLN } @@ -964,7 +964,7 @@ numeric cos(numeric const & x) /** Numeric tangent (trigonometric function). * * @return arbitrary precision numerical tan(x). */ -numeric tan(numeric const & x) +numeric tan(const numeric & x) { return ::tan(*x.value); // -> CLN } @@ -972,7 +972,7 @@ numeric tan(numeric const & x) /** Numeric inverse sine (trigonometric function). * * @return arbitrary precision numerical asin(x). */ -numeric asin(numeric const & x) +numeric asin(const numeric & x) { return ::asin(*x.value); // -> CLN } @@ -980,7 +980,7 @@ numeric asin(numeric const & x) /** Numeric inverse cosine (trigonometric function). * * @return arbitrary precision numerical acos(x). */ -numeric acos(numeric const & x) +numeric acos(const numeric & x) { return ::acos(*x.value); // -> CLN } @@ -990,7 +990,7 @@ numeric acos(numeric const & x) * @param z complex number * @return atan(z) * @exception overflow_error (logarithmic singularity) */ -numeric atan(numeric const & x) +numeric atan(const numeric & x) { if (!x.is_real() && x.real().is_zero() && @@ -1004,7 +1004,7 @@ numeric atan(numeric const & x) * @param x real number * @param y real number * @return atan(y/x) */ -numeric atan(numeric const & y, numeric const & x) +numeric atan(const numeric & y, const numeric & x) { if (x.is_real() && y.is_real()) return ::atan(realpart(*x.value), realpart(*y.value)); // -> CLN @@ -1015,7 +1015,7 @@ numeric atan(numeric const & y, numeric const & x) /** Numeric hyperbolic sine (trigonometric function). * * @return arbitrary precision numerical sinh(x). */ -numeric sinh(numeric const & x) +numeric sinh(const numeric & x) { return ::sinh(*x.value); // -> CLN } @@ -1023,7 +1023,7 @@ numeric sinh(numeric const & x) /** Numeric hyperbolic cosine (trigonometric function). * * @return arbitrary precision numerical cosh(x). */ -numeric cosh(numeric const & x) +numeric cosh(const numeric & x) { return ::cosh(*x.value); // -> CLN } @@ -1031,7 +1031,7 @@ numeric cosh(numeric const & x) /** Numeric hyperbolic tangent (trigonometric function). * * @return arbitrary precision numerical tanh(x). */ -numeric tanh(numeric const & x) +numeric tanh(const numeric & x) { return ::tanh(*x.value); // -> CLN } @@ -1039,7 +1039,7 @@ numeric tanh(numeric const & x) /** Numeric inverse hyperbolic sine (trigonometric function). * * @return arbitrary precision numerical asinh(x). */ -numeric asinh(numeric const & x) +numeric asinh(const numeric & x) { return ::asinh(*x.value); // -> CLN } @@ -1047,7 +1047,7 @@ numeric asinh(numeric const & x) /** Numeric inverse hyperbolic cosine (trigonometric function). * * @return arbitrary precision numerical acosh(x). */ -numeric acosh(numeric const & x) +numeric acosh(const numeric & x) { return ::acosh(*x.value); // -> CLN } @@ -1055,14 +1055,14 @@ numeric acosh(numeric const & x) /** Numeric inverse hyperbolic tangent (trigonometric function). * * @return arbitrary precision numerical atanh(x). */ -numeric atanh(numeric const & x) +numeric atanh(const numeric & x) { return ::atanh(*x.value); // -> CLN } /** Numeric evaluation of Riemann's Zeta function. Currently works only for * integer arguments. */ -numeric zeta(numeric const & x) +numeric zeta(const numeric & x) { // A dirty hack to allow for things like zeta(3.0), since CLN currently // only knows about integer arguments and zeta(3).evalf() automatically @@ -1082,7 +1082,7 @@ numeric zeta(numeric const & x) /** The gamma function. * This is only a stub! */ -numeric gamma(numeric const & x) +numeric gamma(const numeric & x) { clog << "gamma(" << x << "): Does anybody know good way to calculate this numerically?" @@ -1092,7 +1092,7 @@ numeric gamma(numeric const & x) /** The psi function (aka polygamma function). * This is only a stub! */ -numeric psi(numeric const & x) +numeric psi(const numeric & x) { clog << "psi(" << x << "): Does anybody know good way to calculate this numerically?" @@ -1102,7 +1102,7 @@ numeric psi(numeric const & x) /** The psi functions (aka polygamma functions). * This is only a stub! */ -numeric psi(numeric const & n, numeric const & x) +numeric psi(const numeric & n, const numeric & x) { clog << "psi(" << n << "," << x << "): Does anybody know good way to calculate this numerically?" @@ -1113,7 +1113,7 @@ numeric psi(numeric const & n, numeric const & x) /** Factorial combinatorial function. * * @exception range_error (argument must be integer >= 0) */ -numeric factorial(numeric const & nn) +numeric factorial(const numeric & nn) { if (!nn.is_nonneg_integer()) throw (std::range_error("numeric::factorial(): argument must be integer >= 0")); @@ -1126,7 +1126,7 @@ numeric factorial(numeric const & nn) * @param n integer argument >= -1 * @return n!! == n * (n-2) * (n-4) * ... * ({1|2}) with 0!! == 1 == (-1)!! * @exception range_error (argument must be integer >= -1) */ -numeric doublefactorial(numeric const & nn) +numeric doublefactorial(const numeric & nn) { // META-NOTE: The whole shit here will become obsolete and may be moved // out once CLN learns about double factorial, which should be as soon as @@ -1190,7 +1190,7 @@ numeric doublefactorial(numeric const & nn) * integer n and k and positive n this is the number of ways of choosing k * objects from n distinct objects. If n is negative, the formula * binomial(n,k) == (-1)^k*binomial(k-n-1,k) is used to compute the result. */ -numeric binomial(numeric const & n, numeric const & k) +numeric binomial(const numeric & n, const numeric & k) { if (n.is_integer() && k.is_integer()) { if (n.is_nonneg_integer()) { @@ -1212,7 +1212,7 @@ numeric binomial(numeric const & n, numeric const & k) * * @return the nth Bernoulli number (a rational number). * @exception range_error (argument must be integer >= 0) */ -numeric bernoulli(numeric const & nn) +numeric bernoulli(const numeric & nn) { if (!nn.is_integer() || nn.is_negative()) throw (std::range_error("numeric::bernoulli(): argument must be integer >= 0")); @@ -1250,7 +1250,7 @@ numeric bernoulli(numeric const & nn) } /** Absolute value. */ -numeric abs(numeric const & x) +numeric abs(const numeric & x) { return ::abs(*x.value); // -> CLN } @@ -1262,7 +1262,7 @@ numeric abs(numeric const & x) * * @return a mod b in the range [0,abs(b)-1] with sign of b if both are * integer, 0 otherwise. */ -numeric mod(numeric const & a, numeric const & b) +numeric mod(const numeric & a, const numeric & b) { if (a.is_integer() && b.is_integer()) return ::mod(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN @@ -1274,7 +1274,7 @@ numeric mod(numeric const & a, numeric const & b) * Equivalent to Maple's mods. * * @return a mod b in the range [-iquo(abs(m)-1,2), iquo(abs(m),2)]. */ -numeric smod(numeric const & a, numeric const & b) +numeric smod(const numeric & a, const numeric & b) { // FIXME: Should this become a member function? if (a.is_integer() && b.is_integer()) { @@ -1290,7 +1290,7 @@ numeric smod(numeric const & a, numeric const & b) * sign of a or is zero. * * @return remainder of a/b if both are integer, 0 otherwise. */ -numeric irem(numeric const & a, numeric const & b) +numeric irem(const numeric & a, const numeric & b) { if (a.is_integer() && b.is_integer()) return ::rem(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN @@ -1305,7 +1305,7 @@ numeric irem(numeric const & a, numeric const & b) * * @return remainder of a/b and quotient stored in q if both are integer, * 0 otherwise. */ -numeric irem(numeric const & a, numeric const & b, numeric & q) +numeric irem(const numeric & a, const numeric & b, numeric & q) { if (a.is_integer() && b.is_integer()) { // -> CLN cl_I_div_t rem_quo = truncate2(The(cl_I)(*a.value), The(cl_I)(*b.value)); @@ -1322,7 +1322,7 @@ numeric irem(numeric const & a, numeric const & b, numeric & q) * Equivalent to Maple's iquo as far as sign conventions are concerned. * * @return truncated quotient of a/b if both are integer, 0 otherwise. */ -numeric iquo(numeric const & a, numeric const & b) +numeric iquo(const numeric & a, const numeric & b) { if (a.is_integer() && b.is_integer()) return truncate1(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN @@ -1336,7 +1336,7 @@ numeric iquo(numeric const & a, numeric const & b) * * @return truncated quotient of a/b and remainder stored in r if both are * integer, 0 otherwise. */ -numeric iquo(numeric const & a, numeric const & b, numeric & r) +numeric iquo(const numeric & a, const numeric & b, numeric & r) { if (a.is_integer() && b.is_integer()) { // -> CLN cl_I_div_t rem_quo = truncate2(The(cl_I)(*a.value), The(cl_I)(*b.value)); @@ -1356,13 +1356,13 @@ numeric iquo(numeric const & a, numeric const & b, numeric & r) * @return square root of z. Branch cut along negative real axis, the negative * real axis itself where imag(z)==0 and real(z)<0 belongs to the upper part * where imag(z)>0. */ -numeric sqrt(numeric const & z) +numeric sqrt(const numeric & z) { return ::sqrt(*z.value); // -> CLN } /** Integer numeric square root. */ -numeric isqrt(numeric const & x) +numeric isqrt(const numeric & x) { if (x.is_integer()) { cl_I root; @@ -1376,7 +1376,7 @@ numeric isqrt(numeric const & x) * * @return The GCD of two numbers if both are integer, a numerical 1 * if they are not. */ -numeric gcd(numeric const & a, numeric const & b) +numeric gcd(const numeric & a, const numeric & b) { if (a.is_integer() && b.is_integer()) return ::gcd(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN @@ -1388,7 +1388,7 @@ numeric gcd(numeric const & a, numeric const & b) * * @return The LCM of two numbers if both are integer, the product of those * two numbers if they are not. */ -numeric lcm(numeric const & a, numeric const & b) +numeric lcm(const numeric & a, const numeric & b) { if (a.is_integer() && b.is_integer()) return ::lcm(The(cl_I)(*a.value), The(cl_I)(*b.value)); // -> CLN diff --git a/ginac/numeric.h b/ginac/numeric.h index 27153df3..30d4edf3 100644 --- a/ginac/numeric.h +++ b/ginac/numeric.h @@ -64,34 +64,34 @@ private: class numeric : public basic { // friends - friend numeric exp(numeric const & x); - friend numeric log(numeric const & x); - friend numeric sin(numeric const & x); - friend numeric cos(numeric const & x); - friend numeric tan(numeric const & x); - friend numeric asin(numeric const & x); - friend numeric acos(numeric const & x); - friend numeric atan(numeric const & x); - friend numeric atan(numeric const & y, numeric const & x); - friend numeric sinh(numeric const & x); - friend numeric cosh(numeric const & x); - friend numeric tanh(numeric const & x); - 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); - friend numeric smod(numeric const & a, numeric const & b); - friend numeric irem(numeric const & a, numeric const & b); - friend numeric irem(numeric const & a, numeric const & b, numeric & q); - friend numeric iquo(numeric const & a, numeric const & b); - friend numeric iquo(numeric const & a, numeric const & b, numeric & r); - friend numeric sqrt(numeric const & x); - friend numeric isqrt(numeric const & x); - friend numeric gcd(numeric const & a, numeric const & b); - friend numeric lcm(numeric const & a, numeric const & b); + friend numeric exp(const numeric & x); + friend numeric log(const numeric & x); + friend numeric sin(const numeric & x); + friend numeric cos(const numeric & x); + friend numeric tan(const numeric & x); + friend numeric asin(const numeric & x); + friend numeric acos(const numeric & x); + friend numeric atan(const numeric & x); + friend numeric atan(const numeric & y, const numeric & x); + friend numeric sinh(const numeric & x); + friend numeric cosh(const numeric & x); + friend numeric tanh(const numeric & x); + friend numeric asinh(const numeric & x); + friend numeric acosh(const numeric & x); + friend numeric atanh(const numeric & x); + friend numeric zeta(const numeric & x); + friend numeric bernoulli(const numeric & n); + friend numeric abs(const numeric & x); + friend numeric mod(const numeric & a, const numeric & b); + friend numeric smod(const numeric & a, const numeric & b); + friend numeric irem(const numeric & a, const numeric & b); + friend numeric irem(const numeric & a, const numeric & b, numeric & q); + friend numeric iquo(const numeric & a, const numeric & b); + friend numeric iquo(const numeric & a, const numeric & b, numeric & r); + friend numeric sqrt(const numeric & x); + friend numeric isqrt(const numeric & x); + friend numeric gcd(const numeric & a, const numeric & b); + friend numeric lcm(const numeric & a, const numeric & b); // member functions @@ -100,10 +100,10 @@ class numeric : public basic public: numeric(); ~numeric(); - numeric(numeric const & other); - numeric const & operator=(numeric const & other); + numeric(const numeric & other); + const numeric & operator=(const numeric & other); protected: - void copy(numeric const & other); + void copy(const numeric & other); void destroy(bool call_parent); // other constructors @@ -129,7 +129,7 @@ public: ex diff(symbol const & s) const; ex normal(lst &sym_lst, lst &repl_lst, int level=0) const; numeric integer_content(void) const; - ex smod(numeric const &xi) const; + ex smod(const numeric &xi) const; numeric max_coefficient(void) const; protected: int compare_same_type(basic const & other) const; @@ -144,33 +144,26 @@ protected: // non-virtual functions in this class public: - numeric add(numeric const & other) const; - numeric sub(numeric const & other) const; - numeric mul(numeric const & other) const; - numeric div(numeric const & other) const; - numeric power(numeric const & other) const; - numeric const & add_dyn(numeric const & other) const; - numeric const & sub_dyn(numeric const & other) const; - numeric const & mul_dyn(numeric const & other) const; - numeric const & div_dyn(numeric const & other) const; - numeric const & power_dyn(numeric const & other) const; - numeric const & operator=(int i); - numeric const & operator=(unsigned int i); - numeric const & operator=(long i); - numeric const & operator=(unsigned long i); - numeric const & operator=(double d); - numeric const & operator=(char const * s); - /* - numeric add_dyn(numeric const & other) const { return add(other); } - numeric sub_dyn(numeric const & other) const { return sub(other); } - numeric mul_dyn(numeric const & other) const { return mul(other); } - numeric div_dyn(numeric const & other) const { return div(other); } - numeric power_dyn(numeric const & other) const { return power(other); } - */ + numeric add(const numeric & other) const; + numeric sub(const numeric & other) const; + numeric mul(const numeric & other) const; + numeric div(const numeric & other) const; + numeric power(const numeric & other) const; + const numeric & add_dyn(const numeric & other) const; + const numeric & sub_dyn(const numeric & other) const; + const numeric & mul_dyn(const numeric & other) const; + const numeric & div_dyn(const numeric & other) const; + const numeric & power_dyn(const numeric & other) const; + const numeric & operator=(int i); + const numeric & operator=(unsigned int i); + const numeric & operator=(long i); + const numeric & operator=(unsigned long i); + const numeric & operator=(double d); + const numeric & operator=(char const * s); numeric inverse(void) const; int csgn(void) const; - int compare(numeric const & other) const; - bool is_equal(numeric const & other) const; + int compare(const numeric & other) const; + bool is_equal(const numeric & other) const; bool is_zero(void) const; bool is_positive(void) const; bool is_negative(void) const; @@ -184,12 +177,12 @@ public: bool is_real(void) const; bool is_cinteger(void) const; bool is_crational(void) const; - bool operator==(numeric const & other) const; - bool operator!=(numeric const & other) const; - bool operator<(numeric const & other) const; - bool operator<=(numeric const & other) const; - bool operator>(numeric const & other) const; - bool operator>=(numeric const & other) const; + bool operator==(const numeric & other) const; + bool operator!=(const numeric & other) const; + bool operator<(const numeric & other) const; + bool operator<=(const numeric & other) const; + bool operator>(const numeric & other) const; + bool operator>=(const numeric & other) const; int to_int(void) const; double to_double(void) const; numeric real(void) const; @@ -217,42 +210,42 @@ extern _numeric_digits Digits; // global functions -numeric exp(numeric const & x); -numeric log(numeric const & x); -numeric sin(numeric const & x); -numeric cos(numeric const & x); -numeric tan(numeric const & x); -numeric asin(numeric const & x); -numeric acos(numeric const & x); -numeric atan(numeric const & x); -numeric atan(numeric const & y, numeric const & x); -numeric sinh(numeric const & x); -numeric cosh(numeric const & x); -numeric tanh(numeric const & x); -numeric asinh(numeric const & x); -numeric acosh(numeric const & x); -numeric atanh(numeric const & x); -numeric zeta(numeric const & x); -numeric gamma(numeric const & x); -numeric psi(numeric const & x); -numeric psi(numeric const & n, numeric const & x); -numeric factorial(numeric const & n); -numeric doublefactorial(numeric const & n); -numeric binomial(numeric const & n, numeric const & k); -numeric bernoulli(numeric const & n); - -numeric abs(numeric const & x); -numeric mod(numeric const & a, numeric const & b); -numeric smod(numeric const & a, numeric const & b); -numeric irem(numeric const & a, numeric const & b); -numeric irem(numeric const & a, numeric const & b, numeric & q); -numeric iquo(numeric const & a, numeric const & b); -numeric iquo(numeric const & a, numeric const & b, numeric & r); -numeric sqrt(numeric const & x); -numeric isqrt(numeric const & x); - -numeric gcd(numeric const & a, numeric const & b); -numeric lcm(numeric const & a, numeric const & b); +numeric exp(const numeric & x); +numeric log(const numeric & x); +numeric sin(const numeric & x); +numeric cos(const numeric & x); +numeric tan(const numeric & x); +numeric asin(const numeric & x); +numeric acos(const numeric & x); +numeric atan(const numeric & x); +numeric atan(const numeric & y, const numeric & x); +numeric sinh(const numeric & x); +numeric cosh(const numeric & x); +numeric tanh(const numeric & x); +numeric asinh(const numeric & x); +numeric acosh(const numeric & x); +numeric atanh(const numeric & x); +numeric zeta(const numeric & x); +numeric gamma(const numeric & x); +numeric psi(const numeric & x); +numeric psi(const numeric & n, const numeric & x); +numeric factorial(const numeric & n); +numeric doublefactorial(const numeric & n); +numeric binomial(const numeric & n, const numeric & k); +numeric bernoulli(const numeric & n); + +numeric abs(const numeric & x); +numeric mod(const numeric & a, const numeric & b); +numeric smod(const numeric & a, const numeric & b); +numeric irem(const numeric & a, const numeric & b); +numeric irem(const numeric & a, const numeric & b, numeric & q); +numeric iquo(const numeric & a, const numeric & b); +numeric iquo(const numeric & a, const numeric & b, numeric & r); +numeric sqrt(const numeric & x); +numeric isqrt(const numeric & x); + +numeric gcd(const numeric & a, const numeric & b); +numeric lcm(const numeric & a, const numeric & b); /** Exception thrown by numeric members to signal failure */ struct numeric_fail @@ -262,58 +255,61 @@ struct numeric_fail }; // wrapper functions around member functions -inline numeric inverse(numeric const & x) +inline numeric pow(const numeric & x, const numeric & y) +{ return x.power(y); } + +inline numeric inverse(const numeric & x) { return x.inverse(); } -inline bool csgn(numeric const & x) +inline bool csgn(const numeric & x) { return x.csgn(); } -inline bool is_zero(numeric const & x) +inline bool is_zero(const numeric & x) { return x.is_zero(); } -inline bool is_positive(numeric const & x) +inline bool is_positive(const numeric & x) { return x.is_positive(); } -inline bool is_integer(numeric const & x) +inline bool is_integer(const numeric & x) { return x.is_integer(); } -inline bool is_pos_integer(numeric const & x) +inline bool is_pos_integer(const numeric & x) { return x.is_pos_integer(); } -inline bool is_nonneg_integer(numeric const & x) +inline bool is_nonneg_integer(const numeric & x) { return x.is_nonneg_integer(); } -inline bool is_even(numeric const & x) +inline bool is_even(const numeric & x) { return x.is_even(); } -inline bool is_odd(numeric const & x) +inline bool is_odd(const numeric & x) { return x.is_odd(); } -inline bool is_prime(numeric const & x) +inline bool is_prime(const numeric & x) { return x.is_prime(); } -inline bool is_rational(numeric const & x) +inline bool is_rational(const numeric & x) { return x.is_rational(); } -inline bool is_real(numeric const & x) +inline bool is_real(const numeric & x) { return x.is_real(); } -inline bool is_cinteger(numeric const & x) +inline bool is_cinteger(const numeric & x) { return x.is_cinteger(); } -inline bool is_crational(numeric const & x) +inline bool is_crational(const numeric & x) { return x.is_crational(); } -inline numeric real(numeric const & x) +inline numeric real(const numeric & x) { return x.real(); } -inline numeric imag(numeric const & x) +inline numeric imag(const numeric & x) { return x.imag(); } -inline numeric numer(numeric const & x) +inline numeric numer(const numeric & x) { return x.numer(); } -inline numeric denom(numeric const & x) +inline numeric denom(const numeric & x) { return x.denom(); } // numeric evaluation functions for class constant objects: diff --git a/ginac/operators.h b/ginac/operators.h index a60da5d6..f226ff1e 100644 --- a/ginac/operators.h +++ b/ginac/operators.h @@ -41,10 +41,10 @@ ex operator/(ex const & lh, ex const & rh); ex operator%(ex 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); -numeric operator*(numeric const & lh, numeric const & rh); -numeric operator/(numeric const & lh, numeric const & rh); +numeric operator+(const numeric & lh, const numeric & rh); +numeric operator-(const numeric & lh, const numeric & rh); +numeric operator*(const numeric & lh, const numeric & rh); +numeric operator/(const numeric & lh, const numeric & rh); // binary arithmetic assignment operators with ex ex const & operator+=(ex & lh, ex const & rh); @@ -54,17 +54,17 @@ 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 -numeric const & operator+=(numeric & lh, numeric const & rh); -numeric const & operator-=(numeric & lh, numeric const & rh); -numeric const & operator*=(numeric & lh, numeric const & rh); -numeric const & operator/=(numeric & lh, numeric const & rh); +const numeric & operator+=(numeric & lh, const numeric & rh); +const numeric & operator-=(numeric & lh, const numeric & rh); +const numeric & operator*=(numeric & lh, const numeric & rh); +const numeric & operator/=(numeric & lh, const numeric & rh); // unary operators ex operator+(ex const & lh); ex operator-(ex const & lh); -numeric operator+(numeric const & lh); -numeric operator-(numeric const & lh); +numeric operator+(const numeric & lh); +numeric operator-(const numeric & lh); numeric& operator++(numeric & rh); numeric& operator--(numeric & rh); numeric operator++(numeric & lh, int); diff --git a/ginac/power.cpp b/ginac/power.cpp index 51e5bfc5..dfd1a470 100644 --- a/ginac/power.cpp +++ b/ginac/power.cpp @@ -120,11 +120,15 @@ basic * power::duplicate() const 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 << ")"; + if (exponent.is_equal(_ex1_2())) { + os << "sqrt(" << basis << ")"; + } else { + 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 diff --git a/ginac/utils.cpp b/ginac/utils.cpp index e113c04b..1a529363 100644 --- a/ginac/utils.cpp +++ b/ginac/utils.cpp @@ -53,6 +53,118 @@ int compare_pointers(void const * a, void const * b) // `construct on first use' chest of numbers ////////// +// numeric -60 +numeric const & _num_60(void) +{ + const static ex e = ex(numeric(-60)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_60(void) +{ + static ex * e = new ex(_num_60()); + return *e; +} + +// numeric -120 +numeric const & _num_120(void) +{ + const static ex e = ex(numeric(-120)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_120(void) +{ + static ex * e = new ex(_num_120()); + return *e; +} + +// numeric -30 +numeric const & _num_30(void) +{ + const static ex e = ex(numeric(-30)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_30(void) +{ + static ex * e = new ex(_num_30()); + return *e; +} + +// numeric -25 +numeric const & _num_25(void) +{ + const static ex e = ex(numeric(-25)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_25(void) +{ + static ex * e = new ex(_num_25()); + return *e; +} + +// numeric -24 +numeric const & _num_24(void) +{ + const static ex e = ex(numeric(-24)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_24(void) +{ + static ex * e = new ex(_num_24()); + return *e; +} + +// numeric -20 +numeric const & _num_20(void) +{ + const static ex e = ex(numeric(-20)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_20(void) +{ + static ex * e = new ex(_num_20()); + return *e; +} + +// numeric -18 +numeric const & _num_18(void) +{ + const static ex e = ex(numeric(-18)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_18(void) +{ + static ex * e = new ex(_num_18()); + return *e; +} + +// numeric -15 +numeric const & _num_15(void) +{ + const static ex e = ex(numeric(-15)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_15(void) +{ + static ex * e = new ex(_num_15()); + return *e; +} + // numeric -12 numeric const & _num_12(void) { @@ -249,6 +361,20 @@ ex const & _ex_1_3(void) return *e; } +// numeric -1/4 +numeric const & _num_1_4(void) +{ + const static ex e = ex(numeric(-1,4)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex_1_4(void) +{ + static ex * e = new ex(_num_1_4()); + return *e; +} + // numeric 0 numeric const & _num0(void) { @@ -263,6 +389,20 @@ ex const & _ex0(void) return *e; } +// numeric 1/4 +numeric const & _num1_4(void) +{ + const static ex e = ex(numeric(1,4)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex1_4(void) +{ + static ex * e = new ex(_num1_4()); + return *e; +} + // numeric 1/3 numeric const & _num1_3(void) { @@ -459,6 +599,118 @@ ex const & _ex12(void) return *e; } +// numeric 15 +numeric const & _num15(void) +{ + const static ex e = ex(numeric(15)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex15(void) +{ + static ex * e = new ex(_num15()); + return *e; +} + +// numeric 18 +numeric const & _num18(void) +{ + const static ex e = ex(numeric(18)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex18(void) +{ + static ex * e = new ex(_num18()); + return *e; +} + +// numeric 20 +numeric const & _num20(void) +{ + const static ex e = ex(numeric(20)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex20(void) +{ + static ex * e = new ex(_num20()); + return *e; +} + +// numeric 24 +numeric const & _num24(void) +{ + const static ex e = ex(numeric(24)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex24(void) +{ + static ex * e = new ex(_num24()); + return *e; +} + +// numeric 25 +numeric const & _num25(void) +{ + const static ex e = ex(numeric(25)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex25(void) +{ + static ex * e = new ex(_num25()); + return *e; +} + +// numeric 30 +numeric const & _num30(void) +{ + const static ex e = ex(numeric(30)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex30(void) +{ + static ex * e = new ex(_num30()); + return *e; +} + +// numeric 60 +numeric const & _num60(void) +{ + const static ex e = ex(numeric(60)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex60(void) +{ + static ex * e = new ex(_num60()); + return *e; +} + +// numeric 120 +numeric const & _num120(void) +{ + const static ex e = ex(numeric(120)); + const static numeric * n = static_cast(e.bp); + return *n; +} + +ex const & _ex120(void) +{ + static ex * e = new ex(_num120()); + return *e; +} + // comment skeleton for header files diff --git a/ginac/utils.h b/ginac/utils.h index 0f6082a6..133e9f67 100644 --- a/ginac/utils.h +++ b/ginac/utils.h @@ -146,65 +146,101 @@ OutputIterator mymerge3(InputIterator1 first1, InputIterator1 last1, class numeric; class ex; - + +numeric const & _num_120(void); // -120 +ex const & _ex_120(void); +numeric const & _num_60(void); // -60 +ex const & _ex_60(void); +numeric const & _num_30(void); // -30 +ex const & _ex_30(void); +numeric const & _num_25(void); // -25 +ex const & _ex_25(void); +numeric const & _num_24(void); // -24 +ex const & _ex_24(void); +numeric const & _num_20(void); // -20 +ex const & _ex_20(void); +numeric const & _num_18(void); // -18 +ex const & _ex_18(void); +numeric const & _num_15(void); // -15 +ex const & _ex_15(void); numeric const & _num_12(void); // -12 +ex const & _ex_12(void); numeric const & _num_11(void); // -11 +ex const & _ex_11(void); numeric const & _num_10(void); // -10 +ex const & _ex_10(void); numeric const & _num_9(void); // -9 +ex const & _ex_9(void); numeric const & _num_8(void); // -8 +ex const & _ex_8(void); numeric const & _num_7(void); // -7 +ex const & _ex_7(void); numeric const & _num_6(void); // -6 +ex const & _ex_6(void); numeric const & _num_5(void); // -5 +ex const & _ex_5(void); numeric const & _num_4(void); // -4 +ex const & _ex_4(void); numeric const & _num_3(void); // -3 +ex const & _ex_3(void); numeric const & _num_2(void); // -2 +ex const & _ex_2(void); numeric const & _num_1(void); // -1 +ex const & _ex_1(void); numeric const & _num_1_2(void); // -1/2 +ex const & _ex_1_2(void); numeric const & _num_1_3(void); // -1/3 +ex const & _ex_1_3(void); +numeric const & _num_1_4(void); // -1/4 +ex const & _ex_1_4(void); numeric const & _num0(void); // 0 +ex const & _ex0(void); +numeric const & _num1_4(void); // 1/4 +ex const & _ex1_4(void); numeric const & _num1_3(void); // 1/3 +ex const & _ex1_3(void); numeric const & _num1_2(void); // 1/2 +ex const & _ex1_2(void); numeric const & _num1(void); // 1 +ex const & _ex1(void); numeric const & _num2(void); // 2 +ex const & _ex2(void); numeric const & _num3(void); // 3 +ex const & _ex3(void); numeric const & _num4(void); // 4 +ex const & _ex4(void); numeric const & _num5(void); // 5 +ex const & _ex5(void); numeric const & _num6(void); // 6 +ex const & _ex6(void); numeric const & _num7(void); // 7 +ex const & _ex7(void); numeric const & _num8(void); // 8 +ex const & _ex8(void); numeric const & _num9(void); // 9 +ex const & _ex9(void); numeric const & _num10(void); // 10 +ex const & _ex10(void); numeric const & _num11(void); // 11 +ex const & _ex11(void); numeric const & _num12(void); // 12 -ex const & _ex_12(void); // -12 -ex const & _ex_11(void); // -11 -ex const & _ex_10(void); // -10 -ex const & _ex_9(void); // -9 -ex const & _ex_8(void); // -8 -ex const & _ex_7(void); // -7 -ex const & _ex_6(void); // -6 -ex const & _ex_5(void); // -5 -ex const & _ex_4(void); // -4 -ex const & _ex_3(void); // -3 -ex const & _ex_2(void); // -2 -ex const & _ex_1(void); // -1 -ex const & _ex_1_2(void); // -1/2 -ex const & _ex_1_3(void); // -1/3 -ex const & _ex0(void); // 0 -ex const & _ex1_3(void); // 1/3 -ex const & _ex1_2(void); // 1/2 -ex const & _ex1(void); // 1 -ex const & _ex2(void); // 2 -ex const & _ex3(void); // 3 -ex const & _ex4(void); // 4 -ex const & _ex5(void); // 5 -ex const & _ex6(void); // 6 -ex const & _ex7(void); // 7 -ex const & _ex8(void); // 8 -ex const & _ex9(void); // 9 -ex const & _ex10(void); // 10 -ex const & _ex11(void); // 11 -ex const & _ex12(void); // 12 +ex const & _ex12(void); +numeric const & _num15(void); // 15 +ex const & _ex15(void); +numeric const & _num18(void); // 18 +ex const & _ex18(void); +numeric const & _num20(void); // 20 +ex const & _ex20(void); +numeric const & _num24(void); // 24 +ex const & _ex24(void); +numeric const & _num25(void); // 25 +ex const & _ex25(void); +numeric const & _num30(void); // 30 +ex const & _ex30(void); +numeric const & _num60(void); // 60 +ex const & _ex60(void); +numeric const & _num120(void); // 120 +ex const & _ex120(void); #ifndef NO_GINAC_NAMESPACE } // namespace GiNaC -- 2.44.0