X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=blobdiff_plain;f=ginac%2Fnumeric.cpp;h=3bd24a6a8571f8588ed868e531040c6e9999bdfa;hp=c4c420ddb7a0c67cbb1aa846f4577a11d40b247e;hb=094911eb78cacb6f2877a70c9ac74766df58ccea;hpb=252c08252bd3dcab99be3b44c86135a7cb4aa0a6 diff --git a/ginac/numeric.cpp b/ginac/numeric.cpp index c4c420dd..3bd24a6a 100644 --- a/ginac/numeric.cpp +++ b/ginac/numeric.cpp @@ -7,7 +7,7 @@ * of special functions or implement the interface to the bignum package. */ /* - * GiNaC Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany + * GiNaC Copyright (C) 1999-2001 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 @@ -63,14 +63,12 @@ #include #include -#ifndef NO_NAMESPACE_GINAC namespace GiNaC { -#endif // ndef NO_NAMESPACE_GINAC GINAC_IMPLEMENT_REGISTERED_CLASS(numeric, basic) ////////// -// default constructor, destructor, copy constructor assignment +// default ctor, dtor, copy ctor assignment // operator and helpers ////////// @@ -79,58 +77,34 @@ GINAC_IMPLEMENT_REGISTERED_CLASS(numeric, basic) /** default ctor. Numerically it initializes to an integer zero. */ numeric::numeric() : basic(TINFO_numeric) { - debugmsg("numeric default constructor", LOGLEVEL_CONSTRUCT); + debugmsg("numeric default ctor", LOGLEVEL_CONSTRUCT); value = cln::cl_I(0); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); -} - -numeric::~numeric() -{ - debugmsg("numeric destructor" ,LOGLEVEL_DESTRUCT); - destroy(false); -} - -numeric::numeric(const numeric & other) -{ - debugmsg("numeric copy constructor", LOGLEVEL_CONSTRUCT); - copy(other); -} - -const numeric & numeric::operator=(const numeric & other) -{ - debugmsg("numeric operator=", LOGLEVEL_ASSIGNMENT); - if (this != &other) { - destroy(true); - copy(other); - } - return *this; + setflag(status_flags::evaluated | status_flags::expanded); } // protected -void numeric::copy(const numeric & other) +/** For use by copy ctor and assignment operator. */ +void numeric::copy(const numeric &other) { - basic::copy(other); + inherited::copy(other); value = other.value; } void numeric::destroy(bool call_parent) { - if (call_parent) basic::destroy(call_parent); + if (call_parent) inherited::destroy(call_parent); } ////////// -// other constructors +// other ctors ////////// // public numeric::numeric(int i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from int",LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor 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. However, if the integer is small enough, @@ -140,16 +114,13 @@ numeric::numeric(int i) : basic(TINFO_numeric) value = cln::cl_I(i); else value = cln::cl_I((long) i); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } numeric::numeric(unsigned int i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from uint",LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor 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. However, if the integer is small enough, @@ -159,32 +130,23 @@ numeric::numeric(unsigned int i) : basic(TINFO_numeric) value = cln::cl_I(i); else value = cln::cl_I((unsigned long) i); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } numeric::numeric(long i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from long",LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor from long",LOGLEVEL_CONSTRUCT); value = cln::cl_I(i); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } numeric::numeric(unsigned long i) : basic(TINFO_numeric) { - debugmsg("numeric constructor from ulong",LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor from ulong",LOGLEVEL_CONSTRUCT); value = cln::cl_I(i); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } /** Ctor for rational numerics a/b. @@ -192,35 +154,30 @@ 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("numeric ctor from long/long",LOGLEVEL_CONSTRUCT); if (!denom) throw std::overflow_error("division by zero"); value = cln::cl_I(numer) / cln::cl_I(denom); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } numeric::numeric(double d) : basic(TINFO_numeric) { - debugmsg("numeric constructor from double",LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor 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: value = cln::cl_float(d, cln::default_float_format); - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } + /** ctor from C-style string. It also accepts complex numbers in GiNaC * notation like "2+5*I". */ numeric::numeric(const char *s) : basic(TINFO_numeric) { - debugmsg("numeric constructor from string",LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor from string",LOGLEVEL_CONSTRUCT); cln::cl_N ctorval = 0; // parse complex numbers (functional but not completely safe, unfortunately // std::string does not understand regexpese): @@ -292,22 +249,17 @@ numeric::numeric(const char *s) : basic(TINFO_numeric) } } while(delim != std::string::npos); value = ctorval; - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } + /** Ctor from CLN types. This is for the initiated user or internal use * only. */ -numeric::numeric(const cln::cl_N & z) : basic(TINFO_numeric) +numeric::numeric(const cln::cl_N &z) : basic(TINFO_numeric) { - debugmsg("numeric constructor from cl_N", LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor from cl_N", LOGLEVEL_CONSTRUCT); value = z; - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } ////////// @@ -317,7 +269,7 @@ numeric::numeric(const cln::cl_N & z) : basic(TINFO_numeric) /** Construct object from archive_node. */ numeric::numeric(const archive_node &n, const lst &sym_lst) : inherited(n, sym_lst) { - debugmsg("numeric constructor from archive_node", LOGLEVEL_CONSTRUCT); + debugmsg("numeric ctor from archive_node", LOGLEVEL_CONSTRUCT); cln::cl_N ctorval = 0; // Read number as string @@ -349,10 +301,7 @@ numeric::numeric(const archive_node &n, const lst &sym_lst) : inherited(n, sym_l } } value = ctorval; - calchash(); - setflag(status_flags::evaluated | - status_flags::expanded | - status_flags::hash_calculated); + setflag(status_flags::evaluated | status_flags::expanded); } /** Unarchive the object. */ @@ -403,15 +352,6 @@ void numeric::archive(archive_node &n) const // functions overriding virtual functions from bases classes ////////// -// public - -basic * numeric::duplicate() const -{ - debugmsg("numeric duplicate", LOGLEVEL_DUPLICATE); - return new numeric(*this); -} - - /** Helper function to print a real number in a nicer way than is CLN's * default. Instead of printing 42.0L0 this just prints 42.0 to ostream os * and instead of 3.99168L7 it prints 3.99168E7. This is fine in GiNaC as @@ -419,7 +359,7 @@ basic * numeric::duplicate() const * want to visibly distinguish from cl_LF. * * @see numeric::print() */ -static void print_real_number(std::ostream & os, const cln::cl_R & num) +static void print_real_number(std::ostream &os, const cln::cl_R &num) { cln::cl_print_flags ourflags; if (cln::instanceof(num, cln::cl_RA_ring)) { @@ -439,7 +379,7 @@ static void print_real_number(std::ostream & os, const cln::cl_R & num) * with the other routines and produces something compatible to ginsh input. * * @see print_real_number() */ -void numeric::print(std::ostream & os, unsigned upper_precedence) const +void numeric::print(std::ostream &os, unsigned upper_precedence) const { debugmsg("numeric print", LOGLEVEL_PRINT); cln::cl_R r = cln::realpart(cln::the(value)); @@ -504,16 +444,16 @@ void numeric::print(std::ostream & os, unsigned upper_precedence) const } -void numeric::printraw(std::ostream & os) const +void numeric::printraw(std::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(" << cln::the(value) << ")"; + os << class_name() << "(" << cln::the(value) << ")"; } -void numeric::printtree(std::ostream & os, unsigned indent) const +void numeric::printtree(std::ostream &os, unsigned indent) const { debugmsg("numeric printtree", LOGLEVEL_PRINT); os << std::string(indent,' ') << cln::the(value) @@ -524,7 +464,7 @@ void numeric::printtree(std::ostream & os, unsigned indent) const } -void numeric::printcsrc(std::ostream & os, unsigned type, unsigned upper_precedence) const +void numeric::printcsrc(std::ostream &os, unsigned type, unsigned upper_precedence) const { debugmsg("numeric print csrc", LOGLEVEL_PRINT); std::ios::fmtflags oldflags = os.flags(); @@ -610,11 +550,11 @@ bool numeric::info(unsigned inf) const * results: (2+I).has(-2) -> true. But this is consistent, since we also * would like to have (-2+I).has(2) -> true and we want to think about the * sign as a multiplicative factor. */ -bool numeric::has(const ex & other) const +bool numeric::has(const ex &other) const { if (!is_exactly_of_type(*other.bp, numeric)) return false; - const numeric & o = static_cast(const_cast(*other.bp)); + const numeric &o = static_cast(const_cast(*other.bp)); if (this->is_equal(o) || this->is_equal(-o)) return true; if (o.imag().is_zero()) // e.g. scan for 3 in -3*I @@ -651,30 +591,21 @@ ex numeric::evalf(int level) const { // level can safely be discarded for numeric objects. return numeric(cln::cl_float(1.0, cln::default_float_format) * - (cln::the(value))); + (cln::the(value))); } // protected -/** Implementation of ex::diff() for a numeric. It always returns 0. - * - * @see ex::diff */ -ex numeric::derivative(const symbol & s) const -{ - return _ex0(); -} - - -int numeric::compare_same_type(const basic & other) const +int numeric::compare_same_type(const basic &other) const { GINAC_ASSERT(is_exactly_of_type(other, numeric)); - const numeric & o = static_cast(const_cast(other)); + const numeric &o = static_cast(const_cast(other)); return this->compare(o); } -bool numeric::is_equal_same_type(const basic & other) const +bool numeric::is_equal_same_type(const basic &other) const { GINAC_ASSERT(is_exactly_of_type(other,numeric)); const numeric *o = static_cast(&other); @@ -688,6 +619,7 @@ unsigned numeric::calchash(void) const // Use CLN's hashcode. Warning: It depends only on the number's value, not // its type or precision (i.e. a true equivalence relation on numbers). As // a consequence, 3 and 3.0 share the same hashvalue. + setflag(status_flags::hash_calculated); return (hashvalue = cln::equal_hashcode(cln::the(value)) | 0x80000000U); } @@ -705,8 +637,8 @@ unsigned numeric::calchash(void) const // public /** Numerical addition method. Adds argument to *this and returns result as - * a new numeric object. */ -const numeric numeric::add(const numeric & other) const + * a numeric object. */ +const numeric numeric::add(const numeric &other) const { // Efficiency shortcut: trap the neutral element by pointer. static const numeric * _num0p = &_num0(); @@ -720,16 +652,16 @@ const numeric numeric::add(const numeric & other) const /** Numerical subtraction method. Subtracts argument from *this and returns - * result as a new numeric object. */ -const numeric numeric::sub(const numeric & other) const + * result as a numeric object. */ +const numeric numeric::sub(const numeric &other) const { return numeric(cln::the(value)-cln::the(other.value)); } /** Numerical multiplication method. Multiplies *this and argument and returns - * result as a new numeric object. */ -const numeric numeric::mul(const numeric & other) const + * result as a numeric object. */ +const numeric numeric::mul(const numeric &other) const { // Efficiency shortcut: trap the neutral element by pointer. static const numeric * _num1p = &_num1(); @@ -743,10 +675,10 @@ const numeric numeric::mul(const numeric & other) const /** Numerical division method. Divides *this by argument and returns result as - * a new numeric object. + * a numeric object. * * @exception overflow_error (division by zero) */ -const numeric numeric::div(const numeric & other) const +const numeric numeric::div(const numeric &other) const { if (cln::zerop(cln::the(other.value))) throw std::overflow_error("numeric::div(): division by zero"); @@ -754,7 +686,9 @@ const numeric numeric::div(const numeric & other) const } -const numeric numeric::power(const numeric & other) const +/** Numerical exponentiation. Raises *this to the power given as argument and + * returns result as a numeric object. */ +const numeric numeric::power(const numeric &other) const { // Efficiency shortcut: trap the neutral exponent by pointer. static const numeric * _num1p = &_num1(); @@ -775,7 +709,7 @@ const numeric numeric::power(const numeric & other) const } -const numeric & numeric::add_dyn(const numeric & other) const +const numeric &numeric::add_dyn(const numeric &other) const { // Efficiency shortcut: trap the neutral element by pointer. static const numeric * _num0p = &_num0(); @@ -789,14 +723,14 @@ const numeric & numeric::add_dyn(const numeric & other) const } -const numeric & numeric::sub_dyn(const numeric & other) const +const numeric &numeric::sub_dyn(const numeric &other) const { return static_cast((new numeric(cln::the(value)-cln::the(other.value)))-> setflag(status_flags::dynallocated)); } -const numeric & numeric::mul_dyn(const numeric & other) const +const numeric &numeric::mul_dyn(const numeric &other) const { // Efficiency shortcut: trap the neutral element by pointer. static const numeric * _num1p = &_num1(); @@ -810,7 +744,7 @@ const numeric & numeric::mul_dyn(const numeric & other) const } -const numeric & numeric::div_dyn(const numeric & other) const +const numeric &numeric::div_dyn(const numeric &other) const { if (cln::zerop(cln::the(other.value))) throw std::overflow_error("division by zero"); @@ -819,7 +753,7 @@ const numeric & numeric::div_dyn(const numeric & other) const } -const numeric & numeric::power_dyn(const numeric & other) const +const numeric &numeric::power_dyn(const numeric &other) const { // Efficiency shortcut: trap the neutral exponent by pointer. static const numeric * _num1p=&_num1(); @@ -841,37 +775,37 @@ const numeric & numeric::power_dyn(const numeric & other) const } -const numeric & numeric::operator=(int i) +const numeric &numeric::operator=(int i) { return operator=(numeric(i)); } -const numeric & numeric::operator=(unsigned int i) +const numeric &numeric::operator=(unsigned int i) { return operator=(numeric(i)); } -const numeric & numeric::operator=(long i) +const numeric &numeric::operator=(long i) { return operator=(numeric(i)); } -const numeric & numeric::operator=(unsigned long i) +const numeric &numeric::operator=(unsigned long i) { return operator=(numeric(i)); } -const numeric & numeric::operator=(double d) +const numeric &numeric::operator=(double d) { return operator=(numeric(d)); } -const numeric & numeric::operator=(const char * s) +const numeric &numeric::operator=(const char * s) { return operator=(numeric(s)); } @@ -890,7 +824,7 @@ const numeric numeric::inverse(void) const * 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(const numeric & other) */ + * @see numeric::compare(const numeric &other) */ int numeric::csgn(void) const { if (cln::zerop(cln::the(value))) @@ -917,7 +851,7 @@ int numeric::csgn(void) const * * @return csgn(*this-other) * @see numeric::csgn(void) */ -int numeric::compare(const numeric & other) const +int numeric::compare(const numeric &other) const { // Comparing two real numbers? if (cln::instanceof(value, cln::cl_R_ring) && @@ -935,7 +869,7 @@ int numeric::compare(const numeric & other) const } -bool numeric::is_equal(const numeric & other) const +bool numeric::is_equal(const numeric &other) const { return cln::equal(cln::the(value),cln::the(other.value)); } @@ -1025,13 +959,13 @@ bool numeric::is_real(void) const } -bool numeric::operator==(const numeric & other) const +bool numeric::operator==(const numeric &other) const { return equal(cln::the(value), cln::the(other.value)); } -bool numeric::operator!=(const numeric & other) const +bool numeric::operator!=(const numeric &other) const { return !equal(cln::the(value), cln::the(other.value)); } @@ -1070,7 +1004,7 @@ bool numeric::is_crational(void) const /** Numerical comparison: less. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator<(const numeric & other) const +bool numeric::operator<(const numeric &other) const { if (this->is_real() && other.is_real()) return (cln::the(value) < cln::the(other.value)); @@ -1081,7 +1015,7 @@ bool numeric::operator<(const numeric & other) const /** Numerical comparison: less or equal. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator<=(const numeric & other) const +bool numeric::operator<=(const numeric &other) const { if (this->is_real() && other.is_real()) return (cln::the(value) <= cln::the(other.value)); @@ -1092,7 +1026,7 @@ bool numeric::operator<=(const numeric & other) const /** Numerical comparison: greater. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator>(const numeric & other) const +bool numeric::operator>(const numeric &other) const { if (this->is_real() && other.is_real()) return (cln::the(value) > cln::the(other.value)); @@ -1103,7 +1037,7 @@ bool numeric::operator>(const numeric & other) const /** Numerical comparison: greater or equal. * * @exception invalid_argument (complex inequality) */ -bool numeric::operator>=(const numeric & other) const +bool numeric::operator>=(const numeric &other) const { if (this->is_real() && other.is_real()) return (cln::the(value) >= cln::the(other.value)); @@ -1140,6 +1074,15 @@ double numeric::to_double(void) const } +/** Returns a new CLN object of type cl_N, representing the value of *this. + * This method may be used when mixing GiNaC and CLN in one project. + */ +cln::cl_N numeric::to_cl_N(void) const +{ + return cln::cl_N(cln::the(value)); +} + + /** Real part of a number. */ const numeric numeric::real(void) const { @@ -1229,15 +1172,6 @@ int numeric::int_length(void) const } -/** Returns a new CLN object of type cl_N, representing the value of *this. - * This method is useful for casting when mixing GiNaC and CLN in one project. - */ -numeric::operator cln::cl_N() const -{ - return cln::cl_N(cln::the(value)); -} - - ////////// // static member variables ////////// @@ -1250,19 +1184,18 @@ unsigned numeric::precedence = 30; // global constants ////////// -const numeric some_numeric; -const std::type_info & typeid_numeric = typeid(some_numeric); /** Imaginary unit. This is not a constant but a numeric since we are - * natively handing complex numbers anyways. */ + * natively handing complex numbers anyways, so in each expression containing + * an I it is automatically eval'ed away anyhow. */ const numeric I = numeric(cln::complex(cln::cl_I(0),cln::cl_I(1))); /** Exponential function. * * @return arbitrary precision numerical exp(x). */ -const numeric exp(const numeric & x) +const numeric exp(const numeric &x) { - return cln::exp(cln::cl_N(x)); + return cln::exp(x.to_cl_N()); } @@ -1271,56 +1204,56 @@ const numeric exp(const numeric & x) * @param z complex number * @return arbitrary precision numerical log(x). * @exception pole_error("log(): logarithmic pole",0) */ -const numeric log(const numeric & z) +const numeric log(const numeric &z) { if (z.is_zero()) throw pole_error("log(): logarithmic pole",0); - return cln::log(cln::cl_N(z)); + return cln::log(z.to_cl_N()); } /** Numeric sine (trigonometric function). * * @return arbitrary precision numerical sin(x). */ -const numeric sin(const numeric & x) +const numeric sin(const numeric &x) { - return cln::sin(cln::cl_N(x)); + return cln::sin(x.to_cl_N()); } /** Numeric cosine (trigonometric function). * * @return arbitrary precision numerical cos(x). */ -const numeric cos(const numeric & x) +const numeric cos(const numeric &x) { - return cln::cos(cln::cl_N(x)); + return cln::cos(x.to_cl_N()); } /** Numeric tangent (trigonometric function). * * @return arbitrary precision numerical tan(x). */ -const numeric tan(const numeric & x) +const numeric tan(const numeric &x) { - return cln::tan(cln::cl_N(x)); + return cln::tan(x.to_cl_N()); } /** Numeric inverse sine (trigonometric function). * * @return arbitrary precision numerical asin(x). */ -const numeric asin(const numeric & x) +const numeric asin(const numeric &x) { - return cln::asin(cln::cl_N(x)); + return cln::asin(x.to_cl_N()); } /** Numeric inverse cosine (trigonometric function). * * @return arbitrary precision numerical acos(x). */ -const numeric acos(const numeric & x) +const numeric acos(const numeric &x) { - return cln::acos(cln::cl_N(x)); + return cln::acos(x.to_cl_N()); } @@ -1329,13 +1262,13 @@ const numeric acos(const numeric & x) * @param z complex number * @return atan(z) * @exception pole_error("atan(): logarithmic pole",0) */ -const numeric atan(const numeric & x) +const numeric atan(const numeric &x) { if (!x.is_real() && x.real().is_zero() && abs(x.imag()).is_equal(_num1())) throw pole_error("atan(): logarithmic pole",0); - return cln::atan(cln::cl_N(x)); + return cln::atan(x.to_cl_N()); } @@ -1344,11 +1277,11 @@ const numeric atan(const numeric & x) * @param x real number * @param y real number * @return atan(y/x) */ -const numeric atan(const numeric & y, const numeric & x) +const numeric atan(const numeric &y, const numeric &x) { if (x.is_real() && y.is_real()) - return cln::atan(cln::the(cln::cl_N(x)), - cln::the(cln::cl_N(y))); + return cln::atan(cln::the(x.to_cl_N()), + cln::the(y.to_cl_N())); else throw std::invalid_argument("atan(): complex argument"); } @@ -1357,59 +1290,59 @@ const numeric atan(const numeric & y, const numeric & x) /** Numeric hyperbolic sine (trigonometric function). * * @return arbitrary precision numerical sinh(x). */ -const numeric sinh(const numeric & x) +const numeric sinh(const numeric &x) { - return cln::sinh(cln::cl_N(x)); + return cln::sinh(x.to_cl_N()); } /** Numeric hyperbolic cosine (trigonometric function). * * @return arbitrary precision numerical cosh(x). */ -const numeric cosh(const numeric & x) +const numeric cosh(const numeric &x) { - return cln::cosh(cln::cl_N(x)); + return cln::cosh(x.to_cl_N()); } /** Numeric hyperbolic tangent (trigonometric function). * * @return arbitrary precision numerical tanh(x). */ -const numeric tanh(const numeric & x) +const numeric tanh(const numeric &x) { - return cln::tanh(cln::cl_N(x)); + return cln::tanh(x.to_cl_N()); } /** Numeric inverse hyperbolic sine (trigonometric function). * * @return arbitrary precision numerical asinh(x). */ -const numeric asinh(const numeric & x) +const numeric asinh(const numeric &x) { - return cln::asinh(cln::cl_N(x)); + return cln::asinh(x.to_cl_N()); } /** Numeric inverse hyperbolic cosine (trigonometric function). * * @return arbitrary precision numerical acosh(x). */ -const numeric acosh(const numeric & x) +const numeric acosh(const numeric &x) { - return cln::acosh(cln::cl_N(x)); + return cln::acosh(x.to_cl_N()); } /** Numeric inverse hyperbolic tangent (trigonometric function). * * @return arbitrary precision numerical atanh(x). */ -const numeric atanh(const numeric & x) +const numeric atanh(const numeric &x) { - return cln::atanh(cln::cl_N(x)); + return cln::atanh(x.to_cl_N()); } -/*static cln::cl_N Li2_series(const ::cl_N & x, - const ::float_format_t & prec) +/*static cln::cl_N Li2_series(const ::cl_N &x, + const ::float_format_t &prec) { // Note: argument must be in the unit circle // This is very inefficient unless we have fast floating point Bernoulli @@ -1436,8 +1369,8 @@ const numeric atanh(const numeric & x) /** Numeric evaluation of Dilogarithm within circle of convergence (unit * circle) using a power series. */ -static cln::cl_N Li2_series(const cln::cl_N & x, - const cln::float_format_t & prec) +static cln::cl_N Li2_series(const cln::cl_N &x, + const cln::float_format_t &prec) { // Note: argument must be in the unit circle cln::cl_N aug, acc; @@ -1455,8 +1388,8 @@ static cln::cl_N Li2_series(const cln::cl_N & x, } /** Folds Li2's argument inside a small rectangle to enhance convergence. */ -static cln::cl_N Li2_projection(const cln::cl_N & x, - const cln::float_format_t & prec) +static cln::cl_N Li2_projection(const cln::cl_N &x, + const cln::float_format_t &prec) { const cln::cl_R re = cln::realpart(x); const cln::cl_R im = cln::imagpart(x); @@ -1481,7 +1414,7 @@ static cln::cl_N Li2_projection(const cln::cl_N & x, * continuous with quadrant IV. * * @return arbitrary precision numerical Li2(x). */ -const numeric Li2(const numeric & x) +const numeric Li2(const numeric &x) { if (x.is_zero()) return _num0(); @@ -1489,7 +1422,7 @@ const numeric Li2(const numeric & x) // what is the desired float format? // first guess: default format cln::float_format_t prec = cln::default_float_format; - const cln::cl_N value = cln::cl_N(x); + const cln::cl_N value = x.to_cl_N(); // second guess: the argument's format if (!x.real().is_rational()) prec = cln::float_format(cln::the(cln::realpart(value))); @@ -1505,13 +1438,13 @@ const numeric Li2(const numeric & x) - cln::zeta(2, prec) - Li2_projection(cln::recip(value), prec)); else - return Li2_projection(cln::cl_N(x), prec); + return Li2_projection(x.to_cl_N(), prec); } /** Numeric evaluation of Riemann's Zeta function. Currently works only for * integer arguments. */ -const numeric zeta(const numeric & x) +const 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 @@ -1519,12 +1452,12 @@ const numeric zeta(const numeric & x) // 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()) { - const int aux = (int)(cln::double_approx(cln::the(cln::cl_N(x)))); - if (cln::zerop(cln::cl_N(x)-aux)) + const int aux = (int)(cln::double_approx(cln::the(x.to_cl_N()))); + if (cln::zerop(x.to_cl_N()-aux)) return cln::zeta(aux); } std::clog << "zeta(" << x - << "): Does anybody know good way to calculate this numerically?" + << "): Does anybody know a good way to calculate this numerically?" << std::endl; return numeric(0); } @@ -1532,17 +1465,17 @@ const numeric zeta(const numeric & x) /** The Gamma function. * This is only a stub! */ -const numeric lgamma(const numeric & x) +const numeric lgamma(const numeric &x) { std::clog << "lgamma(" << x - << "): Does anybody know good way to calculate this numerically?" + << "): Does anybody know a good way to calculate this numerically?" << std::endl; return numeric(0); } -const numeric tgamma(const numeric & x) +const numeric tgamma(const numeric &x) { std::clog << "tgamma(" << x - << "): Does anybody know good way to calculate this numerically?" + << "): Does anybody know a good way to calculate this numerically?" << std::endl; return numeric(0); } @@ -1550,10 +1483,10 @@ const numeric tgamma(const numeric & x) /** The psi function (aka polygamma function). * This is only a stub! */ -const numeric psi(const numeric & x) +const numeric psi(const numeric &x) { std::clog << "psi(" << x - << "): Does anybody know good way to calculate this numerically?" + << "): Does anybody know a good way to calculate this numerically?" << std::endl; return numeric(0); } @@ -1561,10 +1494,10 @@ const numeric psi(const numeric & x) /** The psi functions (aka polygamma functions). * This is only a stub! */ -const numeric psi(const numeric & n, const numeric & x) +const numeric psi(const numeric &n, const numeric &x) { std::clog << "psi(" << n << "," << x - << "): Does anybody know good way to calculate this numerically?" + << "): Does anybody know a good way to calculate this numerically?" << std::endl; return numeric(0); } @@ -1574,7 +1507,7 @@ const numeric psi(const numeric & n, const numeric & x) * * @param n integer argument >= 0 * @exception range_error (argument must be integer >= 0) */ -const numeric factorial(const numeric & n) +const numeric factorial(const numeric &n) { if (!n.is_nonneg_integer()) throw std::range_error("numeric::factorial(): argument must be integer >= 0"); @@ -1588,9 +1521,9 @@ const numeric factorial(const numeric & n) * @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) */ -const numeric doublefactorial(const numeric & n) +const numeric doublefactorial(const numeric &n) { - if (n == numeric(-1)) + if (n.is_equal(_num_1())) return _num1(); if (!n.is_nonneg_integer()) @@ -1604,7 +1537,7 @@ const numeric doublefactorial(const numeric & n) * 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. */ -const numeric binomial(const numeric & n, const numeric & k) +const numeric binomial(const numeric &n, const numeric &k) { if (n.is_integer() && k.is_integer()) { if (n.is_nonneg_integer()) { @@ -1627,7 +1560,7 @@ const numeric binomial(const numeric & n, const numeric & k) * * @return the nth Bernoulli number (a rational number). * @exception range_error (argument must be integer >= 0) */ -const numeric bernoulli(const numeric & nn) +const numeric bernoulli(const numeric &nn) { if (!nn.is_integer() || nn.is_negative()) throw std::range_error("numeric::bernoulli(): argument must be integer >= 0"); @@ -1704,7 +1637,7 @@ const numeric bernoulli(const numeric & nn) * @param n an integer * @return the nth Fibonacci number F(n) (an integer number) * @exception range_error (argument must be an integer) */ -const numeric fibonacci(const numeric & n) +const numeric fibonacci(const numeric &n) { if (!n.is_integer()) throw std::range_error("numeric::fibonacci(): argument must be integer"); @@ -1734,7 +1667,7 @@ const numeric fibonacci(const numeric & n) cln::cl_I u(0); cln::cl_I v(1); - cln::cl_I m = cln::the(cln::cl_N(n)) >> 1L; // floor(n/2); + cln::cl_I m = cln::the(n.to_cl_N()) >> 1L; // floor(n/2); for (uintL bit=cln::integer_length(m); bit>0; --bit) { // Since a squaring is cheaper than a multiplication, better use // three squarings instead of one multiplication and two squarings. @@ -1760,7 +1693,7 @@ const numeric fibonacci(const numeric & n) /** Absolute value. */ const numeric abs(const numeric& x) { - return cln::abs(cln::cl_N(x)); + return cln::abs(x.to_cl_N()); } @@ -1771,11 +1704,11 @@ const numeric abs(const numeric& x) * * @return a mod b in the range [0,abs(b)-1] with sign of b if both are * integer, 0 otherwise. */ -const numeric mod(const numeric & a, const numeric & b) +const numeric mod(const numeric &a, const numeric &b) { if (a.is_integer() && b.is_integer()) - return cln::mod(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + return cln::mod(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); else return _num0(); } @@ -1785,12 +1718,12 @@ const numeric mod(const numeric & a, const numeric & b) * Equivalent to Maple's mods. * * @return a mod b in the range [-iquo(abs(m)-1,2), iquo(abs(m),2)]. */ -const numeric smod(const numeric & a, const numeric & b) +const numeric smod(const numeric &a, const numeric &b) { if (a.is_integer() && b.is_integer()) { - const cln::cl_I b2 = cln::ceiling1(cln::the(cln::cl_N(b)) >> 1) - 1; - return cln::mod(cln::the(cln::cl_N(a)) + b2, - cln::the(cln::cl_N(b))) - b2; + const cln::cl_I b2 = cln::ceiling1(cln::the(b.to_cl_N()) >> 1) - 1; + return cln::mod(cln::the(a.to_cl_N()) + b2, + cln::the(b.to_cl_N())) - b2; } else return _num0(); } @@ -1802,11 +1735,11 @@ const numeric smod(const numeric & a, const numeric & b) * sign of a or is zero. * * @return remainder of a/b if both are integer, 0 otherwise. */ -const numeric irem(const numeric & a, const numeric & b) +const numeric irem(const numeric &a, const numeric &b) { if (a.is_integer() && b.is_integer()) - return cln::rem(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + return cln::rem(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); else return _num0(); } @@ -1819,11 +1752,11 @@ const numeric irem(const numeric & a, const numeric & b) * * @return remainder of a/b and quotient stored in q if both are integer, * 0 otherwise. */ -const numeric irem(const numeric & a, const numeric & b, numeric & q) +const numeric irem(const numeric &a, const numeric &b, numeric &q) { if (a.is_integer() && b.is_integer()) { - const cln::cl_I_div_t rem_quo = cln::truncate2(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + const cln::cl_I_div_t rem_quo = cln::truncate2(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); q = rem_quo.quotient; return rem_quo.remainder; } else { @@ -1837,11 +1770,11 @@ const numeric irem(const numeric & a, const numeric & 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. */ -const numeric iquo(const numeric & a, const numeric & b) +const numeric iquo(const numeric &a, const numeric &b) { if (a.is_integer() && b.is_integer()) - return truncate1(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + return truncate1(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); else return _num0(); } @@ -1853,11 +1786,11 @@ const numeric iquo(const numeric & a, const numeric & b) * * @return truncated quotient of a/b and remainder stored in r if both are * integer, 0 otherwise. */ -const numeric iquo(const numeric & a, const numeric & b, numeric & r) +const numeric iquo(const numeric &a, const numeric &b, numeric &r) { if (a.is_integer() && b.is_integer()) { - const cln::cl_I_div_t rem_quo = cln::truncate2(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + const cln::cl_I_div_t rem_quo = cln::truncate2(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); r = rem_quo.remainder; return rem_quo.quotient; } else { @@ -1871,11 +1804,11 @@ const numeric iquo(const numeric & a, const numeric & b, numeric & r) * * @return The GCD of two numbers if both are integer, a numerical 1 * if they are not. */ -const numeric gcd(const numeric & a, const numeric & b) +const numeric gcd(const numeric &a, const numeric &b) { if (a.is_integer() && b.is_integer()) - return cln::gcd(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + return cln::gcd(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); else return _num1(); } @@ -1885,11 +1818,11 @@ const numeric gcd(const numeric & a, const numeric & b) * * @return The LCM of two numbers if both are integer, the product of those * two numbers if they are not. */ -const numeric lcm(const numeric & a, const numeric & b) +const numeric lcm(const numeric &a, const numeric &b) { if (a.is_integer() && b.is_integer()) - return cln::lcm(cln::the(cln::cl_N(a)), - cln::the(cln::cl_N(b))); + return cln::lcm(cln::the(a.to_cl_N()), + cln::the(b.to_cl_N())); else return a.mul(b); } @@ -1903,18 +1836,18 @@ const numeric lcm(const numeric & a, const numeric & b) * @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. */ -const numeric sqrt(const numeric & z) +const numeric sqrt(const numeric &z) { - return cln::sqrt(cln::cl_N(z)); + return cln::sqrt(z.to_cl_N()); } /** Integer numeric square root. */ -const numeric isqrt(const numeric & x) +const numeric isqrt(const numeric &x) { if (x.is_integer()) { cln::cl_I root; - cln::isqrt(cln::the(cln::cl_N(x)), &root); + cln::isqrt(cln::the(x.to_cl_N()), &root); return root; } else return _num0(); @@ -1942,13 +1875,15 @@ ex CatalanEvalf(void) } +/** _numeric_digits default ctor, checking for singleton invariance. */ _numeric_digits::_numeric_digits() : digits(17) { // It initializes to 17 digits, because in CLN float_format(17) turns out // to be 61 (<64) while float_format(18)=65. The reason is we want to // have a cl_LF instead of cl_SF, cl_FF or cl_DF. - assert(!too_late); + if (too_late) + throw(std::runtime_error("I told you not to do instantiate me!")); too_late = true; cln::default_float_format = cln::float_format(17); } @@ -1972,14 +1907,14 @@ _numeric_digits::operator long() /** Append global Digits object to ostream. */ -void _numeric_digits::print(std::ostream & os) const +void _numeric_digits::print(std::ostream &os) const { debugmsg("_numeric_digits print", LOGLEVEL_PRINT); os << digits; } -std::ostream& operator<<(std::ostream& os, const _numeric_digits & e) +std::ostream& operator<<(std::ostream &os, const _numeric_digits &e) { e.print(os); return os; @@ -1998,6 +1933,4 @@ bool _numeric_digits::too_late = false; * assignment from C++ unsigned ints and evaluated like any built-in type. */ _numeric_digits Digits; -#ifndef NO_NAMESPACE_GINAC } // namespace GiNaC -#endif // ndef NO_NAMESPACE_GINAC