{VERSION 4 0 "IBM INTEL NT" "4.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{CSTYLE "" -1 256 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 257 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 258 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 259 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 260 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 261 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 262 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 263 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 264 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 265 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 266 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 267 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 268 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 269 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 270 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 271 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 272 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 273 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 274 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 275 "" 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 }{PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 256 1 {CSTYLE "" -1 -1 "" 1 16 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 256 "" 0 "" {TEXT 256 25 "The mathematical pendul um" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 102 "We follow Gerd Baumann's approach in Mathematica to determine the period of the mathematical pendulum." }}{PARA 0 "" 0 "" {TEXT -1 165 "First, we write down the kinetic and potential energies: the kinetic energy \+ is expressed in terms of the angular velocity and the moment of inerti a for a point mass " }{TEXT 274 1 "m" }{TEXT -1 23 " suspended at dist ance " }{TEXT 273 1 "l" }{TEXT -1 1 "." }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 43 "We call the angular velocity d/dt phi=omega" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 24 "Tkin:=1/2*m*l^2*omega^2;" }}}{EXCHG {PARA 0 "" 0 " " {TEXT -1 77 "The potential energy follows from the height of the pen dulum bob (point mass " }{TEXT 275 1 "m" }{TEXT -1 74 "). We choose th e potential energy to be zero if the mass is at the bottom." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "Upot:=m*g*l*(1-cos(phi));" } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 51 "For plotting we have to pick som e parameter values:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "l:=1; " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 6 "g:=10;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "m:=1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 141 "plots[contourplot](Tkin+Upot,phi=-3/2*Pi..3/2*Pi,ome ga=-10..10,contours=[1,5,9,13,17,21,25,29],axes=boxed,filled=true, col oring=[red,green]);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 547 "This phas e space diagram is obviously periodic with period 2 Pi. For small angu lar velocities one has a system approximated by a harmonic oscillator, while for larger velocites (or energies, i.e., higher contours) one h as solutions that simply revolve (winding the pendulum; in the periodi c phase space diagrams one leaves e.g., at the right to re-enter at th e same angular velocity at the left. The motion is left to right in th e top half (for the open curves), and right to the left in the bottom \+ half (look at the sign of the angular velocity)." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 115 "The regularity of the di agram is very striking. The units of energy (chosen for the contour li nes) are in MKSA, if " }{TEXT 259 1 "l" }{TEXT -1 1 "," }{TEXT 258 1 " g" }{TEXT -1 1 "," }{TEXT 257 1 "m" }{TEXT -1 19 " are given in MKSA. " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 358 "Ther e is a critical energy where the separation occurs between oscillating solutions (closed contours, in the somewhat abstract phase space one \+ obtains elliptic or deformed elliptic 'trajectories'; translate those \+ into oscillations using your imagination, by watching the sign of omeg a, passage through phi=0, i.e., the bottom position of the pendulum, e tc.)" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 229 " The locations phi=+Pi, omega=0, or phi=-Pi, omega=0 correspond to the \+ unstable equilibrium (pendulum at the top at rest). One can show that \+ the phase-space 'trajectory' become a cosine curves of phi/2 [Marion-T hornton, eq. 4.31]." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 " " {TEXT -1 385 "For winding solutions we observe how their angular vel ocity is modulated as the go through the (periodic) potential well. Fo r sufficiently high energies (not shown) the phase space 'trajectories ' approach straigh lines, as the pendulum bob's kinetic energy is so h igh that the variation in the potential energy (by being moved up and \+ down in the graviational field) becomes negligible." }}{PARA 0 "" 0 " " {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 298 "Time-dependent perturb ations of the pendulum by an external force can lead to a loss of regu larity in the phase space diagram (CHAOS). The modulation can come fro m a simple coupling to a second pendulum, which leads to a time-varyin g support point (science museum demonstrations + classroom demo!)." }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 102 "The dete rmination of the phase-space 'trajectories' was easy. The actual integ ration of phi(t) is not!" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 145 "The \+ total energy can be expressed through the potential energy at the turn ing point (pendulum bob at rest), i.e., through the initial angle phi0 ." }}{PARA 0 "" 0 "" {TEXT -1 77 "One introduces the energy parameter \+ (initial phase angle phi0) k=sin(phi0/2);" }}{PARA 0 "" 0 "" {TEXT -1 135 "To recover the general expression for the potential energy we una ssign the parameters for length, gravitational acceleration, and mass: " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "l:='l': g:='g': m:='m': Upot;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 46 "Upot:=algsubs(cos (phi)=1-2*sin(phi/2)^2,Upot);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "Etot:=subs(phi=phi0,Upot);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 61 "The total energy for a given trajectory thus depends only on " } {TEXT 261 1 "k" }{TEXT -1 66 "=sin(phi/2). Therefore, the period is ca lculated as a function of " }{TEXT 260 1 "k" }{TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 145 "To calculate t he period for the oscillatory solutions, one uses the energy expressio n which depends on phi and omega=diff(phi,t), isolates omega." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "Tkin+Upot=Etot;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "solve(%,omega);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 41 "DE:=diff(phi(t),t)=subs(phi=phi(t),%[1]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 37 "sol:=dsolve(\{DE,phi(0)=p hi0\},phi(t));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 224 "Maple 5.0 come s up with a formal solution for phi(t). It hints at the possibility of finding the solution through elliptic integrals. In Maple 6.0 no answ er is returned after the program tries to find a solution for a while. " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 36 "For t he determination of the period " }{TEXT 262 1 "T" }{TEXT -1 117 " one \+ can proceed differently: the differential equation DE can be separated so that an expression for dt is obtained:" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 32 "subs(diff(phi(t),t)=dphi/dt,DE);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "solve(%,dt);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 203 "If we integrate the right-hand side over phi from 0 to p hi0, the left-hand side will give us a quarter of the period T (for th e full period phi would have to continue through 0 to -phi0 and back t o phi0." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 71 "The integral to be carried out falls into a category known as elli ptic." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "assume(k>0);" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 61 "A required integral can be express ed as an elliptic function:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "Int(1/sqrt((1-z^2)*(1-k^2*z^2)),z=0..1);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "value(%);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "omega[0]:=sqrt(g/l);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 76 "T o help Maple we restrict the inital phase to be positive, and less tha n Pi." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "assume(phi0>0,phi0 " 0 "" {MPLTEXT 1 0 66 "T:=2/omega[0]*Int( 1/sqrt(sin(phi0/2)^2-sin(phi/2)^2),phi=0..phi0);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "T:=value(T);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 126 "For plotting purposes let us simplify matters by chosing the pendulum length and the gravitational acceleration in MKSA units: " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "g:=10; l:=10;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "T;" }}}{EXCHG {PARA 0 "" 0 " " {TEXT -1 67 "We now explore this answer as a function of the initial angle phi0:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "plot(T,phi0 =0..Pi);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 108 "In the limit of the \+ intial phase going to Pi (unstable equilibrium at the top) we obtain a n infinite period:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "limit (T,phi0=Pi,left);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 " " 0 "" {TEXT -1 108 "What do we get for small initial phase angles (we should recover isochronicity in the small-amplitude limit)" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "plot(T,phi0=0..Pi/16);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "taylor(T,phi0=0,5);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "T1:=algsubs(sin(phi0/2)=k,T) ;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "series(T1,k=0,8);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "T1t:=convert(%,polynom);" }} }{EXCHG {PARA 0 "" 0 "" {TEXT -1 69 "This agrees with eq. 2.25 in Baum ann (or analytical mechanics texts)." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "plot(\{Re(T1),T1t\},k=0..1);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 50 "An alternative would be to introduce the variable " } {TEXT 269 1 "k" }{TEXT -1 26 " directly in the integrand" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 64 "T2:=2/omega[0]*int(1/sqrt(k^2-sin(p hi/2)^2),phi=0..2*arcsin(k));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "simplify(T2);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 30 "This does \+ not work in Maple5! " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "sim plify(T2,symbolic);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 132 "Maple is \+ reluctant to carry out the cancellation of function and its inverse. I n Maple6 T2 is directly calculated as 4*EllipticK(k)." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 75 "A further exploration \+ should look into the calculation of trajectories phi(" }{TEXT 270 1 "t " }{TEXT -1 2 ")!" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 376 "Before we proceed with this step we construct the phase \+ diagram (a plot of the time derivative of the dynamic variable [the ph ase angle phi] against the dynamic variable itself. For this we need t he energy expression (we will look at contour lines of equal height, i .e. lines of constant energy, as the mechanical energy is conserved fo r a frictionless mathematical pendulum)." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "l:='l': g:='g': m:='m':" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "phi:='phi';" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 3 "DE;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 48 "Let us make use of Ge rd Baumann's analytic work:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "omega[0]:=sqrt(g/l);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 " Etot;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "E_s:=Etot/(m*l^2*o mega[0]^2);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 255 "The scaled energy parameter allows to distinguish easily the different types of motion: for E_s<2 we have oscillations, for E_s>2 rotations, and E_s=2 corres ponds to the asymptotic solution where the pendulum bob remains in the unstable equilibrium point." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 55 "We express the scaled energy in terms of the pa rameter " }{TEXT 263 1 "k" }{TEXT -1 13 "=sin(phi0/2):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "k0:=7/8;" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 12 "E_s:=2*k0^2;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "phi:='phi':" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 99 "We define \+ the differential equation such that positive and negative angular velo cities are allowed." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 58 "DE:= (diff(phi(t),t))^2=omega[0]^2*(2*(cos(phi(t))-1+E_s));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 69 "Oscillatory solutions to the equation of \+ motion (E_s<2) are given by:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 49 "theta:=(k,t)->2*arcsin(k*JacobiSN(omega[0]*t,k));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "g:=10; l:=1;" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 28 "subs(phi(t)=theta(k0,t),DE);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "simplify(%);" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 19 "evalf(subs(t=1,%));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 83 "Without further work Maple cannot figure out that the exp ression written for theta(" }{TEXT 271 1 "t" }{TEXT -1 134 ") solves t he differential equation. Nevertheless, we can verify that within nume rical accuracy the differential equation is satisfied." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 124 "plot([theta(0.5,t),theta(0.8,t),th eta(0.9,t),theta(0.999,t)],t=0..2*Pi,color=[red,blue,green,black],view =[0..2*Pi,-Pi..Pi]);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 181 "One obse rves that for moderate energy values the solution still looks like a s ine-curve, although the dependence of the period on the energy is appa rent already. Only for values of " }{TEXT 264 1 "k" }{TEXT -1 164 " cl ose to unity do we see a qualtiative change in the behaviour of the so lution, i.e., a dramatic slowdown of the pendulum bob in the vicinity \+ of the turning point." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 99 "We can follow this up with an animation. To display th e pendulum's motion we need to calculate the " }{TEXT 267 1 "x" } {TEXT -1 5 " and " }{TEXT 266 1 "y" }{TEXT -1 36 " coordinates. We cho ose a length of " }{TEXT 268 1 "l" }{TEXT -1 99 "=1 and define the coo rdinates such that the height is zero at the suspension point of the p endulum." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "x:=(k,t)->evalf (sin(theta(k,t))): y:=(k,t)->-evalf(cos(theta(k,t))):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "x(0.99,1),y(0.99,1);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 50 "We draw lines from the origin to the pend ulum bob:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 151 "PL:=(k,t)->pl ot([[0,0],[x(k,t),y(k,t)]],style=line,color=magenta,title=cat(cat(\"k= \",convert(evalf(k,3),string)),\"time= \",convert(evalf(t,3),string)) ):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "PLs:=seq(PL(0.995,it/ 10),it=1..200):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "with(plo ts):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 60 "display(PLs,inseque nce=true,axes=boxed,scaling=constrained);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 231 "The slowdown of the pendulum bob as the top position is \+ approached is very obvious. Thus the lack of isochronicity, or the dep endence of the period on the maximum amplitude (energy deposited into \+ the pendulum) is made very obvious." }}{PARA 0 "" 0 "" {TEXT -1 298 "I n Maple 6 the title for each frame is actually displayed in the animat ion. However, we found two errors for Maple 6.0 in this worksheet: 1) \+ the contourplot is wrong; 2) the numerical evaluation of the JacobiSN \+ functions is in error. We can only hope that a service pack will remed y these problems." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 215 "While the demonstration of Maple's handling of el liptic integrals has been impressive, we wish to demonstrate also a ge neral approach that can be applied to any differential equation, namel y the numerical solution. " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 63 "DE:=k->(diff(phi(t),t))^2=omega[0]^2*(2*(cos(phi(t))-1+2*k^2));" } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "omega[0]:=1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "IC:=phi(0)=0;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "sol:=dsolve(\{DE(0.999),IC\},phi(t),numeric); " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 100 "Bad news! We need to go back to the original Newton law. Maple 6 is, however, improved in this are a." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 90 "omega:='omega': # it \+ is not good to mix omega(t), and omega[0], which makes omega a table. " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 64 "DE:=diff(phi(t),t)=omeg a(t),diff(omega(t),t)=-(g/l)*sin(phi(t));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 94 "The initial condition is chosen to be the turning point, \+ where all energy is in the potential." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "phifun:=k->2*arcsin(k);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "phi0:=phifun(0.999);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 27 "IC:=omega(0)=0,phi(0)=phi0;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "sol:=dsolve(\{DE,IC\},[omega(t),phi(t)],numeric, output=listprocedure);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "t heta:=subs(sol,phi(t));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 125 "In th e plot command we need to put quotes around theta(t), as the expressio n can only be evaluated for numerical values of t." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 27 "plot('theta(t)',t=0..2*Pi);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "Omega:=subs(sol,omega(t));" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 27 "plot('Omega(t)',t=0..2*Pi); " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 152 "It is evident that there are extended time segments when the pendulum bob moves with low velocitie s. A useful exercise would be to reduce the parameter " }{TEXT 265 1 " k" }{TEXT -1 137 ", and to observe when the relationship between angul ar velocity and angle returns to that of a harmonic oscillator (sine/c osine curves). " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 65 "We can also lo ok at the solutions where the pendulum bob rotates." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "theta:=(t,k)->2*JacobiAM(sqrt(g/l)*t,k); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "g:=10; l:=1;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "plot(Re(theta(t,0.6*I)),t=0. .4);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "evalf(theta(1,5.5*I ));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 24 "The imaginary values of " }{TEXT 272 1 "k" }{TEXT -1 84 " come from the following analytic conti nuation of the energy-initial angle relation:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "E_s:=2.1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 28 "solve(1-E_s=cos(phi1),phi1);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "k:=evalf(sin(%));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 141 "Exercise: Adapt the animation to display the motion in the cas e where the pendulum bob has sufficient energy to rotate over the top \+ position." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{MARK "0 \+ 0 0" 0 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }