I am really horrified! If you integrate elementary function 1/(9x^2+6x+2), you obtain the
wrong result arctan(1+3x) instead of the correct 1/3*arctan(1+3x)! It is only for Maple6.
Maple5 returns the correct result.
It is corrected with Maple 7. (U. Klein)
This bug results from trying to program a great variety of simple cases rather
than just using a handful of algorithms for the general cases. In this problem, the
special case that is handled in procedure `int/ratpoly/ratpoly` is when the
denominator is exactly one more than a perfect square of first-degree binomial. In other
words,
> int(c/((a*x+b)^2+1), x)
will give the wrong answer. Change that 1 to anything else and you’ll get the right answer. This case is so special – why do they even bother with it? It’s hard to imagine that it is worth the bother to treat it separately.
The general case algorithms need to be there anyway. Does it really save that much time to sort out and handle the simple cases? When there are so many special cases, it makes it that much more likely that a trivial bug will go undetected for a long time. Does Maple have an official programming philosophy?
In this case, the bug is truly trivial – essentially a typo – and easy for the user to correct.
The fact that someone did not notice this bug sooner is a sign that it is not worth the
trouble to treat this case separately. I’ve traced it to lines 21 and 23 of procedure
`int/ratpoly/ratpoly`. I show the fix after the code below.
> restart; > showstat(`int/ratpoly/ratpoly`); `int/ratpoly/ratpoly` := proc(f) local ans, den, num, d, n, k, const, q, r, rest, g, const2, A, B, answer; 1 if type(f,polynom(anything,_X)) then 2 RETURN(`int/polynom`(args)) end if; 3 if type(f,`+`) then 4 return map(procname,f) end if; 5 num := 1; 6 const := 1; 7 if type(f,`*`) then 8 rest, const := selectremove(has,f,_X); 9 num, den := numer(rest), denom(rest) else 10 const, num, den := 1, 1, denom(f) end if; 11 den, n := op(`if`(den::`^`,den,[den, 1])); 12 d := degree(den,_X); 13 answer := FAIL; 14 if d = 1 then 15 const := const*num/coeff(den,_X,1); 16 if n = 1 then 17 answer := const*ln(den) else 18 answer := -const/(n-1)/(den^(n-1)) end if elif n = 1 and degree(num,_X) = d-1 and rem(num,diff(den,_X),_X,q) = 0 then 19 answer := const*q*ln(den) elif n = 1 and d = 2+2*degree(num,_X) and type(den,polynom(polynom(rational),_X)) and psqrt(den-1) <> _NOSQRT and rem(diff(den-1,_X),num,_X) = 0 then 20 r := psqrt(den-1); 21 answer := const/content(r,_X)*arctan(r) ^^^^^^^^^^^^^ The bug is that they divide by the common factor (the "content") of the binomial rather than by its leading coefficient. elif n = 1 and d = 2+2*degree(num,_X) and type(den,polynom(polynom(rational),_X)) and psqrt(den+1) <> _NOSQRT and rem(diff(den+1,_X),num,_X) = 0 then 22 r := psqrt(den+1); 23 answer := const/content(r,_X)*arctanh(r) ^^^^^^^^^^^^^
Same bug here.
elif d = 2 then 24 answer := const*`int/ratpoly/quadratic`(num,den,n-1) elif d = 3 then 25 r := `int/ratpoly/cubic`(num,den,n); 26 answer := `if`(r = FAIL,FAIL,const*r) elif 0 < ldegree(num,_X) and `int/ratpoly/subs`(num,den,'k') then 27 num := collect(num/(_X^(k-1)),_X); 28 r := `int/indef1`(subs(_X = _X^(1/k),num/(den^n))); 29 answer := 1/k*const*subs(_X = _X^k,r) elif irem(d,2,'q') = 0 and degree(num,_X) = q-1 and `int/ratpoly/arctan`(num,den^(n+1),'q') then 30 answer := const*q elif d = 4 then 31 r := `int/ratpoly/quartic`(num,den,n); 32 answer := `if`(r = FAIL,FAIL,const*r) elif d = 5 then 33 r := `int/ratpoly/quintic`(num,den,n-1); 34 answer := `if`(r = FAIL,FAIL,const*r) elif d = 6 then 35 r := `int/ratpoly/sextic`(num,den,n-1); 36 answer := `if`(r = FAIL,FAIL,const*r) end if; 37 if answer = FAIL then 38 if type(f,'ratpoly(rational,_X)') then 39 const*`int/risch/ratpoly`(num/(den^n),_X) else 40 g, const2, A, B := `int/ratpoly/horowitz`(num/(den^n),_X); 41 if A = 0 then 42 const*g else 43 const*(g+const2*sum(subs(_X = _R,A/diff(B,_X))*ln(_X-_R),_R = RootOf(B,_X))) end if end if else 44 answer end if end proc
An example of the bug in action:
> int(c/((a*x+b)^2+1), x); c arctan(a x + b)
Before you test my fix, make sure to do a restart to clear the erroneous answers from the remember tables.
> restart;
Since the "content" command is used only on the erroneous lines, you can easily correct this bug:
> `int/ratpoly/ratpoly`:= subs(content= lcoeff, eval(`int/ratpoly/ratpoly`)):
You could also use "diff" instead of "lcoeff". That would be more natural.
Test it:
> int(c/((a*x+b)^2+1), x); c arctan(a x + b) ----------------- a > int(1/(9*x^2+6*x+2), x); 1/3 arctan(1 + 3 x)
In release 5, these special cases are not treated. They are merely passed on to
`int/ratpoly/quadratic` as in line 24 above. Why were these cases singled out in release
6?
I am not absolutely sure that my fix works in all cases. Be wary when doing integrals with a quadratic denominator.