X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=blobdiff_plain;f=ginac%2Fnumeric.cpp;h=ff5d510ca39023a8aa8262ee1cc021bf6bd916ec;hp=5649353f308f68a0c26e715034b84977c5702979;hb=4f0b17af13eb5f7f34fdab171c6d630a77badb3d;hpb=def26469ff96228c66e877bb5594e7d9a24b638f diff --git a/ginac/numeric.cpp b/ginac/numeric.cpp index 5649353f..ff5d510c 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-2008 Johannes Gutenberg University Mainz, Germany + * GiNaC Copyright (C) 1999-2020 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 @@ -24,21 +24,22 @@ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA */ +#ifdef HAVE_CONFIG_H #include "config.h" - -#include -#include -#include -#include -#include +#endif #include "numeric.h" #include "ex.h" #include "operators.h" #include "archive.h" -#include "tostring.h" #include "utils.h" +#include +#include +#include +#include +#include + // CLN should pollute the global namespace as little as possible. Hence, we // include most of it here and include only the part needed for properly // declaring cln::cl_number in numeric.h. This can only be safely done in @@ -88,7 +89,7 @@ numeric::numeric() numeric::numeric(int i) { // 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 + // first. This is due to the behavior of the cl_I-ctor, which // emphasizes efficiency. However, if the integer is small enough // we save space and dereferences by using an immediate type. // (C.f. ) @@ -109,7 +110,7 @@ numeric::numeric(int i) numeric::numeric(unsigned int i) { // 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 + // first. This is due to the behavior of the cl_I-ctor, which // emphasizes efficiency. However, if the integer is small enough // we save space and dereferences by using an immediate type. // (C.f. ) @@ -140,6 +141,17 @@ numeric::numeric(unsigned long i) setflag(status_flags::evaluated | status_flags::expanded); } +numeric::numeric(long long i) +{ + value = cln::cl_I(i); + setflag(status_flags::evaluated | status_flags::expanded); +} + +numeric::numeric(unsigned long long i) +{ + value = cln::cl_I(i); + setflag(status_flags::evaluated | status_flags::expanded); +} /** Constructor for rational numerics a/b. * @@ -223,7 +235,7 @@ numeric::numeric(const char *s) // E to lower case term = term.replace(term.find("E"),1,"e"); // append _ to term - term += "_" + ToString((unsigned)Digits); + term += "_" + std::to_string((unsigned)Digits); // construct float using cln::cl_F(const char *) ctor. if (imaginary) ctorval = ctorval + cln::complex(cln::cl_I(0),cln::cl_F(term.c_str())); @@ -278,8 +290,9 @@ static const cln::cl_F read_real_float(std::istream& s) return x; } -numeric::numeric(const archive_node &n, lst &sym_lst) : inherited(n, sym_lst) +void numeric::read_archive(const archive_node &n, lst &sym_lst) { + inherited::read_archive(n, sym_lst); value = 0; // Read number as string @@ -324,6 +337,7 @@ numeric::numeric(const archive_node &n, lst &sym_lst) : inherited(n, sym_lst) } setflag(status_flags::evaluated | status_flags::expanded); } +GINAC_BIND_UNARCHIVER(numeric); static void write_real_float(std::ostream& s, const cln::cl_R& n) { @@ -373,8 +387,6 @@ void numeric::archive(archive_node &n) const n.add_string("number", s.str()); } -DEFAULT_UNARCHIVE(numeric) - ////////// // functions overriding virtual functions from base classes ////////// @@ -441,7 +453,7 @@ static void print_real_csrc(const print_context & c, const cln::cl_R & x) // Rational number const cln::cl_I numer = cln::numerator(cln::the(x)); const cln::cl_I denom = cln::denominator(cln::the(x)); - if (cln::plusp(x) > 0) { + if (cln::plusp(x)) { c.s << "("; print_integer_csrc(c, numer); } else { @@ -465,7 +477,7 @@ static inline bool coerce(T1& dst, const T2& arg); /** * @brief Check if CLN integer can be converted into int * - * @sa http://www.ginac.de/pipermail/cln-list/2006-October/000248.html + * @sa https://www.ginac.de/pipermail/cln-list/2006-October/000248.html */ template<> inline bool coerce(int& dst, const cln::cl_I& arg) @@ -500,20 +512,20 @@ static void print_real_cl_N(const print_context & c, const cln::cl_R & x) { if (cln::instanceof(x, cln::cl_I_ring)) { - int dst; - // fixnum - if (coerce(dst, cln::the(x))) { - // can be converted to native int - if (dst < 0) - c.s << "(-" << dst << ")"; - else - c.s << dst; - } else { - // bignum - c.s << "cln::cl_I(\""; - print_real_number(c, x); - c.s << "\")"; - } + int dst; + // fixnum + if (coerce(dst, cln::the(x))) { + // can be converted to native int + if (dst < 0) + c.s << '(' << dst << ')'; + else + c.s << dst; + } else { + // bignum + c.s << "cln::cl_I(\""; + print_real_number(c, x); + c.s << "\")"; + } } else if (cln::instanceof(x, cln::cl_RA_ring)) { // Rational number @@ -700,7 +712,7 @@ bool numeric::info(unsigned inf) const case info_flags::negative: return is_negative(); case info_flags::nonnegative: - return !is_negative(); + return is_zero() || is_positive(); case info_flags::posint: return is_pos_integer(); case info_flags::negint: @@ -713,8 +725,6 @@ bool numeric::info(unsigned inf) const return is_odd(); case info_flags::prime: return is_prime(); - case info_flags::algebraic: - return !is_real(); } return false; } @@ -774,10 +784,8 @@ bool numeric::has(const ex &other, unsigned options) const /** Evaluation of numbers doesn't do anything at all. */ -ex numeric::eval(int level) const +ex numeric::eval() const { - // Warning: if this is ever gonna do something, the ex ctors from all kinds - // of numbers should be checking for status_flags::evaluated. return this->hold(); } @@ -787,11 +795,9 @@ ex numeric::eval(int level) const * currently set. In case the object already was a floating point number the * precision is trimmed to match the currently set default. * - * @param level ignored, only needed for overriding basic::evalf. * @return an ex-handle to a numeric. */ -ex numeric::evalf(int level) const +ex numeric::evalf() const { - // level can safely be discarded for numeric objects. return numeric(cln::cl_float(1.0, cln::default_float_format) * value); } @@ -928,9 +934,8 @@ const numeric &numeric::add_dyn(const numeric &other) const return other; else if (&other==_num0_p) return *this; - - return static_cast((new numeric(value + other.value))-> - setflag(status_flags::dynallocated)); + + return dynallocate(value + other.value); } @@ -944,9 +949,8 @@ const numeric &numeric::sub_dyn(const numeric &other) const // hack is supposed to keep the number of distinct numeric objects low. if (&other==_num0_p || cln::zerop(other.value)) return *this; - - return static_cast((new numeric(value - other.value))-> - setflag(status_flags::dynallocated)); + + return dynallocate(value - other.value); } @@ -963,8 +967,7 @@ const numeric &numeric::mul_dyn(const numeric &other) const else if (&other==_num1_p) return *this; - return static_cast((new numeric(value * other.value))-> - setflag(status_flags::dynallocated)); + return dynallocate(value * other.value); } @@ -982,8 +985,8 @@ const numeric &numeric::div_dyn(const numeric &other) const return *this; if (cln::zerop(cln::the(other.value))) throw std::overflow_error("division by zero"); - return static_cast((new numeric(value / other.value))-> - setflag(status_flags::dynallocated)); + + return dynallocate(value / other.value); } @@ -1009,8 +1012,8 @@ const numeric &numeric::power_dyn(const numeric &other) const else return *_num0_p; } - return static_cast((new numeric(cln::expt(value, other.value)))-> - setflag(status_flags::dynallocated)); + + return dynallocate(cln::expt(value, other.value)); } @@ -1374,7 +1377,7 @@ const numeric numeric::numer() const if (cln::instanceof(r, cln::cl_RA_ring) && cln::instanceof(i, cln::cl_RA_ring)) { const cln::cl_I s = cln::lcm(cln::denominator(r), cln::denominator(i)); return numeric(cln::complex(cln::numerator(r)*(cln::exquo(s,cln::denominator(r))), - cln::numerator(i)*(cln::exquo(s,cln::denominator(i))))); + cln::numerator(i)*(cln::exquo(s,cln::denominator(i))))); } } // at least one float encountered @@ -1529,7 +1532,7 @@ const numeric atan(const numeric &y, const numeric &x) return *_num0_p; if (x.is_real() && y.is_real()) return numeric(cln::atan(cln::the(x.to_cl_N()), - cln::the(y.to_cl_N()))); + cln::the(y.to_cl_N()))); // Compute -I*log((x+I*y)/sqrt(x^2+y^2)) // == -I*log((x+I*y)/sqrt((x+I*y)*(x-I*y))) @@ -1747,7 +1750,7 @@ class lanczos_coeffs std::vector *current_vector; }; -std::vector* lanczos_coeffs::coeffs = 0; +std::vector* lanczos_coeffs::coeffs = nullptr; bool lanczos_coeffs::sufficiently_accurate(int digits) { if (digits<=20) { @@ -1773,8 +1776,8 @@ cln::cl_N lanczos_coeffs::calc_lanczos_A(const cln::cl_N &x) const { cln::cl_N A = (*current_vector)[0]; int size = current_vector->size(); - for (int i=1; i(b.to_cl_N()) >> 1) - 1; - return numeric(cln::mod(cln::the(a.to_cl_N()) + b2, - cln::the(b.to_cl_N())) - b2); + * @return a mod b in the range [-iquo(abs(b),2), iquo(abs(b),2)]. */ +const numeric smod(const numeric &a_, const numeric &b_) +{ + if (a_.is_integer() && b_.is_integer()) { + const cln::cl_I a = cln::the(a_.to_cl_N()); + const cln::cl_I b = cln::the(b_.to_cl_N()); + const cln::cl_I b2 = b >> 1; + const cln::cl_I m = cln::mod(a, b); + const cln::cl_I m_b = m - b; + const cln::cl_I ret = m > b2 ? m_b : m; + return numeric(ret); } else return *_num0_p; } @@ -2531,9 +2540,8 @@ _numeric_digits& _numeric_digits::operator=(long prec) cln::default_float_format = cln::float_format(prec); // call registered callbacks - std::vector::const_iterator it = callbacklist.begin(), end = callbacklist.end(); - for (; it != end; ++it) { - (*it)(digitsdiff); + for (auto it : callbacklist) { + (it)(digitsdiff); } return *this;