The Driven Pendulum

Solving the Lagrangian

Setup the problem with Potential Energy (U) and Kinetic Energy (T)
The pendulum suport point moves as a function of time, xd[t].  The coordinates of the pendulum mass are xm and ym.

xm[t]=xd[t]+l Sin[phi[t]];
ym[t]=l Cos[phi[t]];
U=Simplify[-m g ym[t]]
T = Simplify[(1/2) m (D[xm[t],t]^2+D[ym[t],t]^2)]

Solve the Euler-Lagrange Equation

Fphi=Simplify[D[T-U,phi[t]]]
pphi=Simplify[D[T-U,phi'[t]]]
dpphi=Simplify[D[pphi,t]]
eomd=Simplify[Fphi==dpphi,Assumptions->{m>0,l>0}]

Solve the Equations of Motion:
Low frequency drive - superposition

List dependent variables

dvars={phi};

Set initial condisions and drive function

incon={phi[0]==phi0,phi'[0]==v0};
x1[t]=ad Sin[fd Sqrt[g/l] t];
eom=eomd/.{xd''[t]->D[D[x1[t],t],t]};

Define all parameters:
g = gravity
l = length of pendulum
ad = drive amplitude
fd = drive frequency / pendulum resonant frequency
phi0 = phi at t = 0 = initial angle from vertical
v0 = dphi/dt at t = 0 = initial angular velocity

pars={g->π^2,l->1,ad->0.1,fd->0.1,phi0->0,v0->0};
tmax=30 2π/Sqrt[g/l]/.pars;

Solve with NDSolve:
Try different values of fd all less than 1...

sol=NDSolve[{eom,incon}//.pars,dvars,{t,0,tmax},MaxStepSize->N[1/(10 fd)/.pars]]

Plot the solution (phi[t], which doesn’t change much for low frequency drives)

plotpars={PerformanceGoal->"Quality",ImageSize->Full,AspectRatio->1/4};
Plot[Evaluate[phi[t]/.sol],{t,0,tmax},Evaluate[plotpars]]

Solve the Equations of Motion:
Near resonance drive - beats

Try different values of fd near 1...

pars={g->10,l->1,ad->0.01,fd->0.9,phi0->0,v0->0};
tmax=30 2π/Sqrt[g/l]/.pars;
sol=NDSolve[{eom,incon}//.pars,dvars,{t,0,tmax},MaxStepSize->N[1/(10 fd)/.pars]]
Plot[Evaluate[phi[t]/.sol],{t,0,tmax},Evaluate[plotpars]]

Solve the Equations of Motion:
On resonance drive - saturation

Try different values of ad and note how the saturation time and amplitude change

pars={g->10,l->1,ad->0.02,fd->1,phi0->0,v0->0};
tmax=30 2π/Sqrt[g/l]/.pars;
sol=NDSolve[{eom,incon}//.pars,dvars,{t,0,tmax},MaxStepSize->N[1/(10 fd)/.pars]]
Plot[Evaluate[phi[t]/.sol],{t,0,tmax},Evaluate[plotpars]]

Solve the Equations of Motion:
No dirve, big swing - anharmonic

Try different values of phi0 and note how the number of cycles changes (it should be 5 for low amplitude).

pars={g->10,l->1,ad->0.0,fd->0.1,phi0->2 ,v0->0};
tmax=5 2π/Sqrt[g/l]/.pars;
sol=NDSolve[{eom,incon}//.pars,dvars,{t,0,tmax},MaxStepSize->N[1/(10 fd)/.pars]]
Plot[Evaluate[phi[t]/.sol],{t,0,tmax},Evaluate[plotpars]]

Solve the Equations of Motion:
High Frequency Drive - strange

pars={g->10,l->1,ad->0.2,fd->10,phi0->π/2,v0->0};
tmax=5 Sqrt[g/l]/.pars;
sol=NDSolve[{eom,incon}//.pars,dvars,{t,0,tmax},MaxStepSize->N[1/(10 fd)/.pars]]
Plot[Evaluate[phi[t]/.sol],{t,0,tmax},Evaluate[plotpars]]

Solve the Equations of Motion:
Higher Frequency Drive - very strange

Try different values of fd and note how the stability and frequency of the π/2 oscillations change

pars={g->10,l->1,ad->0.2,fd->11,phi0->π/2,v0->0};
tmax=5 Sqrt[g/l]/.pars;
sol=NDSolve[{eom,incon}//.pars,dvars,{t,0,tmax},MaxStepSize->N[1/(10 fd)/.pars]]
Plot[Evaluate[phi[t]/.sol],{t,0,tmax},Evaluate[plotpars]]