Thread Subject:

Efficient integration between +/- infinity

I am trying to calculate the Voigt profile for absorption line-shapes.

The function requires calculating H(a,u) at each value of u (around a central absorption line):

H(a,u) = a/pi. Integral [exp(-y^2)/((u-y)^2+a^2)] .dy with y from -inf to +inf

I see that Matlab has a function for numerical integral - quadl - but how do I do this between +/- infinity?

And how do I do it efficiently?

Appreciate any help, thanks.

"Peter " <pgillies3@gmail.com> wrote in message <im2340$cne$1@fred.mathworks.com>...
> I am trying to calculate the Voigt profile for absorption line-shapes.
>
> The function requires calculating H(a,u) at each value of u (around a central absorption line):
>
> H(a,u) = a/pi. Integral [exp(-y^2)/((u-y)^2+a^2)] .dy with y from -inf to +inf
>
> I see that Matlab has a function for numerical integral - quadl - but how do I do this between +/- infinity?
>
> And how do I do it efficiently?

Change the variable to finite interval.

Bruno

"Peter " <pgillies3@gmail.com> wrote in message <im2340$cne$1@fred.mathworks.com>...
> I am trying to calculate the Voigt profile for absorption line-shapes.
>
> The function requires calculating H(a,u) at each value of u (around a central absorption line):
>
> H(a,u) = a/pi. Integral [exp(-y^2)/((u-y)^2+a^2)] .dy with y from -inf to +inf
>
> I see that Matlab has a function for numerical integral - quadl - but how do I do this between +/- infinity?
>
> And how do I do it efficiently?
>
> Appreciate any help, thanks.

As I recall, a Voigt profile is a convolution between a Gaussian
and a Lorentzian (also known as a Cauchy) distribution.

You should be able to decide how far out to look based on
the parameters of the two components. Beyond that point,
the contributions from each of the component distribution
PDFs will be essentially zero anyway.

John

Bruno & John - thanks for helping out.

I was hoping there was a way that Matlab would "figure it out" but I will run a few different tests to see what limits are appropriate.

"Peter " <pgillies3@gmail.com> wrote in message <im34es$t6$1@fred.mathworks.com>...
> Bruno & John - thanks for helping out.
>
> I was hoping there was a way that Matlab would "figure it out" but I will run a few different tests to see what limits are appropriate.
- - - - - - - - -
  I think Bruno was suggesting you make a change of variable that wouldn't go to infinity at the two integration limits. I would recommend trying

 x = 1/a*atan((u-y)/a)

which is

 y = u - a*tan(a*x)

The integration limits then become finite and the integrand remains finite there, in fact becoming zero, so numerical quadrature should work very well.

Roger Stafford

"Roger Stafford" wrote in message <im3qhq$rh8$1@fred.mathworks.com>...
> "Peter " <pgillies3@gmail.com> wrote in message <im34es$t6$1@fred.mathworks.com>...
> > Bruno & John - thanks for helping out.
> >
> > I was hoping there was a way that Matlab would "figure it out" but I will run a few different tests to see what limits are appropriate.
> - - - - - - - - -
> I think Bruno was suggesting you make a change of variable that wouldn't go to infinity at the two integration limits. I would recommend trying
>
> x = 1/a*atan((u-y)/a)
>
> which is
>
> y = u - a*tan(a*x)
>
> The integration limits then become finite and the integrand remains finite there, in fact becoming zero, so numerical quadrature should work very well.
>
> Roger Stafford
suppose, i want to find the E(x) and E(x^2) for continuous distribution like normal,gamma where the values of random variable ranges from (-inf to inf) and (0 to inf) respectively. how would i do it using matlab.. pls help.. i got stuck with it for many days.. any suggestion is appreciated...

On 22 Jul., 12:38, "simon " <sghosh...@gmail.com> wrote:
> "Roger Stafford" wrote in message <im3qhq$rh...@fred.mathworks.com>...
> > "Peter " <pgilli...@gmail.com> wrote in message <im34es$t...@fred.mathworks.com>...
> > > Bruno & John - thanks for helping out.
>
> > > I was hoping there was a way that Matlab would "figure it out" but I will run a few different tests to see what limits are appropriate.
> > - - - - - - - - -
> > I think Bruno was suggesting you make a change of variable that wouldn't go to infinity at the two integration limits. I would recommend trying
>
> > x = 1/a*atan((u-y)/a)
>
> > which is
>
> > y = u - a*tan(a*x)
>
> > The integration limits then become finite and the integrand remains finite there, in fact becoming zero, so numerical quadrature should work very well.
>
> > Roger Stafford
>
> suppose, i want to find the E(x) and E(x^2) for continuous distribution like normal,gamma where the values of random variable ranges from (-inf to inf) and (0 to inf) respectively. how would i do it using matlab.. pls help.. i got stuck with it for many days.. any suggestion is appreciated...- Zitierten Text ausblenden -
>
> - Zitierten Text anzeigen -


help quadgk

Best wishes
Torsten.