Avoid calling log(1-x) for x=1 in Li_projection.

[ginac.git] / ginac / inifcns_nstdsums.cpp

diff --git a/ginac/inifcns_nstdsums.cpp b/ginac/inifcns_nstdsums.cpp
index 3284ff698d267757e47e6b23b1f79a563e121283..7487b31430de89c01871d8b8911192516e5ecf21 100644 (file)
--- a/ginac/inifcns_nstdsums.cpp
+++ b/ginac/inifcns_nstdsums.cpp
@@ -47,7 +47,7 @@
  */
 
 /*
- *  GiNaC Copyright (C) 1999-2009 Johannes Gutenberg University Mainz, Germany
+ *  GiNaC Copyright (C) 1999-2011 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
@@ -346,7 +346,11 @@ cln::cl_N Li_projection(int n, const cln::cl_N& x, const cln::float_format_t& pr
                } else {
                        // choose the faster algorithm
                        if (cln::abs(cln::realpart(x)) > 0.75) {
-                               return -Li2_do_sum(1-x) - cln::log(x) * cln::log(1-x) + cln::zeta(2);
+                               if ( x == 1 ) {
+                                       return cln::zeta(2);
+                               } else {
+                                       return -Li2_do_sum(1-x) - cln::log(x) * cln::log(1-x) + cln::zeta(2);
+                               }
                        } else {
                                return -Li2_do_sum_Xn(1-x) - cln::log(x) * cln::log(1-x) + cln::zeta(2);
                        }
@@ -368,7 +372,8 @@ cln::cl_N Li_projection(int n, const cln::cl_N& x, const cln::float_format_t& pr
                                return Lin_do_sum_Xn(n, x);
                        }
                } else {
-                       cln::cl_N result = -cln::expt(cln::log(x), n-1) * cln::log(1-x) / cln::factorial(n-1);
+                       cln::cl_N result = 0;
+                       if ( x != 1 ) result = -cln::expt(cln::log(x), n-1) * cln::log(1-x) / cln::factorial(n-1);
                        for (int j=0; j<n-1; j++) {
                                result = result + (S_num(n-j-1, 1, 1) - S_num(1, n-j-1, 1-x))
                                                  * cln::expt(cln::log(x), j) / cln::factorial(j);
@@ -2601,7 +2606,7 @@ ex trafo_H_1tx_prepend_zero(const ex& e, const ex& arg)
                ex addzeta = convert_H_to_zeta(newparameter);
                return e.subs(h == (addzeta-H(newparameter, h.op(1)).hold())).expand();
        } else {
-               return e * (-H(lst(0),1/arg).hold());
+               return e * (-H(lst(ex(0)),1/arg).hold());
        }
 }
 
@@ -2632,7 +2637,7 @@ ex trafo_H_prepend_one(const ex& e, const ex& arg)
                newparameter.prepend(1);
                return e.subs(h == H(newparameter, h.op(1)).hold());
        } else {
-               return e * H(lst(1),1-arg).hold();
+               return e * H(lst(ex(1)),1-arg).hold();
        }
 }
 
@@ -2664,8 +2669,8 @@ ex trafo_H_1tx_prepend_minusone(const ex& e, const ex& arg)
                ex addzeta = convert_H_to_zeta(newparameter);
                return e.subs(h == (addzeta-H(newparameter, h.op(1)).hold())).expand();
        } else {
-               ex addzeta = convert_H_to_zeta(lst(-1));
-               return (e * (addzeta - H(lst(-1),1/arg).hold())).expand();
+               ex addzeta = convert_H_to_zeta(lst(ex(-1)));
+               return (e * (addzeta - H(lst(ex(-1)),1/arg).hold())).expand();
        }
 }
 
@@ -2696,7 +2701,7 @@ ex trafo_H_1mxt1px_prepend_minusone(const ex& e, const ex& arg)
                newparameter.prepend(-1);
                return e.subs(h == H(newparameter, h.op(1)).hold()).expand();
        } else {
-               return (e * H(lst(-1),(1-arg)/(1+arg)).hold()).expand();
+               return (e * H(lst(ex(-1)),(1-arg)/(1+arg)).hold()).expand();
        }
 }
 
@@ -2727,7 +2732,7 @@ ex trafo_H_1mxt1px_prepend_one(const ex& e, const ex& arg)
                newparameter.prepend(1);
                return e.subs(h == H(newparameter, h.op(1)).hold()).expand();
        } else {
-               return (e * H(lst(1),(1-arg)/(1+arg)).hold()).expand();
+               return (e * H(lst(ex(1)),(1-arg)/(1+arg)).hold()).expand();
        }
 }
 
@@ -2806,7 +2811,7 @@ struct map_trafo_H_1mx : public map_function
                                        // leading one
                                        map_trafo_H_1mx recursion;
                                        map_trafo_H_mult unify;
-                                       ex res = H(lst(1), arg).hold() * H(newparameter, arg).hold();
+                                       ex res = H(lst(ex(1)), arg).hold() * H(newparameter, arg).hold();
                                        std::size_t firstzero = 0;
                                        while (parameter.op(firstzero) == 1) {
                                                firstzero++;
@@ -2870,7 +2875,7 @@ struct map_trafo_H_1overx : public map_function
                                        }
                                        if (allthesame) {
                                                map_trafo_H_mult unify;
-                                               return unify((pow(H(lst(-1),1/arg).hold() - H(lst(0),1/arg).hold(), parameter.nops())
+                                               return unify((pow(H(lst(ex(-1)),1/arg).hold() - H(lst(ex(0)),1/arg).hold(), parameter.nops())
                                                       / factorial(parameter.nops())).expand());
                                        }
                                } else {
@@ -2882,7 +2887,7 @@ struct map_trafo_H_1overx : public map_function
                                        }
                                        if (allthesame) {
                                                map_trafo_H_mult unify;
-                                               return unify((pow(H(lst(1),1/arg).hold() + H(lst(0),1/arg).hold() + H_polesign, parameter.nops())
+                                               return unify((pow(H(lst(ex(1)),1/arg).hold() + H(lst(ex(0)),1/arg).hold() + H_polesign, parameter.nops())
                                                       / factorial(parameter.nops())).expand());
                                        }
                                }
@@ -2925,7 +2930,7 @@ struct map_trafo_H_1overx : public map_function
                                        // leading one
                                        map_trafo_H_1overx recursion;
                                        map_trafo_H_mult unify;
-                                       ex res = H(lst(1), arg).hold() * H(newparameter, arg).hold();
+                                       ex res = H(lst(ex(1)), arg).hold() * H(newparameter, arg).hold();
                                        std::size_t firstzero = 0;
                                        while (parameter.op(firstzero) == 1) {
                                                firstzero++;
@@ -2981,7 +2986,7 @@ struct map_trafo_H_1mxt1px : public map_function
                                        }
                                        if (allthesame) {
                                                map_trafo_H_mult unify;
-                                               return unify((pow(-H(lst(1),(1-arg)/(1+arg)).hold() - H(lst(-1),(1-arg)/(1+arg)).hold(), parameter.nops())
+                                               return unify((pow(-H(lst(ex(1)),(1-arg)/(1+arg)).hold() - H(lst(ex(-1)),(1-arg)/(1+arg)).hold(), parameter.nops())
                                                       / factorial(parameter.nops())).expand());
                                        }
                                } else if (parameter.op(0) == -1) {
@@ -2993,7 +2998,7 @@ struct map_trafo_H_1mxt1px : public map_function
                                        }
                                        if (allthesame) {
                                                map_trafo_H_mult unify;
-                                               return unify((pow(log(2) - H(lst(-1),(1-arg)/(1+arg)).hold(), parameter.nops())
+                                               return unify((pow(log(2) - H(lst(ex(-1)),(1-arg)/(1+arg)).hold(), parameter.nops())
                                                       / factorial(parameter.nops())).expand());
                                        }
                                } else {
@@ -3005,7 +3010,7 @@ struct map_trafo_H_1mxt1px : public map_function
                                        }
                                        if (allthesame) {
                                                map_trafo_H_mult unify;
-                                               return unify((pow(-log(2) - H(lst(0),(1-arg)/(1+arg)).hold() + H(lst(-1),(1-arg)/(1+arg)).hold(), parameter.nops())
+                                               return unify((pow(-log(2) - H(lst(ex(0)),(1-arg)/(1+arg)).hold() + H(lst(ex(-1)),(1-arg)/(1+arg)).hold(), parameter.nops())
                                                       / factorial(parameter.nops())).expand());
                                        }
                                }
@@ -3048,7 +3053,7 @@ struct map_trafo_H_1mxt1px : public map_function
                                        // leading one
                                        map_trafo_H_1mxt1px recursion;
                                        map_trafo_H_mult unify;
-                                       ex res = H(lst(1), arg).hold() * H(newparameter, arg).hold();
+                                       ex res = H(lst(ex(1)), arg).hold() * H(newparameter, arg).hold();
                                        std::size_t firstzero = 0;
                                        while (parameter.op(firstzero) == 1) {
                                                firstzero++;