From 0962930400269497db414baba12740496157a63b Mon Sep 17 00:00:00 2001 From: Jens Vollinga Date: Mon, 17 Nov 2003 22:18:02 +0000 Subject: [PATCH] * Removed mZeta * Relaxed convergence check for multiple polylog * Implemented alternating Euler sums as zeta(m,s) --- ginac/inifcns_nstdsums.cpp | 243 ++++++++++++++++++++++++++++--------- ginsh/ginsh_parser.yy | 3 +- 2 files changed, 184 insertions(+), 62 deletions(-) diff --git a/ginac/inifcns_nstdsums.cpp b/ginac/inifcns_nstdsums.cpp index 22deea9e..fbc4a3fa 100644 --- a/ginac/inifcns_nstdsums.cpp +++ b/ginac/inifcns_nstdsums.cpp @@ -4,10 +4,11 @@ * * The functions are: * classical polylogarithm Li(n,x) - * multiple polylogarithm Li(lst(n_1,...,n_k),lst(x_1,...,x_k)) + * multiple polylogarithm Li(lst(m_1,...,m_k),lst(x_1,...,x_k)) * nielsen's generalized polylogarithm S(n,p,x) - * harmonic polylogarithm H(n,x) or H(lst(n_1,...,n_k),x) - * multiple zeta value zeta(n) or zeta(lst(n_1,...,n_k)) + * harmonic polylogarithm H(m,x) or H(lst(m_1,...,m_k),x) + * multiple zeta value zeta(m) or zeta(lst(m_1,...,m_k)) + * alternating Euler sum zeta(m,s) or zeta(lst(m_1,...,m_k),lst(s_1,...,s_k)) * * Some remarks: * @@ -15,6 +16,7 @@ * [Kol] Nielsen's Generalized Polylogarithms, K.S.Kolbig, SIAM J.Math.Anal. 17 (1986), pp. 1232-1258. * [Cra] Fast Evaluation of Multiple Zeta Sums, R.E.Crandall, Math.Comp. 67 (1998), pp. 1163-1172. * [ReV] Harmonic Polylogarithms, E.Remiddi, J.A.M.Vermaseren, Int.J.Mod.Phys. A15 (2000), pp. 725-754 + * [BBB] Special Values of Multiple Polylogarithms, J.Borwein, D.Bradley, D.Broadhurst, P.Lisonek, Trans.Amer.Math.Soc. 353/3 (2001), pp. 907-941 * * - The order of parameters and arguments of H, Li and zeta is defined according to their order in the * nested sums representation. @@ -24,8 +26,8 @@ * one. The parameters for every function (n, p or n_i) must be positive integers. * * - The calculation of classical polylogarithms is speed up by using Bernoulli numbers and - * look-up tables. S uses look-up tables as well. The zeta function applies the algorithm in - * [Cra] for speed up. + * look-up tables. S uses look-up tables as well. The zeta function applies the algorithms in + * [Cra] and [BBB] for speed up. * * - The functions have no series expansion as nested sums. To do it, you have to convert these functions * into the appropriate objects from the nestedsums library, do the expansion and convert the @@ -82,7 +84,7 @@ namespace GiNaC { ////////////////////////////////////////////////////////////////////// // -// Classical polylogarithm Li +// Classical polylogarithm Li(n,x) // // helper functions // @@ -383,7 +385,7 @@ numeric Li_num(int n, const numeric& x) ////////////////////////////////////////////////////////////////////// // -// Multiple polylogarithm Li +// Multiple polylogarithm Li(n,x) // // helper function // @@ -410,7 +412,7 @@ cln::cl_N multipleLi_do_sum(const std::vector& s, const std::vector=0; k--) { t[k] = t[k] + t[k+1] * cln::expt(x[k], q+j-1-k) / cln::expt(cln::cl_I(q+j-1-k), s[k]); } - } while (t[0] != t0buf); + } while ((t[0] != t0buf) || (q<10)); return t[0]; } @@ -421,7 +423,7 @@ cln::cl_N multipleLi_do_sum(const std::vector& s, const std::vector(x1) && is_a(x2)) { + ex conv = 1; for (int i=0; i(x2.op(i))) { return Li(x1,x2).hold(); } - if (x2.op(i) >= 1) { + conv *= x2.op(i); + if ((conv > 1) || ((conv == 1) && (x1.op(0) == 1))) { return Li(x1,x2).hold(); } } @@ -515,7 +519,7 @@ REGISTER_FUNCTION(Li, ////////////////////////////////////////////////////////////////////// // -// Nielsen's generalized polylogarithm S +// Nielsen's generalized polylogarithm S(n,p,x) // // helper functions // @@ -865,7 +869,7 @@ numeric S_num(int n, int p, const numeric& x) ////////////////////////////////////////////////////////////////////// // -// Nielsen's generalized polylogarithm S +// Nielsen's generalized polylogarithm S(n,p,x) // // GiNaC function // @@ -926,7 +930,7 @@ REGISTER_FUNCTION(S, ////////////////////////////////////////////////////////////////////// // -// Harmonic polylogarithm H +// Harmonic polylogarithm H(m,x) // // helper function // @@ -1418,7 +1422,7 @@ cln::cl_N H_do_sum(const std::vector& s, const cln::cl_N& x) ////////////////////////////////////////////////////////////////////// // -// Harmonic polylogarithm H +// Harmonic polylogarithm H(m,x) // // GiNaC function // @@ -1540,7 +1544,7 @@ REGISTER_FUNCTION(H, ////////////////////////////////////////////////////////////////////// // -// Multiple zeta values zeta +// Multiple zeta values zeta(x) and zeta(x,s) // // helper functions // @@ -1797,12 +1801,102 @@ cln::cl_N zeta_do_sum_simple(const std::vector& r) } +// does Hoelder convolution. see [BBB] (7.0) +cln::cl_N zeta_do_Hoelder_convolution(const std::vector& m_, const std::vector& s_) +{ + // prepare parameters + // holds Li arguments in [BBB] notation + std::vector s = s_; + std::vector m_p = m_; + std::vector m_q; + // holds Li arguments in nested sums notation + std::vector s_p(s.size(), cln::cl_N(1)); + s_p[0] = s_p[0] * cln::cl_N("1/2"); + // convert notations + int sig = 1; + for (int i=0; i s_q; + cln::cl_N signum = 1; + + // first term + cln::cl_N res = multipleLi_do_sum(m_p, s_p); + + // middle terms + do { + + // change parameters + if (s.front() > 0) { + if (m_p.front() == 1) { + m_p.erase(m_p.begin()); + s_p.erase(s_p.begin()); + if (s_p.size() > 0) { + s_p.front() = s_p.front() * cln::cl_N("1/2"); + } + s.erase(s.begin()); + m_q.front()++; + } else { + m_p.front()--; + m_q.insert(m_q.begin(), 1); + if (s_q.size() > 0) { + s_q.front() = s_q.front() * 2; + } + s_q.insert(s_q.begin(), cln::cl_N("1/2")); + } + } else { + if (m_p.front() == 1) { + m_p.erase(m_p.begin()); + s_p.erase(s_p.begin()); + if (s_p.size() > 0) { + s_p.front() = s_p.front() * cln::cl_N("1/2"); + } + s.erase(s.begin()); + for (int i=0; i 0) break; + + } + m_q.insert(m_q.begin(), 1); + if (s_q.size() > 0) { + s_q.front() = s_q.front() * 4; + } + s_q.insert(s_q.begin(), cln::cl_N("1/4")); + signum = -signum; + } else { + m_p.front()--; + m_q.insert(m_q.begin(), 1); + if (s_q.size() > 0) { + s_q.front() = s_q.front() * 2; + } + s_q.insert(s_q.begin(), cln::cl_N("1/2")); + } + } + + // exiting the loop + if (m_p.size() == 0) break; + + res = res + signum * multipleLi_do_sum(m_p, s_p) * multipleLi_do_sum(m_q, s_q); + + } while (true); + + // last term + res = res + signum * multipleLi_do_sum(m_q, s_q); + + return res; +} + + } // end of anonymous namespace ////////////////////////////////////////////////////////////////////// // -// Multiple zeta values zeta +// Multiple zeta values zeta(x) // // GiNaC function // @@ -1921,74 +2015,103 @@ unsigned zeta1_SERIAL::serial = ////////////////////////////////////////////////////////////////////// // -// Multiple zeta values mZeta -// -// The use of mZeta is deprecated! This function will be removed -// from GiNaC source soon. Use zeta instead!! +// Alternating Euler sum zeta(x,s) // // GiNaC function // ////////////////////////////////////////////////////////////////////// -static ex mZeta_eval(const ex& x1) -{ - return mZeta(x1).hold(); -} - - -static ex mZeta_evalf(const ex& x1) +static ex zeta2_evalf(const ex& x, const ex& s) { - if (is_a(x1)) { - for (int i=0; i(x)) { - const int j = x1.nops(); + // alternating Euler sum + const int count = x.nops(); + const lst& xlst = ex_to(x); + const lst& slst = ex_to(s); + std::vector xi(count); + std::vector si(count); - std::vector r(j); - for (int i=0; i(x1.op(i)).to_int(); - } + // check parameters and convert them + lst::const_iterator it_xread = xlst.begin(); + lst::const_iterator it_sread = slst.begin(); + std::vector::iterator it_xwrite = xi.begin(); + std::vector::iterator it_swrite = si.begin(); + do { + if (!(*it_xread).info(info_flags::posint)) { + return zeta(x, s).hold(); + } + *it_xwrite = ex_to(*it_xread).to_int(); + if (*it_sread > 0) { + *it_swrite = 1; + } else { + *it_swrite = -1; + } + it_xread++; + it_sread++; + it_xwrite++; + it_swrite++; + } while (it_xwrite != xi.end()); // check for divergence - if (r[0] == 1) { - return mZeta(x1).hold(); + if ((xi[0] == 1) && (si[0] == 1)) { + return zeta(x, s).hold(); } - // if only one argument, use cln::zeta - if (j == 1) { - return numeric(cln::zeta(r[0])); - } + // use Hoelder convolution + return numeric(zeta_do_Hoelder_convolution(xi, si)); + } - // decide on summation algorithm - // this is still a bit clumsy - int limit = (Digits>17) ? 10 : 6; - if (r[0]3 && r[1](s)) { + const lst& l = ex_to(s); + lst::const_iterator it = l.begin(); + while (it != l.end()) { + if ((*it).info(info_flags::negative)) { + return zeta(x, s).hold(); + } + it++; + } + return zeta(x); + } else { + if (s.info(info_flags::positive)) { + return zeta(x); } - } else if (x1.info(info_flags::posint) && (x1 != 1)) { - return numeric(cln::zeta(ex_to(x1).to_int())); } - return mZeta(x1).hold(); + return zeta(x, s).hold(); } -static ex mZeta_deriv(const ex& x, unsigned deriv_param) +static ex zeta2_deriv(const ex& x, const ex& s, unsigned deriv_param) { - return 0; + GINAC_ASSERT(deriv_param==0); + + if (is_exactly_a(x)) { + return _ex0; + } else { + if ((is_exactly_a(s) && (s.op(0) > 0)) || (s > 0)) { + return zeta(_ex1, x); + } + return _ex0; + } } -REGISTER_FUNCTION(mZeta, - eval_func(mZeta_eval). - evalf_func(mZeta_evalf). - do_not_evalf_params(). - derivative_func(mZeta_deriv)); +unsigned zeta2_SERIAL::serial = + function::register_new(function_options("zeta"). + eval_func(zeta2_eval). + evalf_func(zeta2_evalf). + do_not_evalf_params(). + derivative_func(zeta2_deriv). + latex_name("\\zeta"). + overloaded(2)); } // namespace GiNaC diff --git a/ginsh/ginsh_parser.yy b/ginsh/ginsh_parser.yy index 5dc2228b..e1420877 100644 --- a/ginsh/ginsh_parser.yy +++ b/ginsh/ginsh_parser.yy @@ -621,13 +621,12 @@ static const fcn_help_init builtin_help[] = { {"sinh", "hyperbolic sine function"}, {"tan", "tangent function"}, {"tanh", "hyperbolic tangent function"}, - {"zeta", "zeta function\nzeta(x) is Riemann's zeta function, zeta(n,x) its nth derivative"}, + {"zeta", "zeta function\nzeta(x) is Riemann's zeta function, zetaderiv(n,x) its nth derivative.\nIf x is a GiNaC::lst, it is a multiple zeta value\nzeta(x,s) is an alternating Euler sum"}, {"Li2", "dilogarithm"}, {"Li3", "trilogarithm"}, {"Li", "(multiple) polylogarithm"}, {"S", "Nielsen's generalized polylogarithm"}, {"H", "harmonic polylogarithm"}, - {"mZeta", "multiple zeta value"}, {"Order", "order term function (for truncated power series)"}, {"Derivative", "inert differential operator"}, {NULL, NULL} // End marker -- 2.44.0