Message Boards Message Boards

0
|
5296 Views
|
6 Replies
|
0 Total Likes
View groups...
Share
Share this post:
GROUPS:

Integrate gives wrong result

Posted 9 years ago

When I try to analytically solve the Integral

F = Integrate[Sqrt[(x + Cos[?])^2 + (y + Sin[?])^2], {?, 0,2 ?}]

Mathematica gives the solution

F = 0

(so does Maple and other programs).

Since Sin and Cos only give values between -1 and +1 one can clearly see that this can not be true if x or y are high or low numbers. I tried to find a fitting curve and found that

F = Power[(2 ?)^E + (2 ? Sqrt[x^2 + y^2])^E, (E)^-1]

I have no mathematical proof for that, but it seems elegant enough to be true. Of course I know that Integrate has some "possible issues" and I always recommend to test the result with NIntegrate. But anyway, since I seem to have found a better analytical solution than zero I thought it might help if I share it. Maybe this can be fixed in future versions of Mathematica.

POSTED BY: tyran simon
6 Replies

As pointed out above the phase in the integral is of no importance because a and b in the integrand Sqrt [ a + b Cos [ w ] ] are constant and the Integration is around the whole circle from 0 to 2 Pi. Therefore we can write

In[1]:= j1a =  Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[phi - t]] /. t -> 0

Out[1]= Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[phi]]

Avoiding lots of Messages by

In[2]:= SetOptions[Integrate, GenerateConditions -> False]

Out[2]= {Assumptions :> $Assumptions, GenerateConditions -> False, 
 PrincipalValue -> False}

we get

In[3]:= res = Integrate[j1a, {phi, 0, 2 Pi}]

Out[3]= 4 Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2]]
  EllipticE[(4 Sqrt[x^2 + y^2])/(1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])]

and as above

In[4]:= res /. x -> 0 /. y -> 0

Out[4]= 2 \[Pi]

In[5]:= res /. x -> 5. /. y -> 3.8

Out[5]= 39.7097
POSTED BY: Hans Dolhaine

Hello,

first of all, the Integrand as a (positive) square-root of a sum of squares is greater or equal Zero. As is pointed out later, I think it can be Zero at only one Point, and that means, you integrate over a positive Integrand, so the integral should be positive and non-Zero.

Look at it from a geometrical Point of view. Your Integrand Looks like the norm of a vector, in fact the norm of the sum of two vectors,

XX = {x, y};
VV = {Cos[phi], Sin[phi]};

And this is a circle of radius 1 centered a the Point { - x, - y }. So your Integrand is in fact the length of a vector from the origin to a Point of that circle. That length is always positive, only if the circle passes through the origin it can be Zero, at exactly one Point. So, the Integrand is (mostly) positive ant the integral non-Zero.

Let's check it. Define the norm or length of a vector by

norm[a_] := Sqrt[a.a]

and we get

 In[11]:= norm[XX + VV]

    Out[11]= Sqrt[(x + Cos[phi])^2 + (y + Sin[phi])^2]

exactly your Integrand.

In other words you integrate over the norm[ XX + YY ].

The vector nearest to the origin is XX = { 0, 0 } and as norm[ YY ] is 1 the integral should then become 2 Pi.

Taking into account the inequality

norm[a + b ] <= norm[a] + norm[ b] 

and norm [ a ] = constant we arrive at (at least in my opinion)

2 \[Pi] <=  integral <= 2 \[Pi] ( 1 + norm [XX] )

Now let us rotate the Frame of reference. Define a Rotation Matrix dd

dd = ( {
    {Cos[t], Sin[t]},
    {-Sin[t], Cos[t]}
   } );

and apply it to the vector in question. Note that this is an orthogonal Transformation and lengths are not changed

  In[21]:= vecnew = dd.(XX + VV) // Expand // TrigReduce

  Out[21]= {Cos[phi - t] + x Cos[t] + y Sin[t], 
   y Cos[t] + Sin[phi - t] - x Sin[t]}

We Chose the angle t in a way, that the Point { x, y } is rotated to the x-axes of the new System

In[22]:= vecnew1 = 
 vecnew /. Cos[t] -> x/Sqrt[x^2 + y^2] /. Sin[t] -> y/Sqrt[x^2 + y^2]

Out[22]= {x^2/Sqrt[x^2 + y^2] + y^2/Sqrt[x^2 + y^2] + Cos[phi - t], 
 Sin[phi - t]}

Note that this is the same Situation as before, only the Point in question is tranferred to the x-axes of a new coordinate - System. And now we have

   j1 = norm[vecnew1] // FullSimplify

    Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[phi - t]]

But this is of the form Sqrt [ a + b Cos [ w ] ], as the Phase shift doesn't matter when we integrate over phi from 0 to 2 Pi.

Indeed we get

In[28]:= J1 = \[Integral]j1 \[DifferentialD]phi

Out[28]= (2 Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[phi - t]]
  EllipticE[(phi - t)/2, (4 Sqrt[x^2 + y^2])/(
  1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])])/Sqrt[(1 + x^2 + y^2 + 
 2 Sqrt[x^2 + y^2] Cos[phi - t])/(1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])]

(hmm, interstingly the rotation angle is maintained in the Elliptical thing...) and

In[30]:= T2 = J1 /. phi -> 2 Pi

Out[30]= (2 Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[t]]
  EllipticE[1/2 (2 \[Pi] - t), (4 Sqrt[x^2 + y^2])/(
  1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])])/Sqrt[(1 + x^2 + y^2 + 
 2 Sqrt[x^2 + y^2] Cos[t])/(1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])]

In[31]:= T1 = J1 /. phi -> 0

Out[31]= -((
 2 Sqrt[1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[t]]
   EllipticE[t/2, (4 Sqrt[x^2 + y^2])/(
   1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])])/Sqrt[((
 1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[t])/(
 1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2]))])

and finally

In[34]:= res = T2 - T1 // FullSimplify

Out[34]= (2 Sqrt[
   1 + x^2 + y^2 + 
    2 Sqrt[x^2 + y^2]
      Cos[t]] (EllipticE[\[Pi] - t/2, (4 Sqrt[x^2 + y^2])/(
      1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])] + 
     EllipticE[t/2, (4 Sqrt[x^2 + y^2])/(
      1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])]))/(Sqrt[(
  1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2] Cos[t])/(
  1 + x^2 + y^2 + 2 Sqrt[x^2 + y^2])])

This should be the result of the Integration. If Mathematica says 0 that is wrong indeed.

Two final checks

In[35]:= res /. x -> 0 /. y -> 0

Out[35]= 2 \[Pi]

ok, And

In[40]:= res /. t -> ArcSin[y/Sqrt[x^2 + y^2]] /. x -> 5. /. 
  y -> 3.8 // N

Out[40]= 39.7097

compared to

In[41]:= NIntegrate[((Sqrt[(x - Cos[a])^2 + (y - Sin[a])^2]) /. 
    x -> 5. /. y -> 3.8), {a, 0, 2 Pi}]

Out[41]= 39.7097

Regards Hans

POSTED BY: Hans Dolhaine

If you first expand and simplify the argument of the square root, the integration works:

In[1]:= 
Simplify[Expand[(x + Cos[\[CurlyPhi]])^2 + (y + Sin[\[CurlyPhi]])^2]]

Out[1]= 1 + x^2 + y^2 + 2 x Cos[\[CurlyPhi]] + 2 y Sin[\[CurlyPhi]]

In[2]:= Integrate[
 Sqrt[1 + x^2 + y^2 + 2 x Cos[\[CurlyPhi]] + 
   2 y Sin[\[CurlyPhi]]], {\[CurlyPhi], 0, 2 \[Pi]}]

Out[2]= ConditionalExpression[
 4 Sqrt[1 + x^2 + y^2 - 2 Sqrt[x^2 + y^2]]
   EllipticE[-((4 Sqrt[x^2 + y^2])/(
    1 + x^2 + y^2 - 2 Sqrt[x^2 + y^2]))], 
 2 Sqrt[x^2 + y^2] < 1 + x^2 + y^2 && (1 + x^2 + y^2)/Sqrt[
   x^2 + y^2] >= 2]
POSTED BY: S M Blinder

The difference between F and G is of the order of .2 or less. It is too small to perceive it in the two plots you have made. Try plotting the difference:

With[{n = E, y = 0}, 
 Plot[NIntegrate[
    Sqrt[(x + Cos[\[CurlyPhi]])^2 + (y + 
         Sin[\[CurlyPhi]])^2], {\[CurlyPhi], 0, 2 \[Pi]}] - 
   Power[(2 \[Pi])^E + (2 \[Pi] Sqrt[x^2 + y^2])^E, (E)^-1], {x, 0, 
   2 Pi}]]
POSTED BY: Gianluca Gorni

I tested your formula for F numerically, and it doesn't seem true:

With[{x = 1, 
   y = 2}, {Integrate[
    Sqrt[(x + Cos[\[CurlyPhi]])^2 + (y + 
         Sin[\[CurlyPhi]])^2], {\[CurlyPhi], 0, 2 \[Pi]}],
   Power[(2 \[Pi])^E + (2 \[Pi] Sqrt[x^2 + y^2])^E, (E)^-1]}] // N
With[{y = 2}, 
 Plot[NIntegrate[
    Sqrt[(x + Cos[\[CurlyPhi]])^2 + (y + 
         Sin[\[CurlyPhi]])^2], {\[CurlyPhi], 0, 2 \[Pi]}] -
   Power[(2 \[Pi])^E + (2 \[Pi] Sqrt[x^2 + y^2])^E, (E)^-1], {x, -20, 
   20}]]
POSTED BY: Gianluca Gorni
Posted 9 years ago

On the Plot it fits:

x,y

r

POSTED BY: tyran simon
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract