Exact solution for pendulum

1. Jan 25, 2006

gulsen

Yup. How can I solve:
$$d^2\theta / dt^2 = -(g/L)sin(\theta)$$

2. Jan 25, 2006

arildno

There doesn't exist a nice, exact solution of this diff. eq for general initial values.

Last edited: Jan 25, 2006
3. Jan 25, 2006

benorin

Look here

The definite integral that gives so little as the period of the pendulum involves a complete elliptic integral of the first kind! and indeed this is not nice either.

4. Jan 25, 2006

saltydog

I disagree. It is very nice. Beautiful in fact. I'll start it for you:

We obtain the first integral:

$$\frac{1}{2}\left(\frac{d\theta}{dt}\right)^2-\frac{g}{L}Cos(\theta)=C$$

You know, integrate the second order ODE directly to get this. Try it.

If $\theta=\omega$ is the maximum displacement of the pendulum from its equlibrium position, then $\theta^{'}=0$ for this value and we can evaluate C, which is:

$$C=-\frac{gCos(\omega)}{L}$$

Now we solve for $\frac{d\theta}{dt}$ and obtain:

$$\frac{d\theta}{dt}=\sqrt{\frac{2g}{L}}\sqrt{Cos(\theta)-Cos(\omega)}$$

or:

$$dt=\sqrt{\frac{L}{2g}}\frac{d\theta}{\sqrt{Cos(\theta)-Cos(\omega)}}$$

Letting:

$$Cos(\theta)=1-2k^2 Sin^2(\phi)$$

$$k=Sin(\frac{\omega}{2})$$

we can then reduce the last dif. eq. to an elliptic integral and solve it accordingly . . .

. . . equal rights for special functions . . .

Last edited: Jan 25, 2006
5. Jan 26, 2006

arildno

Equal rights for special functions?
As long as you don't count in the Haenkel function, I'll accept that.

6. Jan 26, 2006

dextercioby

The solution angle(t) is proportional to sinus amplitudinis, one of Jacobi's elliptical functions.

Daniel.

7. Jan 26, 2006

saltydog

Well I' d like to finish this since I think the analysis is as I saidit was above.

First some clarification:

1. I obtained the first integral by writing the ODE as:

$$d\left[\theta^{'}\right]=-\frac{g}{L}Sin[\theta]$$

multiplying both sides by $\theta^{'}[/tex] and writing it as: $$\theta^{'}d\left[\theta^{'}\right]=-\frac{g}{L}Sin[\theta]d\theta$$ The variables have been separated and can now be directly integrated to give the first integral (pretty sure that's poke-a-poke, kinda' rough for me too). 2. I glossed over the fact that I chose the initial velocity to be zero. This makes the algebra easier although we can't model round-and-round motion with this solution. I left the problem above as: $$dt=\sqrt{\frac{L}{2g}}\frac{d\theta}{\sqrt{Cos(\theta)-Cos(\omega)}}$$ and seek a transformation that will convert the RHS to an elliptic integral: $$F[\phi,k]=\int_0^{\phi} \frac{dv}{\sqrt{1-k^2 Sin^2(v)}}$$ Well, it turns out that if: $$Cos(\theta)=1-2k^2Sin^2(\phi)\quad \text{and}\quad k=Sin(\omega/2)$$ then the following can be shown to be true: $$Cos(\theta)-Cos(\omega)=2k^2Cos^2(\phi)$$ $$Sin(\theta)=2k Sin(\phi)\sqrt{1-k^2Sin^2(\phi)}$$ $$Sin(\theta)d\theta=4k^2 Sin(\phi) Cos(\phi)d\phi$$ Making this substitution, which I think is quite a challenge in itself, we obtain: $$dt=\sqrt{\frac{L}{g}}\frac{d\phi}{\sqrt{1-k^2 Sin^2(\phi)}}$$ Integrating from t=0 to some time t: $$t=\sqrt{\frac{L}{g}}\int_{0}^{\phi}\frac{d\phi}{\sqrt{1-k^2 Sin^2(\phi)}}$$ Now, this last integral is an elliptic integral. It has an inverse. We can then take the inverse elliptic integral of both sides in order to isolate [itex]\phi$. A special function called the JacobiSN function expresses the sin of the inverse. That is if:

$$F[\phi,k]=\int_0^{\phi} \frac{dv}{\sqrt{1-k^2 Sin^2(v)}}$$

then:

$$\text{JacobiSN}\left\{F[\phi,k]\right\}=Sin(\phi)$$

Doing this to the equation above yields:

$$\text{JacobiSN}\left(t\sqrt{\frac{g}{L}},k\right)=Sin(\phi)=\frac{1}{k}Sin(\theta/2)$$

Isolating theta:

$$\theta(t)=2\text{ArcSin}\left[k\text{JacobiSN}\left(t\sqrt{g/L},k\right)\right],\;k=sin(\omega/2)$$

See, told you it was beautiful

A sample plot of this function (for a 3 foot pendulum raised 1 radian) is attached.

Attached Files:

• non-linear pendulum.JPG
File size:
8.6 KB
Views:
309
Last edited: Jan 26, 2006
8. Feb 4, 2006

Hootenanny

Staff Emeritus
Saltydog,

$$\theta(t)=2\text{ArcSin}\left[k\text{JacobiSN}\left(t\sqrt{g/L},k\right)\right],\;k=sin(\omega/2)$$
but I was wondering if there was anyway you could find an exact solution for the period of a pendulum when $\sin\theta \not\approx \theta$

Thanks

9. Feb 4, 2006

HallsofIvy

Staff Emeritus
No there is not exact formula- just as there is no exact formula for circumference of an ellipse.

10. Feb 4, 2006

Hootenanny

Staff Emeritus
So how would I solve it? Could I do an expansion?

11. Feb 4, 2006

saltydog

Well hey Hootenanny. Do I need to look up what that means so I know who I'm talkin' to? Anyway, not sure what Hall is referring to but my understanding is that there is a way:

Using the three transformation rules above with the expression:

$$\frac{d\theta}{dt}=\sqrt{\frac{2g}{L}}\sqrt{Cos(\theta)-Cos(\omega)}$$

we obtain:

$$\frac{d\theta}{dt}=2k\sqrt{\frac{g}{L}} Cos(\phi)$$

Since the desired value of theta (the positions of maximum displacement in order to calculate the period) is that for which $\frac{d\theta}{dt}=0$, this then corresponds to:

$$\phi=\pi/2$$

Letting P(k) be the period of the pendulum, we get:

$$P(k)=4\sqrt{\frac{L}{g}}\int_0^{\pi/2} \frac{d\phi}{\sqrt{1-k^2 Sin^2(\phi)}}$$

or:

$$P(k)=4\sqrt{\frac{L}{g}}K(k)$$

where K(k) is the complete elliptic integral of the first kind . . . equal rights . . . ok, you guys are getting annoyed by that.

Anyway, when k=0, this reduces to:

$$P=2\pi\sqrt{\frac{L}{g}}$$

Last edited: Feb 4, 2006
12. Feb 4, 2006

Hootenanny

Staff Emeritus
Sorry, I'm sure this is gonna sound really stupid but which three transformation rules?

13. Feb 4, 2006

saltydog

These:

$$Cos(\theta)=1-2k^2Sin^2(\phi)\quad \text{and}\quad k=Sin(\omega/2)$$

$$Cos(\theta)-Cos(\omega)=2k^2Cos^2(\phi)$$

$$Sin(\theta)=2k Sin(\phi)\sqrt{1-k^2Sin^2(\phi)}$$

$$Sin(\theta)d\theta=4k^2 Sin(\phi) Cos(\phi)d\phi$$

Tell you what though, it would be a good idea if you went through the change of variables using these transformations that lead to the elliptic form of the integral as I described above.

14. Feb 4, 2006

saltydog

Alright, I don't understand this. What, I can say that can't I?

So I plot the period function:

$$P(k)=4\sqrt{\frac{L}{g}}K(k)$$

with:

$$k=sin[\omega/2]$$

with omega being the maximum displacement in radians. The plot is attached. Well, anyway, when k=0, that corresponds to the pendulum at rest. Why then is the period:

$$P(0)=2\pi \sqrt{L/g}\approx 1.92 \text{sec}$$

(for a 3 foot pendulum I mean)

when it's not moving?

Tell you what though, if I had easy access to a college physics lab you can bet I'd be in there with a rigid pendulum checking all this.

Attached Files:

• period of pendulum.JPG
File size:
3.2 KB
Views:
186
Last edited: Feb 4, 2006
15. Feb 4, 2006

Hootenanny

Staff Emeritus
I'm just thinking that because just because k = 0 doesn't have to mean that the pendulum is at rest:

$$k= \sin [\omega/2]$$
$\sin = 0$ when $\omega/2 = 0$.
$$\omega = \frac{v}{r}$$

So when k = 0 could it just indicate that the pendulum is at its maximum displacement because v = 0??

16. Feb 4, 2006

HallsofIvy

Staff Emeritus
No, $\omega$ is NOT the angle of motion, it is the angle of maximum displacement and is a constant. If that maximum displacement if 0, so that k= 0, then there is no motion.

17. Feb 4, 2006

saltydog

The analysis I gave above is taken from "Introduction to Non-Linear Differential and Integral Equations", by H. Davis., Ch. 7. What, you didn't think I did all that from scratch did you?

For example, Davis gives this problem at the end of the section:

A pendulum is displaced through an angle of 45 degrees. Compute the peroid. Answer: $6.53 \sqrt{L/g}$ seconds.

Ok, so I use:

$$\text{maxdisp}=(45/360) 2\pi$$

$$k[\omega]=Sin[\omega/2]$$

I then calculate:

$$4 \text{EllipticK[(k[maxdisp])}^2]}$$

Mathematica returns:

$$6.534$$

so this gives me some confidence the formula is correct.

I tell you what, if we can't get it resolved here I may risk stepping into the General Physics Forum and presenting this apparent discrepancy to them to explain. It's not just a "limiting" issue is is? I mean even if the max displacement is close to zero, the period will be around 1.9 seconds and approaches the limiting case P(0) as omega goes to zero?

Last edited: Feb 4, 2006
18. Feb 5, 2006

Hootenanny

Staff Emeritus
So your saying that there may be a minimum displacement before a significant error is induced?

19. Feb 5, 2006

mathwonk

instead of, or perhaps in addition to, these complicated formulas, it is enlightening to draw the phase curves.

i.e. consider it as a system with y = dx/dt, and then dy/dx = - Csin(x). Then find and plot all the equilibria, i.e. points where both dy/dt = 0 = dx/dt.

Then linearize the system near each equilibrium, i.e. solve the usual approximate system d^2x/dt^2 = -x, and get an idea of the possible local behavior near the equilibria.

in fact these linearizations do give a good approximation to that behavior.
then try to draw in the full phase curves guided by these local pictures.

this gives a rather enlightening picture of the whole situation, although some further analysis is needed to be sure the linear models do reflect the non linear solutions in this case; i.e. since the linear equilibria are "centers", it is not immediate that the non - linear equilibria are also centers, but here they are.

20. Feb 6, 2006

saltydog

I agree with Mathwonk regarding the importance of qualitative analysis of differential equations and am a big advocator of such. I only presented the analytical approach for it's sheer beauty. Although some here have criticized Bob Devaney and his approach to Differential Equations (emphasis on qualitative analysis), I find such absolutely essential in obtaining an intutitive understanding of dynamics. Tell you what though, I say yes, do the phase-portrait for that one but really one should then proceed to open the throttle wide:

$$\frac{d^2\theta}{dt^2}+\frac{1}{q}\frac{d\theta}{dt}+Sin(\theta)=gCos(\omega t)$$

Although I've not studied this one in detail, I understand it exhibits chaos and allows for the existence of strange attractors. Thus one should be able to draw both a Feigenbaum diagram and a strange attractor for it. Yep, yep, I think this is a good one for Hootenanny to work on. Please post a complete report of your findings, plots too.

Edit: Oh yea, just set it up in Mathematica and solve it numerically then plot the attractor parametrically and the Feigenbaum plot as a Poincare' section.

Last edited: Feb 6, 2006