# Double Pendulum

 In[7]:= XAnimate[Graphics[{{PointSize[.025], {Red, Point[{x1[t], y1[t]}]}, {Blue, Point[{x2[t], y2[t]}]}, Line[{{0, 0}, {x1[t], y1[t]}, {x2[t], y2[t]}}]} /. soldp, {Gray, Line[Map[Function[Evaluate[{x2[#], y2[#]} /. soldp]], Range[0, t, 0.025]]]}}, PlotRange -> {{-2, 2}, {-2, 0}}, Axes -> True, Ticks -> False, ImageSize -> 500], {t, 0, 10, .025}, SaveDefinitions -> True]
 In[8]:= XWith[{t = 3.7}, Graphics[{{Thick, Line[{{0, 0}, {x1[t], y1[t]}, {x2[t], y2[t]}}], PointSize[.025], {Red, Point[{x1[t], y1[t]}]}, {Blue, Point[{x2[t], y2[t]}]},} /. soldp, {Gray, Line[Map[Function[Evaluate[{x2[#], y2[#]} /. soldp]], Range[0, t, 0.025]]]}}, PlotRange -> {{-1.6, 1.6}, {-2, 0}}, Axes -> True, Ticks -> False, ImageSize -> 500]]
 Play Animation »Stop Animation »

#### Derive the governing equations using Newton's second law of motion, and .

 In[1]:= Xdeqns = {Subscript[m, 1] x1''[t] == (\[Lambda]1[t]/Subscript[l, 1]) x1[ t] - (\[Lambda]2[t]/Subscript[l, 2]) (x2[t] - x1[t]), Subscript[m, 1] y1''[t] == (\[Lambda]1[t]/Subscript[l, 1]) y1[ t] - (\[Lambda]2[t]/Subscript[l, 2]) (y2[t] - y1[t]) - Subscript[m, 1] g, Subscript[m, 2] x2''[t] == (\[Lambda]2[t]/Subscript[l, 2]) (x2[t] - x1[t]), Subscript[m, 2] y2''[t] == (\[Lambda]2[t]/Subscript[l, 2]) (y2[t] - y1[t]) - Subscript[m, 2] g};

#### The lengths of the pendulum rods are fixed. These are expressed as algebraic constraints.

 In[2]:= Xaeqns = {x1[t]^2 + y1[t]^2 == Subscript[l, 1]^2, (x2[t] - x1[t])^2 + (y2[t] - y1[t])^2 == Subscript[l, 2]^2};

#### Specify the initial state of the system as initial conditions.

 In[3]:= Xics = {x1[0] == 1, y1[0] == 0, x1'[0] == 0, y1'[0] == 0, x2[0] == 1, y2[0] == -1, x2'[0] == 0, y2'[0] == 0};

#### Specify the physical parameters for the pendulum system.

 In[4]:= Xparams = {g -> 9.81, Subscript[m, 1] -> 1, Subscript[m, 2] -> 1, Subscript[l, 1] -> 1, Subscript[l, 2] -> 1};

#### Solve and visualize the system.

 In[5]:= Xsoldp = First[ NDSolve[{deqns, aeqns, ics} /. params, {x1, y1, x2, y2, \[Lambda]1, \[Lambda]2}, {t, 0, 15}, Method -> {"IndexReduction" -> {Automatic, "IndexGoal" -> 0}}]];
 In[6]:= XParametricPlot[ Evaluate[{{x1[t], y1[t]}, {x2[t], y2[t]}} /. soldp], {t, 0, 15}, PlotStyle -> {Red, Blue}, ImageSize -> Medium, PlotLegends -> {"Trajectory of pendulum 1", "Trajectory of pendulum 2"}]
 Out[6]=