/************************************************************************/
/*
opt-powell: optimize policy using Numerical Recipes Powell
*/
/************************************************************************/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#include "main.h"
#include "main2.h"

/*****************************************************************************/
/*****************************************************************************/
/*****************************************************************************/

extern SIM sim;

/*****************************************************************************/
/* Powell stuff */

#define FTOL 1.0e-7
static float p[MAX_N_PARAMETERS+1];
static float p_min[MAX_N_PARAMETERS+1];
static float p_max[MAX_N_PARAMETERS+1];
static float **xi;

/*****************************************************************************/
/************************************************************************/
/************************************************************************/
/* Main: initialize and then call powell optimizer */

int main( int argc, char **argv )
{
  int i, j, k;
  float **matrix();
  float fret;
  double step = 0.01;
  double cost2;
  PARAMETER *read_parameter_file();
  
  init_default_parameters( &sim );

  if ( argc < 2 )
    {
      fprintf( stderr, "Optimize using which parameter file?\n" );
      fprintf( stderr, "%s parameter-file\n", argv[0] );
      exit( -1 );
    }
  sim.params = read_parameter_file( argv[1] );
  sim.n_parameters = process_parameters( sim.params, &sim, 1 );
  if ( sim.n_parameters > MAX_N_PARAMETERS )
    {
      fprintf( stderr, "Too many parameters %d > %d\n",
	       sim.n_parameters, MAX_N_PARAMETERS );
      exit( -1 );
    }
  sprintf( sim.output_file, "%s.new", argv[1] );

  init_sim( &sim );

  sim.n_func_calls_per_eval = 1;
  sim.all_time_low_cost = 1e20;
  sim.debug_criterion = 1;

  for( i = 1; i <= sim.n_parameters; i++ )
    {
      p_min[i] = -1e10;
      p_max[i] = 1e10;
    }
  xi = matrix( 1, sim.n_parameters, 1, sim.n_parameters );

  for ( i = 1; i <= sim.n_parameters; i++ )
    for ( j = 1; j <= sim.n_parameters; j++ )
      xi[i][j] = (i == j ? step : 0.0 );
  parameters_to_vector( sim.params, &p[1] );
  call_many_func( p );

  printf( "Press return to continue.\n" );
  getchar();

  /* return() here if all you want to see is one policy run. */
  powell( p, p_min, p_max, xi, sim.n_parameters, FTOL, &(sim.iter), 
	  &fret, call_many_func, sim.debug_criterion );

  cost2 = call_many_func( p );

  printf( "Iterations: %d (%d)\n\n", sim.iter, sim.func_calls );
  printf( "Minimum found at: \n" );
  for ( i = 1; i <= sim.n_parameters; i++ )
    printf( "p[%d] = %18.12f;\n", i-1, p[i] );
  printf( "\n\nMinimum function value = %12.6f, should equal %12.6f \n\n",
	  fret, cost2 );
  printf( "All time low cost: %g\n", sim.all_time_low_cost );
}

/************************************************************************/
