This week only. Pastebin PRO Accounts Christmas Special! Don't miss out!Want more features on Pastebin? Sign Up, it's FREE!
Guest

Untitled

By: a guest on Dec 7th, 2013  |  syntax: None  |  size: 0.82 KB  |  views: 40  |  expires: Never
download  |  raw  |  embed  |  report abuse  |  print
Text below is selected. Please press Ctrl+C to copy to your clipboard. (⌘+C on Mac)
  1. (*Initial Conditions*)
  2.  
  3. m = 0.050;
  4. M = 2;
  5. r = 0.25;
  6. g = 9.8;
  7. \[Epsilon] = 0.1;
  8. \[Theta]0 = \[Pi]/3;
  9. tf = 100;
  10.  
  11. IR = m r^2;
  12.  
  13. (*Code*)
  14.  
  15. nPend = 3;
  16.  
  17. p[i_][t_] := \[Theta][i]''[
  18.    t] + ((m r g)/
  19.     IR) Sin[\[Theta][i][t]] + \[Epsilon] ((\[Theta][i][t]/\[Theta]0)^2 -
  20.      1) \[Theta][i]'[t] + ((r m Cos[\[Theta][i][t]])/IR) y''[t]
  21. y[t_] := -(m/(M + nPend m)) r (
  22.    Total[ Table[Sin[\[Theta][j][t]], {j, 1, nPend}]])
  23.  
  24. eqns = Thread[Table[p[i][t], {i, 1, nPend}] == 0]; (*Equations of motion*)
  25. initpos = Thread[Table[\[Theta][i][t], {i, 1, nPend}] == 1]; (*Initial position*)
  26. initvel = Thread[Table[\[Theta][i]'[t], {i, 1, nPend}] == 0]; (*Initial velocity*)
  27.  
  28. alleqns = Join[eqns, initpos, initvel];
  29.  
  30. slist = Table[\[Theta][i], {i, 1, nPend}]; (*What NDSolve solves for*)
  31.  
  32. NDSolve[alleqns, slist, {t, 1, tf}]
clone this paste RAW Paste Data