• API
• FAQ
• Tools
• Archive
SHARE
TWEET

# Double Pendulum

Matthen Jul 23rd, 2011 1,846 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. (* Drawing of diagram *)
2. DrawPendula[\[Theta]1_, \[Theta]2_, m1_, m2_, l1_, l2_] :=
4.     radius1 = 0.2 (m1/Max[m1, m2])^(1/3);
5.     radius2 = 0.2 (m2/Max[m1, m2])^(1/3);
6.     Graphics[{
7.       PointSize[Large],
8.       Blue, Point[{0, 0}],
9.       Orange,
10.       Line[{{0,
11.          0}, {l1 Sin[\[Theta]1], -l1 Cos[\[Theta]1]}, {l1 \
12. Sin[\[Theta]1] + l2 Sin[\[Theta]2], -l1 Cos[\[Theta]1] -
13.           l2 Cos[\[Theta]2]}}],
14.       Red, Disk[{l1 Sin[\[Theta]1], -l1 Cos[\[Theta]1]}, radius1],
15.       Disk[{l1 Sin[\[Theta]1] +
16.          l2 Sin[\[Theta]2], -l1 Cos[\[Theta]1] - l2 Cos[\[Theta]2]},
18.       }, PlotRange -> Table[{-1 - l1 - l2, 1 + l1 + l2}, {2}],
19.      ImageSize -> {500}]
20.     )];
21. (* Interface *)
22. doublepend[\[Theta]1_, \[Theta]2_, d\[Theta]1_, d\[Theta]2_, m1_, m2_,
23.     l1_, l2_, g_, \[Alpha]_, ta_] :=
24.   Module[{eqns, soln, A, B, plot, \[Alpha]2},
25.    \[Alpha]2 = \[Alpha] B'[t];
26.    eqns = {m2 l2 B''[t] + m2 l1 A''[t] Cos[A[t] - B[t]] -
27.        m2 l1 (A'[t]^2) Sin[A[t] - B[t]] +
28.        m2 g Sin[B[t]] + \[Alpha]2 m2 l2 ==
29.       0, (m1 + m2) l1 A''[t] + m2 l2 B''[t] Cos[A[t] - B[t]] +
30.        m2 l2 (B'[t]^2) Sin[A[t] - B[t]] +
31.        g (m1 + m2) Sin[A[t]] + \[Alpha] A'[t] (m1 + m2) == 0,
32.      A[0] == \[Theta]1, A'[0] == d\[Theta]1, B[0] == \[Theta]2,
33.      B'[0] == d\[Theta]2};
34.    soln =
35.     NDSolve[eqns, {A, B}, {t, ta - ta/5 - 0.6, ta},
36.      MaxSteps -> Infinity, PrecisionGoal -> 10];
37.    If[ta < 25,
38.     plot =
39.      ParametricPlot[{{l1 Sin[A[t]], -l1 Cos[A[t]]}, {l1 Sin[A[t]] +
40.           l2 Sin[B[t]], -l1 Cos[A[t]] - l2 Cos[B[t]]}} /. soln, {t,
41.        ta - ta/5 - 0.5, ta}];
42.     Show[First[DrawPendula[A[ta], B[ta], m1, m2, l1, l2] /. soln],
43.      plot],
44.     Graphics[Text[Style["(Solution now unstable, probs)", Large]],
45.      ImageSize -> {500}]
46.     ]
47.    ];
48.
49. Manipulate[
50.  doublepend[\[Theta]1, \[Theta]2, d\[Theta]1, d\[Theta]2, m1, m2, l1,
51.   l2, g, \[Alpha], ta]
52.  ,
53.  {{\[Theta]1, 7 Pi/8,
54.    "Initial \!\(\*SubscriptBox[\"\[Theta]\", \"1\"]\)"}, 0, 2 Pi},
55.  {{\[Theta]2, Pi/3,
56.    "Initial \!\(\*SubscriptBox[\"\[Theta]\", \"2\"]\)"}, 0, 2 Pi},
57.  {{d\[Theta]1, 2,
58.    "Initial \!\(\*SubscriptBox[\"\[Theta]\", \"1\"]\)'"}, -5, 5},
59.  {{d\[Theta]2, 1,
60.    "Initial \!\(\*SubscriptBox[\"\[Theta]\", \"2\"]\)'"}, -5, 5},
61.  {{m1, 1, "\!\(\*SubscriptBox[\"m\", \"1\"]\)"}, 0.1, 5},
62.  {{m2, 2, "\!\(\*SubscriptBox[\"m\", \"2\"]\)"}, 0.1, 5},
63.  {{l1, 1.5, "\!\(\*SubscriptBox[\"l\", \"1\"]\)"}, 0.1, 5},
64.  {{l2, 1, "\!\(\*SubscriptBox[\"l\", \"2\"]\)"}, 0.1, 5},
65.  {{g, 10, "gravity"}, 0, 20},
66.  {{\[Alpha], 0.1, "friction"}, 0, 0.5},
67.  {{ta, 0.001, "animate"}, 0.001, Infinity, 0.05,
68.   ControlType -> Trigger}
69.  ]
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.

Top