Implementing multiplication operator for mathematical functions C++

Question

I have the following abstract base class:

class Function {
 virtual double Eval(double x) const = 0;
};     

I want to be able to use expressions like f * g or f->operator *(g), where f and g are concrete objects of the class Function, in my main file, say for example when I want to calculate a definite integral so that I can write:

AnIntegrationMethod(f*g);

A rather unsophisticated method I came up with consists of declaring a class Product (only header file shown, the implementation is obvious):

class Product: public Function {
 public: Product(Function *g, Function *f); 
  ~Product(); 
  virtual double Eval(double x) const; //return _g->Eval(x)*_f->Eval(x)
 private: Function *_g; Function *_f;
 };

and then in any of my functions

#include "Product.h"

class ConcreteFunction: public Function {
 public: 
   ConcreteFunction(); 
  ~ConcreteFunction(); 
  virtual double Eval(double x) const;
 Function* operator*(Function *f) {return new Product(this, f);}
 };

This actually works for simple stuff but the problem is that the operator * is only defined within single derived classes of the base class instead of being defined for every possible derived class. This means, for instance, that if I have a concrete object f representing a mathematical function I can call f->operator *g but if I want to call again the operator * to get the object (f->operator * g)->operator * f I am not going to be able to because the function f* g does not have the * operator defined as f.

I suppose I should define the operator * directly in my base class but then I can't figure out how to implement the operator because I don't really know how to get the proper Eval function for the product since I cannot use the class Product now, it wouldn't make sense to use the class Product derived from the class Function in the class Function itself. I think I'm also quite confused over whether in general is correct to write something like the following:

 Function* Function::operator*(Function *f) {
 Function *g;
 ...
 //operations that allow g to be have the right Eval function
 //which should of the sort this->Eval(x) * f->Eval(x)
 ...
 return g;
 }

Any hints or suggestions on how to proceed are appreciated. My level is very basic, I've been programming two month now.

Answer 1

Just a sketch, you might do something like this:

#include <memory>

// Base Function: f(x) = x
class Function
{
    protected:
    struct Implementation
    {
        virtual ~Implementation() {}
        virtual double evaluate(double x) const { return x; }
    };

    public:
    Function()
    :   self(std::make_shared<Implementation>())
    {}

    double operator () (double x) const { return self->evaluate(x); }

    protected:
    Function(std::shared_ptr<Implementation> self)
    :   self(self)
    {}

    private:
    std::shared_ptr<Implementation> self;
};
typedef Function Identity;


// Unary Function: u(-f(x))
class UnaryMinus : public Function
{
    protected:
    struct Implementation : Function::Implementation
    {
        Function f;
        Implementation(Function f)
        :   f(f)
        {};

        virtual double evaluate(double x) const override { return -f(x); }
    };

    public:
    UnaryMinus(Function f)
    :   Function(std::make_shared<Implementation>(f))
    {}
};

// Binary Function: u(f(x) + g(x))
class BinaryAdd : public Function
{
    protected:
    struct Implementation : Function::Implementation
    {
        Function f;
        Function g;
        Implementation(Function f, Function g)
        :   f(f), g(g)
        {};

        virtual double evaluate(double x) const override { return f(x) + g(x); }
    };

    public:
    BinaryAdd(Function f, Function g)
    :   Function(std::make_shared<Implementation>(f, g))
    {}
};

// Binary Function: u(f(x) * g(x))
class BinaryMultiply : public Function
{
    protected:
    struct Implementation : Function::Implementation
    {
        Function f;
        Function g;
        Implementation(Function f, Function g)
        :   f(f), g(g)
        {};

        virtual double evaluate(double x) const override { return f(x) * g(x); }
    };

    public:
    BinaryMultiply(Function f, Function g)
    :   Function(std::make_shared<Implementation>(f, g))
    {}
};

inline UnaryMinus operator - (Function f) { return UnaryMinus(f); }
inline BinaryAdd operator + (Function f, Function g) { return BinaryAdd(f, g); }
inline BinaryMultiply operator * (Function f, Function g) { return BinaryMultiply(f, g); }

#include <iostream>
int main() {
    Identity x;
    Function result = -x * (x + x) + x;
    std::cout << result(2) << '\n';
}

Answer 2

You could overload operator* as a free-standing function. Put it in Product.h/cpp even if it is not a member of Product, since it is tightly related to it.

In Product.h:

Function* operator*(Function *f, Function *g);

In Product.cpp:

Function* operator*(Function *f, Function *g) {
    return new Product(f, g);
}

Same with addition, subtraction, etc.

Alternatively, you can implement them as member functions, but put the definition in Function.cpp and include Product.h, etc. there.

Note that your design has a huge flaw. You create new Function objects on the heap and pass around pointers. You need to delete them somewhere in your code, I assume in your deconstructors. But then you also need to take care about copying your objects. Usually, it is a nightmare to care about proper deletion manually, and there are automatic ways to do it (called "memory management"). You could consider using smart pointers, for example. Have a look at std::shared_ptr. While not being most efficient in every case, it is a good thing to use in general when you first want to learn the language and not too many details about memory management. To apply shared_ptr to your code, replace every Function* with shared_ptr<Function>, and replace new Function(...) with make_shared<Function>(...) (and the same with other types).

Also note that * for math functions is ambiguous: In some contexts / literature, f*g means multiplying the result, while in others it means function convolution.

Answer 3

Here's a C++11 generic programming solution that doesn't rely on inherited polymorphism (i.e. - virtual functions) and doesn't require dynamic allocation.

I'm not a C++ expert and this can likely be improved upon significantly, but it works and gets the idea across. In particular, the code below only works for functions of double's. You can probably make the operands and return types be a template type too, so that this can generically work on different types (e.g. - complex) too. I don't know the proper way to scope the template operators so that you can use the shorthand operator notation and not have them accidentally invoked on (or make ambiguous) other types that have operator()(double x). If anyone has any suggestions to improve upon this answer, then please chime in and I'll edit my answer.

#include <iostream>

using namespace std;

struct Identity
{
  double operator() (double x) const { return x; }
};

struct Constant
{
  template<typename T1>
  Constant(const T1 &x) : _x(x) {}

  double operator()(double x) const { return _x; }

private:
  double _x;
};


template<typename T1>
struct Negate
{
  Negate(const T1 &f) : _f(f) {}

  double operator() (double x) const { return -_f(x); }

private:
  T1 _f;
};

template<typename T1>
struct Reciprocal
{
  Reciprocal(const T1 &f) : _f(f) {}

  double operator() (double x) const { return 1 / _f(x); }

private:
  T1 _f;
};

template<typename T1, typename T2>
struct Sum
{
  Sum(const T1 &f, const T2 &g) : _f(f), _g(g) {}

  double operator() (double x) const { return _f(x) + _g(x); }

private:
  T1 _f;
  T2 _g;
};

template<typename T1, typename T2>
struct Product
{
  Product(const T1 &f, const T2 &g) : _f(f), _g(g) {}

  double operator() (double x) const { return _f(x) * _g(x); }

private:
  T1 _f;
  T2 _g;
};

template<typename T1> 
Negate<T1> operator-(const T1 &f) 
{ return Negate<T1>(f); }

template<typename T1, typename T2> 
Sum<T1, T2> operator+(const T1 &f, const T2 &g) 
{ return Sum<T1, T2>(f, g); }

template<typename T1, typename T2> 
Sum<T1, Negate<T2> > operator-(const T1 &f, const T2 &g) 
{ return Sum<T1, Negate<T2> >(f, Negate<T2>(g)); }

template<typename T1, typename T2> 
Product<T1, T2> operator*(const T1 &f, const T2 &g) 
{ return Product<T1, T2>(f, g); }

template<typename T1, typename T2> 
Product<T1, Reciprocal<T2> > operator/(const T1 &f, const T2 &g) 
{ return Product<T1, Reciprocal<T2> >(f, Reciprocal<T2>(g)); }


int main()
{
  auto f = (Identity() * Constant(4.0) + Constant(5)) / Identity();  // f(x) = (x * 4 + 5) / x; f(2) = 6.5
  auto g = f * f;                                                    // g(x) = f(x) * f(x);     g(2) = 42.25

  cout << f(2) << " " << g(2) << " " << (g / f)(2) << endl;  // prints 6.5 42.25 6.5

  return 0;
}

EDIT: The main drawback of this approach is that the type of a "formula" must be fully known at compile time and embedded in the template generated classes. That means that very complicated formulas will generate lots of different classes and code. So, this approach could lead to nasty code bloat. Also, you can't do something like:

for (i = 1; i < j; ++i)  // raise f to the jth power (other than 0)
  f *= f;  

Since the type of f must be fully known at compile time and the multiplication is invoking new types iteratively. The other approaches that use a class hierarchy, dynamic allocation (with auto-cleanup) and polymorphism can do this though and don't have the problem of code bloat. Still, it was interesting to try.