Re: Integrate (a difficult integral)
- To: mathgroup at smc.vnet.net
- Subject: [mg74593] Re: Integrate (a difficult integral)
- From: "dimitris" <dimmechan at yahoo.com>
- Date: Tue, 27 Mar 2007 04:06:18 -0500 (EST)
- References: <esls78$q0v$1@smc.vnet.net><eu1rm1$l6j$1@smc.vnet.net>
>Since Mathematica return unevaluated both Integrate[F[z], {z, -
>Infinity, Infinity}] and
>Integrate[F[z], {z, -Infinity, 4 - I, Infinity}] we will use the
>previous setting with Limits.
I should have written...
"Since Mathematica return unevaluated both Integrate[f[z], {z, -
Infinity, Infinity}]
and Integrate[f[z], {z, -Infinity, 4 - I, Infinity}] we will use the
previous setting with Limits."
I apologize for any inconvience!
Dimitris
=CF/=C7 dimitris =DD=E3=F1=E1=F8=E5:
> As I said in another working with CAS don't mean that we... stop
> thinking with our mind.
>
> To obtain the best results (which are possible of course!)-not only in
> integration but in general doing serious Mathematics with Mathematica-
> you need an a combination of art, experience, reading, thinking (and
> guessing sometimes!).
>
> Let demonstrate one way of working for one integral that appeared in a
> recent thread in another forum.
> (The question was by Vladimir Bondarenko and I do consider that the
> integral should be appeared also in MathGroup.)
>
> Note that I don't say that it is the best one for the present
> situation but rather than it is a way that helped
> Mathematica to overcome its shothcomings in the present situation and
> return results in reasonable time (or someone would prefer to say: "Oh
> man, don't say many things; it's...just what you think! No more, no
> less!")
>
> Anyway...
>
> $Version
> "5.2 for Microsoft Windows (June 20, 2005)"
>
> Consider the function
>
> f[z_] := 1/Sqrt[z^3 - I*z^2 + z - I]
>
> The following integral returns unevaluated
>
> Integrate[F[z], {z, -Infinity, Infinity}]
>
> It is exactly this integral that we are looking for!
>
> The indefinite integral is
>
> F[z_] = Integrate[1/Sqrt[z^3 - I*z^2 + z - I], z]
> -(((1/4 + I/4)*(-I + z)*Sqrt[I + z]*(2*ArcTan[Sqrt[I + z]/(-2 + Sqrt[I
> + z])] - 2*ArcTan[Sqrt[I + z]/(2 + Sqrt[I + z])] +
> I*(Log[(2 + I) + z - 2*Sqrt[I + z]] - Log[(2 + I) + z + 2*Sqrt[I
> + z]])))/Sqrt[(-I + z)^2*(I + z)])
>
> Of course "blinded" application of the Nweton-Leibniz formula leads
> to incorrect result
>
> Limit[F[z], z -> Infinity] - Limit[F[z], z -> -Infinity]
> {N[%], NIntegrate[f[x], {x, -Infinity, Infinity}]}
>
> (-(1/2) - I/2)*Pi
> {-1.5707963267948966 - 1.5707963267948966*I, 3.1415926536059606 +
> 3=2E14159265360596*I}
>
> Here are some plots
>
> (Plot[#1[f[z]], {z, -10, 10}] & ) /@ {Re, Im};
> (Plot[#1[F[z]], {z, -20, 20}, Ticks -> {Range[-20, 20, 4], Automatic},
> Axes -> False, Frame -> True] & ) /@
> {Re, Im};
>
> (ContourPlot[#1[F[x + I*y]], {x, -5, 5}, {y, -5, 5}, PlotPoints ->
> 100, Contours -> 100, ContourShading -> False] & ) /@ {Re, Im};
> (ContourPlot[#1[f[x + I*y]], {x, -2, 2}, {y, -2, 2}, PlotPoints ->
> 100, Contours -> 100, ContourShading -> False] & ) /@ {Re, Im}
>
> As we see both the real and imaginary part has jump discontinuities in
> the real axis.
>
> The following plots are very hepful
>
> (Plot[#1[ArcTan[Sqrt[I + z]/(-2 + Sqrt[I + z])]], {z, -10, 10}] & ) /@
> {Re, Im};
> (Plot[#1[ArcTan[Sqrt[I + z]/(2 + Sqrt[I + z])]], {z, -10, 10}] & ) /@
> {Re, Im};
> (Plot[#1[Log[(2 + I) + z - 2*Sqrt[I + z]]], {z, -10, 10}] & ) /@ {Re,
> Im};
> (Plot[#1[Log[(2 + I) + z + 2*Sqrt[I + z]]], {z, -10, 10}] & ) /@ {Re,
> Im};
>
> These plots indicate that the first function introduce the jump
> discontinuities in the real axis
>
> Here is the branch point of this function
>
> Solve[-2 + Sqrt[I + z] == 0, z]
> {{z -> 4 - I}}
>
> Before obtaing Vladimir's integral let's see a simpler example
>
> g = HoldForm[Integrate[f[z], {z - 100, 100}]]
>
> The result returned by Mathematica is wrong as we see below
>
> FullSimplify[ReleaseHold[g]]
> {N[%], ReleaseHold[g /. Integrate[x___] :> NIntegrate[x, MaxRecursion
> -
> > 12]]}
>
> (1/2 + I/2)*Pi + Log[((101620289 - 320288*Sqrt[10001] +
> 8*I*Sqrt[-21072330472601 + 1017117472601*Sqrt[10001]])/100020001)^(-
> (1/4) + I/4)]
> {1.3704665093461683 + 1.370466509346168*I, 2.941262836141067 +
> 2=2E941262836141067*I}
>
> However using the UNDOCUNTATED (for Integrate; it is well documentated
> for NIntegrate in the Help Browser. BTW, I believe that its first
> appearance-which I saw
> myself of course!-is an old post by Ted Ersek) setting of the third
> element in the iterator
> of Integrate as follows
>
> forces Mathematica to consider the above mentioned discontinuity,
> giving in this way the correct result!
>
> FullSimplify[Integrate[f[z], {z, -100, 4 - I, 100}]]
> N[%]
>
> (1/4 + I/4)*(4*Pi + Log[((101620289 - 320288*Sqrt[10001] +
> 8*I*Sqrt[-21072330472601 + 1017117472601*Sqrt[10001]])/100020001)^I])
> 2=2E941262836141065 + 2.941262836141065*I
>
>
> Previous setting is equivalent to
>
> FullSimplify[Limit[F[z], z -> 100, Direction -> 1] - Limit[F[z], z ->
> 4 - I, Direction -> -1] +
> Limit[F[z], z -> 4 - I, Direction -> 1] - Limit[F[z], z -> -100,
> Direction -> -1]]
> N[%]
>
> (1/4 + I/4)*(4*Pi + Log[((101620289 - 320288*Sqrt[10001] +
> 8*I*Sqrt[-21072330472601 + 1017117472601*Sqrt[10001]])/100020001)^I])
> 2=2E941262836141065 + 2.941262836141065*I
>
> Since Mathematica return unevaluated both Integrate[F[z], {z, -
> Infinity, Infinity}] and
> Integrate[F[z], {z, -Infinity, 4 - I, Infinity}] we will use the
> previous setting with Limits
> for our needs.
>
> We have
>
> Limit[F[z], z -> Infinity] (*correct*)
> 0
>
> Limit[F[z], z -> -Infinity] (*incorrect; it misses a minus sign*)
> (1/2 + I/2)*Pi
>
> (Limit[F[z], z -> 4 - I, Direction -> #1] & ) /@ {-1, 1} (*correct*)
> {(-(1/4) - I/4)*(Pi - 2*ArcTan[1/2] - I*Log[5]), (1/4 + I/4)*(Pi +
> 2*ArcTan[1/2] + I*Log[5])}
>
> So, (at last!) we finally have
>
> (0 - (-(1/4) - I/4)*(Pi - 2*ArcTan[1/2] - I*Log[5])) + ((1/4 + I/
> 4)*(Pi + 2*ArcTan[1/2] + I*Log[5]) - (-(1/2 + I/2))*Pi)
> Simplify[%]
> {(N[#1, 22] & )[%], NIntegrate[f[x], {x, -Infinity, Infinity},
> WorkingPrecision -> 50, PrecisionGoal -> 20]}
>
> (1/2 + I/2)*Pi + (1/4 + I/4)*(Pi - 2*ArcTan[1/2] - I*Log[5]) + (1/4 +
> I/4)*(Pi + 2*ArcTan[1/2] + I*Log[5])
>
> (1 + I)*Pi
>
> {3.1415926535897932384626433832795028842`22. +
> 3=2E1415926535897932384626433832795028842`22.*I,
> 3.14159265358979323846264338327950287583`21.64106815259695 +
> 3=2E14159265358979323846264303665222933418`21.64106815259695*I}
>
> I=2Ee the desired integral is equal to (1 + I)*Pi.
>
>
> Best Regards
> Dimitris
>
>
>
>
>
>
>
> =CF/=C7 Michael Weyrauch =DD=E3=F1=E1=F8=E5:
> > Hello,
> >
> > Dimitris, this is a marvelous solution to my problem. I really apprec=
ia=
> te
> > your help. I will now see if I can solve all my other (similar) integra=
ls=
> using the same trick.
> > Timing is not really the big issue if I get results in a reasonable amo=
un=
> t of time.
> >
> > Also the references you cited are quite interesting, because they give =
so=
> me insight
> > what might go on inside Mathematica concerning integration.
> >
> > I also agree with you that it is my task as a programmer to "help" Math=
em=
> atica
> > and guide it into the direction I want it to go. Sometimes this is an "=
ar=
> t"...
> >
> > Nevertheless, in the present situation I do not really understand why M=
at=
> hematica wants
> > me to do that rather trivial variable transformation, which is at the h=
ea=
> rt of your solution.
> > The integrand is still a rather complicated rational function of the sa=
me=
> order. The form
> > of the integrand did not really change
> > substantially as it is the case with some other ingenious substitutions=
o=
> ne uses in order to
> > do some complicated looking integrals "by hand".
> >
> > I think the fact that we are forced to such tricks shows that the Mathe=
ma=
> tica integrator
> > is still a bit "immature" in special cases, as also the very interestin=
g =
> article by D. Lichtblau,
> > which you cite, seems to indicate. So all this is probably perfectly kn=
ow=
> n to the
> > Mathematica devellopers. And I hope the next Mathematica version has al=
l =
> this "ironed out"??
> >
> > Many thanks again, Michael