List Archives >
Maple User Group List Archive >
Archive by date >
This Month By Date >
This Month By Topic
[MUG] integral - correct for evalf(Int(f)), incorrect for int(f)
| [MUG] integral - correct for evalf(Int(f)), incorrect for int(f) |
|
Author: Robert J Gooding
Posted: 29/01/2001 17:53:54 GMT
|
>> From: "Robert J. Gooding"
The following is a pair of integrals that are clearly convergent
for y < -3, and one of them gives a complete elliptic integral,
as it should. However, the second fails analytically, but
is correct numerically - I am using 6.01.
I have not seen this specific bad result before, although I am becoming used
to being exceedingly cautious with elliptic integrals around singularities.
> restart:assume(y<-3):assume(qy,real):assume(qx,real);
Firstly, I do one integral analytically, and it works.
> (1/(2*Pi)^2)*int(int(1/(y+cos(qx)+cos(qy)+cos(qx)*cos(qy)),
qy=-Pi..Pi),qx=-Pi..Pi);
2*EllipticK(2/sqrt(-y+1))/pi/(y-1)
It also works for a particular value of y, and evaluates to the correct
expression.
> evalf((1/(2*Pi)^2)*int(int(1/(s-3.1+cos(qx)+cos(qy)+cos(qx)*cos(qy)),
qy=-Pi..Pi),qx=-Pi..Pi)); # NOTE the use of int and not Int
-0.8805
For the second integral, which by comparison to the first, must
also be convergent, is however not done properly by maple.
>(1/(2*Pi)^2)*int(int((cos(qx)*cos(qy))/(y+cos(qx)+cos(qy)+cos(qx)*cos(qy)),
qy=-Pi..Pi),qx=-Pi..Pi); # if you do not assume qx,qy real, maple will not
complete this integral.
infinity
WHICH IS INCORRECT.
Try this with a fixed value of y (for which the integral will be convergent).
>(1/(2*Pi)^2)*int(int((cos(qx)*cos(qy))/(-4+cos(qx)+cos(qy)+cos(qx)*cos(qy)),
qy=-Pi..Pi),qx=-Pi..Pi);
infinity
WHICH IS ALSO INCORRECT.
Now do the integration one independent variable at a time.
>
f:=unapply(int((cos(qx)*cos(qy))/(y+cos(qx)+cos(qy)+cos(qx)*cos(qy)),qx=-Pi..Pi)
,qy,y);
f(qy,y):= some function that I'm too lazy to
type out. there can be zeroes in the denominator
(see below)
To be sure that you're not going crazy, plot this function, say at y=-4,
and see that the integrand is a boring well-behaved function.
>plot(f(x,-4),x=-Pi..Pi);
NOW DO THE INTEGRATION ANALYTICALLY AND NUMERICALLY AND COMPARE.
>int(f(qy,-4),qy=-Pi..Pi);
infinity
which is again incorrect; but
>evalf(Int(f(qy,-4),qy=-Pi..Pi);
-1.486
which is presumably correct - it's at least reasonable!
So, this is not a 2D integral problem - it arises in the 1D integrals
as well.
One can examine the integrand of the 1D integral and discover that it's
denominator goes as (1+cos(qy)). However, the series expansion for this
integrand about qy=Pi is Pi/50*(qy-Pi)^2 + ..., which is very well
behaved.
So, like previous versions of maple, we still seem to be faced with
always tracking down every possible singularity in our integrands
before we can trust our results, although the numerical integrations
are sometimes OK.
--
Bob Gooding,
|
[View Complete Thread]
Previous by date: [MUG] logical operations in maple, Guy-Armand Kamendje
Next by date: [MUG] Re: Does the rand() function really give random numbers?, Maple Group
Previous thread: [MUG] step by step OR steps are too simple for Maple, Boris Alexeev
Next thread: [MUG] Does the rand() function really give random numbers?, Nathan Sokalski
|
Datasheet