Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | File List | Namespace Members | Class Members | Related Pages | Search for

# FunDef.cpp

    00001 /*
    00002  *             Copyright (C) 2005 Maarten Keijzer
    00003  *
    00004  *          This program is free software; you can redistribute it and/or modify
    00005  *          it under the terms of version 2 of the GNU General Public License as
    00006  *          published by the Free Software Foundation.
    00007  *
    00008  *          This program is distributed in the hope that it will be useful,
    00009  *          but WITHOUT ANY WARRANTY; without even the implied warranty of
    00010  *          MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    00011  *          GNU General Public License for more details.
    00012  *
    00013  *          You should have received a copy of the GNU General Public License
    00014  *          along with this program; if not, write to the Free Software
    00015  *          Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
    00016  */
    00017
    00018
    00019 #include <sstream>
    00020 #include "Sym.h"
    00021 #include "FunDef.h"
    00022 #include <LanguageTable.h>
    00023
    00024 using namespace std;
    00025 using namespace boost::numeric;
    00026
    00027 vector<const FunDef*> language;
    00028
    00029 token_t add_function(FunDef* function) {
    00030     language.push_back(function);
    00031     return token_t(language.size()-1);
    00032 }
    00033
    00034 const FunDef& get_element(token_t token) { return *language[token]; }
    00035
    00036 /* Printing */
    00037
    00038 string c_print(const Sym& sym) {
    00039     return c_print(sym, vector<string>());
    00040 }
    00041
    00042 string c_print(const Sym& sym, const vector<string>& vars) {
    00043     const SymVec& args = sym.args();
    00044     vector<string> names(args.size());
    00045     for (unsigned i = 0; i < args.size(); ++i) {
    00046         names[i] = c_print(args[i], vars);
    00047     }
    00048     return language[sym.token()]->c_print(names, vars);
    00049 }
    00050
    00051 /* Evaluation */
    00052
    00053
    00054 double eval(const Sym& sym, const std::vector<double>& inputs) {
    00055     return language[sym.token()]->eval(sym.args(), inputs);
    00056 }
    00057
    00058
    00059 /* Interval Logic */
    00060 Interval eval(const Sym& sym, const vector<Interval>& inputs) {
    00061     const SymVec& args = sym.args();
    00062     vector<Interval> interv(args.size());
    00063     for (unsigned i = 0; i < args.size(); ++i) {
    00064         interv[i] = eval(args[i], inputs);
    00065
    00066         if (!valid(interv[i])) throw interval_error();
    00067     }
    00068     return language[sym.token()]->eval(interv, inputs);
    00069 }
    00070
    00071 /*  */
    00072 void add_function_to_table(LanguageTable& table, token_t token) {
    00073     const FunDef& fundef = *language[token];
    00074
    00075     if (fundef.has_varargs() == false) {
    00076         table.add_function(token, fundef.min_arity());
    00077     } else { // sum or prod (or min or max)
    00078         table.add_function(token, 2);
    00079     }
    00080 }
    00081
    00082
    00083 // by default it is eager
    00084 double FunDef::eval(const SymVec& args, const vector<double>& inputs) const {
    00085     vector<double> values(args.size());
    00086     for (unsigned i = 0; i < args.size(); ++i) {
    00087         values[i] = ::eval(args[i], inputs);
    00088     }
    00089
    00090     return eval(values, inputs);
    00091 }
    00092
    00093 /* Variable Handling */
    00094 FunDef* make_var(int idx); // defined in FunDefs.h
    00095 static vector<token_t> var_token;
    00096
    00097 Sym SymVar(unsigned idx) {
    00098     if (var_token.size() <= idx) {
    00099         // it is new
    00100         var_token.resize(idx+1, token_t(-1));
    00101         var_token[idx] = add_function( make_var(idx) );
    00102     } else if (var_token[idx] == token_t(-1)) {
    00103         var_token[idx] = add_function( make_var(idx) );
    00104     }
    00105     return Sym(var_token[idx]);
    00106 }
    00107
    00108
    00109 /* Constant Handling */
    00110
    00111 struct HashDouble{
    00112     size_t operator()(double val) const {
    00113         unsigned long h = 0;
    00114         char* s = (char*)&val;
    00115         for (unsigned i=0 ; i<sizeof(double); ++i)
    00116             h = 5*h + s[i];
    00117         return size_t(h);
    00118     }
    00119 };
    00120
    00121 #if USE_TR1
    00122 typedef std::tr1::unordered_map<double, token_t> DoubleSet;
    00123 typedef std::tr1::unordered_map<Sym, token_t>    LambdaSet;
    00124 #else
    00125 typedef hash_map<double, token_t, HashDouble> DoubleSet;
    00126 typedef hash_map<Sym, token_t, HashSym> LambdaSet;
    00127 #endif
    00128
    00129 static DoubleSet doubleSet; // for quick checking if a constant already exists
    00130 static vector<double> token_value;
    00131
    00132 static LambdaSet lambdaSet;
    00133 static vector<Sym> token_lambda;
    00134
    00135 static std::vector<token_t> free_list;
    00136
    00137 void delete_val(token_t token) { // clean up the information about this value
    00138
    00139     if (is_constant(token)) {
    00140         //cout << "Deleting constant token " << token << endl;
    00141         double value = token_value[token];
    00142         doubleSet.erase(value);
    00143
    00144         delete language[token];
    00145         language[token] = 0;
    00146         free_list.push_back(token);
    00147     }
    00148     else if (is_lambda(token)) {
    00149         //cout << "Deleting lambda token " << token << endl;
    00150
    00151         Sym expression = token_lambda[token];
    00152         lambdaSet.erase(expression);
    00153
    00154         delete language[token];
    00155         language[token] = 0;
    00156         free_list.push_back(token);
    00157     }
    00158 }
    00159
    00160
    00161 FunDef* make_const(double value);
    00162
    00163 void extend_free_list();
    00164
    00165 Sym SymConst(double value) {
    00166
    00167     DoubleSet::iterator it = doubleSet.find(value);
    00168
    00169     if (it != doubleSet.end()) {
    00170         return Sym(it->second); // already exists
    00171     }
    00172
    00173
    00174     if (free_list.empty()) { // make space for tokens;
    00175         extend_free_list();
    00176     }
    00177
    00178     token_t token = free_list.back();
    00179     free_list.pop_back();
    00180     //cout << "Creating constant with token " << token << endl;
    00181     assert(language[token] == 0);
    00182
    00183     language[token] = make_const(value);
    00184
    00185     doubleSet[value] = token;
    00186     if (token_value.size() < token) token_value.resize(token+1);
    00187     token_value[token] = value;
    00188
    00189     return Sym(token);
    00190 }
    00191
    00192 /* LanguageTable depends on this one, XXX move somewhere safe.*/
    00193 #include <utils/eoRNG.h>
    00194 extern Sym default_const() { return SymConst(rng.normal()); }
    00195
    00196 /* The functions */
    00197 namespace {
    00198
    00199 class Var : public FunDef {
    00200     public :
    00201         unsigned idx;
    00202         string default_str;
    00203
    00204         Var(unsigned _idx) : idx(_idx) {
    00205             ostringstream os;
    00206             os << "x[" << idx << ']'; // CompiledCode expects this form
    00207             default_str = os.str();
    00208         }
    00209
    00210         double eval(const vector<double>& _, const vector<double>& inputs) const { return inputs[idx]; }
    00211         double eval(const SymVec& _, const vector<double>& inputs) const { return inputs[idx]; }
    00212         string c_print(const vector<string>& _, const vector<string>& names) const {
    00213             if (names.empty()) {
    00214                 return default_str;
    00215             }
    00216             return names[idx];
    00217         }
    00218
    00219         Interval eval(const vector<Interval>& _, const vector<Interval>& inputs) const {
    00220             return inputs[idx];
    00221         }
    00222
    00223         unsigned min_arity() const { return 0; }
    00224
    00225         string name() const { return "var"; }
    00226
    00227 };
    00228
    00229 class Const : public FunDef {
    00230     public:
    00231         double value;
    00232         string value_str;
    00233
    00234         Const(double _value) : value(_value) {
    00235             ostringstream os;
    00236             os.precision(17);
    00237             os.setf(ios::showpoint);
    00238             os << '(' << value << ')';
    00239             value_str = os.str();
    00240         }
    00241
    00242
    00243         double eval(const vector<double>& _, const vector<double>& inputs) const { return value; }
    00244         double eval(const SymVec& _, const vector<double>& inputs) const { return value; }
    00245         string c_print(const vector<string>& _, const vector<string>& names) const {
    00246             return value_str;
    00247         }
    00248
    00249         Interval eval(const vector<Interval>& _, const vector<Interval>& inputs) const {
    00250             // Profil/Bias seems to have a problem with 0 * inf when the Interval is exact zero (fpe)
    00251             //if (value == 0.0) return Interval(-BiasEpsilon,BiasEpsilon);
    00252             return Interval(value);
    00253         }
    00254
    00255         unsigned min_arity() const { return 0; }
    00256
    00257         string name() const { return "parameter"; }
    00258 };
    00259
    00260 } // namespace
    00261
    00262 void get_constants(Sym sym, vector<double>& ret) {
    00263     token_t token = sym.token();
    00264     if (is_constant(token)) {
    00265         double val = static_cast<const Const*>(language[token])->value;
    00266         ret.push_back(val);
    00267     }
    00268
    00269     const SymVec& args = sym.args();
    00270     for (unsigned i = 0; i < args.size(); ++i) {
    00271          get_constants(args[i], ret);
    00272     }
    00273
    00274 }
    00275
    00276 double get_constant_value(token_t token) {
    00277     return static_cast<const Const*>(language[token])->value;
    00278 }
    00279
    00281 vector<double> get_constants(Sym sym) {
    00282     vector<double> retval;
    00283     get_constants(sym, retval);
    00284     return retval;
    00285 }
    00286
    00289 Sym set_constants(Sym sym, vector<double>::const_iterator& it) {
    00290
    00291     token_t token = sym.token();
    00292     if (is_constant(token)) {
    00293         return SymConst(*it++);
    00294     }
    00295
    00296     SymVec args = sym.args();
    00297     for (unsigned i = 0; i < args.size(); ++i) {
    00298          args[i] = set_constants(args[i], it);
    00299     }
    00300
    00301     return Sym(token, args);
    00302 }
    00303
    00304 Sym set_constants(Sym sym, const vector<double>& constants) {
    00305     vector<double>::const_iterator it = constants.begin();
    00306     return set_constants(sym, it);
    00307 }
    00308
    00309 // Get functions out, excluding Const and Var
    00310 vector<const FunDef*> get_defined_functions() {
    00311     vector<const FunDef*> res;
    00312     for (unsigned i = 0; i < language.size(); ++i) {
    00313         res.push_back(language[i]);
    00314
    00315         if (is_constant(i) || is_variable(i)) {
    00316             res.back() = 0; // erase
    00317         }
    00318     }
    00319
    00320     return res;
    00321 }
    00322
    00323 FunDef* make_var(int idx) { return new Var(idx); }
    00324 FunDef* make_const(double value) { return new Const(value); }
    00325
    00326 bool is_constant(token_t token) {
    00327     const Const* cnst = dynamic_cast<const Const*>( language[token] );
    00328     return cnst != 0;
    00329 }
    00330
    00331 bool is_variable(token_t token) {
    00332     const Var* var = dynamic_cast<const Var*>( language[token] );
    00333     return var != 0;
    00334 }
    00335
    00336 unsigned get_variable_index(token_t token) {
    00337     const Var* var = static_cast<const Var*>( language[token] );
    00338     return var->idx;
    00339 }
    00340
    00341 namespace {
    00342 class Lambda : public FunDef {
    00343     public:
    00344         Sym expression;
    00345         int arity;
    00346
    00347         Lambda(Sym expr, int arity_) : expression(expr), arity(arity_) {}
    00348
    00349         double eval(const vector<double>& vals, const vector<double>& _) const {
    00350             return ::eval(expression, vals);
    00351         }
    00352
    00353         string c_print(const vector<string>& args, const vector<string>& _) const {
    00354             return string("/*f*/") + ::c_print(expression, args) + string("/*eof*/");
    00355         }
    00356
    00357         Interval eval(const vector<Interval>& args, const vector<Interval>& _) const {
    00358             return ::eval(expression, args);
    00359         }
    00360
    00361         unsigned min_arity() const { return arity; }
    00362         string name() const { return "F"; }
    00363
    00364 };
    00365     Sym normalize(Sym sym, SymVec& args) {
    00366         // check if it's a variable
    00367         token_t token = sym.token();
    00368         const Var* var = dynamic_cast< const Var*>(language[token]);
    00369
    00370         if (var != 0) {
    00371             for (unsigned i = 0; i < args.size(); ++i) {
    00372                 if (sym == args[i]) {
    00373                     return SymVar(i); // replace with reference to arg
    00374                 }
    00375             }
    00376
    00377             // not replaced, add it
    00378             args.push_back(sym);
    00379             return SymVar(args.size()-1);
    00380
    00381         }
    00382
    00383         SymVec a = sym.args();
    00384         for (unsigned i = 0; i < a.size(); ++i) {
    00385             a[i] = normalize(a[i], args);
    00386         }
    00387
    00388         return Sym(token, a);
    00389     }
    00390 }
    00391
    00392 bool is_lambda(token_t token) {
    00393     const Lambda* lambda = dynamic_cast<const Lambda*>( language[token]);
    00394     return lambda != 0;
    00395 }
    00396
    00397 std::ostream& print_list(Sym sym, ostream& os) {
    00398     os << sym.token() << ' ';
    00399
    00400     const SymVec& args = sym.args();
    00401     for (unsigned i = 0; i < args.size(); ++i) {
    00402         print_list(args[i], os);
    00403     }
    00404     return os;
    00405 }
    00406
    00407 token_t new_lambda(Sym sym, int arity) {
    00408     // check if already present
    00409
    00410     LambdaSet::iterator it = lambdaSet.find(sym);
    00411     if (it != lambdaSet.end()) {
    00412         return it->second;
    00413     }
    00414
    00415
    00416     // new, insert
    00417     Lambda* lambda = new Lambda(sym, arity);
    00418
    00419     if (free_list.empty()) {
    00420         extend_free_list();
    00421     }
    00422
    00423     token_t token = free_list.back();
    00424     free_list.pop_back();
    00425     language[token] = lambda;
    00426
    00427     lambdaSet[sym] = token;
    00428     if (token_lambda.size() <= token) token_lambda.resize(token+1);
    00429     token_lambda[token] = sym;
    00430
    00431     return token;
    00432 }
    00433
    00434 /* Compression */
    00435 typedef hash_map<Sym, unsigned, HashSym> OccMap;
    00436
    00437 void count_occurances(Sym sym, OccMap& occ) {
    00438     occ[sym]++;
    00439     const SymVec& args = sym.args();
    00440     for (unsigned i = 0; i < args.size(); ++i) {
    00441         count_occurances(args[i], occ);
    00442     }
    00443 }
    00444
    00445 Sym create_lambda(Sym sym, OccMap& occ, unsigned nvars, vector<Sym>& args) {
    00446     unsigned o = occ[sym];
    00447     unsigned sz = sym.size();
    00448
    00449     if (o * sz > o + sz + nvars || is_variable(sym.token()) ) {
    00450         // check if it's already present
    00451         for (unsigned i = 0; i < args.size(); ++i) {
    00452             if (args[i] == sym) {
    00453                 return SymVar(i);
    00454             }
    00455         }
    00456         // push_back
    00457         args.push_back(sym);
    00458         return SymVar(args.size()-1);
    00459     }
    00460
    00461     SymVec sym_args = sym.args();
    00462     for (unsigned i = 0; i < sym_args.size(); ++i) {
    00463         sym_args[i] = create_lambda(sym_args[i], occ, nvars, args);
    00464     }
    00465
    00466     return Sym(sym.token(), sym_args);
    00467
    00468 }
    00469
    00470 Sym compress(Sym sym) {
    00471     OccMap occ(sym.size());
    00472     count_occurances(sym, occ);
    00473
    00474     unsigned nvars = 0;
    00475     for (OccMap::iterator it = occ.begin(); it != occ.end(); ++it) {
    00476         if (is_variable(it->first.token())) nvars++;
    00477     }
    00478
    00479     SymVec args;
    00480     Sym body = create_lambda(sym, occ, nvars, args);
    00481
    00482
    00483     if (body.size() < sym.size()) {
    00484         // see if the body can be compressed some more
    00485         body = compress(body);
    00486
    00487         token_t token = new_lambda(body, args.size());
    00488
    00489         for (unsigned i = 0; i < args.size(); ++i) {
    00490             args[i] = compress(args[i]);
    00491         }
    00492
    00493         Sym result = Sym(token, args);
    00494         return compress(result); // see if it can be compressed some more
    00495     }
    00496
    00497     return sym;
    00498 }
    00499
    00500 Sym SymLambda(Sym expr) { return compress(expr); }
    00501
    00502 Sym expand(Sym expr, const SymVec& args) {
    00503
    00504     const Var* var = dynamic_cast<const Var*>( language[expr.token()] );
    00505     if (var != 0) {
    00506         return args[var->idx];
    00507     }
    00508
    00509     SymVec expr_args = expr.args();
    00510     for (unsigned i = 0; i < expr_args.size(); ++i) {
    00511         expr_args[i] = expand(expr_args[i], args);
    00512     }
    00513
    00514     return Sym(expr.token(), expr_args);
    00515 }
    00516
    00517 Sym SymUnlambda(Sym sym) {
    00518     Sym retval = sym;
    00519     const Lambda* lambda = dynamic_cast<const Lambda*>( language[sym.token()] );
    00520     if (lambda != 0) {
    00521         retval = expand(lambda->expression, sym.args());
    00522     }
    00523
    00524     return retval;
    00525 }
    00526
    00527 Sym expand_all(Sym sym) {
    00528     SymVec args = sym.args();
    00529     for (unsigned i = 0; i < args.size(); ++i) {
    00530         args[i] = expand_all(args[i]);
    00531     }
    00532
    00533     Sym nw = SymUnlambda( Sym(sym.token(), args) );
    00534
    00535     if (nw != sym) {
    00536         nw = expand_all(nw);
    00537     }
    00538
    00539     return nw;
    00540 }
    00541
    00542 namespace {
    00543
    00544 class Sum : public FunDef {
    00545
    00546     public :
    00547
    00548         double eval(const vector<double>& vals, const vector<double>& _) const {
    00549             double res = 0;
    00550             for (unsigned i = 0; i < vals.size(); ++i) res += vals[i];
    00551             return res;
    00552         }
    00553
    00554         string c_print(const vector<string>& args, const vector<string>& _) const {
    00555             if (args.empty()) { return "0.0"; }
    00556
    00557             ostringstream os;
    00558             os << "(" << args[0];
    00559             for (unsigned i = 1; i < args.size(); ++i) {
    00560                 os << "+" << args[i];
    00561             }
    00562             os << ")";
    00563             return os.str();
    00564         }
    00565
    00566         Interval eval(const vector<Interval>& args, const vector<Interval>& inputs) const {
    00567             Interval interv(0.0);
    00568             for (unsigned i = 0; i < args.size(); ++i) {
    00569                 interv += args[i];
    00570             }
    00571             return interv;
    00572         }
    00573
    00574         unsigned min_arity() const { return 0; }
    00575         bool has_varargs() const { return true; }
    00576
    00577         string name() const { return "sum"; }
    00578 };
    00579
    00580
    00581 class Prod : public FunDef {
    00582
    00583     public :
    00584
    00585         double eval(const vector<double>& vals, const vector<double>& _) const {
    00586             double res = 1;
    00587             for (unsigned i = 0; i < vals.size(); ++i) res *= vals[i];
    00588             return res;
    00589         }
    00590
    00591         string c_print(const vector<string>& args, const vector<string>& _) const {
    00592             if (args.empty()) { return "1.0"; }
    00593
    00594             ostringstream os;
    00595             os << "(" << args[0];
    00596             for (unsigned i = 1; i < args.size(); ++i) {
    00597                 os << "*" << args[i];
    00598             }
    00599             os << ")";
    00600
    00601             return os.str();
    00602         }
    00603
    00604         Interval eval(const vector<Interval>& args, const vector<Interval>& inputs) const {
    00605             Interval interv(1.0);
    00606             for (unsigned i = 0; i < args.size(); ++i) {
    00607                 interv *= args[i];
    00608             }
    00609             return interv;
    00610         }
    00611
    00612         unsigned min_arity() const { return 0; }
    00613         bool has_varargs() const { return true; }
    00614
    00615         string name() const { return "prod"; }
    00616 };
    00617
    00618
    00619 class Power : public FunDef {
    00620     public :
    00621         double eval(const vector<double>& vals, const vector<double>& _) const {
    00622             return pow(vals[0], vals[1]);
    00623         }
    00624
    00625         string c_print(const vector<string>& args, const vector<string>& _) const {
    00626             return "pow(" + args[0] + ',' + args[1] + ')';
    00627         }
    00628
    00629         Interval eval(const vector<Interval>& args, const vector<Interval>& _) const {
    00630             Interval first = args[0];
    00631             Interval second = args[1];
    00632             Interval lg = log(first);
    00633             if (!valid(lg)) throw interval_error();
    00634             return exp(second * lg);
    00635         }
    00636
    00637         unsigned min_arity() const { return 2; }
    00638
    00639         string name() const { return "pow"; }
    00640 };
    00641
    00642 class IsNeg : public FunDef {
    00643
    00644     public:
    00645         double eval(const vector<double>& vals, const vector<double>& _) const {
    00646             if (vals[0] < 0.0) return vals[1];
    00647             return vals[2];
    00648         }
    00649
    00650         double eval(const Sym& sym, const vector<double>& inputs) const {
    00651             const SymVec& args = sym.args();
    00652             double arg0 = ::eval(args[0], inputs);
    00653             if (arg0 < 0.0) {
    00654                 return ::eval(args[1], inputs);
    00655             }
    00656             return ::eval(args[2], inputs);
    00657         }
    00658
    00659         string c_print(const vector<string>& args, const vector<string>& _) const {
    00660             return "((" + args[0] + "<0.0)?" + args[1] + ":" + args[2]+")";
    00661         }
    00662
    00663         Interval eval(const vector<Interval>& args, const vector<Interval>& _) const {
    00664             Interval a0 = args[0];
    00665             if (a0.upper() < 0.0) return args[1];
    00666             if (a0.lower() >= 0.0) return args[2];
    00667
    00668             return Interval( std::min(args[1].lower(), args[2].lower()), std::max(args[1].upper(), args[2].upper()));
    00669         }
    00670
    00671         unsigned min_arity() const { return 3; }
    00672
    00673         string name() const { return "ifltz"; }
    00674 };
    00675
    00676 template <typename Func>
    00677 class Unary : public FunDef {
    00678
    00679     Func un;
    00680
    00681     double eval(const vector<double>& vals, const vector<double>& _) const {
    00682         return un(vals[0]);
    00683     }
    00684
    00685     string c_print(const vector<string>& args, const vector<string>& _) const {
    00686         return un(args[0]);
    00687     }
    00688
    00689     Interval eval(const vector<Interval>& args, const vector<Interval>& _) const {
    00690         return un(args[0]);
    00691     }
    00692
    00693     unsigned min_arity() const { return 1; }
    00694
    00695     string name() const { return un.name(); }
    00696
    00697 };
    00698
    00699 struct Inv {
    00700     double operator()(double val) const { return 1.0 / val; }
    00701     string operator()(string v)   const { return "(1./" + v + ")"; }
    00702     Interval operator()(Interval v) const { return 1.0 / v; }
    00703
    00704     string name() const { return "inv"; }
    00705 };
    00706
    00707 struct Min {
    00708     double operator()(double val) const { return -val; }
    00709     string operator()(string v) const   { return "(-" + v + ")"; }
    00710     Interval operator()(Interval v) const { return -v; }
    00711
    00712     string name() const { return "min"; }
    00713 };
    00714
    00715 } // namespace
    00716
    00717 string prototypes = "double pow(double, double);";
    00718 string get_prototypes() { return prototypes; }
    00719 unsigned add_prototype(string str) { prototypes += string("double ") + str + "(double);"; return prototypes.size(); }
    00720
    00721 token_t add_function(FunDef* function, token_t where) {
    00722     if (language.size() <= where) language.resize(where+1);
    00723     language[where] = function;
    00724     return 0;
    00725 }
    00726
    00727 namespace {
    00728
    00729 #define FUNCDEF(funcname) struct funcname##_struct { \
    00730     double operator()(double val) const { return funcname(val); }\
    00731     string operator()(string val) const  { return string(#funcname) + '(' + val + ')'; }\
    00732     Interval operator()(Interval val) const { return funcname(val); }\
    00733     string name() const { return string(#funcname); }\
    00734 };\
    00735 static const token_t funcname##_token_static = add_function( new Unary<funcname##_struct>, funcname##_token);\
    00736 unsigned funcname##_size = add_prototype(#funcname);
    00737
    00738 static token_t ssum_token  = add_function( new Sum , sum_token);
    00739 static token_t sprod_token = add_function( new Prod, prod_token);
    00740 static token_t sinv_token  = add_function( new Unary<Inv>, inv_token);
    00741 static token_t smin_token  = add_function( new Unary<Min>, min_token);
    00742 static token_t spow_token  = add_function( new Power, pow_token);
    00743 static token_t sifltz_token = add_function( new IsNeg, ifltz_token);
    00744
    00745 FUNCDEF(sin);
    00746 FUNCDEF(cos);
    00747 FUNCDEF(tan);
    00748 FUNCDEF(asin);
    00749 FUNCDEF(acos);
    00750 FUNCDEF(atan);
    00751
    00752 FUNCDEF(sinh);
    00753 FUNCDEF(cosh);
    00754 FUNCDEF(tanh);
    00755 FUNCDEF(asinh);
    00756 FUNCDEF(acosh);
    00757 FUNCDEF(atanh);
    00758
    00759 FUNCDEF(exp);
    00760 FUNCDEF(log);
    00761 } // namespace
    00762
    00763 double sqr(double x) { return x*x; }
    00764
    00765 namespace {
    00766 FUNCDEF(sqr);
    00767 FUNCDEF(sqrt);
    00768
    00769 const int buildInFunctionOffset = language.size();
    00770 } // namespace
    00771
    00772 void add_tokens() {
    00773     unsigned sz = language.size();
    00774     language.resize(sz + sz+1); // double
    00775
    00776     for (unsigned i = sz; i < language.size(); ++i) {
    00777        free_list.push_back(i);
    00778     }
    00779 }
    00780
    00781 void extend_free_list() {
    00782     // first check if we can clean up unused tokens;
    00783     const vector<unsigned>& refcount = Sym::token_refcount();
    00784     for (unsigned i = buildInFunctionOffset; i < refcount.size(); ++i) {
    00785         if (language[i] == 0) continue;
    00786
    00787         bool c = is_constant(i);
    00788         bool l = is_lambda(i);
    00789
    00790         if (refcount[i] == 0 && (c || l)) {
    00791
    00792             if (c) {
    00793                 doubleSet.erase(token_value[i]);
    00794             }
    00795
    00796             if (l) {
    00797                 lambdaSet.erase(token_lambda[i]);
    00798             }
    00799
    00800             delete language[i];
    00801             language[i] = 0;
    00802             free_list.push_back(i);
    00803         }
    00804     }
    00805
    00806     // if still empty, add new tokens
    00807     if (free_list.empty()) {
    00808         add_tokens();
    00809     }
    00810 }
    00811
    00812
    00813 /* Serialization */
    00814 void write_raw(ostream& os, const Sym& sym) {
    00815     token_t token = sym.token();
    00816     const SymVec& args = sym.args();
    00817
    00818     if (is_constant(token)) {
    00819         os << "c" << language[token]->c_print(vector<string>(), vector<string>());
    00820     } else {
    00821
    00822         const Var* var = dynamic_cast<const Var*>( language[token] );
    00823
    00824         if (var != 0) {
    00825             os << "v" << var->idx;
    00826         } else {
    00827             os << "f" << token << ' ' << args.size();
    00828         }
    00829     }
    00830
    00831     for (unsigned i = 0; i < args.size(); ++i) {
    00832         write_raw(os, args[i]);
    00833     }
    00834 }
    00835
    00836 string write_raw(const Sym& sym) {
    00837
    00838     ostringstream os;
    00839     write_raw(os, sym);
    00840
    00841     return os.str();
    00842 }
    00843
    00844 Sym read_raw(istream& is) {
    00845     char id = is.get();
    00846
    00847     switch (id) {
    00848         case 'c' :
    00849             {
    00850                 double val;
    00851                 is.get(); // skip '('
    00852                 is >> val;
    00853                 is.get(); // skip ')'
    00854                 return SymConst(val);
    00855             }
    00856         case 'v' :
    00857             {
    00858                 unsigned idx;
    00859                 is >> idx;
    00860                 return SymVar(idx);
    00861             }
    00862         case 'f' :
    00863             {
    00864                 token_t token;
    00865                 unsigned arity;
    00866                 is >> token;
    00867                 is >> arity;
    00868                 SymVec args(arity);
    00869                 for (unsigned i = 0; i < arity; ++i) {
    00870                     args[i] = read_raw(is);
    00871                 }
    00872
    00873                 return Sym(token, args);
    00874             }
    00875         default : {
    00876                       cerr << "Character = " << id << " Could not read formula from stream" << endl;
    00877                       exit(1);
    00878                   }
    00879
    00880     }
    00881
    00882     return Sym();
    00883 }
    00884
    00885 Sym read_raw(string str) {
    00886     istringstream is(str);
    00887     return read_raw(is);
    00888 }
    00889
    00890 void read_raw(istream& is, Sym& sym) {
    00891     sym = read_raw(is);
    00892 }
    00893