From cda7fd319d75f95fbdac21213a10a04de724d228 Mon Sep 17 00:00:00 2001 From: Stefan Weinzierl Date: Sun, 12 Oct 2014 18:30:41 +0000 Subject: [PATCH] Fix log() error in multiple polylog G. Tobias Huber reported the following bug: > G({-0.18784426860496276*I,-0.1579006477353911*I,+0.18784426860496276*I,-0.18784426860496276*I},1); log_eval(): log(0) The problem is the sorting operation in G_do_trafo. This operation should 1) put the elements in increasing order of |x| 2) put equal elements next to each other (which then will avoid the log(0) problem). The current code fails for goal 2) if we have complex numbers of equal absolute value, one number occuring more than once, in an initial order like in the example above: Entries 1,3,4 of the list all have the same absolute value, entries 1 and 4 are equal. In the sorting operation 1 and 4 should be put next to each other. Previously, the sorting operation would give 2,1,3,4. What we would like to have is either 2,1,4,3 or 2,3,1,4. It is clear that the key for the sorting operation cannot be |x| alone. If |x1|=|x2| we have to use the phase as well. This patch takes the phase into account in the sorting. --- ginac/inifcns_nstdsums.cpp | 25 ++++++++++++++++++++++--- 1 file changed, 22 insertions(+), 3 deletions(-) diff --git a/ginac/inifcns_nstdsums.cpp b/ginac/inifcns_nstdsums.cpp index 866472f4..a42c1bd3 100644 --- a/ginac/inifcns_nstdsums.cpp +++ b/ginac/inifcns_nstdsums.cpp @@ -1043,6 +1043,25 @@ G_do_hoelder(std::vector x, /* yes, it's passed by value */ return result; } +class less_object_for_cl_N +{ +public: + bool operator() (const cln::cl_N & a, const cln::cl_N & b) const + { + // absolute value? + if (abs(a) != abs(b)) + return (abs(a) < abs(b)) ? true : false; + + // complex phase? + if (phase(a) != phase(b)) + return (phase(a) < phase(b)) ? true : false; + + // equal, therefore "less" is not true + return false; + } +}; + + // convergence transformation, used for numerical evaluation of G function. // the parameter x, s and y must only contain numerics static cln::cl_N @@ -1050,17 +1069,17 @@ G_do_trafo(const std::vector& x, const std::vector& s, const cln::cl_N& y, bool flag_trailing_zeros_only) { // sort (|x|<->position) to determine indices - typedef std::multimap sortmap_t; + typedef std::multimap sortmap_t; sortmap_t sortmap; std::size_t size = 0; for (std::size_t i = 0; i < x.size(); ++i) { if (!zerop(x[i])) { - sortmap.insert(std::make_pair(abs(x[i]), i)); + sortmap.insert(std::make_pair(x[i], i)); ++size; } } // include upper limit (scale) - sortmap.insert(std::make_pair(abs(y), x.size())); + sortmap.insert(std::make_pair(y, x.size())); // generate missing dummy-symbols int i = 1; -- 2.37.1