7.73 Bug in int in Maple 6 (28.6.00)

7.73.1 Martin Hirschmanner
7.73.2 Robert Israel (29.6.00)
7.73.3 Helmut Kahovec (18.7.00)
7.73.4 Willard, Daniel Dr (21.7.00)

7.73.1 Martin Hirschmanner

There seems to be a bug in Maple 6. When I do the following computation I get different results for the numerical Intergration. In the first case I took as limits q=0..2*Pi and in the second q=0..2*evalf(Pi). In the first case the result is wrong. In Maple V both computations give the same result.

Is this a bug ?

> restart: 
> with(linalg): 
> h:=matrix(9,1,[a1*cos(q),a1*sin(q),q,r*cos(q)+a2*sqrt(1-lambda^2*(sin(q))^2), 
                 (r-lambda*a2)*sin(q),-arcsin(lambda*sin(q)), 
                 r*cos(q)+l*sqrt(1-lambda^2*(sin(q))^2),0,0]); 
> hst:=map(diff,h,q); 
> M:=diag(m1,m1,I1,m2,m2,I2,m3,m3,I3); 
> Theta:=evalm(transpose(hst)&* M &* hst); 
> l:=0.25; 
> r:=0.03; 
> lambda:=r/l; 
> a1:=21.877e-3; 
> m1:=50.32; 
> I1:=1.0267; 
> a2:=0.034345; 
> m2:=16.645; 
> I2:=0.23251; 
> m3:=44.55; 
 
> Theta_m:=evalf(int(1/2/Pi*Theta[1,1],q=0..2*Pi)); 
 
                        Theta_m := 1.058110253 
 
> Theta_m:=evalf(int(1/2/evalf(Pi)*Theta[1,1],q=0..2*evalf(Pi))); 
 
                        Theta_m := 1.085648004
 

7.73.2 Robert Israel (29.6.00)

es, this is a bug. But it’s not numerical integration, it’s symbolic integration. If you want numerical integration, use evalf(Int(...)) not evalf(int(...)).

Here is this bug in a simpler form:

> h:= (1-cos(t)/sqrt(1-a*sin(t)^2))^2; 
> int(expand(h),t) - int(h,t); 
 
                        t
 

The integral of the expanded form of \(h\) would be correct (at least when \(a\) is small and positive). Somehow it’s as if int(h,t) expands \(h\) but forgets about the 1. In your example, if you use expand(Theta[1,1]) instead of Theta[1,1] you should get the correct answer 1.085648004 both times.

7.73.3 Helmut Kahovec (18.7.00)

All I found out is that if the (trivial) changes described below are made in `int/indefinite`() and `int/definite`() then the example input will give the correct results.

Unfortunately, that statement is not correct, it is a bug, really. I am sorry and congratulate Robert Israel!

Well, I did not at all expect the symbolic integration to be wrong. Note that releases 4 and 5 of Maple V cannot evaluate that integral symbolically and thus numerically integrate definite integrals like

> lambda:=1/2: 
> j:=(1+lambda^2*cos(q)/(sqrt(1-lambda^2*sin(q)^2)))^2; 
> evalf(int(j,q=0..2*Pi));
 

maple6, however, can indeed integrate

> int(j,q);
 

and

> int(j,q=0..2*Pi);
 

symbolically but, as Robert has shown, will give the wrong result. `int/indefinite`() and `int/definite`() should expand their first argument in this case. I cannot pinpoint the bug more closely than that. If we make the following trivial changes to `int/indefinite`() and `int/definite`()

`int/indefinite` := proc(FF, zz) 
local ff, opts, answer; 
global StandardFunctions; 
   1   ff := expand(FF); 
   2   if nargs = 2 then 
   ... 
  19   if answer = FAIL then 
  20     answer := ('int')(ff,zz) 
       end if; 
  21   answer 
end proc 
 
`int/definite` := proc(FF, zz, aa, bb) 
local ff, answer, opts; 
   1   ff := expand(FF); 
   2   if nargs = 5 then 
   ... 
  33   if answer = FAIL then 
  34     answer := ('int')( 
           ff,zz = aa .. bb,op(subs(('formula') = NULL,opts)) 
         ) 
       end if; 
  35   RETURN(answer) 
end proc
 

then the following integrals will give the correct results:

> int(j,q); 
 
  q + arcsin(1/2 sin(q)) 
 
         - 1/8 sqrt(3) arctan(1/3 (2 tan(1/2 q) - 1) sqrt(3)) 
 
         - 1/8 sqrt(3) arctan(1/3 (2 tan(1/2 q) + 1) sqrt(3)) 
 
         + 1/2 arctan(tan(1/2 q)) 
 
> int(j,q=0..2*Pi); 
 
                       5/2 Pi - 1/4 sqrt(3) Pi
 

Note, too, that maple6 integrates

> int(j,q=0..2*evalf(Pi));
 

numerically right from the beginning, cf. lines 1-4 of the procedure int().

7.73.4 Willard, Daniel Dr (21.7.00)

Please note that V 5.1 will integrate it symbolically with the answer given in terms of elliptic integrals.