Re: Integration qpint [message #81060 is a reply to message #81059] |
Mon, 13 August 2012 02:00  |
Craig Markwardt
Messages: 1869 Registered: November 1996
|
Senior Member |
|
|
On Sunday, August 12, 2012 1:42:57 PM UTC-4, Gompie wrote:
> Hi
>
> Can anyone help. I am using qpint1d to integrate but i get
>
> "Program caused arithmetic error: Floating overflow:"
>
>
>
> The code is
>
> C1 = 1.1910427E-12
>
> c2=1.43883
>
> T=300
>
> c=c1*(c2^(-4))*T^4
>
> llimit=c2*2380/T
>
> ulimit=c2*2940/T
>
> print,"X=",llimit,ulimit
>
> lin= c*qpint1d('x^3/(EXP(X)-1)', llimit, +inf, /expr)
>
> rin=c*qpint1d('x^3/(EXP(X)-1)', ulimit, +inf, /expr)
>
> print,abs(lin-rin)
>
> end
>
>
>
> This outputs a number ( i am not sure if it is correct) but it also says floating point error.
Check the status variable. It's not an error. You can verify this yourself by
IDL reports *any* overflow or underflow that occurs. That's not necessarily bad.
What do you think happens when you try to evaluate EXP(+inf)? Or EXP(+3000)? Your user function generates an overflow, that's what happens. IDL is warning you that it is treating EXP(+3000) equivalent to +infinity. That is what happens when IDL cannot evaluate the full range of your user function with full precision. EXP(+3000) simply cannot be represented on a modern computer.
If you don't want overflow errors, then rewrite your expression to use EXP(-X) instead. Of course, then you will get underflow errors. Life sucks, doesn't it? :-)
If you really don't want overflow errors for this function, then you must avoid using +inf as your upper bound.
Craig Markwardt
|
|
|