--- /dev/null
+
+///////////////////////////////////////////////////////////////////////////
+// This program is free software: you can redistribute it and/or modify //
+// it under the terms of the version 3 of the GNU General Public License //
+// as published by the Free Software Foundation. //
+// //
+// This program is distributed in the hope that it will be useful, but //
+// WITHOUT ANY WARRANTY; without even the implied warranty of //
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU //
+// General Public License for more details. //
+// //
+// You should have received a copy of the GNU General Public License //
+// along with this program. If not, see <http://www.gnu.org/licenses/>. //
+// //
+// Written and (C) by Francois Fleuret //
+// Contact <francois.fleuret@idiap.ch> for comments & bug reports //
+///////////////////////////////////////////////////////////////////////////
+
+#include <math.h>
+#include "mappings.h"
+
+#ifdef DEBUG
+int nb_refs = 0;
+#endif
+
+// I looked a bit at string.h of the g++ library to do this ref
+// counting things in the grab() and release(). This strategy allows
+// to avoid useless copies.
+
+// Each time you build an object which has a reference to a Map, you
+// have to grab() the Map, and each time you destroy an object with a
+// reference to a Map, you have to release() the Map. You should not
+// have grab() or release() out of constructors / destructors and
+// assignement operator.
+
+// With such a politic, you should not need to worry about references
+// anywhere else
+
+class Map {
+ int ref;
+public:
+ // If DEBUG is defined, we keep a count of the total number of
+ // constructed / deleted Maps. "NO REF!" is printed each time this
+ // counter goes to 0.
+
+ Map() { ref = 0;
+#ifdef DEBUG
+ nb_refs++;
+#endif
+}
+ virtual ~Map() {
+#ifdef DEBUG
+ nb_refs--; if(nb_refs == 0) cout << "NO REF!" << endl;
+#endif
+ }
+ virtual float eval(float x) = 0;
+ virtual Map *derivative() = 0;
+ virtual Map *compose(Map *g) = 0;
+ void release() { if(--ref == 0) delete this; }
+ Map *grab() { ref++; return this; }
+ virtual void print(ostream &o) = 0;
+};
+
+class Cst : public Map {
+ float value;
+public:
+ Cst(float v) : value(v) { }
+ virtual float eval(float x) { return value; }
+ virtual Map *derivative() { return new Cst(0.0); }
+ virtual Map *compose(Map *g) { return this; }
+ virtual void print(ostream &o) { o << value; }
+};
+
+class Var : public Map {
+public:
+ Var() {}
+ virtual float eval(float x) { return x; }
+ virtual Map *derivative() { return new Cst(1.0); }
+ virtual Map *compose(Map *g) { return g; }
+ virtual void print(ostream &o) { o << "X"; }
+};
+
+class Sum : public Map {
+ Map *f1, *f2;
+public:
+ Sum(Map *g1, Map *g2) : f1(g1->grab()), f2(g2->grab()) { }
+ ~Sum() { f1->release(); f2->release(); }
+ virtual float eval(float x) { return f1->eval(x) + f2->eval(x); }
+ virtual Map *derivative() { return new Sum(f1->derivative(), f2->derivative()); }
+ virtual Map *compose(Map *g) { return new Sum(f1->compose(g), f2->compose(g)); }
+ virtual void print(ostream &o) { o << "( "; f1->print(o); o<< " ) + ( " ; f2->print(o); o << " )"; }
+};
+
+class Prd : public Map {
+ Map *f1, *f2;
+public:
+ Prd(Map *g1, Map *g2) : f1(g1->grab()), f2(g2->grab()) { }
+ ~Prd() { f1->release(); f2->release(); }
+ virtual float eval(float x) { return f1->eval(x) * f2->eval(x); }
+ virtual Map *derivative() { return new Sum(new Prd(f1->derivative(), f2), new Prd(f1, f2->derivative())); }
+ virtual Map *compose(Map *g) { return new Prd(f1->compose(g), f2->compose(g)); }
+ virtual void print(ostream &o) { o << "( "; f1->print(o); o<< " ) * ( " ; f2->print(o); o << " )"; }
+};
+
+class Sin : public Map {
+ Map *f;
+public:
+ Sin(Map *g) : f(g->grab()) { }
+ ~Sin() { f->release(); }
+ virtual float eval(float x) { return sin(f->eval(x)); }
+ virtual Map *derivative() { return new Prd(f->derivative(), new Sin(new Sum(f, new Cst(M_PI/2)))); }
+ virtual Map *compose(Map *g) { return new Sin(f->compose(g)); }
+ virtual void print(ostream &o) { o << "sin( "; f->print(o); o<< " )"; }
+};
+
+class Inv : public Map {
+ Map *f;
+public:
+ Inv(Map *g) : f(g->grab()) {}
+ ~Inv() { f->release(); }
+ virtual float eval(float x) { return 1.0/f->eval(x); }
+ virtual Map *derivative() { return new Prd(new Prd(new Cst(-1.0), f->derivative()), new Inv(new Prd(f, f))); }
+ virtual Map *compose(Map *g) { return new Inv(f->compose(g)); }
+ virtual void print(ostream &o) { o << "1 / ( "; f->print(o); o<< " )"; }
+};
+
+class Exp : public Map {
+ Map *f;
+public:
+ Exp(Map *g) : f(g->grab()) {}
+ ~Exp() { f->release(); }
+ virtual float eval(float x) { return exp(f->eval(x)); }
+ virtual Map *derivative() { return new Prd(f->derivative(), new Exp(f)); }
+ virtual Map *compose(Map *g) { return new Exp(f->compose(g)); }
+ virtual void print(ostream &o) { o << "exp( "; f->print(o); o<< " )"; }
+};
+
+class Log : public Map {
+ Map *f;
+public:
+ Log(Map *g) : f(g->grab()) {}
+ ~Log() { f->release(); }
+ virtual float eval(float x) { return log(f->eval(x)); }
+ virtual Map *derivative() { return new Prd(f->derivative(), new Prd(f->derivative(), new Inv(f))); }
+ virtual Map *compose(Map *g) { return new Log(f->compose(g)); }
+ virtual void print(ostream &o) { o << "log( "; f->print(o); o<< " )"; }
+};
+
+class Pow : public Map {
+ Map *f;
+ int exponent;
+public:
+ Pow(Map *g, int e) : f(g->grab()), exponent(e) { }
+ ~Pow() { f->release(); }
+ virtual float eval(float x) { return pow(f->eval(x), exponent); }
+ virtual Map *derivative() { return new Prd(new Prd(new Cst(float(exponent)), f->derivative()), new Pow(f, exponent-1)); }
+ virtual Map *compose(Map *g) { return new Pow(f->compose(g), exponent); }
+ virtual void print(ostream &o) { o << "( "; f->print(o); o<< " )^" << exponent; }
+};
+
+const Mapping Mapping::X(new Var());
+
+Mapping::Mapping() : f(0) { }
+Mapping::Mapping(Map *g) : f(g ? g->grab(): 0) { }
+Mapping::Mapping(const Mapping &s) : f(s.f ? s.f->grab() : 0) { }
+Mapping::Mapping(float x) : f((new Cst(x))->grab()) { }
+Mapping::~Mapping() { if(f) f->release(); }
+
+// This is the only place where grab() and release() are not in a
+// constructor / destructor
+Mapping &Mapping::operator = (const Mapping &m) {
+ if(&m != this) {
+ if(f) f->release();
+ if(m.f) f = m.f->grab(); else f = 0;
+ }
+ return *this;
+}
+
+float Mapping::operator () (float x) const { return f->eval(x); }
+Mapping Mapping::derivative() const { return Mapping(f->derivative()); };
+Mapping Mapping::compose(const Mapping &m) const { return Mapping(f->compose(m.f)); }
+
+Mapping Mapping::add(const Mapping &m) const { return Mapping(new Sum(f, m.f)); }
+Mapping Mapping::sub(const Mapping &m) const { return Mapping(new Sum(f, new Prd(new Cst(-1.0), m.f))); }
+Mapping Mapping::mul(const Mapping &m) const { return Mapping(new Prd(f, m.f)); }
+Mapping Mapping::div(const Mapping &m) const { return Mapping(new Prd(f, new Inv(m.f))); }
+Mapping Mapping::pow(int k) const { return Mapping(new Pow(f, k)); }
+Mapping Mapping::neg() const { return Mapping(new Prd(new Cst(-1), f)); }
+
+Mapping Mapping::sin() const { return Mapping(new Sin(f)); }
+Mapping Mapping::cos() const { return Mapping(new Sin(new Sum(f, new Cst(M_PI/2)))); }
+Mapping Mapping::log() const { return Mapping(new Log(f)); }
+Mapping Mapping::exp() const { return Mapping(new Exp(f)); }
+
+void Mapping::print(ostream &o) const { if(f) f->print(o); else o << "UNDEFINED"; }
+
+ostream &operator << (ostream &o, const Mapping &m) { m.print(o); return o; }
+
+Mapping operator + (const Mapping &mr) { return mr; }
+Mapping operator + (const Mapping &ml, const Mapping &mr) { return ml.add(mr); }
+Mapping operator - (const Mapping &mr) { return mr.neg(); }
+Mapping operator - (const Mapping &ml, const Mapping &mr) { return ml.sub(mr); }
+Mapping operator * (const Mapping &ml, const Mapping &mr) { return ml.mul(mr); }
+Mapping operator / (const Mapping &ml, const Mapping &mr) { return ml.div(mr); }
+Mapping operator ^ (const Mapping &ml, int k) { return ml.pow(k); }
+
+Mapping sin(const Mapping &m) { return m.sin(); }
+Mapping cos(const Mapping &m) { return m.cos(); }
+Mapping log(const Mapping &m) { return m.log(); }
+Mapping exp(const Mapping &m) { return m.exp(); }
+
--- /dev/null
+
+///////////////////////////////////////////////////////////////////////////
+// This program is free software: you can redistribute it and/or modify //
+// it under the terms of the version 3 of the GNU General Public License //
+// as published by the Free Software Foundation. //
+// //
+// This program is distributed in the hope that it will be useful, but //
+// WITHOUT ANY WARRANTY; without even the implied warranty of //
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU //
+// General Public License for more details. //
+// //
+// You should have received a copy of the GNU General Public License //
+// along with this program. If not, see <http://www.gnu.org/licenses/>. //
+// //
+// Written and (C) by Francois Fleuret //
+// Contact <francois.fleuret@idiap.ch> for comments & bug reports //
+///////////////////////////////////////////////////////////////////////////
+
+#ifndef MAPPINGS_H
+#define MAPPINGS_H
+
+#include <iostream>
+
+using namespace std;
+
+// We have a bunch of such classes inside
+class Map;
+
+class Mapping {
+ Map *f;
+
+ // We'll need this one from time to time, but nobody is supposed to
+ // use it from "outside"
+ Mapping(Map *g);
+
+public:
+ // This is the variable
+ const static Mapping X;
+
+ // This build a non-defined Mapping. Everything is illegal on such a
+ // Mapping except delete, = and <<. All other operations will *crash*
+ Mapping();
+
+ // Cast from float to allow operations with constant values
+ Mapping(float x);
+
+ // Standard copy constructor and destructor
+ Mapping(const Mapping &s);
+ ~Mapping();
+
+ // Assignment
+ Mapping &operator = (const Mapping &m);
+
+ // Evaluation
+ float operator () (float x) const;
+
+ // Derivation
+ Mapping derivative() const;
+
+ // Composition of Mappings
+ Mapping compose(const Mapping &m) const;
+
+ // Standard operations
+ Mapping pow(int k) const;
+ Mapping mul(const Mapping &m) const;
+ Mapping div(const Mapping &m) const;
+ Mapping add(const Mapping &m) const;
+ Mapping sub(const Mapping &m) const;
+ Mapping neg() const;
+ Mapping sin() const;
+ Mapping cos() const;
+ Mapping log() const;
+ Mapping exp() const;
+
+ // Print the expression
+ void print(ostream &s) const;
+};
+
+// I guess we have all the ones we need?
+Mapping operator + (const Mapping &m);
+Mapping operator + (const Mapping &ml, const Mapping &mr);
+Mapping operator - (const Mapping &m);
+Mapping operator - (const Mapping &ml, const Mapping &mr);
+Mapping operator * (const Mapping &ml, const Mapping &mr);
+Mapping operator / (const Mapping &ml, const Mapping &mr);
+
+// Be very careful with this operator, it has a lower precedence than
+// +-*/
+Mapping operator ^ (const Mapping &ml, int k);
+
+// stream stuff, I dont have the courage to do the >>
+ostream &operator << (ostream &o, const Mapping &m);
+
+Mapping sin(const Mapping &m);
+Mapping cos(const Mapping &m);
+Mapping log(const Mapping &m);
+Mapping exp(const Mapping &m);
+
+#endif