yanera_mdmin_simplex.c

Go to the documentation of this file.

/*============================================================================
Copyright 2006 Thad Harroun, Kevin Yager

This file is part of "yanera 2.0".

"yanera 2.0" is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.

"yanera 2.0" is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with "yanera 2.0"; if not, write to the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA 
=============================================================================*/
#include "yanera.h"
#include "yanera_mdmin_simplex.h"
#include "yanera_mdmin.h"
#include "yanera_util.h"
/*============================================================================*/
int yanera_mdmin_simplex_main(yanera_container *yanera)
{
  int                         i, j, count;
  int                         status;
  gsl_vector                 *guess, *step_size;
  gsl_multimin_function       mdmin_func;
  gsl_multimin_fminimizer    *solvr;
  
  /* Inital guess is the paramters as they were passed into here */
  count = 0;
  for (i=0; i<yanera->parameters.n; i++) 
  {
    if (yanera->parameters.f[i] == NO)
    { count += 1; }
  }
  guess = gsl_vector_alloc(count);
  step_size = gsl_vector_alloc(count);
  for (i=0, j=0; i<yanera->parameters.n; i++) 
  {
    if (yanera->parameters.f[i] == NO)
    {
      gsl_vector_set(guess, j, yanera->parameters.p[i]);
      gsl_vector_set(step_size, j, yanera->parameters.p[i]/100.0);
      j++;
    }
  }

  mdmin_func.f      = &mdmin_simplex_f;
  mdmin_func.n      = count;
  mdmin_func.params = (void *)yanera;
  
  solvr = gsl_multimin_fminimizer_alloc( \
                  yanera->misc.fit_method_mdmin_simplex, count);
  
  /*==================================================*/
  count = 0;
  gsl_multimin_fminimizer_set(solvr, &mdmin_func, guess, step_size);
  yanera_report(yanera, count, solvr);

  do
  {
    count++;
    
    status = gsl_multimin_fminimizer_iterate(solvr);
    if (status != GSL_SUCCESS) 
    {
      fprintf(stderr,"STATUS: %s\n", gsl_strerror(status));
      break;
    }
    
    yanera_report(yanera, count, solvr);
    
    status = gsl_multimin_test_size(solvr->size, 1.0e-5);
    
    if (status == GSL_SUCCESS) 
    {
      fprintf(stderr,"STATUS: %s\n", gsl_strerror(status));
    }
    
  } while ((status == GSL_CONTINUE) && (count < yanera->misc.fit_iterations));

  yanera_report(yanera, -1, solvr);

  gsl_multimin_fminimizer_free(solvr);
  gsl_vector_free(guess);
  gsl_vector_free(step_size);
  return GSL_SUCCESS;
}
/*============================================================================*/
double mdmin_simplex_f(const gsl_vector *parameters, \
                       void *yanera_pointer)
{
  int     i, j;
  double  tmp;
  yanera_container  *yanera;
  /*
  ** r =>> residual
  */
  yanera = (yanera_container *)yanera_pointer;
  /*
  **  Load the params back into the parameter array.
  */
  for(i=0, j=0; i<yanera->parameters.n; i++)
  {
    if (yanera->parameters.f[i] == NO)
    {
      yanera->parameters.p[i] = gsl_vector_get(parameters, j);
      j++;
    }
  }

  tmp = mdmin_chisq(yanera);
  
  return tmp;
}

---