Hello all,
I had a question about the gsl simulated annealing code and I was hoping that
maybe someone here might have a solution. I have been using the simulated
annealing functions from the gsl library in an attempt to perform basic linear
regressions. I have successfully utilized it for a simple one-variable case,
however I am now in the process of trying to fit more variables. (i.e) y = b0 +
b1x1 + b2x2. I am confused on how to pass multiple variables to the
gsl_siman_solve function. I have tried a simple array, and also tried to use
the gsl_vector structure. I am either having problems retrieving the current
condition in the energy, step, and metric evaluations, or am passing condition
to the gsl_siman_solve function improperly.Ultimately I would like to pass in 2
variables (B0 and B1, or perhaps even more, and minimize them using
experimental data. I have tried to alter the step function and metric function
for multiple variables to no success. I have attached a sample code I thought
might have worked, but every time I try to run it, I get a segmentation fault.
If you have the time, maybe you could take a peek and see what I am doing wrong.
Thanks,
Jesse Coyle
Nuclear & Radiological Engineering
Georgia Institute of Technology
#include <math.h>
#include <stdlib.h>
#include <string.h>
#include <gsl/gsl_siman.h>
#include <gsl/gsl_vector.h>
/* set up parameters for this simulated annealing run */
/* how many points do we try before stepping */
#define N_TRIES 200
/* how many iterations for each T? */
#define ITERS_FIXED_T 10
/* max step size in random walk */
#define STEP_SIZE 1.0
/* Boltzmann constant */
#define K 1.0
/* initial temperature */
#define T_INITIAL 0.008
/* damping factor for temperature */
#define MU_T 1.003
#define T_MIN 2.0e-10
gsl_siman_params_t params
= {N_TRIES, ITERS_FIXED_T, STEP_SIZE,
K, T_INITIAL, MU_T, T_MIN};
/* now some functions to test in one dimension */
double yvals[10]={70,65,90,95,110,115,120,140,155,150};
double xvals[10]={80,100,120,140,160,180,200,220,240,260};
gsl_vector * initialconfiguration;
double xavg = 170.0;
double yavg = 111.0;
double E1(void *xp){
double beta0 = gsl_vector_get(xp,0);
double beta1 = gsl_vector_get(xp,1);
double yscore=0.0;
int j=0;
for (j=0; j<10; j++){
yscore = pow((yvals[j]-(beta0)-(beta1*xvals[j])),2.0)+yscore;
}
return yscore;
}
double M1(void *xp, void *yp)
{
double x0 = gsl_vector_get(xp,0);
double x1 = gsl_vector_get(xp,1);
double y0 = gsl_vector_get(yp,0);
double y1 = gsl_vector_get(yp,1);
return sqrt(pow((y0-x0),2.0) + pow((y1-x1),2.0));
}
void S1(const gsl_rng * r, void *xp, double step_size)
{
double old_x0 = gsl_vector_get(xp,0);
double new_x0;
double u = gsl_rng_uniform(r);
new_x0 = u * 2 * step_size - step_size + old_x0;
u = gsl_rng_uniform(r);
double old_x1 = gsl_vector_get(xp,1);
double new_x1;
new_x1 = u * 2 * step_size - step_size + old_x1;
/*memcpy(gsl_vector_ptr(xp,0), &new_x0, sizeof(new_x0));
memcpy(gsl_vector_ptr(xp,1), &new_x1, sizeof(new_x1));*/
gsl_vector_set(xp,0,new_x0);
gsl_vector_set(xp,1,new_x1);
}
void P1(void *xp)
{
printf ("%12g %12g", gsl_vector_get(xp,0), gsl_vector_get(xp,1));
}
int
main(int argc, char *argv[])
{
initialconfiguration = gsl_vector_alloc(2);
const gsl_rng_type * T;
gsl_rng * r;
double x_initial = 25.0;
gsl_vector_set_all(initialconfiguration, 25.0);
gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc(T);
gsl_siman_solve(r, &initialconfiguration, E1, S1, M1, P1,
NULL, NULL, NULL,
sizeof(initialconfiguration), params);
gsl_rng_free (r);
return 0;
}
_______________________________________________
Help-gsl mailing list
[email protected]
http://lists.gnu.org/mailman/listinfo/help-gsl