// ---------------------------------------
//  funtest.cpp
// ---------------------------------------

#include "w.h"
/*
	C++ is not perfectly suited to functional programming.
	However, the virtual base class fun provides a generic
	base type for functions: pointers to instances of
	derived classes, such as gfun (below), are sufficiently
	'equivalent' that tools (like ftab) can be designed once,
	for all of them.  Thus we can go on to design other
	global tools such as root finders, integrators,
	differential-equation solvers, and graph plotters.
	Moreover, new classes such as dif (below), construct
   (for example) from a gfun object its derivative: all the
	tools designed will again apply to this new 'function'.
	As a matter of policy, we shall usually overload our
	tools so that they also apply to top-level functions,
	for example the fptype1 functions defined in w.h

												Richard Hall

	*/
 // ----------------------------------------------
 class gfun : public fun
	{ private:
		  double a;
	  public:
		  gfun(double ai = 1): a(ai) {}
		  void seta(double ai) { a = ai; }
		  double geta() { return a; }
		  double f(double x)
			 { return exp(-a*x*x); } };
 // ----------------------------------------------

	class dif: public fun         // derivative
		  { private: fun *f1;
				double h,h2;
			 public:
				dif(fun *f1i, double hi = 0.01)
				  { f1 = f1i; h = hi; h2 = 0.5/h; }
				double f(double x)
					{ double y =  (f1->f)(x+h);
					  y -= (f1->f)(x-h);
					  return y*h2; } };
// -----------------------------------------------


 void main()

	{
	  nl(0);
	  banner("funtest.cpp");
	  gfun g1; // a = 1 (default)
	  nl();
	  pfl("g1(x) = exp(-ax^2),  a = %.6f   g1(1) = %0.6f",
				 g1.geta(), g1.f(1));
	  nl();
	  gfun *g = new gfun(0.25);
	  pfl("g(x) = exp(-ax^2),   a = %.6f    g(2) = %.6f"
				,g->geta(),g->f(2));
	  nl();
	  g->seta(1);
	  pfl("Now a in g = %.6f",g->geta());
	  ftab(g,"g",-2,2,4);   // g is gfun pointer
	  dif *dg = new dif(g);
	  ftab(dg,"dg",-2,2,4); // dg is dif pointer
	  wait();
	  delete g,dg;
		}

 /* output

   -------------
    funtest.cpp
   -------------

   g1(x) = exp(-ax^2),  a = 1.000000   g1(1) = 0.367879

   g(x) = exp(-ax^2),   a = 0.250000    g(2) = 0.367879

   Now a in g = 1.000000


   g( -2.000000) =     0.018316
   g( -1.000000) =     0.367879
   g(  0.000000) =     1.000000
   g(  1.000000) =     0.367879
   g(  2.000000) =     0.018316


   dg( -2.000000) =     0.073275
   dg( -1.000000) =     0.735734
   dg(  0.000000) =     0.000000
   dg(  1.000000) =    -0.735734
   dg(  2.000000) =    -0.073275


  */