{"id":1204,"date":"2021-11-27T15:40:12","date_gmt":"2021-11-27T23:40:12","guid":{"rendered":"https:\/\/gantovnik.com\/bio-tips\/?p=1204"},"modified":"2024-05-30T19:06:09","modified_gmt":"2024-05-31T02:06:09","slug":"210-parametric-curve-in-3d-2-2","status":"publish","type":"post","link":"https:\/\/gantovnik.com\/bio-tips\/2021\/11\/210-parametric-curve-in-3d-2-2\/","title":{"rendered":"#212 The double pendulum simulation using python"},"content":{"rendered":"<p><img data-recalc-dims=\"1\" decoding=\"async\" src=\"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2021\/11\/ex212.png?resize=500%2C400&#038;ssl=1\" alt=\"\" width=\"500\" height=\"400\" class=\"alignnone size-full wp-image-1205\" srcset=\"https:\/\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2021\/11\/ex212.png 500w, https:\/\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2021\/11\/ex212-480x384.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 500px, 100vw\" \/><\/p>\n<pre class=\"brush: python; title: ; notranslate\" title=\"\">\r\nfrom numpy import sin, cos\r\nimport numpy as np\r\nimport matplotlib.pyplot as plt\r\nimport scipy.integrate as integrate\r\nimport matplotlib.animation as animation\r\nfrom collections import deque\r\nG = 9.8  # acceleration due to gravity, in m\/s^2\r\nL1 = 1.0  # length of pendulum 1 in m\r\nL2 = 1.0  # length of pendulum 2 in m\r\nL = L1 + L2  # maximal length of the combined pendulum\r\nM1 = 1.0  # mass of pendulum 1 in kg\r\nM2 = 1.0  # mass of pendulum 2 in kg\r\nt_stop = 25  # how many seconds to simulate\r\nhistory_len = 10000  # how many trajectory points to display\r\n\r\ndef derivs(state, t):\r\n    dydx = np.zeros_like(state)\r\n    dydx&#x5B;0] = state&#x5B;1]\r\n    delta = state&#x5B;2] - state&#x5B;0]\r\n    den1 = (M1+M2) * L1 - M2 * L1 * cos(delta) * cos(delta)\r\n    dydx&#x5B;1] = ((M2 * L1 * state&#x5B;1] * state&#x5B;1] * sin(delta) * cos(delta)\r\n                + M2 * G * sin(state&#x5B;2]) * cos(delta)\r\n                + M2 * L2 * state&#x5B;3] * state&#x5B;3] * sin(delta)\r\n                - (M1+M2) * G * sin(state&#x5B;0])) \/ den1)\r\n    dydx&#x5B;2] = state&#x5B;3]\r\n    den2 = (L2\/L1) * den1\r\n    dydx&#x5B;3] = ((- M2 * L2 * state&#x5B;3] * state&#x5B;3] * sin(delta) * cos(delta)\r\n                + (M1+M2) * G * sin(state&#x5B;0]) * cos(delta)\r\n                - (M1+M2) * L1 * state&#x5B;1] * state&#x5B;1] * sin(delta)\r\n                - (M1+M2) * G * sin(state&#x5B;2])) \/ den2)\r\n    return dydx\r\n\r\n# create a time array from 0..t_stop sampled at 0.02 second steps\r\ndt = 0.02\r\nt = np.arange(0, t_stop, dt)\r\n# th1 and th2 are the initial angles (degrees)\r\n# w10 and w20 are the initial angular velocities (degrees per second)\r\nth1 = 120.0\r\nw1 = 0.0\r\nth2 = -10.0\r\nw2 = 0.0\r\n# initial state\r\nstate = np.radians(&#x5B;th1, w1, th2, w2])\r\n# integrate your ODE using scipy.integrate.\r\ny = integrate.odeint(derivs, state, t)\r\nx1 = L1*sin(y&#x5B;:, 0])\r\ny1 = -L1*cos(y&#x5B;:, 0])\r\nx2 = L2*sin(y&#x5B;:, 2]) + x1\r\ny2 = -L2*cos(y&#x5B;:, 2]) + y1\r\nfig = plt.figure(figsize=(5, 4))\r\nax = fig.add_subplot(autoscale_on=False, xlim=(-L, L), ylim=(-L, 1.))\r\nax.set_aspect(&#039;equal&#039;)\r\nax.grid()\r\nline, = ax.plot(&#x5B;], &#x5B;], &#039;o-&#039;, lw=1,c=&quot;black&quot;)\r\ntrace, = ax.plot(&#x5B;], &#x5B;], &#039;o-&#039;, lw=1, ms=2)\r\ntime_template = &#039;time = %.1fs&#039;\r\ntime_text = ax.text(0.05, 0.9, &#039;&#039;, transform=ax.transAxes)\r\nhistory_x, history_y = deque(maxlen=history_len), deque(maxlen=history_len)\r\n\r\ndef animate(i):\r\n    thisx = &#x5B;0, x1&#x5B;i], x2&#x5B;i]]\r\n    thisy = &#x5B;0, y1&#x5B;i], y2&#x5B;i]]\r\n    if i == 0:\r\n        history_x.clear()\r\n        history_y.clear()\r\n\r\n    history_x.appendleft(thisx&#x5B;2])\r\n    history_y.appendleft(thisy&#x5B;2])\r\n    line.set_data(thisx, thisy)\r\n    trace.set_data(history_x, history_y)\r\n    time_text.set_text(time_template % (i*dt))\r\n    return line, trace, time_text\r\n\r\nani = animation.FuncAnimation(fig, animate, len(y), interval=dt*1000, blit=True,repeat=False)\r\nplt.savefig(&#039;ex212.png&#039;, dpi=72)\r\nplt.show()\r\n<\/pre>\n","protected":false},"excerpt":{"rendered":"<p>from numpy import sin, cos import numpy as np import matplotlib.pyplot as plt import scipy.integrate as integrate import matplotlib.animation as animation from collections import deque G = 9.8 # acceleration due to gravity, in m\/s^2 L1 = 1.0 # length of pendulum 1 in m L2 = 1.0 # length of pendulum 2 in m [&hellip;]<\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"nf_dc_page":"","_et_pb_use_builder":"","_et_pb_old_content":"","_et_gb_content_width":"","_lmt_disableupdate":"yes","_lmt_disable":"","_jetpack_newsletter_access":"","_jetpack_dont_email_post_to_subs":false,"_jetpack_newsletter_tier_id":0,"_jetpack_memberships_contains_paywalled_content":false,"_jetpack_memberships_contains_paid_content":false,"footnotes":"","jetpack_post_was_ever_published":false},"categories":[89,2],"tags":[],"class_list":["post-1204","post","type-post","status-publish","format-standard","hentry","category-animation","category-python"],"modified_by":"gantovnik","jetpack_featured_media_url":"","jetpack_sharing_enabled":true,"jetpack_shortlink":"https:\/\/wp.me\/p8bH0k-jq","jetpack_likes_enabled":true,"jetpack-related-posts":[{"id":193,"url":"https:\/\/gantovnik.com\/bio-tips\/2019\/01\/double-pendulum\/","url_meta":{"origin":1204,"position":0},"title":"#39 Double pendulum using python","author":"gantovnik","date":"2019-01-10","format":false,"excerpt":"import os import numpy as np import matplotlib.pyplot as plt from scipy import integrate import sympy os.chdir(r'D:\\projects\\wordpress\\ex39') os.getcwd() t, g, m1, l1, m2, l2 = sympy.symbols(\"t, g, m_1, l_1, m_2, l_2\") theta1, theta2 = sympy.symbols(\"theta_1, theta_2\", cls=sympy.Function) ode1 = sympy.Eq((m1+m2)*l1 * theta1(t).diff(t,t) + m2*l2 * theta2(t).diff(t,t) + m2*l2 * theta2(t).diff(t)**2\u2026","rel":"","context":"In &quot;python&quot;","block_context":{"text":"python","link":"https:\/\/gantovnik.com\/bio-tips\/category\/python\/"},"img":{"alt_text":"example39","src":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example39.png?resize=350%2C200","width":350,"height":200,"srcset":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example39.png?resize=350%2C200 1x, https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example39.png?resize=525%2C300 1.5x"},"classes":[]},{"id":190,"url":"https:\/\/gantovnik.com\/bio-tips\/2019\/01\/coupled-damped-springs-jacobian-is-available\/","url_meta":{"origin":1204,"position":1},"title":"Coupled damped springs (Jacobian is available)","author":"gantovnik","date":"2019-01-10","format":false,"excerpt":"\u00a0 import os import numpy as np import matplotlib.pyplot as plt from scipy import integrate os.chdir(r'D:\\projects\\wordpress\\ex38') os.getcwd() def f(t, y, args): m1, k1, g1, m2, k2, g2 = args return [y[1], - k1\/m1 * y[0] + k2\/m1 * (y[2] - y[0]) - g1\/m1 * y[1], y[3], - k2\/m2 * (y[2]\u2026","rel":"","context":"In &quot;python&quot;","block_context":{"text":"python","link":"https:\/\/gantovnik.com\/bio-tips\/category\/python\/"},"img":{"alt_text":"example38","src":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example38.png?resize=350%2C200","width":350,"height":200,"srcset":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example38.png?resize=350%2C200 1x, https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example38.png?resize=525%2C300 1.5x"},"classes":[]},{"id":187,"url":"https:\/\/gantovnik.com\/bio-tips\/2019\/01\/coupled-damped-springs\/","url_meta":{"origin":1204,"position":2},"title":"Coupled damped springs","author":"gantovnik","date":"2019-01-09","format":false,"excerpt":"\u00a0 import os import numpy as np import matplotlib.pyplot as plt from scipy import integrate os.chdir(r'D:\\projects\\wordpress\\ex37') os.getcwd() def f(t, y, args): m1, k1, g1, m2, k2, g2 = args return [y[1], - k1\/m1 * y[0] + k2\/m1 * (y[2] - y[0]) - g1\/m1 * y[1], y[3], - k2\/m2 * (y[2]\u2026","rel":"","context":"In &quot;python&quot;","block_context":{"text":"python","link":"https:\/\/gantovnik.com\/bio-tips\/category\/python\/"},"img":{"alt_text":"example37","src":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example37.png?resize=350%2C200","width":350,"height":200,"srcset":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example37.png?resize=350%2C200 1x, https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2019\/01\/example37.png?resize=525%2C300 1.5x"},"classes":[]},{"id":2048,"url":"https:\/\/gantovnik.com\/bio-tips\/2024\/01\/405-cardioid-in-polar-coordinates-by-matplot-in-python\/","url_meta":{"origin":1204,"position":3},"title":"#405 Cardioid in polar coordinates by matplotlib in python","author":"gantovnik","date":"2024-01-14","format":false,"excerpt":"A cardioid is a two-dimensional plane figure that has a heart-shaped curve described in polar coordinates by the equation r = 2a(1 + cos \u03b8) for 0 \u2264 \u03b8 \u2264 2\u03c0: [code language=\"python\"] import numpy as np import matplotlib.pyplot as plt theta = np.linspace(0, 2.*np.pi, 1000) a = 1.0 r\u2026","rel":"","context":"In &quot;math&quot;","block_context":{"text":"math","link":"https:\/\/gantovnik.com\/bio-tips\/category\/math\/"},"img":{"alt_text":"","src":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2024\/01\/ex405.png?resize=350%2C200&ssl=1","width":350,"height":200,"srcset":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2024\/01\/ex405.png?resize=350%2C200&ssl=1 1x, https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2024\/01\/ex405.png?resize=525%2C300&ssl=1 1.5x"},"classes":[]},{"id":1758,"url":"https:\/\/gantovnik.com\/bio-tips\/2023\/01\/335-solution-of-nonlinear-equation-using-brent-method-in-python\/","url_meta":{"origin":1204,"position":4},"title":"#335 Solution of nonlinear equation by Brent&#8217;s method in python","author":"gantovnik","date":"2023-01-06","format":false,"excerpt":"Brent's method is a hybrid root-finding algorithm combining the bisection method, the secant method, and inverse quadratic interpolation. ex335.py [code language=\"python\"] import numpy as np import matplotlib.pyplot as plt import scipy.optimize as optimize def f(x): return 2*x - 1 + 2*np.cos(np.pi*x) x = np.linspace(0.0,2.0,201) y=f(x) plt.plot(x,y,label='$2x - 1 + 2\\\\cos(\\\\pi\u2026","rel":"","context":"In &quot;optimization&quot;","block_context":{"text":"optimization","link":"https:\/\/gantovnik.com\/bio-tips\/category\/optimization\/"},"img":{"alt_text":"","src":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2023\/01\/ex335.png?resize=350%2C200&ssl=1","width":350,"height":200},"classes":[]},{"id":1242,"url":"https:\/\/gantovnik.com\/bio-tips\/2021\/11\/210-parametric-curve-in-3d-2-2-2-2-2-2-2-2-2-2-2-2-2\/","url_meta":{"origin":1204,"position":5},"title":"#223 Initial value problem in python","author":"gantovnik","date":"2021-11-29","format":false,"excerpt":"[code language=\"python\"] import numpy as np import matplotlib.pyplot as plt class IV_Problem: \"\"\" Initial value problem (IVP) class \"\"\" def __init__(self, rhs, y0, interval, name='IVP'): \"\"\" rhs 'right hand side' function of the ordinary differential equation f(t,y) y0 array with initial values interval start and end value of the interval\u2026","rel":"","context":"In &quot;python&quot;","block_context":{"text":"python","link":"https:\/\/gantovnik.com\/bio-tips\/category\/python\/"},"img":{"alt_text":"","src":"https:\/\/i0.wp.com\/gantovnik.com\/bio-tips\/wp-content\/uploads\/2021\/11\/ex223.png?resize=350%2C200&ssl=1","width":350,"height":200},"classes":[]}],"_links":{"self":[{"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/posts\/1204","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/comments?post=1204"}],"version-history":[{"count":3,"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/posts\/1204\/revisions"}],"predecessor-version":[{"id":2207,"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/posts\/1204\/revisions\/2207"}],"wp:attachment":[{"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/media?parent=1204"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/categories?post=1204"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/gantovnik.com\/bio-tips\/wp-json\/wp\/v2\/tags?post=1204"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}