# Help using gsl_multiroot_fsolver in C

1. Nov 14, 2012

### Yuu Suzumi

Dear all,

Unfortunately I am not skilled in GSL. I am using gsl_multiroot_fsolver_hybrid in C to solve a system of 3 equations with three unknowns. I have tried very hard to find my
error and repeatedly consulted the GSL manual, but without success. I hope
there is someone who can offer some insight! I will first post the part of the
code where the error occurs. Please tell me if it is necessary to send a
simpler but standalone code displaying the problem. Sofar I post only the two
functions to keep it short.

The problem occurs in the second (iter =1) iteration, inside int quasinormal_f
, where the vector x (independent variable of the problem, which the solver
should vary) takes values (nan, nan, nan). I cannot see why.

The functions I defined are:
/* --- --- */
int quasinormal_f (const gsl_vector * x, void *params,
gsl_vector * f)
{

const double x0 = gsl_vector_get (x, 0);
const double x1 = gsl_vector_get (x, 1);
const double x2 = gsl_vector_get (x, 2);

qparams->N =x0;
qparams->A =x1;
qparams->S =x2;

printf("INNER POINT 1 -- Attention!\n");
printf("N,A,S = %e %e %e \n",
x0, x1, x2);
const double y0 = norm_quasinormal(qparams) - 1.0;
//const double y0 = 0;

const double y1 = mean_quasinormal(qparams) - 0.0;
//const double y1 =0;

const double y2 = var_quasinormal(qparams) - kappa_var(gcosmo-> theta);

gsl_vector_set (f, 0, y0);
gsl_vector_set (f, 1, y1);
gsl_vector_set (f, 2, y2);

return GSL_SUCCESS;
}

int COMPUTE_QUASINORMAL_PARAMETERS()
{

const gsl_multiroot_fsolver_type *T;
gsl_multiroot_fsolver *s;

int status;
size_t i, iter = 0;

const size_t ndim = 3;

gsl_multiroot_function F;

F.f = &quasinormal_f;
F.params = NULL;
F.n = 3;

double x_init[3] = {1.5, 0.1,0.5*gcosmo->ausgelagertes_lognormal_sigma};
gsl_vector *x = gsl_vector_calloc (ndim);

printf("Check point 2: x = %e %e %e \n",
gsl_vector_get(x,0), gsl_vector_get(x,1), gsl_vector_get(x,2));

gsl_vector_set (x, 0, x_init[0]);
gsl_vector_set (x, 1, x_init[1]);
gsl_vector_set (x, 2, x_init[2]);

printf("Check point 3: x = %e %e %e \n",
gsl_vector_get(x,0), gsl_vector_get(x,1), gsl_vector_get(x,2));

T = gsl_multiroot_fsolver_hybrid;
s = gsl_multiroot_fsolver_alloc (T, 3);

gsl_multiroot_fsolver_set (s, &F, x);
printf("THE SOLVER WAS SET!\n");
do
{
iter++;
printf("iter = %i\n",
iter);

printf("x_init = %e %e %e \n",
x_init[0], x_init[1], x_init[2]);

status = gsl_multiroot_fsolver_iterate (s);

printf("THE SOLVER WAS ITERATED\n");
printf("iter = %i\n",
iter);

printf("x_init = %e %e %e \n",
x_init[0], x_init[1], x_init[2]);

if (status) /* check if solver is stuck */
break;

status =
gsl_multiroot_test_residual (s->f, 1e-7);
}
while (status == GSL_CONTINUE && iter < 1000);

printf ("status = %s\n", gsl_strerror (status));

gsl_vector_free (x);

return 0;
}
/* --- --- */

The output of the programme is:

Check point 2: x = 0.000000e+00 0.000000e+00 0.000000e+00
Check point 3: x = 1.500000e+00 1.000000e-01 1.035659e-01
INNER POINT 1 -- Attention!
N,A,S = 1.500000e+00 1.000000e-01 1.035659e-01
INNER POINT 1 -- Attention!
N,A,S = 1.500000e+00 1.000000e-01 1.035659e-01
INNER POINT 1 -- Attention!
N,A,S = 1.500000e+00 1.000000e-01 1.035659e-01
INNER POINT 1 -- Attention!
N,A,S = 1.500000e+00 1.000000e-01 1.035659e-01
THE SOLVER WAS SET!