Definite Integration tests

1. Definite Integration tests tests suggested by Dave Rusin 5/7/2002 sci.math.symbolic

3. The story A colleague teaching a graduate complex analysis course mentioned that one of his better students reported that he liked computing definite (real) integrals with the standard line-integral-and-residue techniques, but that, in the end, these were nothing that Mathematica could not do. So the instructor asked me for some challenge problems for the student to feed to Mma. I found some in old notes (taken from this newsgroup and elsewhere). Since software evolves, it is reasonable to expect that old challenge problems are challenges no longer; indeed some progress seems to be evident. (I tested with Mathematica version 4.1 and Maple version 7 on a Sun Solaris.) The case of definite integrals is particularly interesting to me since the great hype which accompanied early releases of Mathematica induced people to gleefully report its errors, many of which were indeed definite integrals computed in "closed form". So this became an afternoon project of tracking down some definite integrals which show software struggling. I am interested in (i) other examples of definite integrals which can be done by a grad student in complex analysis, but not by these packages; (Examples from the residue calculus particularly welcome.) (ii) more examples of definite integrals which the packages (still) do wrong; (iii) lessons about the science -- as opposed to the art -- of computing closed-form expressions for definite integrals; what's the algorithm? I am aware of the Risch "algorithm" for computing _indefinite_ integrals, so my interest is in definite integrals which can be computed only because the limits of integration are "special". I don't know, for example, how to determine whether the integrals of rational functions of coordinates over _closed loops_ on a Riemann surface can be expressed in some closed form. (For example, it was recently discussed in this group that the length of the curve y=x2(4-x) from x=0 to x=4 can be given in a "closed form" which is quite lengthy.) For all but the third, fourth, and fifth of the examples below, the antiderivatives actually _may_ be written using the primitives known to the software, but -- particularly for some important choices of the parameters a, b -- the definite integrals are clearly simpler than the antiderivatives in form. (I'm not wholly sure how to quantify "simpler".) Also, sometimes the software knows the general antiderivative but cannot seem to compute the definite integral correctly! Attached below are the integrals I considered, showing input and output for first Maple and then Mathematica. Here they are in TeX form in case you're writing a final exam right about now... $$\int_01 (1-x^a)^{1/a} \, dx$$ $$\int_{-\infty}^\infty {{e^{i a x}}\over{x2+b2}} \, dx$$ $$\int_0^\infty x^{-b} (\log(x))^a \log(1+x+x2) \, dx$$ $$\int_{-\pi}^\pi (1-a \cos(x)-b \cos(2 x))^{-1/2} \, dx$$ $$\int_0^{\pi/2} {1\over{1+(\tan(x))^a}} \, dx$$ $$\int_0^{2\pi} \log(a+\sin(x)) \, dx$$ $$\int_0^{\pi/2} {1\over{(a+\sqrt{\cos(x)})2}} \, dx$$ dave ================ int((1-x^a)^(1/a), x=0 .. 1); Definite integration: Can't determine if the integral is convergent. Need to know the sign of --> -1/a Will now try indefinite integration and then take limits. (1 - 2/a) 1/2 2 Pi GAMMA(1/a + 1) ------------------------------- (2/a - 1) GAMMA(1/a - 1/2) #Hmm: if it needs to know the sign of -1/a, shouldn't that be #in the answer somewhere? ................ Integrate[(1-x^a)^(1/a), {x, 0, 1}] 1 Sqrt[Pi] Gamma[1 + -] 1 a If[Re[-] > -1 && Re[a] > 0, ---------------------, a 2/a 1 1 2 Gamma[- + -] 2 a (otherwise no answer) ================ int(exp(I*a*x)/(x2+b2), x= -infinity .. infinity); # Gives the answer 0 ; but if I add the assumptions a>0 and b>0 I get exp(-a~ b~) Pi -------------- b~ ................ Integrate[Exp[I a x]/(x2+b2), {x, -Infinity, Infinity}] -2 2 Sqrt[b ] Pi (* otherwise no answer *) ================ Note: the integrals converge, and closed-form solutions exist, for int( x^(-b)*log(x)^a*log(1+x+x2), x=0 .. infinity); Integrate[ x^(-b) Log[x]^a Log[1+x+x2], {x,0,Infinity}] for all b in (1, 2) and all counting numbers a ; this was in a sci.math thread last November. (Hint: differentiation with respect to b allows us to assume a=0. The other logarithm term is ln(1-x3) - ln(1-x). ) But both systems leave these integrals unevaluated. So I restricted to the case b=3/2. First Maple: int( x^(-3/2)*log(x)^a*log(1+x+x2), x=0 .. infinity); infinity #Uh oh... It even gives this response when a = 2, and yet it decides evalf( Int( x^(-3/2)*log(x)2*log(1+x+x2), x=0 .. infinity) ); 208.8558862 ................ Mathematica declined to evaluate the integral for general a, and even for a=3 (and a=4 and a=5, and a=-1). On the other hand, it _was_ willing to supply an answer for a=0, 1, and 2. The last is Integrate[ x^(-3/2) Log[x]^2 Log[1+x+x2], {x,0,Infinity}] (* long, slowly-computed, complex-valued answer deleted. Simplify[%] shortens it; FullSimplify[%] gets it to one line, still complex-valued *) ComplexExpand[%] 2 3 8 Pi 16 Pi 16 Sqrt[3] Pi + ----- + --------- 3 3 Sqrt[3] ================ int(1/(1-a*cos(x)-b*cos(2*x))^(1/2),x=-Pi .. Pi); # Lengthy answer is slighly shortened by simplify(%); but not pretty! # It is expressed using lim_{x->0+} and lim_{x->0-} , and # runs about 16 lines using lprint . # Here is the simplified result of running the computation with a=1/2, b=2/5: / 1/2 1/2 1/2 1/2 | 2 (11 + 473 ) 8 5 |2 EllipticK(%1) - EllipticF(1/2 ---------------------, %1) | 1/2 1/2 \ (-5 + 473 ) 1/2 1/2 1/2 \ 2 (21 + 473 ) | / 1/2 1/2 - EllipticF(1/2 ---------------------, %1)| / (44 + 2 473 ) 1/2 1/2 | / (5 + 473 ) / 1/4 473 %1 := 8 --------------------------------- 1/2 1/2 1/2 1/2 (21 + 473 ) (11 + 473 ) ................ Integrate[1/(1 - a Cos[x] - b Cos[2 x])^(1/2), {x, -Pi, Pi}] -2 Power::infy: Infinite expression 0 encountered. -2 Power::infy: Infinite expression 0 encountered. Out[17]= Indeterminate (* The special case mentioned above triggers additional warnings before these, and then still the answer "Indeterminate". *) ================ int(1/(1+(tan(x))^a), x= 0.. Pi/2); Integrate[1/(1+(Tan[x])^a), {x, 0, Pi/2} ] No answer from either of them. (Hint: replace x by pi/2 - x .) (In a Jan. 1999 thread it was shown that MuPad can do this one without help, and the others can be coaxed into giving an answer.) Note that antiderivates _can_ be given for many values of a but that this is neither necessary nor particularly helpful, and as far as I can tell the antiderivative is not expressible in closed form for all a. ================ int(log(a+sin(x)), x=0..2*Pi); #The answer, after using simplify(%) , is 2 1/2 2 1/2 a + I (-a + 1) + I a - I (-a + 1) + I 2 Pi ln(a) - 2 Pi ln(----------------------) - 2 Pi ln(----------------------) 2 1/2 2 1/2 a + I (-a + 1) a - I (-a + 1) #But note that branches of the logarithms must be chosen carefully! #The answer must be real if a is real and greater than 1. In that case, #the integeral= 2 pi (cosh^{-1}(a) - log 2). This is zero if a=5/4, so try it: int(log(5/4+sin(x)), x=0..2*Pi); 2 2 Pi ln(5) - 2 Pi ln(1 + 2 I) - 2 Pi ln(2 + I) + I Pi #simplify(%) does give zero. ................ Integrate[Log[a+Sin[x]], {x, 0, 2*Pi}] 2 2 -I (Pi + 2 I Pi Log[a] - 2 I Pi Log[1 + I a - I Sqrt[-1 + a ]] - 2 2 I Pi Log[1 + I a + I Sqrt[-1 + a ]]) (* I was unable to coerce it present a real answer for real a. I tried the same special case as above: *) Integrate[Log[5/4+Sin[x]], {x, 0, 2*Pi}] (* A very long answer involving PolyLog, etc. has been deleted. With FullSimplify[%] we eventually get the (complex!) answer *) 2 4 I Pi (* The problem lies with the integration, not the simplification: the numerical value of the first long answer is also imaginary. *) ================ int(1/(a+sqrt(cos(x)))2, x=0..Pi/2); #Here is the answer in compressed form (obtained very quickly): 2*((a4-1)^(1/2)+arctan((a2+1)/((a2+1)*(a2-1))^(1/2))+a4*arctan( (a2+1)/((a2+1)*(a2-1))^(1/2)))/(a2+1)/(a2-1)/(a4-1)^(1/2)+a/(a2+1)* 2^(1/2)*EllipticK(1/2*2^(1/2))+2*a/(a2+1)/(a2-1)*2^(1/2)* EllipticK(1/2*2^(1/2))-2*a/(a2+1)/(a2-1)*2^(1/2)*EllipticE(1/2*2^(1/2)) -a*(a4+1)/(a2+1)/(a2-1)2*2^(1/2)*EllipticPi(-1/(a2-1),1/2*2^(1/2)) #This function has a singularity at a=1; the numerical value when #a=.99999 is about .5406924657909949240642166 - 24836532.758693 i, and #that imaginary part grows proportionally to 1/(1-a) for a < 1. #On the other hand, the integrals are all real and nicely behaved at a=1; #the real part of the previous expression seems to get the integral. #Maple will do the case a=1 separately if asked; numerically it's OK 1/2 1/2 1/2 -1 - 2 I EllipticF(1/2 2 , 2 ) + 2 I EllipticK(2 ) - 2/3 -------------------------------------------------------- 1/2 1/2 (2 + 2 ) (2 - 2) ................ Integrate[1/(a+Sqrt[Cos[x]])2, {x, 0, Pi/2}] (* Mathematica has a similar-looking answer, with the same behaviour as a -> 1- . However, taking the integral with a=1 separately gives a patently wrong answer: *) Integrate[1/(1+Sqrt[Cos[x]])2, {x, 0, Pi/2}] Pi 2 EllipticF[--, 2] 1 4 -(-) - ------------------ 3 3 (* Unfortunately Mathematica numerically evaluates this to -1.20735 (plus some imaginary fuzz) which is not reassuring, the integrand being positive. *) ................ Incidentally, in the case a=1, changing the exponent in the last integrand from a 2 to a 1 or a 3 makes Mathematica provide promising-looking answers using elliptic functions, but these evaluate numerically to a values which are again noticeably different from the numerical value of the integral. After some time Maple found the value for exponent 1 to be 1/2 1/2 1/2 1/2 -1 - I EllipticF(1/2 2 , 2 ) + I EllipticE(1/2 2 , 2 ) 1/2 1/2 + I EllipticK(2 ) - I EllipticE(2 ) I let Maple try the case with exponent 3 but gave up after about an hour.