#ifndef FUNCTION_H
#define FUNCTION_H

#include <functional>
#include <vector>
#include <iostream>
#include <time.h>
#include <math.h>

using namespace std;

//
// The following group of classes were an attempt to develop a way to call functions on functions
//
class Function : public unary_function<double,double>
{
public:
	Function() {};
	virtual ~Function() 
	{
		//do garbage collection on tempFunctions //crashes program?
		for (vector<Function*>::iterator  iter = tempFunctions.begin(); iter != tempFunctions.end(); iter++)
			delete *iter;
	};
	virtual double operator()(double x){return x;};
	virtual Function& operator^(Function& f);
	virtual Function& operator*(Function& f);
	virtual Function& operator+(Function& f);
	virtual Function& operator-(Function& f);
	virtual Function& operator()(Function& f);

	virtual Function& operator^(double x);
	virtual Function& operator*(double x);
	virtual Function& operator+(double x);
	virtual Function& operator-(double x);

	virtual Function& operator-();
	virtual Function& getDerivative() { return *this;};
protected:
	 vector<Function*> tempFunctions;
};

	

class PowFunc : public Function
{
public:
	PowFunc(Function&f1, Function& f2) : a(f1), b(f2) {}
	virtual double operator() ( double x) { return pow(a(x),b(x));};
	virtual Function& getDerivative();
	
protected:
	Function& a;
	Function& b;
};

Function& Function::operator^(Function& f)
{
	PowFunc* a;
	tempFunctions.insert(tempFunctions.end(), (a = new PowFunc(*this, f)));
	return *a;
};

Function& PowFunc::getDerivative() 
{ 
	cerr << "can not evaluate f(x)^g(x).  Need to develop numerical approximation" <<endl;
	return a*b;
};


class AdditionFunc : public Function
{
public:
	AdditionFunc(Function&f1, Function& f2) : a(f1), b(f2) {}
	virtual double operator() ( double x) { return (a(x)+b(x));};
	virtual Function& getDerivative();
protected:
	Function& a;
	Function& b;
};

Function& Function::operator+(Function& f)
{
	AdditionFunc* a;
	tempFunctions.insert(tempFunctions.end(), (a = new AdditionFunc(*this, f)));
	return *a;
};

Function& AdditionFunc::getDerivative() 
{ 
	return (a.getDerivative() + b.getDerivative());
};

class SubtractionFunc : public Function
{
public:
	SubtractionFunc(Function&f1, Function& f2) : a(f1), b(f2) {}
	virtual double operator() ( double x) { return (a(x)-b(x));};
	virtual Function& getDerivative();
	
protected:
	Function& a;
	Function& b;
};

Function& Function::operator-(Function& f)
{
	SubtractionFunc* s;
	tempFunctions.insert(tempFunctions.end(), (s = new SubtractionFunc(*this, f)));
	return *s;
};

Function& SubtractionFunc::getDerivative() 
{ 
	return (a.getDerivative() - b.getDerivative());
};

class MultiplyFunc : public Function
{
public:
	MultiplyFunc(Function& f1, Function& f2) :a(f1), b(f2) {};
	virtual double operator() ( double x) { return a(x)*b(x);};
	virtual Function& getDerivative();
	
protected:
	Function& a;
	Function& b;
};

Function& Function::operator*(Function& f)
{
	MultiplyFunc * m;
	tempFunctions.insert(tempFunctions.end(), (m = new MultiplyFunc(*this, f)));
	return *m;
};

Function& MultiplyFunc::getDerivative() 
{ 
	return (a.getDerivative()*b + a*b.getDerivative());
};

class CompositionFunc : public Function
{
public:
	CompositionFunc(Function& f1, Function& f2) :a(f1), b(f2) {};
	virtual double operator() ( double x) { return a(b(x));};
	virtual Function& getDerivative();
protected:
	Function& a;
	Function& b;
};

Function& Function::operator()(Function& f)
{
	CompositionFunc* c;
	tempFunctions.insert(tempFunctions.end(), (c = new CompositionFunc(*this, f)));
	return *c;
};

Function& CompositionFunc::getDerivative() 
{ 
	return (a.getDerivative()(b)*b.getDerivative());
};

class ConstFunc : public Function
{
public:
	ConstFunc(double x) : _x(x) {};
	virtual double operator()(double x){return _x;};
protected:
	double _x;
};

class IdentityFunction : public ConstFunc
{
public:
	IdentityFunction() : ConstFunc(1) {};
};

class VectorFunc : public Function
{
public:
	VectorFunc(vector<double> &list) : v(list) {};
	virtual double operator()(double n) { return v[(int)n]; };
protected:
	vector<double> & v;
};

Function& Function::operator^(double x)
{
	ConstFunc* c;
	tempFunctions.insert(tempFunctions.end(), (c = new ConstFunc(x)));
	return ((*this)^(*c));
}
Function& Function::operator*(double x)
{
	ConstFunc* c;
	tempFunctions.insert(tempFunctions.end(), (c = new ConstFunc(x)));
	return ((*this)*(*c));
};
Function& Function::operator+(double x)
{
	ConstFunc* c;
	tempFunctions.insert(tempFunctions.end(), (c = new ConstFunc(x)));
	return ((*this)+(*c));
}
;
Function& Function::operator-(double x)
{
	ConstFunc* c;
	tempFunctions.insert(tempFunctions.end(), (c = new ConstFunc(x)));
	return ((*this)-(*c));
}
;
Function& Function::operator-()
{
	ConstFunc* c;
	tempFunctions.insert(tempFunctions.end(), (c = new ConstFunc(0)));
	return ((*c) - (*this));
};

#define functor(NAME, FUNC) \
class NAME : public Function \
{\
public: \
	virtual Function&  operator()(Function& x) { return Function::operator ()(x);}; \
	virtual double operator()(double x) { return FUNC(x);}; \
}; 

#define FunctionObject(FUNCTION) \
	functor(FunctionObject##FUNCTION, FUNCTION) \
	FunctionObject##FUNCTION FUNCTION##F;

#define FunctionObjectSpec(NAME, FUNCTION, OBJECT) \
	functor(NAME,FUNCTION) \
		NAME OBJECT;\

#define FunctionObjectSpecPtr(NAME, FUNCTION, OBJECT) \
	functor(NAME,FUNCTION) \
		NAME* OBJECT;\

 
#define FunctionObject_(NAME, FUNCTION) \
	functor(NAME,FUNCTION) \
	NAME FUNCTION##F;

FunctionObjectSpec(ExponentialFunc, exp, exponentialFunc);

FunctionObjectSpec(CosFunc, cos, cosFunc)

FunctionObjectSpec(SinFunc, sin, sinFunc);

class SumFunc : public Function
{
public:
	enum Type {Int, Dub};
	SumFunc(int lower, int upper, int increment = 1) : a(lower), b(upper), inc(increment), type(Int) {};
	SumFunc(double lower, double upper, double increment = 1) : a(lower), b(upper), inc(increment), type(Dub){};
	virtual double operator()(double x) { return ((b-a)/inc)*x; };
	virtual Function& operator()(Function& f) 
	{ 
		double val = 0;
		for (double n = a;n<=b;n+=inc)
			val += f(n);
		ConstFunc *c = new ConstFunc(val);
		tempFunctions.insert(tempFunctions.end(), c);
		return *c;
	};


protected:
	double a, b, inc;
	Type type;
};

#endif