Week 7

Coupled Simple Harmonic Oscillators - The double pendulum

Let us consider the coupled pendula of the class discussion.  In this case, the potential energy of the system is given by U = mgh_1+mgh_2 .  The heights can be approximated using the small angle formula:  h_1= x_1^2/2/l and h_2= h_1 +(x1-x2)^2/2/l.

U[x1_,x2_]=m g (x1^2/2/l)+m g (x1^2/2/l+(x1-x2)^2/2/l);

We can now figure out the equations of motion by calculating the forces:

F1=-D[U[x1,x2],x1]//Simplify

F2=-D[U[x1,x2],x2]//Simplify

We now apply Newton’s 2nd law, and take the assumption that we are exciting a normal mode of oscillation - each of EOM1 and EOM2 are expressions that are equal to zero.

EOM1=m D[A1 Cos[ω t],{t,2}]-(F1/.{x1->A1 Cos[ω t],x2->A2 Cos[ω t]})

EOM2=m D[A2 Cos[ω t],{t,2}]-(F2/.{x1->A1 Cos[ω t],x2->A2 Cos[ω t]})

These equations become algebraic equations after dividing by Cos(ω t).  We can also eliminate A1, and focus only on the ratio A1/A2=r by setting A1=r and A2=1.  Note we also express ω as the square root of ω^2.  This allows us to avoid extra solutions with negative ω that would otherwise show up:

EOM1Alg=Simplify[(EOM1/Cos[ω t])/.{A1->r,A2->1,ω->Sqrt[ωsqr]}]

EOM2Alg=Simplify[(EOM2/Cos[ω t])/.{A1->r,A2->1,ω->Sqrt[ωsqr]}]

Now we can solve for both ω^2 and r:

NMinfo=Simplify[Solve[{EOM1Alg==0,EOM2Alg==0},{r,ωsqr}]]

The set of equations has two solutions that we are calling NMinfo for “Normal Mode info.” Each solution has a given frequency, and each a given ratio of amplitudes of oscillation for each degree of freedom.  Let us try to animate these, asuming g=10m/s^2 and l=10cm=1/10 m:

ωNM1=Sqrt[ωsqr/.NMinfo[[1]]]/.{l->1/10,g->10}

ωNM2=Sqrt[ωsqr/.NMinfo[[2]]]/.{l->1/10,g->10}

rNM1=r/.NMinfo[[1]]

rNM2=r/.NMinfo[[2]]

Animate[ListPlot[{{0,0},{1/20rNM1 Cos[ωNM1 t],-1},{1/20 Cos[ωNM1 t],-2}},Joined->True,PlotMarkers->{Automatic,Large},PlotRange->{{-1,1},{-2.1,0}}],{t,0,10 2π/ωNM1},AnimationRate->1/(10*2π/ωNM1),RefreshRate->100]

Animate[ListPlot[{{0,0},{1/10rNM2 Cos[ωNM2 t],-1},{1/10 Cos[ωNM2 t],-2}},Joined->True,PlotMarkers->{Automatic,Large},PlotRange->{{-1,1},{-2.1,0}}],{t,0,10 2π/ωNM2},AnimationRate->1/(10*2π/ωNM2),RefreshRate->100]

Now let us finally consider superposition of the two normal modes.  You can adjust the amplitudes of each:

Manipulate[Animate[ListPlot[{{0,0},{ANM1 rNM1 Cos[ωNM1 t]+ANM2 rNM2 Cos[ωNM2 t],-1},{ANM1 Cos[ωNM1 t]+ANM2 Cos[ωNM2 t],-2}},Joined->True,PlotMarkers->{Automatic,Large},PlotRange->{{-1,1},{-2.1,0}}],{t,0,10 2π/ωNM2},AnimationRate->1/(10*2π/ωNM2),RefreshRate->100],{{ANM1,1/10},0,1/2},{{ANM2,1/10},0,1/2}]