How to obtain a function from gsl library integration?

In summary, the conversation discusses a problem with using the gsl library for integration via Monte Carlo. The person has obtained the correct result, but needs to know how to use it when there is a dependence on another variable. An example is given where the integration is only on one variable, but the function also has a dependence on another variable. The code used to solve this problem is provided, along with a solution that involves defining a struct and using the gsl_monte_vegas_integrate function.
  • #1
Fabio Kopp
16
0
I am having some problems to use the gsl library. I do the integration via monte carlo and I have obtained the right result, but a need to know how to use this result when I have a dependence on another variable. For example, f(x,y)=x*y, but the integration is only on x. Because I want to evaluate this function on y after the integration.
For the simple case, only one integration without the dependence in another variables, I obtained using this code:
Mod note: Edited to add [ code ] and [ /code ] tags.
C:
#include <iostream>
#include <gsl/gsl_math.h>
#include <gsl/gsl_monte.h>
#include <gsl/gsl_monte_vegas.h>
#include <cmath>

using namespace std;

double g (double *k, size_t dim, void *params)
{
  double A;
   A=1*k[0];
    return A;}int main(void)
{
  double res, err;  double xl[10] = { 0.0, 0.0,0.0};
  double xu[10] = { 1.0, 1.0,1.0};

  const gsl_rng_type *T;
  gsl_rng *r;

  gsl_monte_function G = { &g, 3, 0 };
  gsl_rng_env_setup ();
  T = gsl_rng_default;
  r = gsl_rng_alloc (T);

{

gsl_monte_vegas_state *s = gsl_monte_vegas_alloc (3);
gsl_monte_vegas_integrate (&G, xl, xu, 3, 10000, r, s,
&res, &err);

do
{
size_t calls = 50000;

gsl_monte_vegas_integrate (&G, xl, xu, 3, calls, r, s,
&res, &err);
}
while (fabs (gsl_monte_vegas_chisq (s) - 1.0) > 0.1);

cout<<"The result is "<< res <<"."<<endl;
gsl_monte_vegas_free (s);
}
  gsl_rng_free (r);

  return 0;
}

To compile the above code, I used:

c++ -Wall -I/usr/local/include -c teste00.cc
c++ -L/usr/local/lib teste00.o -lgsl -lgslcblas -lmI have a function that has the following form
A=x*k[0].
As I said before, I do not know how to implement this dependence in x to be evaluated after the integration on GSL.
Does anybody have any idea how to do this?
 
Last edited by a moderator:
Technology news on Phys.org
  • #2
Fabio Kopp said:
I am having some problems to use the gsl library. I do the integration via monte carlo and I have obtained the right result, but a need to know how to use this result when I have a dependence on another variable. For example, f(x,y)=x*y, but the integration is only on x. Because I want to evaluate this function on y after the integration.
For the simple case, only one integration without the dependence in another variables, I obtained using this code:
Mod note: Edited to add [ code ] and [ /code ] tags.
C:
#include <iostream>
#include <gsl/gsl_math.h>
#include <gsl/gsl_monte.h>
#include <gsl/gsl_monte_vegas.h>
#include <cmath>

using namespace std;

double g (double *k, size_t dim, void *params)
{
  double A;
   A=1*k[0];
    return A;}int main(void)
{
  double res, err;  double xl[10] = { 0.0, 0.0,0.0};
  double xu[10] = { 1.0, 1.0,1.0};

  const gsl_rng_type *T;
  gsl_rng *r;

  gsl_monte_function G = { &g, 3, 0 };
  gsl_rng_env_setup ();
  T = gsl_rng_default;
  r = gsl_rng_alloc (T);

{

gsl_monte_vegas_state *s = gsl_monte_vegas_alloc (3);
gsl_monte_vegas_integrate (&G, xl, xu, 3, 10000, r, s,
&res, &err);

do
{
size_t calls = 50000;

gsl_monte_vegas_integrate (&G, xl, xu, 3, calls, r, s,
&res, &err);
}
while (fabs (gsl_monte_vegas_chisq (s) - 1.0) > 0.1);

cout<<"The result is "<< res <<"."<<endl;
gsl_monte_vegas_free (s);
}
  gsl_rng_free (r);

  return 0;
}

To compile the above code, I used:

c++ -Wall -I/usr/local/include -c teste00.cc
c++ -L/usr/local/lib teste00.o -lgsl -lgslcblas -lmI have a function that has the following form
A=x*k[0].
As I said before, I do not know how to implement this dependence in x to be evaluated after the integration on GSL.
Does anybody have any idea how to do this?
I have found the solution for this problem and the code is the following.
#include<iostream>
#include<gsl/gsl_math.h>
#include <gsl/gsl_monte_vegas.h>
#include <gsl/gsl_monte.h>

using namespace std;

struct my_f_params { double a; double b; double c; };

double my_f (double x[], size_t dim, void * p) {
struct my_f_params * fp = (struct my_f_params *)p;
if (dim != 2)
{ fprintf (stderr, "error: dim != 2");
abort ();
}
return fp->a * x[0] * x[0]
+ fp->b * x[0] * x[1]
+ fp->c * x[1] * x[1];
}int main(){

gsl_monte_function F;
for (double i=1.0; i<10;++i){
struct my_f_params params = { i, 2.0, 1.0 };
F.f = &my_f;
F.dim = 2;
F.params = &params;double res, err;
double xl[10] = { 0.0, 0.0,0.0};
double xu[10] = { 1.0, 1.0,1.0};

const gsl_rng_type *T;
gsl_rng *r;
gsl_rng_env_setup ();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

{

gsl_monte_vegas_state *s = gsl_monte_vegas_alloc (3);
gsl_monte_vegas_integrate (&F, xl, xu, 3, 10000, r, s,
&res, &err);

do
{
size_t calls = 50000;gsl_monte_vegas_integrate (&F, xl, xu, 3, calls, r, s,
&res, &err);
}
while (fabs (gsl_monte_vegas_chisq (s) - 1.0) > 0.1);

cout<<" For i = "<< i << " the result is " << res << "."<<endl;gsl_monte_vegas_free (s);

}
gsl_rng_free (r);
}
return 0;
}
 

1. How do I use the gsl library for integration?

To use the gsl library for integration, you will need to include the necessary header files and link the library to your code. The gsl library provides various integration functions, such as gsl_integration_qag() and gsl_integration_qng(), which can be used to obtain a function from the library.

2. What is the difference between gsl_integration_qag() and gsl_integration_qng()?

The main difference between these two functions is the algorithm used for integration. gsl_integration_qag() uses an adaptive integration algorithm, while gsl_integration_qng() uses a non-adaptive algorithm. This means that gsl_integration_qag() is more accurate but can take longer to compute compared to gsl_integration_qng().

3. How do I pass my own function to the gsl integration function?

You can pass your own function to the gsl integration function by defining a C function that takes a double as an input and returns a double as an output. This function will represent the integrand. You can then pass this function as an argument to gsl_integration_qag() or gsl_integration_qng() along with the other necessary parameters.

4. Can I integrate a function with multiple variables using the gsl library?

Yes, the gsl library allows for integration of functions with multiple variables. You will need to use the gsl_integration_ntuple() function and provide the appropriate number of integration variables and limits.

5. How can I check the accuracy of the integration using the gsl library?

The gsl library provides an error estimate for the computed integral. You can use this error estimate to check the accuracy of the integration. Additionally, the library also allows for setting a desired absolute or relative error limit for the integration, which can help in achieving a more accurate result.

Similar threads

  • Programming and Computer Science
Replies
11
Views
1K
  • Programming and Computer Science
Replies
7
Views
1K
  • Programming and Computer Science
Replies
6
Views
916
  • Programming and Computer Science
Replies
21
Views
3K
  • Programming and Computer Science
Replies
1
Views
646
  • Programming and Computer Science
Replies
3
Views
1K
  • Programming and Computer Science
Replies
1
Views
1K
  • Programming and Computer Science
Replies
15
Views
5K
  • Programming and Computer Science
Replies
4
Views
1K
  • Programming and Computer Science
Replies
3
Views
3K
Back
Top