ginac/function.cppy

   1 /** @file function.cpp
   2  *
   3  *  Implementation of class of symbolic functions. */
   4
   5 /*
   6  *  This file was generated automatically by function.py.
   7  *  Please do not modify it directly, edit function.cppy instead!
   8  *  function.py options: maxargs=@maxargs@
   9  *
  10  *  GiNaC Copyright (C) 1999-2010 Johannes Gutenberg University Mainz, Germany
  11  *
  12  *  This program is free software; you can redistribute it and/or modify
  13  *  it under the terms of the GNU General Public License as published by
  14  *  the Free Software Foundation; either version 2 of the License, or
  15  *  (at your option) any later version.
  16  *
  17  *  This program is distributed in the hope that it will be useful,
  18  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20  *  GNU General Public License for more details.
  21  *
  22  *  You should have received a copy of the GNU General Public License
  23  *  along with this program; if not, write to the Free Software
  24  *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
  25  */
  26
  27 #include "function.h"
  28 #include "operators.h"
  29 #include "fderivative.h"
  30 #include "ex.h"
  31 #include "lst.h"
  32 #include "symmetry.h"
  33 #include "print.h"
  34 #include "power.h"
  35 #include "archive.h"
  36 #include "inifcns.h"
  37 #include "tostring.h"
  38 #include "utils.h"
  39 #include "hash_seed.h"
  40 #include "remember.h"
  41
  42 #include <iostream>
  43 #include <limits>
  44 #include <list>
  45 #include <stdexcept>
  46 #include <string>
  47
  48 namespace GiNaC {
  49
  50 //////////
  51 // helper class function_options
  52 //////////
  53
  54 function_options::function_options()
  55 {
  56         initialize();
  57 }
  58
  59 function_options::function_options(std::string const & n, std::string const & tn)
  60 {
  61         initialize();
  62         set_name(n, tn);
  63 }
  64
  65 function_options::function_options(std::string const & n, unsigned np)
  66 {
  67         initialize();
  68         set_name(n, std::string());
  69         nparams = np;
  70 }
  71
  72 function_options::~function_options()
  73 {
  74         // nothing to clean up at the moment
  75 }
  76
  77 void function_options::initialize()
  78 {
  79         set_name("unnamed_function", "\\\\mbox{unnamed}");
  80         nparams = 0;
  81         eval_f = evalf_f = real_part_f = imag_part_f = conjugate_f = expand_f
  82                 = derivative_f = power_f = series_f = 0;
  83         evalf_params_first = true;
  84         use_return_type = false;
  85         eval_use_exvector_args = false;
  86         evalf_use_exvector_args = false;
  87         conjugate_use_exvector_args = false;
  88         real_part_use_exvector_args = false;
  89         imag_part_use_exvector_args = false;
  90         expand_use_exvector_args = false;
  91         derivative_use_exvector_args = false;
  92         power_use_exvector_args = false;
  93         series_use_exvector_args = false;
  94         print_use_exvector_args = false;
  95         use_remember = false;
  96         functions_with_same_name = 1;
  97         symtree = 0;
  98 }
  99
 100 function_options & function_options::set_name(std::string const & n,
 101                                               std::string const & tn)
 102 {
 103         name = n;
 104         if (tn==std::string())
 105                 TeX_name = "\\\\mbox{"+name+"}";
 106         else
 107                 TeX_name = tn;
 108         return *this;
 109 }
 110
 111 function_options & function_options::latex_name(std::string const & tn)
 112 {
 113         TeX_name = tn;
 114         return *this;
 115 }
 116
 117 // the following lines have been generated for max. @maxargs@ parameters
 118 +++ for method, N in [ (f, N) for f in methods[0:-1] for N in range(1, maxargs + 1)]:
 119 function_options & function_options::@method@_func(@method@_funcp_@N@ e)
 120 {
 121         test_and_set_nparams(@N@);
 122         @method@_f = @method@_funcp(e);
 123         return *this;
 124 }
 125 ---
 126 // end of generated lines
 127
 128 +++ for method in methods[0:-1]:
 129 function_options & function_options::@method@_func(@method@_funcp_exvector e)
 130 {
 131         @method@_use_exvector_args = true;
 132         @method@_f = @method@_funcp(e);
 133         return *this;
 134 }
 135 ---
 136
 137 // end of generated lines
 138
 139 function_options & function_options::set_return_type(unsigned rt, const return_type_t* rtt)
 140 {
 141         use_return_type = true;
 142         return_type = rt;
 143         if (rtt != 0)
 144                 return_type_tinfo = *rtt;
 145         else
 146                 return_type_tinfo = make_return_type_t<function>();
 147         return *this;
 148 }
 149
 150 function_options & function_options::do_not_evalf_params()
 151 {
 152         evalf_params_first = false;
 153         return *this;
 154 }
 155
 156 function_options & function_options::remember(unsigned size,
 157                                               unsigned assoc_size,
 158                                               unsigned strategy)
 159 {
 160         use_remember = true;
 161         remember_size = size;
 162         remember_assoc_size = assoc_size;
 163         remember_strategy = strategy;
 164         return *this;
 165 }
 166
 167 function_options & function_options::overloaded(unsigned o)
 168 {
 169         functions_with_same_name = o;
 170         return *this;
 171 }
 172
 173 function_options & function_options::set_symmetry(const symmetry & s)
 174 {
 175         symtree = s;
 176         return *this;
 177 }
 178
 179 void function_options::test_and_set_nparams(unsigned n)
 180 {
 181         if (nparams==0) {
 182                 nparams = n;
 183         } else if (nparams!=n) {
 184                 // we do not throw an exception here because this code is
 185                 // usually executed before main(), so the exception could not
 186                 // be caught anyhow
 187                 std::cerr << "WARNING: " << name << "(): number of parameters ("
 188                           << n << ") differs from number set before ("
 189                           << nparams << ")" << std::endl;
 190         }
 191 }
 192
 193 void function_options::set_print_func(unsigned id, print_funcp f)
 194 {
 195         if (id >= print_dispatch_table.size())
 196                 print_dispatch_table.resize(id + 1);
 197         print_dispatch_table[id] = f;
 198 }
 199
 200 /** This can be used as a hook for external applications. */
 201 unsigned function::current_serial = 0;
 202
 203
 204 GINAC_IMPLEMENT_REGISTERED_CLASS(function, exprseq)
 205
 206 //////////
 207 // default constructor
 208 //////////
 209
 210 // public
 211
 212 function::function() : serial(0)
 213 {
 214 }
 215
 216 //////////
 217 // other constructors
 218 //////////
 219
 220 // public
 221
 222 function::function(unsigned ser) : serial(ser)
 223 {
 224 }
 225
 226 // the following lines have been generated for max. @maxargs@ parameters
 227 +++ for N in range(1, maxargs + 1):
 228 function::function(unsigned ser, @seq('const ex & param%(n)d', N)@)
 229         : exprseq(@seq('param%(n)d', N)@), serial(ser)
 230 {
 231 }
 232 ---
 233
 234 function::function(unsigned ser, const exprseq & es) : exprseq(es), serial(ser)
 235 {
 236
 237         // Force re-evaluation even if the exprseq was already evaluated
 238         // (the exprseq copy constructor copies the flags)
 239         clearflag(status_flags::evaluated);
 240 }
 241
 242 function::function(unsigned ser, const exvector & v, bool discardable)
 243   : exprseq(v,discardable), serial(ser)
 244 {
 245 }
 246
 247 function::function(unsigned ser, std::auto_ptr<exvector> vp)
 248   : exprseq(vp), serial(ser)
 249 {
 250 }
 251
 252 //////////
 253 // archiving
 254 //////////
 255
 256 /** Construct object from archive_node. */
 257 void function::read_archive(const archive_node& n, lst& sym_lst)
 258 {
 259         inherited::read_archive(n, sym_lst);
 260         // Find serial number by function name
 261         std::string s;
 262         if (n.find_string("name", s)) {
 263                 unsigned int ser = 0;
 264                 std::vector<function_options>::const_iterator i = registered_functions().begin(), iend = registered_functions().end();
 265                 while (i != iend) {
 266                         if (s == i->name) {
 267                                 serial = ser;
 268                                 return;
 269                         }
 270                         ++i; ++ser;
 271                 }
 272                 throw (std::runtime_error("unknown function '" + s + "' in archive"));
 273         } else
 274                 throw (std::runtime_error("unnamed function in archive"));
 275 }
 276
 277 /** Archive the object. */
 278 void function::archive(archive_node &n) const
 279 {
 280         inherited::archive(n);
 281         GINAC_ASSERT(serial < registered_functions().size());
 282         n.add_string("name", registered_functions()[serial].name);
 283 }
 284
 285 GINAC_BIND_UNARCHIVER(function);
 286
 287 //////////
 288 // functions overriding virtual functions from base classes
 289 //////////
 290
 291 // public
 292
 293 void function::print(const print_context & c, unsigned level) const
 294 {
 295         GINAC_ASSERT(serial<registered_functions().size());
 296         const function_options &opt = registered_functions()[serial];
 297         const std::vector<print_funcp> &pdt = opt.print_dispatch_table;
 298
 299         // Dynamically dispatch on print_context type
 300         const print_context_class_info *pc_info = &c.get_class_info();
 301
 302 next_context:
 303         unsigned id = pc_info->options.get_id();
 304         if (id >= pdt.size() || pdt[id] == NULL) {
 305
 306                 // Method not found, try parent print_context class
 307                 const print_context_class_info *parent_pc_info = pc_info->get_parent();
 308                 if (parent_pc_info) {
 309                         pc_info = parent_pc_info;
 310                         goto next_context;
 311                 }
 312
 313                 // Method still not found, use default output
 314                 if (is_a<print_tree>(c)) {
 315
 316                         c.s << std::string(level, ' ') << class_name() << " "
 317                             << opt.name << " @" << this
 318                             << std::hex << ", hash=0x" << hashvalue << ", flags=0x" << flags << std::dec
 319                             << ", nops=" << nops()
 320                             << std::endl;
 321                         unsigned delta_indent = static_cast<const print_tree &>(c).delta_indent;
 322                         for (size_t i=0; i<seq.size(); ++i)
 323                                 seq[i].print(c, level + delta_indent);
 324                         c.s << std::string(level + delta_indent, ' ') << "=====" << std::endl;
 325
 326                 } else if (is_a<print_csrc>(c)) {
 327
 328                         // Print function name in lowercase
 329                         std::string lname = opt.name;
 330                         size_t num = lname.size();
 331                         for (size_t i=0; i<num; i++)
 332                                 lname[i] = tolower(lname[i]);
 333                         c.s << lname;
 334                         printseq(c, '(', ',', ')', exprseq::precedence(), function::precedence());
 335
 336                 } else if (is_a<print_latex>(c)) {
 337                         c.s << opt.TeX_name;
 338                         printseq(c, '(', ',', ')', exprseq::precedence(), function::precedence());
 339                 } else {
 340                         c.s << opt.name;
 341                         printseq(c, '(', ',', ')', exprseq::precedence(), function::precedence());
 342                 }
 343
 344         } else {
 345
 346                 // Method found, call it
 347                 current_serial = serial;
 348                 if (opt.print_use_exvector_args)
 349                         ((print_funcp_exvector)pdt[id])(seq, c);
 350                 else switch (opt.nparams) {
 351                         // the following lines have been generated for max. @maxargs@ parameters
 352 +++ for N in range(1, maxargs + 1):
 353                         case @N@:
 354                                 ((print_funcp_@N@)(pdt[id]))(@seq('seq[%(n)d]', N, 0)@, c);
 355                                 break;
 356 ---
 357                         // end of generated lines
 358                 default:
 359                         throw(std::logic_error("function::print(): invalid nparams"));
 360                 }
 361         }
 362 }
 363
 364 ex function::eval(int level) const
 365 {
 366         if (level>1) {
 367                 // first evaluate children, then we will end up here again
 368                 return function(serial,evalchildren(level));
 369         }
 370
 371         GINAC_ASSERT(serial<registered_functions().size());
 372         const function_options &opt = registered_functions()[serial];
 373
 374         // Canonicalize argument order according to the symmetry properties
 375         if (seq.size() > 1 && !(opt.symtree.is_zero())) {
 376                 exvector v = seq;
 377                 GINAC_ASSERT(is_a<symmetry>(opt.symtree));
 378                 int sig = canonicalize(v.begin(), ex_to<symmetry>(opt.symtree));
 379                 if (sig != std::numeric_limits<int>::max()) {
 380                         // Something has changed while sorting arguments, more evaluations later
 381                         if (sig == 0)
 382                                 return _ex0;
 383                         return ex(sig) * thiscontainer(v);
 384                 }
 385         }
 386
 387         if (opt.eval_f==0) {
 388                 return this->hold();
 389         }
 390
 391         bool use_remember = opt.use_remember;
 392         ex eval_result;
 393         if (use_remember && lookup_remember_table(eval_result)) {
 394                 return eval_result;
 395         }
 396         current_serial = serial;
 397         if (opt.eval_use_exvector_args)
 398                 eval_result = ((eval_funcp_exvector)(opt.eval_f))(seq);
 399         else
 400         switch (opt.nparams) {
 401                 // the following lines have been generated for max. @maxargs@ parameters
 402 +++ for N in range(1, maxargs + 1):
 403                 case @N@:
 404                         eval_result = ((eval_funcp_@N@)(opt.eval_f))(@seq('seq[%(n)d]', N, 0)@);
 405                         break;
 406 ---
 407                 // end of generated lines
 408         default:
 409                 throw(std::logic_error("function::eval(): invalid nparams"));
 410         }
 411         if (use_remember) {
 412                 store_remember_table(eval_result);
 413         }
 414         return eval_result;
 415 }
 416
 417 ex function::evalf(int level) const
 418 {
 419         GINAC_ASSERT(serial<registered_functions().size());
 420         const function_options &opt = registered_functions()[serial];
 421
 422         // Evaluate children first
 423         exvector eseq;
 424         if (level == 1 || !(opt.evalf_params_first))
 425                 eseq = seq;
 426         else if (level == -max_recursion_level)
 427                 throw(std::runtime_error("max recursion level reached"));
 428         else {
 429                 eseq.reserve(seq.size());
 430                 --level;
 431                 exvector::const_iterator it = seq.begin(), itend = seq.end();
 432                 while (it != itend) {
 433                         eseq.push_back(it->evalf(level));
 434                         ++it;
 435                 }
 436         }
 437
 438         if (opt.evalf_f==0) {
 439                 return function(serial,eseq).hold();
 440         }
 441         current_serial = serial;
 442         if (opt.evalf_use_exvector_args)
 443                 return ((evalf_funcp_exvector)(opt.evalf_f))(seq);
 444         switch (opt.nparams) {
 445                 // the following lines have been generated for max. @maxargs@ parameters
 446 +++ for N in range(1, maxargs + 1):
 447                 case @N@:
 448                         return ((evalf_funcp_@N@)(opt.evalf_f))(@seq('eseq[%(n)d]', N, 0)@);
 449 ---
 450                 // end of generated lines
 451         }
 452         throw(std::logic_error("function::evalf(): invalid nparams"));
 453 }
 454
 455 /**
 456  *  This method is defined to be in line with behaviour of function::return_type()
 457  */
 458 ex function::eval_ncmul(const exvector & v) const
 459 {
 460         // If this function is called then the list of arguments is non-empty
 461         // and the first argument is non-commutative, see  function::return_type()
 462         return seq.begin()->eval_ncmul(v);
 463 }
 464
 465 unsigned function::calchash() const
 466 {
 467         unsigned v = golden_ratio_hash(make_hash_seed(typeid(*this)) ^ serial);
 468         for (size_t i=0; i<nops(); i++) {
 469                 v = rotate_left(v);
 470                 v ^= this->op(i).gethash();
 471         }
 472
 473         if (flags & status_flags::evaluated) {
 474                 setflag(status_flags::hash_calculated);
 475                 hashvalue = v;
 476         }
 477         return v;
 478 }
 479
 480 ex function::thiscontainer(const exvector & v) const
 481 {
 482         return function(serial, v);
 483 }
 484
 485 ex function::thiscontainer(std::auto_ptr<exvector> vp) const
 486 {
 487         return function(serial, vp);
 488 }
 489
 490 /** Implementation of ex::series for functions.
 491  *  \@see ex::series */
 492 ex function::series(const relational & r, int order, unsigned options) const
 493 {
 494         GINAC_ASSERT(serial<registered_functions().size());
 495         const function_options &opt = registered_functions()[serial];
 496
 497         if (opt.series_f==0) {
 498                 return basic::series(r, order);
 499         }
 500         ex res;
 501         current_serial = serial;
 502         if (opt.series_use_exvector_args) {
 503                 try {
 504                         res = ((series_funcp_exvector)(opt.series_f))(seq, r, order, options);
 505                 } catch (do_taylor) {
 506                         res = basic::series(r, order, options);
 507                 }
 508                 return res;
 509         }
 510         switch (opt.nparams) {
 511                 // the following lines have been generated for max. @maxargs@ parameters
 512 +++ for N in range(1, maxargs + 1):
 513                 case @N@:
 514                         try {
 515                                 res = ((series_funcp_@N@)(opt.series_f))(@seq('seq[%(n)d]', N, 0)@, r, order, options);
 516                         } catch (do_taylor) {
 517                                 res = basic::series(r, order, options);
 518                         }
 519                         return res;
 520 ---
 521                 // end of generated lines
 522         }
 523         throw(std::logic_error("function::series(): invalid nparams"));
 524 }
 525
 526 /** Implementation of ex::conjugate for functions. */
 527 ex function::conjugate() const
 528 {
 529         GINAC_ASSERT(serial<registered_functions().size());
 530         const function_options & opt = registered_functions()[serial];
 531
 532         if (opt.conjugate_f==0) {
 533                 return conjugate_function(*this).hold();
 534         }
 535
 536         if (opt.conjugate_use_exvector_args) {
 537                 return ((conjugate_funcp_exvector)(opt.conjugate_f))(seq);
 538         }
 539
 540         switch (opt.nparams) {
 541                 // the following lines have been generated for max. @maxargs@ parameters
 542 +++ for N in range(1, maxargs + 1):
 543                 case @N@:
 544                         return ((conjugate_funcp_@N@)(opt.conjugate_f))(@seq('seq[%(n)d]', N, 0)@);
 545 ---
 546                 // end of generated lines
 547         }
 548         throw(std::logic_error("function::conjugate(): invalid nparams"));
 549 }
 550
 551 /** Implementation of ex::real_part for functions. */
 552 ex function::real_part() const
 553 {
 554         GINAC_ASSERT(serial<registered_functions().size());
 555         const function_options & opt = registered_functions()[serial];
 556
 557         if (opt.real_part_f==0)
 558                 return basic::real_part();
 559
 560         if (opt.real_part_use_exvector_args)
 561                 return ((real_part_funcp_exvector)(opt.real_part_f))(seq);
 562
 563         switch (opt.nparams) {
 564                 // the following lines have been generated for max. @maxargs@ parameters
 565 +++ for N in range(1, maxargs + 1):
 566                 case @N@:
 567                         return ((real_part_funcp_@N@)(opt.real_part_f))(@seq('seq[%(n)d]', N, 0)@);
 568 ---
 569                 // end of generated lines
 570         }
 571         throw(std::logic_error("function::real_part(): invalid nparams"));
 572 }
 573
 574 /** Implementation of ex::imag_part for functions. */
 575 ex function::imag_part() const
 576 {
 577         GINAC_ASSERT(serial<registered_functions().size());
 578         const function_options & opt = registered_functions()[serial];
 579
 580         if (opt.imag_part_f==0)
 581                 return basic::imag_part();
 582
 583         if (opt.imag_part_use_exvector_args)
 584                 return ((imag_part_funcp_exvector)(opt.imag_part_f))(seq);
 585
 586         switch (opt.nparams) {
 587                 // the following lines have been generated for max. @maxargs@ parameters
 588 +++ for N in range(1, maxargs + 1):
 589                 case @N@:
 590                         return ((imag_part_funcp_@N@)(opt.imag_part_f))(@seq('seq[%(n)d]', N, 0)@);
 591 ---
 592                 // end of generated lines
 593         }
 594         throw(std::logic_error("function::imag_part(): invalid nparams"));
 595 }
 596
 597 // protected
 598
 599 /** Implementation of ex::diff() for functions. It applies the chain rule,
 600  *  except for the Order term function.
 601  *  \@see ex::diff */
 602 ex function::derivative(const symbol & s) const
 603 {
 604         ex result;
 605
 606         if (serial == Order_SERIAL::serial) {
 607                 // Order Term function only differentiates the argument
 608                 return Order(seq[0].diff(s));
 609         } else {
 610                 // Chain rule
 611                 ex arg_diff;
 612                 size_t num = seq.size();
 613                 for (size_t i=0; i<num; i++) {
 614                         arg_diff = seq[i].diff(s);
 615                         // We apply the chain rule only when it makes sense.  This is not
 616                         // just for performance reasons but also to allow functions to
 617                         // throw when differentiated with respect to one of its arguments
 618                         // without running into trouble with our automatic full
 619                         // differentiation:
 620                         if (!arg_diff.is_zero())
 621                                 result += pderivative(i)*arg_diff;
 622                 }
 623         }
 624         return result;
 625 }
 626
 627 int function::compare_same_type(const basic & other) const
 628 {
 629         GINAC_ASSERT(is_a<function>(other));
 630         const function & o = static_cast<const function &>(other);
 631
 632         if (serial != o.serial)
 633                 return serial < o.serial ? -1 : 1;
 634         else
 635                 return exprseq::compare_same_type(o);
 636 }
 637
 638 bool function::is_equal_same_type(const basic & other) const
 639 {
 640         GINAC_ASSERT(is_a<function>(other));
 641         const function & o = static_cast<const function &>(other);
 642
 643         if (serial != o.serial)
 644                 return false;
 645         else
 646                 return exprseq::is_equal_same_type(o);
 647 }
 648
 649 bool function::match_same_type(const basic & other) const
 650 {
 651         GINAC_ASSERT(is_a<function>(other));
 652         const function & o = static_cast<const function &>(other);
 653
 654         return serial == o.serial;
 655 }
 656
 657 unsigned function::return_type() const
 658 {
 659         GINAC_ASSERT(serial<registered_functions().size());
 660         const function_options &opt = registered_functions()[serial];
 661
 662         if (opt.use_return_type) {
 663                 // Return type was explicitly specified
 664                 return opt.return_type;
 665         } else {
 666                 // Default behavior is to use the return type of the first
 667                 // argument. Thus, exp() of a matrix behaves like a matrix, etc.
 668                 if (seq.empty())
 669                         return return_types::commutative;
 670                 else
 671                         return seq.begin()->return_type();
 672         }
 673 }
 674
 675 return_type_t function::return_type_tinfo() const
 676 {
 677         GINAC_ASSERT(serial<registered_functions().size());
 678         const function_options &opt = registered_functions()[serial];
 679
 680         if (opt.use_return_type) {
 681                 // Return type was explicitly specified
 682                 return opt.return_type_tinfo;
 683         } else {
 684                 // Default behavior is to use the return type of the first
 685                 // argument. Thus, exp() of a matrix behaves like a matrix, etc.
 686                 if (seq.empty())
 687                         return make_return_type_t<function>();
 688                 else
 689                         return seq.begin()->return_type_tinfo();
 690         }
 691 }
 692
 693 //////////
 694 // new virtual functions which can be overridden by derived classes
 695 //////////
 696
 697 // none
 698
 699 //////////
 700 // non-virtual functions in this class
 701 //////////
 702
 703 // protected
 704
 705 ex function::pderivative(unsigned diff_param) const // partial differentiation
 706 {
 707         GINAC_ASSERT(serial<registered_functions().size());
 708         const function_options &opt = registered_functions()[serial];
 709
 710         // No derivative defined? Then return abstract derivative object
 711         if (opt.derivative_f == NULL)
 712                 return fderivative(serial, diff_param, seq);
 713
 714         current_serial = serial;
 715         if (opt.derivative_use_exvector_args)
 716                 return ((derivative_funcp_exvector)(opt.derivative_f))(seq, diff_param);
 717         switch (opt.nparams) {
 718                 // the following lines have been generated for max. @maxargs@ parameters
 719 +++ for N in range(1, maxargs + 1):
 720                 case @N@:
 721                         return ((derivative_funcp_@N@)(opt.derivative_f))(@seq('seq[%(n)d]', N, 0)@, diff_param);
 722 ---
 723                 // end of generated lines
 724         }
 725         throw(std::logic_error("function::pderivative(): no diff function defined"));
 726 }
 727
 728 ex function::power(const ex & power_param) const // power of function
 729 {
 730         GINAC_ASSERT(serial<registered_functions().size());
 731         const function_options &opt = registered_functions()[serial];
 732
 733         if (opt.power_f == NULL)
 734                 return (new GiNaC::power(*this, power_param))->setflag(status_flags::dynallocated |
 735                                                                        status_flags::evaluated);
 736
 737         current_serial = serial;
 738         if (opt.power_use_exvector_args)
 739                 return ((power_funcp_exvector)(opt.power_f))(seq,  power_param);
 740         switch (opt.nparams) {
 741                 // the following lines have been generated for max. @maxargs@ parameters
 742 +++ for N in range(1, maxargs + 1):
 743                 case @N@:
 744                         return ((power_funcp_@N@)(opt.power_f))(@seq('seq[%(n)d]', N, 0)@, power_param);
 745 ---
 746                 // end of generated lines
 747         }
 748         throw(std::logic_error("function::power(): no power function defined"));
 749 }
 750
 751 ex function::expand(unsigned options) const
 752 {
 753         GINAC_ASSERT(serial<registered_functions().size());
 754         const function_options &opt = registered_functions()[serial];
 755
 756         // No expand defined? Then return the same function with expanded arguments (if required)
 757         if (opt.expand_f == NULL) {
 758                 // Only expand arguments when asked to do so
 759                 if (options & expand_options::expand_function_args)
 760                         return inherited::expand(options);
 761                 else
 762                         return (options == 0) ? setflag(status_flags::expanded) : *this;
 763         }
 764
 765         current_serial = serial;
 766         if (opt.expand_use_exvector_args)
 767                 return ((expand_funcp_exvector)(opt.expand_f))(seq,  options);
 768         switch (opt.nparams) {
 769                 // the following lines have been generated for max. @maxargs@ parameters
 770 +++ for N in range(1, maxargs + 1):
 771                 case @N@:
 772                         return ((expand_funcp_@N@)(opt.expand_f))(@seq('seq[%(n)d]', N, 0)@, options);
 773 ---
 774                 // end of generated lines
 775         }
 776         throw(std::logic_error("function::expand(): no expand of function defined"));
 777 }
 778
 779 std::vector<function_options> & function::registered_functions()
 780 {
 781         static std::vector<function_options> rf = std::vector<function_options>();
 782         return rf;
 783 }
 784
 785 bool function::lookup_remember_table(ex & result) const
 786 {
 787         return remember_table::remember_tables()[this->serial].lookup_entry(*this,result);
 788 }
 789
 790 void function::store_remember_table(ex const & result) const
 791 {
 792         remember_table::remember_tables()[this->serial].add_entry(*this,result);
 793 }
 794
 795 // public
 796
 797 unsigned function::register_new(function_options const & opt)
 798 {
 799         size_t same_name = 0;
 800         for (size_t i=0; i<registered_functions().size(); ++i) {
 801                 if (registered_functions()[i].name==opt.name) {
 802                         ++same_name;
 803                 }
 804         }
 805         if (same_name>=opt.functions_with_same_name) {
 806                 // we do not throw an exception here because this code is
 807                 // usually executed before main(), so the exception could not
 808                 // caught anyhow
 809                 std::cerr << "WARNING: function name " << opt.name
 810                           << " already in use!" << std::endl;
 811         }
 812         registered_functions().push_back(opt);
 813         if (opt.use_remember) {
 814                 remember_table::remember_tables().
 815                         push_back(remember_table(opt.remember_size,
 816                                                  opt.remember_assoc_size,
 817                                                  opt.remember_strategy));
 818         } else {
 819                 remember_table::remember_tables().push_back(remember_table());
 820         }
 821         return registered_functions().size()-1;
 822 }
 823
 824 /** Find serial number of function by name and number of parameters.
 825  *  Throws exception if function was not found. */
 826 unsigned function::find_function(const std::string &name, unsigned nparams)
 827 {
 828         std::vector<function_options>::const_iterator i = function::registered_functions().begin(), end = function::registered_functions().end();
 829         unsigned serial = 0;
 830         while (i != end) {
 831                 if (i->get_name() == name && i->get_nparams() == nparams)
 832                         return serial;
 833                 ++i;
 834                 ++serial;
 835         }
 836         throw (std::runtime_error("no function '" + name + "' with " + ToString(nparams) + " parameters defined"));
 837 }
 838
 839 /** Return the print name of the function. */
 840 std::string function::get_name() const
 841 {
 842         GINAC_ASSERT(serial<registered_functions().size());
 843         return registered_functions()[serial].name;
 844 }
 845
 846 } // namespace GiNaC
 847