ginac/inifcns_nstdsums.cpp

   1 /** @file inifcns_nstdsums.cpp
   2  *
   3  *  Implementation of some special functions that have a representation as nested sums.
   4  *  The functions are:
   5  *    classical polylogarithm              Li(n,x)
   6  *    multiple polylogarithm               Li(lst(n_1,...,n_k),lst(x_1,...,x_k)
   7  *    nielsen's generalized polylogarithm  S(n,p,x)
   8  *    harmonic polylogarithm               H(lst(m_1,...,m_k),x)
   9  *    multiple zeta value                  mZeta(lst(m_1,...,m_k))
  10  *
  11  *  Some remarks:
  12  *    - All formulae used can be looked up in the following publication:
  13  *      Nielsen's Generalized Polylogarithms, K.S.Kolbig, SIAM J.Math.Anal. 17 (1986), pp. 1232-1258.
  14  *      This document will be referenced as [Kol] throughout this source code.
  15  *    - Classical polylogarithms (Li) and nielsen's generalized polylogarithms (S) can be numerically
  16  *      evaluated in the whole complex plane except for S(n,p,-1) when p is not unit (no formula yet
  17  *      to tackle these points). And of course, there is still room for speed optimizations ;-).
  18  *    - The remaining functions can only be numerically evaluated with arguments lying in the unit sphere
  19  *      at the moment. Sorry. The evaluation especially for mZeta is very slow ... better not use it
  20  *      right now.
  21  *    - The functions have no series expansion. To do it, you have to convert these functions
  22  *      into the appropriate objects from the nestedsums library, do the expansion and convert the
  23  *      result back.
  24  *    - Numerical testing of this implementation has been performed by doing a comparison of results
  25  *      between this software and the commercial M.......... 4.1.
  26  *
  27  */
  28
  29 /*
  30  *  GiNaC Copyright (C) 1999-2003 Johannes Gutenberg University Mainz, Germany
  31  *
  32  *  This program is free software; you can redistribute it and/or modify
  33  *  it under the terms of the GNU General Public License as published by
  34  *  the Free Software Foundation; either version 2 of the License, or
  35  *  (at your option) any later version.
  36  *
  37  *  This program is distributed in the hope that it will be useful,
  38  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
  39  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  40  *  GNU General Public License for more details.
  41  *
  42  *  You should have received a copy of the GNU General Public License
  43  *  along with this program; if not, write to the Free Software
  44  *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  45  */
  46
  47 #include <stdexcept>
  48 #include <vector>
  49 #include <cln/cln.h>
  50
  51 #include "inifcns.h"
  52 #include "lst.h"
  53 #include "numeric.h"
  54 #include "operators.h"
  55 #include "relational.h"
  56
  57
  58 namespace GiNaC {
  59
  60
  61 //////////////////////
  62 // helper functions //
  63 //////////////////////
  64
  65
  66 // helper function for classical polylog Li
  67 static cln::cl_N Li_series(int n, const cln::cl_N& x, const cln::float_format_t& prec)
  68 {
  69         // Note: argument must be in the unit circle
  70         cln::cl_N aug, acc;
  71         cln::cl_N num = cln::complex(cln::cl_float(1, prec), 0);
  72         cln::cl_N den = 0;
  73         int i = 1;
  74         do {
  75                 num = num * x;
  76                 cln::cl_R ii = i;
  77                 den = cln::expt(ii, n);
  78                 i++;
  79                 aug = num / den;
  80                 acc = acc + aug;
  81         } while (acc != acc+aug);
  82         return acc;
  83 }
  84
  85
  86 // helper function for classical polylog Li
  87 static cln::cl_N Li_projection(int n, const cln::cl_N& x, const cln::float_format_t& prec)
  88 {
  89         return Li_series(n, x, prec);
  90 }
  91
  92
  93 // helper function for classical polylog Li
  94 static numeric Li_num(int n, const numeric& x)
  95 {
  96         if (n == 1) {
  97                 // just a log
  98                 return -cln::log(1-x.to_cl_N());
  99         }
 100         if (x.is_zero()) {
 101                 return 0;
 102         }
 103         if (x == 1) {
 104                 // [Kol] (2.22)
 105                 return cln::zeta(n);
 106         }
 107         else if (x == -1) {
 108                 // [Kol] (2.22)
 109                 return -(1-cln::expt(cln::cl_I(2),1-n)) * cln::zeta(n);
 110         }
 111
 112         // what is the desired float format?
 113         // first guess: default format
 114         cln::float_format_t prec = cln::default_float_format;
 115         const cln::cl_N value = x.to_cl_N();
 116         // second guess: the argument's format
 117         if (!x.real().is_rational())
 118                 prec = cln::float_format(cln::the<cln::cl_F>(cln::realpart(value)));
 119         else if (!x.imag().is_rational())
 120                 prec = cln::float_format(cln::the<cln::cl_F>(cln::imagpart(value)));
 121
 122         // [Kol] (5.15)
 123         if (cln::abs(value) > 1) {
 124                 cln::cl_N result = -cln::expt(cln::log(-value),n) / cln::factorial(n);
 125                 // check if argument is complex. if it is real, the new polylog has to be conjugated.
 126                 if (cln::zerop(cln::imagpart(value))) {
 127                         if (n & 1) {
 128                                 result = result + conjugate(Li_projection(n, cln::recip(value), prec));
 129                         }
 130                         else {
 131                                 result = result - conjugate(Li_projection(n, cln::recip(value), prec));
 132                         }
 133                 }
 134                 else {
 135                         if (n & 1) {
 136                                 result = result + Li_projection(n, cln::recip(value), prec);
 137                         }
 138                         else {
 139                                 result = result - Li_projection(n, cln::recip(value), prec);
 140                         }
 141                 }
 142                 cln::cl_N add;
 143                 for (int j=0; j<n-1; j++) {
 144                         add = add + (1+cln::expt(cln::cl_I(-1),n-j)) * (1-cln::expt(cln::cl_I(2),1-n+j))
 145                                         * Li_num(n-j,1).to_cl_N() * cln::expt(cln::log(-value),j) / cln::factorial(j);
 146                 }
 147                 result = result - add;
 148                 return result;
 149         }
 150         else {
 151                 return Li_projection(n, value, prec);
 152         }
 153 }
 154
 155
 156 // helper function for S(n,p,x)
 157 static cln::cl_N numeric_nielsen(int n, int step)
 158 {
 159         if (step) {
 160                 cln::cl_N res;
 161                 for (int i=1; i<n; i++) {
 162                         res = res + numeric_nielsen(i, step-1) / cln::cl_I(i);
 163                 }
 164                 return res;
 165         }
 166         else {
 167                 return 1;
 168         }
 169 }
 170
 171
 172 // forward declaration needed by function C below
 173 static numeric S_num(int n, int p, const numeric& x);
 174
 175
 176 // helper function for S(n,p,x)
 177 // [Kol] (7.2)
 178 static cln::cl_N C(int n, int p)
 179 {
 180         cln::cl_N result;
 181
 182         for (int k=0; k<p; k++) {
 183                 for (int j=0; j<=(n+k-1)/2; j++) {
 184                         if (k == 0) {
 185                                 if (n & 1) {
 186                                         if (j & 1) {
 187                                                 result = result - 2 * cln::expt(cln::pi(),2*j) * S_num(n-2*j,p,1).to_cl_N() / cln::factorial(2*j);
 188                                         }
 189                                         else {
 190                                                 result = result + 2 * cln::expt(cln::pi(),2*j) * S_num(n-2*j,p,1).to_cl_N() / cln::factorial(2*j);
 191                                         }
 192                                 }
 193                         }
 194                         else {
 195                                 if (k & 1) {
 196                                         if (j & 1) {
 197                                                 result = result + cln::factorial(n+k-1)
 198                                                         * cln::expt(cln::pi(),2*j) * S_num(n+k-2*j,p-k,1).to_cl_N()
 199                                                         / (cln::factorial(k) * cln::factorial(n-1) * cln::factorial(2*j));
 200                                         }
 201                                         else {
 202                                                 result = result - cln::factorial(n+k-1)
 203                                                         * cln::expt(cln::pi(),2*j) * S_num(n+k-2*j,p-k,1).to_cl_N()
 204                                                         / (cln::factorial(k) * cln::factorial(n-1) * cln::factorial(2*j));
 205                                         }
 206                                 }
 207                                 else {
 208                                         if (j & 1) {
 209                                                 result = result - cln::factorial(n+k-1) * cln::expt(cln::pi(),2*j) * S_num(n+k-2*j,p-k,1).to_cl_N()
 210                                                         / (cln::factorial(k) * cln::factorial(n-1) * cln::factorial(2*j));
 211                                         }
 212                                         else {
 213                                                 result = result + cln::factorial(n+k-1)
 214                                                         * cln::expt(cln::pi(),2*j) * S_num(n+k-2*j,p-k,1).to_cl_N()
 215                                                         / (cln::factorial(k) * cln::factorial(n-1) * cln::factorial(2*j));
 216                                         }
 217                                 }
 218                         }
 219                 }
 220         }
 221         int np = n+p;
 222         if ((np-1) & 1) {
 223                 if (((np)/2+n) & 1) {
 224                         result = -result - cln::expt(cln::pi(),np) / (np * cln::factorial(n-1) * cln::factorial(p));
 225                 }
 226                 else {
 227                         result = -result + cln::expt(cln::pi(),np) / (np * cln::factorial(n-1) * cln::factorial(p));
 228                 }
 229         }
 230
 231         return result;
 232 }
 233
 234
 235 // helper function for S(n,p,x)
 236 // [Kol] remark to (9.1)
 237 static cln::cl_N a_k(int k)
 238 {
 239         cln::cl_N result;
 240
 241         if (k == 0) {
 242                 return 1;
 243         }
 244
 245         result = result;
 246         for (int m=2; m<=k; m++) {
 247                 result = result + cln::expt(cln::cl_N(-1),m) * cln::zeta(m) * a_k(k-m);
 248         }
 249
 250         return -result / k;
 251 }
 252
 253
 254 // helper function for S(n,p,x)
 255 // [Kol] remark to (9.1)
 256 static cln::cl_N b_k(int k)
 257 {
 258         cln::cl_N result;
 259
 260         if (k == 0) {
 261                 return 1;
 262         }
 263
 264         result = result;
 265         for (int m=2; m<=k; m++) {
 266                 result = result + cln::expt(cln::cl_N(-1),m) * cln::zeta(m) * b_k(k-m);
 267         }
 268
 269         return result / k;
 270 }
 271
 272
 273 // helper function for S(n,p,x)
 274 static cln::cl_N S_series(int n, int p, const cln::cl_N& x, const cln::float_format_t& prec)
 275 {
 276         n++;
 277         int i = p;
 278         p--;
 279         cln::cl_N aug, acc;
 280         cln::cl_N num = cln::expt(x,p);
 281         cln::cl_N converter = cln::complex(cln::cl_float(1, prec), 0);
 282         cln::cl_N den = 0;
 283         do {
 284                 num = num * x;
 285                 den = cln::expt(cln::cl_I(i), n);
 286                 aug = num / den * numeric_nielsen(i, p);
 287                 i++;
 288                 acc = acc + aug;
 289         } while (acc != acc+aug);
 290
 291         return acc;
 292 }
 293
 294
 295 // helper function for S(n,p,x)
 296 static cln::cl_N S_projection(int n, int p, const cln::cl_N& x, const cln::float_format_t& prec)
 297 {
 298         // [Kol] (5.3)
 299         if (cln::abs(cln::realpart(x)) > cln::cl_F("0.5")) {
 300
 301                 cln::cl_N result = cln::expt(cln::cl_I(-1),p) * cln::expt(cln::log(x),n)
 302                         * cln::expt(cln::log(1-x),p) / cln::factorial(n) / cln::factorial(p);
 303
 304                 for (int s=0; s<n; s++) {
 305                         cln::cl_N res2;
 306                         for (int r=0; r<p; r++) {
 307                                 res2 = res2 + cln::expt(cln::cl_I(-1),r) * cln::expt(cln::log(1-x),r)
 308                                         * S_series(p-r,n-s,1-x,prec) / cln::factorial(r);
 309                         }
 310                         result = result + cln::expt(cln::log(x),s) * (S_num(n-s,p,1).to_cl_N() - res2) / cln::factorial(s);
 311                 }
 312
 313                 return result;
 314         }
 315
 316         return S_series(n, p, x, prec);
 317 }
 318
 319
 320 // helper function for S(n,p,x)
 321 static numeric S_num(int n, int p, const numeric& x)
 322 {
 323         if (x == 1) {
 324                 if (n == 1) {
 325                     // [Kol] (2.22) with (2.21)
 326                         return cln::zeta(p+1);
 327                 }
 328
 329                 if (p == 1) {
 330                     // [Kol] (2.22)
 331                         return cln::zeta(n+1);
 332                 }
 333
 334                 // [Kol] (9.1)
 335                 cln::cl_N result;
 336                 for (int nu=0; nu<n; nu++) {
 337                         for (int rho=0; rho<=p; rho++) {
 338                                 result = result + b_k(n-nu-1) * b_k(p-rho) * a_k(nu+rho+1)
 339                                         * cln::factorial(nu+rho+1) / cln::factorial(rho) / cln::factorial(nu+1);
 340                         }
 341                 }
 342                 result = result * cln::expt(cln::cl_I(-1),n+p-1);
 343
 344                 return result;
 345         }
 346         else if (x == -1) {
 347                 // [Kol] (2.22)
 348                 if (p == 1) {
 349                         return -(1-cln::expt(cln::cl_I(2),-n)) * cln::zeta(n+1);
 350                 }
 351                 throw std::runtime_error("don't know how to evaluate this function!");
 352         }
 353
 354         // what is the desired float format?
 355         // first guess: default format
 356         cln::float_format_t prec = cln::default_float_format;
 357         const cln::cl_N value = x.to_cl_N();
 358         // second guess: the argument's format
 359         if (!x.real().is_rational())
 360                 prec = cln::float_format(cln::the<cln::cl_F>(cln::realpart(value)));
 361         else if (!x.imag().is_rational())
 362                 prec = cln::float_format(cln::the<cln::cl_F>(cln::imagpart(value)));
 363
 364
 365         // [Kol] (5.3)
 366         if (cln::realpart(value) < -0.5) {
 367
 368                 cln::cl_N result = cln::expt(cln::cl_I(-1),p) * cln::expt(cln::log(value),n)
 369                         * cln::expt(cln::log(1-value),p) / cln::factorial(n) / cln::factorial(p);
 370
 371                 for (int s=0; s<n; s++) {
 372                         cln::cl_N res2;
 373                         for (int r=0; r<p; r++) {
 374                                 res2 = res2 + cln::expt(cln::cl_I(-1),r) * cln::expt(cln::log(1-value),r)
 375                                         * S_num(p-r,n-s,1-value).to_cl_N() / cln::factorial(r);
 376                         }
 377                         result = result + cln::expt(cln::log(value),s) * (S_num(n-s,p,1).to_cl_N() - res2) / cln::factorial(s);
 378                 }
 379
 380                 return result;
 381
 382         }
 383         // [Kol] (5.12)
 384         else if (cln::abs(value) > 1) {
 385
 386                 cln::cl_N result;
 387
 388                 for (int s=0; s<p; s++) {
 389                         for (int r=0; r<=s; r++) {
 390                                 result = result + cln::expt(cln::cl_I(-1),s) * cln::expt(cln::log(-value),r) * cln::factorial(n+s-r-1)
 391                                         / cln::factorial(r) / cln::factorial(s-r) / cln::factorial(n-1)
 392                                         * S_num(n+s-r,p-s,cln::recip(value)).to_cl_N();
 393                         }
 394                 }
 395                 result = result * cln::expt(cln::cl_I(-1),n);
 396
 397                 cln::cl_N res2;
 398                 for (int r=0; r<n; r++) {
 399                         res2 = res2 + cln::expt(cln::log(-value),r) * C(n-r,p) / cln::factorial(r);
 400                 }
 401                 res2 = res2 + cln::expt(cln::log(-value),n+p) / cln::factorial(n+p);
 402
 403                 result = result + cln::expt(cln::cl_I(-1),p) * res2;
 404
 405                 return result;
 406         }
 407         else {
 408                 return S_projection(n, p, value, prec);
 409         }
 410 }
 411
 412
 413 // helper function for multiple polylogarithm
 414 static cln::cl_N numeric_zsum(int n, std::vector<cln::cl_N>& x, std::vector<cln::cl_N>& m)
 415 {
 416         cln::cl_N res;
 417         if (x.empty()) {
 418                 return 1;
 419         }
 420         for (int i=1; i<n; i++) {
 421                 std::vector<cln::cl_N>::iterator be;
 422                 std::vector<cln::cl_N>::iterator en;
 423                 be = x.begin();
 424                 be++;
 425                 en = x.end();
 426                 std::vector<cln::cl_N> xbuf(be, en);
 427                 be = m.begin();
 428                 be++;
 429                 en = m.end();
 430                 std::vector<cln::cl_N> mbuf(be, en);
 431                 res = res + cln::expt(x[0],i) / cln::expt(i,m[0]) * numeric_zsum(i, xbuf, mbuf);
 432         }
 433         return res;
 434 }
 435
 436
 437 // helper function for harmonic polylogarithm
 438 static cln::cl_N numeric_harmonic(int n, std::vector<cln::cl_N>& m)
 439 {
 440         cln::cl_N res;
 441         if (m.empty()) {
 442                 return 1;
 443         }
 444         for (int i=1; i<n; i++) {
 445                 std::vector<cln::cl_N>::iterator be;
 446                 std::vector<cln::cl_N>::iterator en;
 447                 be = m.begin();
 448                 be++;
 449                 en = m.end();
 450                 std::vector<cln::cl_N> mbuf(be, en);
 451                 res = res + cln::recip(cln::expt(i,m[0])) * numeric_harmonic(i, mbuf);
 452         }
 453         return res;
 454 }
 455
 456
 457 /////////////////////////////
 458 // end of helper functions //
 459 /////////////////////////////
 460
 461
 462 // Polylogarithm and multiple polylogarithm
 463
 464 static ex Li_eval(const ex& x1, const ex& x2)
 465 {
 466         if (x2.is_zero()) {
 467                 return 0;
 468         }
 469         else {
 470                 return Li(x1,x2).hold();
 471         }
 472 }
 473
 474 static ex Li_evalf(const ex& x1, const ex& x2)
 475 {
 476         // classical polylogs
 477         if (is_a<numeric>(x1) && is_a<numeric>(x2)) {
 478                 return Li_num(ex_to<numeric>(x1).to_int(), ex_to<numeric>(x2));
 479         }
 480         // multiple polylogs
 481         else if (is_a<lst>(x1) && is_a<lst>(x2)) {
 482                 for (int i=0; i<x1.nops(); i++) {
 483                         if (!is_a<numeric>(x1.op(i)))
 484                                 return Li(x1,x2).hold();
 485                         if (!is_a<numeric>(x2.op(i)))
 486                                 return Li(x1,x2).hold();
 487                         if (x2 >= 1)
 488                                 return Li(x1,x2).hold();
 489                 }
 490
 491                 cln::cl_N m_1 = ex_to<numeric>(x1.op(x1.nops()-1)).to_cl_N();
 492                 cln::cl_N x_1 = ex_to<numeric>(x2.op(x2.nops()-1)).to_cl_N();
 493                 std::vector<cln::cl_N> x;
 494                 std::vector<cln::cl_N> m;
 495                 const int nops = ex_to<numeric>(x1.nops()).to_int();
 496                 for (int i=nops-2; i>=0; i--) {
 497                         m.push_back(ex_to<numeric>(x1.op(i)).to_cl_N());
 498                         x.push_back(ex_to<numeric>(x2.op(i)).to_cl_N());
 499                 }
 500
 501                 cln::cl_N res;
 502                 cln::cl_N resbuf;
 503                 for (int i=nops; true; i++) {
 504                         resbuf = res;
 505                         res = res + cln::expt(x_1,i) / cln::expt(i,m_1) * numeric_zsum(i, x, m);
 506                         if (cln::zerop(res-resbuf))
 507                                 break;
 508                 }
 509
 510                 return numeric(res);
 511
 512         }
 513
 514         return Li(x1,x2).hold();
 515 }
 516
 517 REGISTER_FUNCTION(Li, eval_func(Li_eval).evalf_func(Li_evalf).do_not_evalf_params());
 518
 519
 520 // Nielsen's generalized polylogarithm
 521
 522 static ex S_eval(const ex& x1, const ex& x2, const ex& x3)
 523 {
 524         if (x2 == 1) {
 525                 return Li(x1+1,x3);
 526         }
 527         return S(x1,x2,x3).hold();
 528 }
 529
 530 static ex S_evalf(const ex& x1, const ex& x2, const ex& x3)
 531 {
 532         if (is_a<numeric>(x1) && is_a<numeric>(x2) && is_a<numeric>(x3)) {
 533                 if ((x3 == -1) && (x2 != 1)) {
 534                         // no formula to evaluate this ... sorry
 535                         return S(x1,x2,x3).hold();
 536                 }
 537                 return S_num(ex_to<numeric>(x1).to_int(), ex_to<numeric>(x2).to_int(), ex_to<numeric>(x3));
 538         }
 539         return S(x1,x2,x3).hold();
 540 }
 541
 542 REGISTER_FUNCTION(S, eval_func(S_eval).evalf_func(S_evalf).do_not_evalf_params());
 543
 544
 545 // Harmonic polylogarithm
 546
 547 static ex H_eval(const ex& x1, const ex& x2)
 548 {
 549         return H(x1,x2).hold();
 550 }
 551
 552 static ex H_evalf(const ex& x1, const ex& x2)
 553 {
 554         if (is_a<lst>(x1) && is_a<numeric>(x2)) {
 555                 for (int i=0; i<x1.nops(); i++) {
 556                         if (!is_a<numeric>(x1.op(i)))
 557                                 return H(x1,x2).hold();
 558                 }
 559
 560                 cln::cl_N m_1 = ex_to<numeric>(x1.op(x1.nops()-1)).to_cl_N();
 561                 cln::cl_N x_1 = ex_to<numeric>(x2).to_cl_N();
 562                 std::vector<cln::cl_N> m;
 563                 const int nops = ex_to<numeric>(x1.nops()).to_int();
 564                 for (int i=nops-2; i>=0; i--) {
 565                         m.push_back(ex_to<numeric>(x1.op(i)).to_cl_N());
 566                 }
 567
 568                 cln::cl_N res;
 569                 cln::cl_N resbuf;
 570                 for (int i=nops; true; i++) {
 571                         resbuf = res;
 572                         res = res + cln::expt(x_1,i) / cln::expt(i,m_1) * numeric_harmonic(i, m);
 573                         if (cln::zerop(res-resbuf))
 574                                 break;
 575                 }
 576
 577                 return numeric(res);
 578
 579         }
 580
 581         return H(x1,x2).hold();
 582 }
 583
 584 REGISTER_FUNCTION(H, eval_func(H_eval).evalf_func(H_evalf).do_not_evalf_params());
 585
 586
 587 // Multiple zeta value
 588
 589 static ex mZeta_eval(const ex& x1)
 590 {
 591         return mZeta(x1).hold();
 592 }
 593
 594 static ex mZeta_evalf(const ex& x1)
 595 {
 596         if (is_a<lst>(x1)) {
 597                 for (int i=0; i<x1.nops(); i++) {
 598                         if (!is_a<numeric>(x1.op(i)))
 599                                 return mZeta(x1).hold();
 600                 }
 601
 602                 cln::cl_N m_1 = ex_to<numeric>(x1.op(x1.nops()-1)).to_cl_N();
 603                 std::vector<cln::cl_N> m;
 604                 const int nops = ex_to<numeric>(x1.nops()).to_int();
 605                 for (int i=nops-2; i>=0; i--) {
 606                         m.push_back(ex_to<numeric>(x1.op(i)).to_cl_N());
 607                 }
 608
 609                 cln::float_format_t prec = cln::default_float_format;
 610                 cln::cl_N res = cln::complex(cln::cl_float(0, prec), 0);
 611                 cln::cl_N resbuf;
 612                 for (int i=nops; true; i++) {
 613                         // to infinity and beyond ... timewise
 614                         resbuf = res;
 615                         res = res + cln::recip(cln::expt(i,m_1)) * numeric_harmonic(i, m);
 616                         if (cln::zerop(res-resbuf))
 617                                 break;
 618                 }
 619
 620                 return numeric(res);
 621
 622         }
 623
 624         return mZeta(x1).hold();
 625 }
 626
 627 REGISTER_FUNCTION(mZeta, eval_func(mZeta_eval).evalf_func(mZeta_evalf).do_not_evalf_params());
 628
 629
 630 } // namespace GiNaC
 631