bug-gsl
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [Bug-gsl] bug in GSL numerical integration


From: Pedro Gonnet
Subject: Re: [Bug-gsl] bug in GSL numerical integration
Date: Mon, 12 Aug 2013 16:37:04 +0100

Hi Nikhil,

The problem is that your shifted function, once transformed as described
in [1], just looks like a flat line with a spike in it, and I'm guessing
that the nodes of the quadrature rule just miss it.

You can probably make it work by increasing the relative tolerance.

Cheers,
Pedro

[1] 
http://www.gnu.org/software/gsl/manual/html_node/QAGI-adaptive-integration-on-infinite-intervals.html


On Mon, 2013-08-12 at 14:14 +0200, Nikhil Jayant Joshi wrote:
> Hi,
> 
> I am doing integration over (-infty, +infty) using the gsl_integration_qagi
> routine. My expectation is that upon translation along the x-axis the
> result of integration (area under the curve) should not change. However
> that is not, what I observer. Am I making a mistake somewhere? The code is
> attached below:
> 
> The variable offset basically creates the translation. The area remains
> same for offset values of 0, 10.0, 20.0, but suddenly drop to zero after
> offset ~ 40.0
> 
> 
> double offset=200.0;
> 
> 
> double f (double x, void * params) {
> 
>     double alpha = *(double *) params;
> 
>     x += offset;
> 
>     double f = exp(-x*x);
> 
>     return f;
> 
> }
> 
> 
> int main(int argc, const char * argv[])
> 
> {
> 
> 
> 
>     gsl_integration_workspace * w
> 
>     = gsl_integration_workspace_alloc (1000);
> 
> 
> 
>     double result, error;
> 
>     double expected = -4.0;
> 
>     double alpha = 1.0;
> 
> 
> 
>     gsl_function F;
> 
>     F.function = &f;
> 
>     F.params = α
> 
> 
> 
>     gsl_integration_qagi (&F, 0, 0.001, 1000,
> 
>                           w, &result, &error);
> 
> 
> 
>     printf ("result = % .18f\n", result);
> 
> 
> 
>     return 0;
> 
> }
> 
> Thanks in anticipation,
> Nikhil
> --
> My problem:
> I have my brain divided into two parts: The left part has nothing right in
> it, while the right part has nothing left in it.





reply via email to

[Prev in Thread] Current Thread [Next in Thread]