Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

GSL root finding

24 views
Skip to first unread message

the...@hotmail.com

unread,
Mar 6, 2017, 10:37:33 AM3/6/17
to
Hi all,
I am having some troubles with the following code:


#include <stdlib.h>
#include <stdio.h>
#include <gsl/gsl_vector.h>
#include <gsl/gsl_multiroots.h>

using namespace std;

// create a data structure which will contain the parameters for the main function,
// f(x, y, z) = k1*x^2*y - k2*z^2, and the secondary reactions needed to calculate
// x, y and z of the main function, f(x, y) = 2*x - y - 2*x0 + y0 and
// f(y, z) = x + z - x0 -z0
struct rootFind_Params{
// k1, k2 rate constants
// x0, y0, z0 initial conditions
double k1, k2, x0, y0, z0;
};

// create a function that will store the vector result f(v, p), in f for argument
// v and parameters p
int rootFind(const gsl_vector* v, void* p, gsl_vector* f){
struct rootFind_Params* params = (struct rootFind_Params*)p;
const double k1 = (params->k1);
const double k2 = (params->k2);
const double x0 = (params->x0);
const double y0 = (params->y0);
const double z0 = (params->z0);
// return the i-th element of the vector v
const double x = gsl_vector_get(v,0);
const double y = gsl_vector_get(v,1);
const double z = gsl_vector_get(v,2);
// set the value of the i-th element of a vector f to f(x, y, z), f(x, y) or
// f(x, z)
gsl_vector_set(f, 0, (k1*x*x*y - k2*z*z));
gsl_vector_set(f, 1, (2*x - y - 2*x0 + y0));
gsl_vector_set(f, 2, (x + z - x0 -z0));
// return GSL_SUCCESS if calculation was completed successfully
return GSL_SUCCESS;
};

// do main computation
int main(){
// initialize solver:
// create a hybrid solver of type T
const gsl_multiroot_fsolver_type* T = gsl_multiroot_fsolver_hybrids;
// return a pointer to a newly allocated instance of a solver of type T
// for a system of three dimensions
gsl_multiroot_fsolver* s = gsl_multiroot_fsolver_alloc(T, 3);
// provide the function to solve:
// create a variable status which will report an error if the root-finding function
// is unsuccessful
int status;
// specify initial values of indices
size_t i, iter = 0;
// specify dimension of the system, ie. the number of components of the vectors
// v and f
const size_t n = 3;
// assign numerical values to the parameters in the data structure
struct rootFind_Params params = {1, 0.7, 0.5, 1, 0};
// create a function object f which can be passed to a solver
gsl_multiroot_function f = {&rootFind, n, &params};
// create an initial guess for x, y and z
double v_init[3] = {-2.0, -2.0, -2.0};
gsl_vector* v = gsl_vector_alloc(n);
gsl_vector_set(v, 0, v_init[0]);
gsl_vector_set(v, 1, v_init[1]);
gsl_vector_set(v, 2, v_init[2]);
// set the already allocated solver
gsl_multiroot_fsolver_set(s, &f, v);

print_state(iter, s);

do
{
// increment the index
iter++;
// iterate the solver
status = gsl_multiroot_fsolver_iterate(s);
print_state(iter, s);
// check if solver is stuck
if (status)
break;
// tests the residual value of f against the absolute error bound, 1e-7;
// returns GSL_SUCCESS if the condition sum_i |f_i| < 1e-7, else returns
// GSL_CONTINUE
status = gsl_multiroot_test_residual(s->f, 1e-7);
}
while (status == GSL_CONTINUE && iter < 1000);

// report error to user
printf("status = %s\n", gsl_strerror(status));

// free solver
gsl_multiroot_fsolver_free(s);
// free vector
gsl_vector_free(v);

return 0;
}

// display intermediate state of the solution
int print_state(size_t iter, gsl_multiroot_fsolver* s){
printf ("iter = %3u v = % .3f % .3f " "f(v) = % .3e % .3e\n", iter,
gsl_vector_get(s->v, 0),
gsl_vector_get(s->v, 1),
gsl_vector_get(s->v, 2),
gsl_vector_get(s->f, 0),
gsl_vector_get(s->f, 1),
gsl_vector_get(s->f, 2));
}







I am getting two errors:
rootFinding.cc:70:24: error: 'print_state' was not declared in this scope
print_state(iter, s);
and also this one,
rootFinding.cc:103:27: error: 'struct gsl_multiroot_fsolver' has no member named 'v'
gsl_vector_get(s->v, 0),

Why am I having these problems? Any help would be much appreciated!

Louis Krupp

unread,
Mar 6, 2017, 2:58:05 PM3/6/17
to
On Mon, 6 Mar 2017 07:37:06 -0800 (PST), the...@hotmail.com wrote:

>Hi all,
>I am having some troubles with the following code:
>
>
>#include <stdlib.h>
>#include <stdio.h>
>#include <gsl/gsl_vector.h>
>#include <gsl/gsl_multiroots.h>
>
<snip>

> print_state(iter, s);

<snip>

>int print_state(size_t iter, gsl_multiroot_fsolver* s){

<snip>
>
>I am getting two errors:
>rootFinding.cc:70:24: error: 'print_state' was not declared in this scope
> print_state(iter, s);
>and also this one,
>rootFinding.cc:103:27: error: 'struct gsl_multiroot_fsolver' has no member named 'v'
> gsl_vector_get(s->v, 0),
>
>Why am I having these problems? Any help would be much appreciated!

You need to declare print_state() before you use it. See:

http://www.cplusplus.com/articles/yAqpX9L8/

I did a search for gsl_multiroot_fsolver, and I found this
declaration:

typedef struct
71 {
72 const gsl_multiroot_fsolver_type * type;
73 gsl_multiroot_function * function ;
74 gsl_vector * x ;
75 gsl_vector * f ;
76 gsl_vector * dx ;
77 void *state;
78 }
79 gsl_multiroot_fsolver;

in this page:

https://fossies.org/dox/gsl-2.3/gsl__multiroots_8h_source.html

so it's possible that gsl_multiroot_fsolver really has no member named
'v'.

Wild guess: Try passing 'v' as an argument to print_state and using
'v' instead of 's->v'.

Louis
0 new messages