]> www.ginac.de Git - ginac.git/blobdiff - ginac/inifcns_nstdsums.cpp
Removed exits from numeric.cpp and use array indices instead of push_back.
[ginac.git] / ginac / inifcns_nstdsums.cpp
index 1fd80aefe2e453b269fcf34e403d9149669e2014..e449570d8359c211fbc1a1bdf0236355e3e220f4 100644 (file)
@@ -47,7 +47,7 @@
  */
 
 /*
- *  GiNaC Copyright (C) 1999-2005 Johannes Gutenberg University Mainz, Germany
+ *  GiNaC Copyright (C) 1999-2007 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
@@ -471,7 +471,13 @@ namespace {
 // performs the actual series summation for multiple polylogarithms
 cln::cl_N multipleLi_do_sum(const std::vector<int>& s, const std::vector<cln::cl_N>& x)
 {
+       // ensure all x <> 0.
+       for (std::vector<cln::cl_N>::const_iterator it = x.begin(); it != x.end(); ++it) {
+               if ( *it == 0 ) return cln::cl_float(0, cln::float_format(Digits));
+       }
+
        const int j = s.size();
+       bool flag_accidental_zero = false;
 
        std::vector<cln::cl_N> t(j);
        cln::cl_F one = cln::cl_float(1, cln::float_format(Digits));
@@ -480,19 +486,13 @@ cln::cl_N multipleLi_do_sum(const std::vector<int>& s, const std::vector<cln::cl
        int q = 0;
        do {
                t0buf = t[0];
-               // do it once ...
                q++;
                t[j-1] = t[j-1] + cln::expt(x[j-1], q) / cln::expt(cln::cl_I(q),s[j-1]) * one;
                for (int k=j-2; k>=0; k--) {
+                       flag_accidental_zero = cln::zerop(t[k+1]);
                        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]);
                }
-               // ... and do it again (to avoid premature drop out due to special arguments)
-               q++;
-               t[j-1] = t[j-1] + cln::expt(x[j-1], q) / cln::expt(cln::cl_I(q),s[j-1]) * one;
-               for (int k=j-2; k>=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) || flag_accidental_zero );
 
        return t[0];
 }
@@ -1911,7 +1911,7 @@ numeric S_num(int n, int p, const numeric& x)
                prec = cln::float_format(cln::the<cln::cl_F>(cln::imagpart(value)));
 
        // [Kol] (5.3)
-       if ((cln::realpart(value) < -0.5) || (n == 0)) {
+       if ((cln::realpart(value) < -0.5) || (n == 0) || ((cln::abs(value) <= 1) && (cln::abs(value) > 0.95))) {
 
                cln::cl_N result = cln::expt(cln::cl_I(-1),p) * cln::expt(cln::log(value),n)
                                   * cln::expt(cln::log(1-value),p) / cln::factorial(n) / cln::factorial(p);
@@ -2635,7 +2635,7 @@ struct map_trafo_H_1mx : public map_function
                                        if (allthesame) {
                                                lst newparameter;
                                                for (int i=parameter.nops(); i>0; i--) {
-                                                       newparameter.append(0);
+                                                       newparameter.append(1);
                                                }
                                                return pow(-1, parameter.nops()) * H(newparameter, 1-arg).hold();
                                        }
@@ -2651,7 +2651,7 @@ struct map_trafo_H_1mx : public map_function
                                        if (allthesame) {
                                                lst newparameter;
                                                for (int i=parameter.nops(); i>0; i--) {
-                                                       newparameter.append(1);
+                                                       newparameter.append(0);
                                                }
                                                return pow(-1, parameter.nops()) * H(newparameter, 1-arg).hold();
                                        }
@@ -2681,7 +2681,7 @@ struct map_trafo_H_1mx : public map_function
                                        // leading one
                                        map_trafo_H_1mx recursion;
                                        map_trafo_H_mult unify;
-                                       ex res;
+                                       ex res = H(lst(1), arg).hold() * H(newparameter, arg).hold();
                                        int firstzero = 0;
                                        while (parameter.op(firstzero) == 1) {
                                                firstzero++;
@@ -2698,11 +2698,9 @@ struct map_trafo_H_1mx : public map_function
                                                }
                                                res -= H(newparameter, arg).hold();
                                        }
-                                       return (unify((-H(lst(0), 1-arg).hold() * recursion(H(newparameter, arg).hold())).expand()) +
-                                                       recursion(res)) / firstzero;
-
+                                       res = recursion(res).expand() / firstzero;
+                                       return unify(res);
                                }
-
                        }
                }
                return e;