#ifndef _Minimizers_h_
#define _Minimizers_h_
/* David Weenink, 19 aug 1996 */

#ifndef _Thing_h_
	#include "Thing.h"
#endif

/*********** deferred class Minimizer **********************************/

#define Minimizer_members Thing_members							\
    long nParameters;	/* the number of parameters */			\
    double *p;			/* the parameters */					\
    double minimum;		/* current minimum */					\
    double *history;	/* previous minima */					\
    double tolerance;	/* stop criterium */					\
    Any object;		/* reference to the object that uses this Minimizer */	\
    long funcCalls;		/* the number of times 'func' has been called */	\
    long success;		/* indicates whether I'm done */					\
    long start;			/* start iteration series */						\
    long maxNumOfIterations; /* the current maximum number of iterations */	\
    long iteration;     /* the total number of iterations */				\
    int (*after) (I, Any aclosure); /* to be called after each iteration */	\
    Any aclosure;															\
    Any gmonitor;		/* graphics to monitor the minimization process */
#define Minimizer_methods Thing_methods				\
    int (*minimize) (I);  /* does the work */		\
    void (*reset) (I); /* reset the minimizer */	\
    void (*setParameters) (I, Any parameters);

class_create (Minimizer, Thing)
 
int Minimizer_init (I, long nParameters, Any object);
/*
	Preconditions:
		nParameters > 0;
	Postconditions:
		if (gmonitor) gmonitor != NULL
*/

void Minimizer_reset (I, const double guess[]);
/* reset the start values for the minimizer
 * Preconditions:
 *    guess != NULL;
 * Post conditions:
 *    p[] = guess[];
 *    my minimum = 1.0e30;
 *    success = maxNumOfIterations = iteration = funcCalls = 0;
 *    reset (me);
 */

void Minimizer_setAfterEachIteration (I, int (*after) (I, Any aclosure), 
	Any aclosure);
/* set the procedure that is executed after each iteration. */

void Minimizer_setParameters (I, Any parameters); /* for inheritors */

int Minimizer_minimize (I, long maxNumOfIterations, double tolerance, 
	int monitor);
/* Minimizes during maximally maxNumOfIterations. The gmonitor is initialized
 * before minimization and cleared afterwards.
 * Preconditions:
 *    maxNumOfIterations >= 0;
 *    tolerance > 0;
 * Postconditions:
 *    if (reset) Minimizer_reset called with xopt as initial guess.
 *    after each function call: funcCalls++
 *    after each iteration: iteration++
 */

int Minimizer_minimizeManyTimes (I, long numberOfTimes, long maxIterationsPerTime,
	double tolerance);
   
void Minimizer_drawHistory (I, Any graphics, long itmin, long itmax,
    double minimum, double maximum, int garnish);
    
/********** deferred class LineMinimizer ************************************/

#define LineMinimizer_members Minimizer_members						\
	/* the function to be minimized */		\
    double (*func) (Any object, const double p[]);	\
	double maxLineStep;	/*maximum step in line search direction */	\
    double *direction;	/* search direction vector */				\
    double *ptry;		/* point in search direction */
#define LineMinimizer_methods Minimizer_methods \
    void (*linmin) (I, double p[], double fp, double direction[], double *fret);	/* line minimization */
class_create (LineMinimizer, Minimizer)

int LineMinimizer_init (I, long nParameters, Any object, double (*func) 
	(Any object, const double p[]));

/*************** class PowellMinimizer **************************************/

#define PowellMinimizer_members LineMinimizer_members \
	double **directions;
#define PowellMinimizer_methods LineMinimizer_methods
class_create (PowellMinimizer, LineMinimizer)

Any PowellMinimizer_create (long nParameters, Any object, double (*func) 
	(Any object, const double p[]));

/******************  class PolakRibiereMinimizer ****************************/

#define PolakRibiereMinimizer_members LineMinimizer_members	\
	double *dp;	/* gradient */								\
    void  (*dfunc) (Any object, const double p[], double dp[]); 
	/* calculates gradient at position p */
#define PolakRibiereMinimizer_methods LineMinimizer_methods
class_create (PolakRibiereMinimizer, LineMinimizer)

Any PolakRibiereMinimizer_create (long nParameters, Any object, double (*func) 
	(Any object, const double p[]), void (*dfunc) (Any object, const double p[],
	double dp[]));

/******************  class SteepestDescentMinimizer**************************/

typedef struct structSteepestDescentMinimizer_parameters {
	double eta, momentum;
} *structSteepestDescentMinimizer_parameters;

#define SteepestDescentMinimizer_members Minimizer_members	\
	double eta, momentum;									\
    double (*func) (Any object, const double p[]);			\
    void  (*dfunc) (Any object, const double p[], double dp[]);	
	/* calculates gradient at position p */
#define SteepestDescentMinimizer_methods Minimizer_methods
class_create (SteepestDescentMinimizer, Minimizer)

Any SteepestDescentMinimizer_create (long nParameters, Any object, double (*func) 
	(Any object, const double p[]), void (*dfunc) (Any object, const double p[],
	double dp[]));
	 
/*****************  class LevenbergMarquardtMinimiz ************************/

#define LevenbergMarquardtMinimizer_members Minimizer_members \
    double (*hfunc) (Any object, const double p[], double dp[], double **ddp);
#define LevenbergMarquardtMinimizer_methods Minimizer_methods
class_create (LevenbergMarquardtMinimizer, Minimizer)

Any LevenbergMarquardtMinimizer_create (long nParameters, Any object,
	 double (*hfunc) (Any object, const double p[], double dp[], double **ddp));

/******************  class SimplexMinimizer ********************************/

/*
	Down-hill simplex method due to Nelder, J.A. and Mead, R. (1965),
	Computer Journal, vol. 7, pp. 308-313.
	Numerical Recipes, chpater 10.
*/

#define SimplexMinimizer_members Minimizer_members 			\
    double (*func) (Any object, const double p[]);			\
	double **simplex;										\
	double *y;	/* function value at simplex vertices */	\
	double scale;
#define SimplexMinimizer_methods Minimizer_methods
class_create (SimplexMinimizer, Minimizer)

int SimplexMinimizer_init (I, long nParameters, Any object, double (*func) 
	(Any object, const double p[]), double scale);
Any SimplexMinimizer_create (long nParameters, Any object, double (*func) 
	(Any object, const double p[]), double scale);

void SimplexMinimizer_initSimplex (I);

/***************  class SimulatedAnnealingMinimizer ************************/

/*
	Multidimensional minimization by simulated annealing combined with the
	downhill simplex method of Nelder & Mead.
	Source: Numerical Recipes chapter 10.9
*/

typedef struct structSimulatedAnnealingMinimizer_step {
	double t0, fraction; long step;
} *structSimulatedAnnealingMinimizer_step;

#define SimulatedAnnealingMinimizer_members SimplexMinimizer_members 		\
	double t0, tc; /* start and current temperature */						\
    double (*temperature) (I, Any tclosure);	/* the annealing scheme */	\
    Any tclosure;
#define SimulatedAnnealingMinimizer_methods SimplexMinimizer_methods
class_create (SimulatedAnnealingMinimizer, SimplexMinimizer)

/*
	The simmulated annealing minimization consists of the following elements:
	1. A description of possible system configurations.
	2. A generator of random changes in the configuration.
	3. An objective function 'func' whose minimization is the goal.
	4. A control parameter 'temperature' and an annealing schedule which
	   tells how it is lowered from high to low values.
	
	The implementation is with a modified downhill simplex method.
	We add a positive, logarithmically distributed random variable, proportional
	to the temperature (T),  to the stored function value associated with every
	vertex of the simplex, and we subtract a similar random variable from the
	function value of every new point that is tried as a replacement point. This
	method always accepts a true downhill step but sometimes accepts an uphill
	one.
	In the limit T-->0,  this algorithm reduces exactly to the downhill simplex
	and converges to a local minimum.
	At a finite value of T, the simplex expands to a scale that approximates the
	the size of the region that can be reached at this temperature.

	simplex	: matrix p[1...nParameters+1][1...nParameters]. Each row is a 
	  vertex of the starting simplex.
	y	: vector[1...nParameters+1]. Function values evaluated at each vertex.
	The nParameters+1 points of the simplex are initialized according to:
	  P(i) = P(0) + scale*E(i); i=1..N
	  where P(0) is your initial guess and the E(i) are  N unit vectors.
	temperature:
	
	The annealing scheme can be given in the function temperature. Its return 
	value must be >= 0.
	As an example take a look at the function 
	'_SimulatedAnnealingMinimizer_temperatureStep'
	that reduces the temperature by a factor T(1-eps) after every m iterations.

	PRECONDITIONS:
	  nParameters   >= 2;
 */


Any SimulatedAnnealingMinimizer_create (long nParameters, Any object, double 
	(*func) (Any object, const double p[]), double scale, double (*temperature)
	(I, Any tclosure), Any tclosure);

void SimulatedAnnealingMinimizer_setTemperatureProc (I, double (*temperature) 
	(I, Any tclosure), Any tclosure);

double _SimulatedAnnealingMinimizer_temperatureStep (I, Any tclosure);

/*************  class FletcherPowellMinimizer **************************/

#define FletcherPowellMinimizer_members LineMinimizer_members	\
	double *dp;													\
    void  (*dfunc) (Any object, const double p[], double dp[]);	
	/* calculates gradient at position p */
#define FletcherPowellMinimizer_methods LineMinimizer_methods
class_create (FletcherPowellMinimizer, LineMinimizer)

Any FletcherPowellMinimizer_create (long nParameters, Any object, double (*func) 
	(Any object, const double p[]), void (*dfunc) (Any object, const double p[],
	double dp[]));

/**********  class VDSmagtMinimizer ********************************/

typedef struct structVDSmagtMinimizer_parameters {
	double lineSearchGradient;
	long lineSearchMaxNumOfIterations;
} *structVDSmagtMinimizer_parameters;

#define VDSmagtMinimizer_members Minimizer_members                          \
    double (*func) (Any object, const double p[]);							\
	void  (*dfunc) (Any object, const double p[], double dp[]);             \
	double *dp;																\
    double lineSearchGradient;      	                                    \
    long lineSearchMaxNumOfIterations;                                      \
    double gr0, gropt, df, alplim, alpha, dalpha, alphamin;					\
    double *pc;	/* position of current point */								\
    double *gc;	/* gradient of current point */								\
    double *g0;	/* gradient at beginning of line search */					\
    double *s;	/* search direction for line search */						\
    double *srst;/* search direction for first iteration after restart */	\
    double *grst; /* gradient for first iteration after restart */			\
    double fc, grc, fch, gr2s, temp, grs, beta, gcg0;                       \
    double gamma, gamma_in, f0, gsq, gopt_sq;                               \
    long lineSearch_iteration, flag, again, one_up, restart;		    	\
    long restart_flag;
#define VDSmagtMinimizer_methods Minimizer_methods
class_create (VDSmagtMinimizer, Minimizer)

Any VDSmagtMinimizer_create (long dimension, void *object, double (*func) 
	(Any object, const double p[]), void (*dfunc) (Any object, const double p[],
	double dp[]));

#endif /* _Minimizer_h_ */
