[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.