- added prefix and postfix increment and decrement operators for class numeric

- made all function_index const
- added overloaded psi(x) and psi(n,x) with some sensible behaviour

diff --git a/ginac/function.pl b/ginac/function.pl
index 06bf2d674c3c92d6aa4b615ad92b2b7ecf49f0cd..0e6cacfca222c8e7ff6a212ae332d232a1890655 100755 (executable)
--- a/ginac/function.pl
+++ b/ginac/function.pl
@@ -34,7 +34,7 @@ sub generate {
 $declare_function_macro=generate(
     <<'END_OF_DECLARE_FUNCTION_MACRO','GiNaC::ex const & p${N}','p${N}');
 #define DECLARE_FUNCTION_${N}P(NAME) \\
-extern unsigned function_index_##NAME; \\
+extern const unsigned function_index_##NAME; \\
 inline GiNaC::function NAME(${SEQ1}) { \\
     return GiNaC::function(function_index_##NAME, ${SEQ2}); \\
 }
@@ -155,7 +155,7 @@ $declare_function_macro
 // end of generated lines
 
 #define REGISTER_FUNCTION(NAME,E,EF,D,S) \\
-unsigned function_index_##NAME=GiNaC::function::register_new(#NAME,E,EF,D,S);
+const unsigned function_index_##NAME=GiNaC::function::register_new(#NAME,E,EF,D,S);
 
 #define BEGIN_TYPECHECK \\
 bool automatic_typecheck=true;

diff --git a/ginac/inifcns.h b/ginac/inifcns.h
index 725c1ee72e0eb817bcc0b3834e084110448030c3..8be87d64cf15a9e18082bff5ff07c1a49ef2204a 100644 (file)
--- a/ginac/inifcns.h
+++ b/ginac/inifcns.h
@@ -87,11 +87,11 @@ DECLARE_FUNCTION_1P(zeta)
 DECLARE_FUNCTION_1P(gamma)
 
 /** Psi-function (aka polygamma-function). */
-extern unsigned function_index_psi1;
+extern const unsigned function_index_psi1;
 inline function psi(ex const & p1) {
     return function(function_index_psi1, p1);
 }
-extern unsigned function_index_psi2;
+extern const unsigned function_index_psi2;
 inline function psi(ex const & p1, ex const & p2) {
     return function(function_index_psi2, p1, p2);
 }

diff --git a/ginac/inifcns_gamma.cpp b/ginac/inifcns_gamma.cpp
index a208592f82eab47de23bcfde1d83aeaea4095297..5610d48eb14c010a159bbb537cd0a189b80a9717 100644 (file)
--- a/ginac/inifcns_gamma.cpp
+++ b/ginac/inifcns_gamma.cpp
@@ -113,10 +113,30 @@ REGISTER_FUNCTION(gamma, gamma_eval, gamma_evalf, gamma_diff, gamma_series);
 static ex psi1_eval(ex const & x)
 {
     if (x.info(info_flags::numeric)) {
-        // do some stuff...
+        if (x.info(info_flags::integer) && !x.info(info_flags::positive))
+            throw (std::domain_error("psi_eval(): simple pole"));
+        if (x.info(info_flags::positive)) {
+            // psi(n) -> 1 + 1/2 +...+ 1/(n-1) - EulerGamma
+            if (x.info(info_flags::integer)) {
+                numeric rat(0);
+                for (numeric i(ex_to_numeric(x)-numONE()); i.is_positive(); --i)
+                    rat += i.inverse();
+                return rat-EulerGamma;
+            }
+            // psi((2m+1)/2) -> 2/(2m+1) + 2/2m +...+ 2/1 - EulerGamma - 2log(2)
+            if ((exTWO()*x).info(info_flags::integer)) {
+                numeric rat(0);
+                for (numeric i((ex_to_numeric(x)-numONE())*numTWO()); i.is_positive(); i-=numTWO())
+                    rat += numTWO()*i.inverse();
+                return rat-EulerGamma-exTWO()*log(exTWO());
+            }
+            if (x.compare(exONE())==1) {
+                // should call numeric, since >1
+            }
+        }
     }
     return psi(x).hold();
-}    
+}
 
 static ex psi1_evalf(ex const & x)
 {
@@ -134,12 +154,7 @@ static ex psi1_diff(ex const & x, unsigned diff_param)
     return psi(exONE(), x);
 }
 
-static ex psi1_series(ex const & x, symbol const & s, ex const & point, int order)
-{
-    throw(std::logic_error("Nobody told me how to series expand the psi function. :-("));
-}
-
-unsigned function_index_psi1 = function::register_new("psi", psi1_eval, psi1_evalf, psi1_diff, psi1_series);
+const unsigned function_index_psi1 = function::register_new("psi", psi1_eval, psi1_evalf, psi1_diff, NULL);
 
 //////////
 // Psi-functions (aka polygamma-functions)  psi(0,x)==psi(x)
@@ -151,7 +166,7 @@ static ex psi2_eval(ex const & n, ex const & x)
 {
     // psi(0,x) -> psi(x)
     if (n.is_zero())
-        return psi(x).hold();
+        return psi(x);
     if (n.info(info_flags::numeric) && x.info(info_flags::numeric)) {
         // do some stuff...
     }
@@ -180,11 +195,6 @@ static ex psi2_diff(ex const & n, ex const & x, unsigned diff_param)
     return psi(n+1, x);
 }
 
-static ex psi2_series(ex const & n, ex const & x, symbol const & s, ex const & point, int order)
-{
-    throw(std::logic_error("Nobody told me how to series expand the psi functions. :-("));
-}
-
-unsigned function_index_psi2 = function::register_new("psi", psi2_eval, psi2_evalf, psi2_diff, psi2_series);
+const unsigned function_index_psi2 = function::register_new("psi", psi2_eval, psi2_evalf, psi2_diff, NULL);
 
 } // namespace GiNaC

diff --git a/ginac/inifcns_zeta.cpp b/ginac/inifcns_zeta.cpp
index 1a850ac43ebfe7f9dac2583eee6fac7b2afebf96..cd63876506fdb566fdaa681339f4308cfb84442c 100644 (file)
--- a/ginac/inifcns_zeta.cpp
+++ b/ginac/inifcns_zeta.cpp
@@ -44,7 +44,7 @@ static ex zeta_eval(ex const & x)
         if (y.is_integer()) {
             if (y.is_zero())
                 return -exHALF();
-            if (!x.compare(exONE()))
+            if (x.is_equal(exONE()))
                 throw(std::domain_error("zeta(1): infinity"));
             if (x.info(info_flags::posint)) {
                 if (x.info(info_flags::odd))

diff --git a/ginac/numeric.h b/ginac/numeric.h
index 13c790fc6a754e2f939d8e4d6e2c47efc57fc1aa..eb3565907148ffa34537f7e9b6d65afc841f26ea 100644 (file)
--- a/ginac/numeric.h
+++ b/ginac/numeric.h
@@ -317,7 +317,8 @@ inline numeric numer(numeric const & x)
 inline numeric denom(numeric const & x)
 { return x.denom(); }
 
-ex IEvalf(void);
+// numeric evaluation functions for class constant objects:
+
 ex PiEvalf(void);
 ex EulerGammaEvalf(void);
 ex CatalanEvalf(void);

diff --git a/ginac/operators.cpp b/ginac/operators.cpp
index d9882c0e235f5ba5c36e22461c7dc4e1c27e2c38..4d762b305c2424ccb623ad5a6e11a61bfa65681e 100644 (file)
--- a/ginac/operators.cpp
+++ b/ginac/operators.cpp
@@ -272,7 +272,39 @@ numeric operator+(numeric const & lh)
 
 numeric operator-(numeric const & lh)
 {
-    return (numeric(-1)*lh);
+    return numMINUSONE()*lh;
+}
+
+/** Numeric prefix increment.  Adds 1 and returns incremented number. */
+numeric& operator++(numeric & rh)
+{
+    rh = rh+numONE();
+    return rh;
+}
+
+/** Numeric prefix decrement.  Subtracts 1 and returns decremented number. */
+numeric& operator--(numeric & rh)
+{
+    rh = rh-numONE();
+    return rh;
+}
+
+/** Numeric postfix increment.  Returns the number and leaves the original
+ *  incremented by 1. */
+numeric operator++(numeric & lh, int)
+{
+    numeric tmp = lh;
+    lh = lh+numONE();
+    return tmp;
+}
+
+/** Numeric Postfix decrement.  Returns the number and leaves the original
+ *  decremented by 1. */
+numeric operator--(numeric & lh, int)
+{
+    numeric tmp = lh;
+    lh = lh-numONE();
+    return tmp;
 }
 
 // binary relational operators ex with ex

diff --git a/ginac/operators.h b/ginac/operators.h
index 4c95a55164a81adf4b7c67325a009b3361b74346..fab2870b4fff53152dd838342f5bdd174962879c 100644 (file)
--- a/ginac/operators.h
+++ b/ginac/operators.h
@@ -92,6 +92,10 @@ ex operator-(ex const & lh);
 
 numeric operator+(numeric const & lh);
 numeric operator-(numeric const & lh);
+numeric& operator++(numeric & rh);
+numeric& operator--(numeric & rh);
+numeric operator++(numeric & lh, int);
+numeric operator--(numeric & lh, int);
 
 // binary relational operators ex with ex
 relational operator==(ex const & lh, ex const & rh);