Driven Quadruple Pendulum

import%20marimo%0A%0A__generated_with%20%3D%20%220.23.9%22%0Aapp%20%3D%20marimo.App(width%3D%22medium%22%2C%20app_title%3D%22Driven%20Quadruple%20Pendulum%22)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_()%3A%0A%20%20%20%20import%20json%0A%20%20%20%20import%20marimo%20as%20mo%0A%20%20%20%20import%20numpy%20as%20np%0A%0A%20%20%20%20return%20json%2C%20mo%2C%20np%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20Driven%20damped%20quadruple%20pendulum%0A%0A%20%20%20%20Four%20equal%20masses%20on%20rigid%20rods%20with%20a%20periodic%20torque%20applied%20at%20the%20first%20pivot%20and%20viscous%0A%20%20%20%20damping%20on%20each%20joint.%20PC3236%20at%20NUS.%0A%0A%20%20%20%20The%20mass%20matrix%20%24M(%5Ctheta)%24%20couples%20all%20four%20angles%20through%20cosine%20terms%3B%20the%20RHS%20carries%20centrifugal%0A%20%20%20%20couplings%2C%20gravity%2C%20and%20the%20drive%20torque%20%24%5Ctau_0%5Ccos(%5Comega_d%20t)%24%20at%20pivot%201%20only.%20Even%20mild%20drive%0A%20%20%20%20amplitudes%20push%20the%20outer%20joints%20into%20irregular%2C%20non-repeating%20motion.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20tau0_sl%20%3D%20mo.ui.slider(0.0%2C%2012.0%2C%20value%3D8.0%2C%20step%3D0.25%2C%20label%3D%22drive%20%CF%84%E2%82%80%20(N%C2%B7m)%22%2C%20show_value%3DTrue)%0A%20%20%20%20omega_sl%20%3D%20mo.ui.slider(1.0%2C%2012.0%2C%20value%3D5.03%2C%20step%3D0.1%2C%20label%3D%22%CF%89_d%20(rad%2Fs)%22%2C%20show_value%3DTrue)%0A%20%20%20%20b_sl%20%3D%20mo.ui.slider(0.0%2C%200.3%2C%20value%3D0.10%2C%20step%3D0.01%2C%20label%3D%22damping%20b%22%2C%20show_value%3DTrue)%0A%20%20%20%20mo.md(f%22%22%22%7Bmo.hstack(%5Btau0_sl%2C%20omega_sl%2C%20b_sl%5D%2C%20justify%3D%22start%22%2C%20gap%3D2)%7D%22%22%22)%0A%20%20%20%20return%20b_sl%2C%20omega_sl%2C%20tau0_sl%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(b_sl%2C%20np%2C%20omega_sl%2C%20tau0_sl)%3A%0A%20%20%20%20_g%20%3D%209.81%0A%20%20%20%20_L%20%3D%20%5B0.40%2C%200.35%2C%200.30%2C%200.25%5D%0A%20%20%20%20_m%20%3D%20%5B1.0%2C%201.0%2C%201.0%2C%201.0%5D%0A%20%20%20%20_tau0%20%3D%20tau0_sl.value%0A%20%20%20%20_wd%20%20%20%3D%20omega_sl.value%0A%20%20%20%20_b%20%20%20%20%3D%20b_sl.value%0A%0A%20%20%20%20def%20_M(th)%3A%0A%20%20%20%20%20%20%20%20L%2C%20m%20%3D%20_L%2C%20_m%0A%20%20%20%20%20%20%20%20Mm%20%3D%20np.zeros((4%2C%204))%0A%20%20%20%20%20%20%20%20cmsum%20%3D%20%5Bsum(m%5Bk%3A%5D)%20for%20k%20in%20range(4)%5D%0A%20%20%20%20%20%20%20%20for%20i%20in%20range(4)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Mm%5Bi%2C%20i%5D%20%3D%20cmsum%5Bi%5D%20*%20L%5Bi%5D%20**%202%0A%20%20%20%20%20%20%20%20for%20i%20in%20range(4)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20j%20in%20range(i%20%2B%201%2C%204)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20v%20%3D%20sum(m%5Bk%5D%20for%20k%20in%20range(j%2C%204))%20*%20L%5Bi%5D%20*%20L%5Bj%5D%20*%20np.cos(th%5Bi%5D%20-%20th%5Bj%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20Mm%5Bi%2C%20j%5D%20%3D%20Mm%5Bj%2C%20i%5D%20%3D%20v%0A%20%20%20%20%20%20%20%20return%20Mm%0A%0A%20%20%20%20def%20_F(th%2C%20om%2C%20t)%3A%0A%20%20%20%20%20%20%20%20L%2C%20m%2C%20g%20%3D%20_L%2C%20_m%2C%20_g%0A%20%20%20%20%20%20%20%20F%20%3D%20np.zeros(4)%0A%20%20%20%20%20%20%20%20cmsum%20%3D%20%5Bsum(m%5Bk%3A%5D)%20for%20k%20in%20range(4)%5D%0A%20%20%20%20%20%20%20%20for%20i%20in%20range(4)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20F%5Bi%5D%20-%3D%20cmsum%5Bi%5D%20*%20g%20*%20L%5Bi%5D%20*%20np.sin(th%5Bi%5D)%0A%20%20%20%20%20%20%20%20for%20i%20in%20range(4)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20j%20in%20range(4)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20if%20i%20%3D%3D%20j%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20continue%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20c%20%3D%20sum(m%5Bk%5D%20for%20k%20in%20range(max(i%2C%20j)%2C%204))%20*%20L%5Bi%5D%20*%20L%5Bj%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20if%20j%20%3C%20i%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20F%5Bi%5D%20%2B%3D%20c%20*%20om%5Bj%5D%20**%202%20*%20np.sin(th%5Bj%5D%20-%20th%5Bi%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20F%5Bi%5D%20-%3D%20c%20*%20om%5Bj%5D%20**%202%20*%20np.sin(th%5Bi%5D%20-%20th%5Bj%5D)%0A%20%20%20%20%20%20%20%20F%5B0%5D%20%2B%3D%20_tau0%20*%20np.cos(_wd%20*%20t)%0A%20%20%20%20%20%20%20%20F%20-%3D%20_b%20*%20om%0A%20%20%20%20%20%20%20%20return%20F%0A%0A%20%20%20%20def%20_derivs(t%2C%20y)%3A%0A%20%20%20%20%20%20%20%20th%2C%20om%20%3D%20y%5B%3A4%5D%2C%20y%5B4%3A%5D%0A%20%20%20%20%20%20%20%20return%20np.concatenate(%5Bom%2C%20np.linalg.solve(_M(th)%2C%20_F(th%2C%20om%2C%20t))%5D)%0A%0A%20%20%20%20def%20_rk4(t%2C%20y%2C%20dt)%3A%0A%20%20%20%20%20%20%20%20k1%20%3D%20_derivs(t%2C%20y)%0A%20%20%20%20%20%20%20%20k2%20%3D%20_derivs(t%20%2B%20dt%20%2F%202%2C%20y%20%2B%20dt%20%2F%202%20*%20k1)%0A%20%20%20%20%20%20%20%20k3%20%3D%20_derivs(t%20%2B%20dt%20%2F%202%2C%20y%20%2B%20dt%20%2F%202%20*%20k2)%0A%20%20%20%20%20%20%20%20k4%20%3D%20_derivs(t%20%2B%20dt%2C%20y%20%2B%20dt%20*%20k3)%0A%20%20%20%20%20%20%20%20return%20y%20%2B%20(dt%20%2F%206)%20*%20(k1%20%2B%202%20*%20k2%20%2B%202%20*%20k3%20%2B%20k4)%0A%0A%20%20%20%20DT%20%3D%200.01%3B%20T_END%20%3D%2012.0%3B%20MAX_FRAMES%20%3D%20800%0A%20%20%20%20_n%20%3D%20int(T_END%20%2F%20DT)%0A%20%20%20%20_step%20%3D%20max(1%2C%20_n%20%2F%2F%20MAX_FRAMES)%0A%20%20%20%20_y%20%3D%20np.zeros(8)%3B%20_y%5B0%5D%20%3D%200.3%3B%20_y%5B1%5D%20%3D%200.1%0A%20%20%20%20_t%20%3D%200.0%3B%20_hist%20%3D%20%5B%5D%0A%20%20%20%20for%20_i%20in%20range(_n)%3A%0A%20%20%20%20%20%20%20%20if%20_i%20%25%20_step%20%3D%3D%200%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20_hist.append(_y%5B%3A4%5D.tolist())%0A%20%20%20%20%20%20%20%20_y%20%3D%20_rk4(_t%2C%20_y%2C%20DT)%3B%20_t%20%2B%3D%20DT%0A%0A%20%20%20%20sim%20%3D%20%7B%22th%22%3A%20_hist%2C%20%22n%22%3A%20len(_hist)%2C%20%22L%22%3A%20_L%7D%0A%20%20%20%20return%20(sim%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(json%2C%20mo%2C%20sim)%3A%0A%20%20%20%20_JS%20%3D%20r%22%22%22%3C!doctype%20html%3E%3Chtml%3E%3Chead%3E%3Cmeta%20charset%3D%22utf-8%22%3E%3Cstyle%3E%0A%20%20%20%20html%2Cbody%7Bmargin%3A0%3Bbackground%3A%230f1117%3Boverflow%3Ahidden%7Dcanvas%7Bdisplay%3Ablock%3Bwidth%3A100vw%3Bheight%3A100vh%7D%0A%20%20%20%20%3C%2Fstyle%3E%3C%2Fhead%3E%3Cbody%3E%3Ccanvas%20id%3D%22c%22%3E%3C%2Fcanvas%3E%3Cscript%3E%0A%20%20%20%20const%20S%3D__SPEC__%2Ccv%3Ddocument.getElementById(%22c%22)%2Cctx%3Dcv.getContext(%222d%22)%3B%0A%20%20%20%20const%20DPR%3DMath.min(window.devicePixelRatio%7C%7C1%2C2)%3B%0A%20%20%20%20function%20fit()%7Bconst%20w%3Dcv.offsetWidth%7C%7Cwindow.innerWidth%7C%7C800%2Ch%3Dcv.offsetHeight%7C%7Cwindow.innerHeight%7C%7C520%3Bif(w%3E0%26%26h%3E0)%7Bcv.width%3Dw*DPR%3Bcv.height%3Dh*DPR%3B%7D%7D%0A%20%20%20%20addEventListener(%22resize%22%2Cfit)%3Bfit()%3B%0A%20%20%20%20const%20VP%3DS.panels.map(()%3D%3E(%7Bzoom%3A1%2Cpx%3A0%2Cpy%3A0%7D))%3B%0A%20%20%20%20let%20frame%3D0%2Cplaying%3Dtrue%2Clast%3D0%3B%0A%20%20%20%20function%20legend(lx%2Cly%2Clw%2Clh%2Centries)%7B%0A%20%20%20%20%20%20const%20y0%3Dly%2Blh-14*DPR%3B%0A%20%20%20%20%20%20entries.forEach((%5Bl%2Cc%5D%2Ci)%3D%3E%7Bctx.fillStyle%3Dc%3Bctx.fillRect(lx%2B5*DPR%2Bi*72*DPR%2Cy0%2C8*DPR%2C8*DPR)%3Bctx.fillStyle%3D%22%238a93a6%22%3Bctx.font%3D(9*DPR)%2B%22px%20monospace%22%3Bctx.textAlign%3D%22left%22%3Bctx.fillText(l%2Clx%2B15*DPR%2Bi*72*DPR%2Cy0%2B7*DPR)%3B%7D)%3B%0A%20%20%20%20%7D%0A%20%20%20%20function%20axes2D(title%2Cxlabel%2Clx%2Cly%2Clw%2Clh%2Cxd%2Cyd)%7B%0A%20%20%20%20%20%20const%20pL%3D34*DPR%2CpR%3D8*DPR%2CpT%3D18*DPR%2CpB%3D22*DPR%3B%0A%20%20%20%20%20%20const%20mx%3Dlx%2BpL%2Cmy%3Dly%2BpT%2Cmw%3Dlw-pL-pR%2Cmh%3Dlh-pT-pB%3B%0A%20%20%20%20%20%20function%20toX(v)%7Breturn%20mx%2B(v-xd.min)%2F(xd.max-xd.min%7C%7C1)*mw%3B%7D%0A%20%20%20%20%20%20function%20toY(v)%7Breturn%20my%2Bmh-(v-yd.min)%2F(yd.max-yd.min%7C%7C1)*mh%3B%7D%0A%20%20%20%20%20%20ctx.strokeStyle%3D%22%231e2535%22%3Bctx.lineWidth%3D1%3B%0A%20%20%20%20%20%20for(let%20i%3D0%3Bi%3C%3D5%3Bi%2B%2B)%7Bconst%20x%3Dmx%2Bi*mw%2F5%3Bctx.beginPath()%3Bctx.moveTo(x%2Cmy)%3Bctx.lineTo(x%2Cmy%2Bmh)%3Bctx.stroke()%3B%7D%0A%20%20%20%20%20%20for(let%20i%3D0%3Bi%3C%3D4%3Bi%2B%2B)%7Bconst%20y%3Dmy%2Bi*mh%2F4%3Bctx.beginPath()%3Bctx.moveTo(mx%2Cy)%3Bctx.lineTo(mx%2Bmw%2Cy)%3Bctx.stroke()%3B%7D%0A%20%20%20%20%20%20ctx.strokeStyle%3D%22%232a3448%22%3Bctx.lineWidth%3D1%3Bctx.beginPath()%3Bctx.rect(mx%2Cmy%2Cmw%2Cmh)%3Bctx.stroke()%3B%0A%20%20%20%20%20%20ctx.fillStyle%3D%22%234a5568%22%3Bctx.font%3D(8*DPR)%2B%22px%20monospace%22%3Bctx.textAlign%3D%22center%22%3B%0A%20%20%20%20%20%20for(let%20i%3D0%3Bi%3C%3D5%3Bi%2B%3D2)%7Bconst%20v%3Dxd.min%2Bi*(xd.max-xd.min)%2F5%3Bctx.fillText(v.toFixed(0)%2Cmx%2Bi*mw%2F5%2Cmy%2Bmh%2B10*DPR)%3B%7D%0A%20%20%20%20%20%20ctx.textAlign%3D%22right%22%3B%0A%20%20%20%20%20%20for(let%20i%3D0%3Bi%3C%3D4%3Bi%2B%3D2)%7Bconst%20v%3Dyd.min%2Bi*(yd.max-yd.min)%2F4%3Bctx.fillText(v.toFixed(1)%2Cmx-3*DPR%2Cmy%2Bmh-i*mh%2F4%2B4*DPR)%3B%7D%0A%20%20%20%20%20%20if(title)%7Bctx.fillStyle%3D%22%239aa3b0%22%3Bctx.font%3D(10*DPR)%2B%22px%20monospace%22%3Bctx.textAlign%3D%22center%22%3Bctx.fillText(title%2Clx%2Blw%2F2%2Cly%2B12*DPR)%3B%7D%0A%20%20%20%20%20%20if(xlabel)%7Bctx.fillStyle%3D%22%23626b7c%22%3Bctx.font%3D(9*DPR)%2B%22px%20monospace%22%3Bctx.textAlign%3D%22center%22%3Bctx.fillText(xlabel%2Clx%2Blw%2F2%2Cly%2Blh-4*DPR)%3B%7D%0A%20%20%20%20%20%20return%7Bmx%2Cmy%2Cmw%2Cmh%2CtoX%2CtoY%7D%3B%0A%20%20%20%20%7D%0A%20%20%20%20function%20drawPanel(p%2Ci%2Clx%2Cly%2Clw%2Clh)%7B%0A%20%20%20%20%20%20const%20f%3Dframe%3B%0A%20%20%20%20%20%20ctx.save()%3Bctx.beginPath()%3Bctx.rect(lx%2Cly%2Clw%2Clh)%3Bctx.clip()%3B%0A%20%20%20%20%20%20if(p.kind%3D%3D%3D%22pendulum_chain%22)%7B%0A%20%20%20%20const%20Ls%3Dp.L%2CthF%3Dp.th%5Bf%5D%3B%0A%20%20%20%20const%20cx%3Dlx%2Blw*0.5%2Ccy%3Dly%2Blh*0.10%3B%0A%20%20%20%20const%20totalL%3DLs.reduce((a%2Cb)%3D%3Ea%2Bb%2C0)%3B%0A%20%20%20%20const%20sc%3DMath.min(lw%2Clh)*0.76%2FtotalL%3B%0A%20%20%20%20let%20bx%3Dcx%2Cby%3Dcy%3Bconst%20bobs%3D%5B%5D%3B%0A%20%20%20%20for(let%20j%3D0%3Bj%3CLs.length%3Bj%2B%2B)%7B%0A%20%20%20%20%20%20const%20nx%3Dbx%2Bsc*Ls%5Bj%5D*Math.sin(thF%5Bj%5D)%2Cny%3Dby%2Bsc*Ls%5Bj%5D*Math.cos(thF%5Bj%5D)%3B%0A%20%20%20%20%20%20bobs.push(%5Bnx%2Cny%5D)%3Bbx%3Dnx%3Bby%3Dny%3B%0A%20%20%20%20%7D%0A%20%20%20%20let%20rx%3Dcx%2Cry%3Dcy%3B%0A%20%20%20%20for(let%20j%3D0%3Bj%3Cbobs.length%3Bj%2B%2B)%7B%0A%20%20%20%20%20%20ctx.strokeStyle%3Dp.colors%5Bj%5D%7C%7C%22%237dd3fc%22%3Bctx.lineWidth%3D2.5*DPR%3B%0A%20%20%20%20%20%20ctx.beginPath()%3Bctx.moveTo(rx%2Cry)%3Bctx.lineTo(bobs%5Bj%5D%5B0%5D%2Cbobs%5Bj%5D%5B1%5D)%3Bctx.stroke()%3B%0A%20%20%20%20%20%20rx%3Dbobs%5Bj%5D%5B0%5D%3Bry%3Dbobs%5Bj%5D%5B1%5D%3B%0A%20%20%20%20%7D%0A%20%20%20%20const%20rs%3D%5B9%2C7%2C6%2C5%5D%3B%0A%20%20%20%20for(let%20j%3D0%3Bj%3Cbobs.length%3Bj%2B%2B)%7B%0A%20%20%20%20%20%20ctx.fillStyle%3Dp.colors%5Bj%5D%7C%7C%22%237dd3fc%22%3B%0A%20%20%20%20%20%20ctx.beginPath()%3Bctx.arc(bobs%5Bj%5D%5B0%5D%2Cbobs%5Bj%5D%5B1%5D%2Crs%5Bj%5D*DPR%2C0%2C7)%3Bctx.fill()%3B%0A%20%20%20%20%7D%0A%20%20%20%20ctx.fillStyle%3D%22%238a93a6%22%3Bctx.beginPath()%3Bctx.arc(cx%2Ccy%2C3*DPR%2C0%2C7)%3Bctx.fill()%3B%0A%20%20%20%20legend(lx%2Cly%2Clw%2Clh%2C(p.colors%7C%7C%5B%5D).map((c%2Cj)%3D%3E%5B%22%CE%B8%22%2B(j%2B1)%2Cc%5D))%3B%0A%20%20%20%20%20%20%7D%0A%20%20%20%20%20%20else%20if(p.kind%3D%3D%3D%22timeseries%22)%7B%0A%20%20%20%20const%20ys_arr%3Dp.ys%7C%7C%5B%5D%2Cn%3Dys_arr%5B0%5D%3Fys_arr%5B0%5D.length%3A0%3B%0A%20%20%20%20const%20ts%3Dp.t%7C%7CArray.from(%7Blength%3An%7D%2C(_%2Ck)%3D%3Ek)%3B%0A%20%20%20%20let%20ymin%3DInfinity%2Cymax%3D-Infinity%3B%0A%20%20%20%20for(const%20ya%20of%20ys_arr)for(const%20v%20of%20ya)%7Bif(v%3Cymin)ymin%3Dv%3Bif(v%3Eymax)ymax%3Dv%3B%7D%0A%20%20%20%20ymin-%3DMath.abs(ymax-ymin)*0.07%3Bymax%2B%3DMath.abs(ymax-ymin)*0.07%3B%0A%20%20%20%20const%20xd%3D%7Bmin%3Ats%5B0%5D%7C%7C0%2Cmax%3Ats%5Bn-1%5D%7C%7C1%7D%2Cyd%3D%7Bmin%3Aymin%7C%7C0%2Cmax%3Aymax%7C%7C1%7D%3B%0A%20%20%20%20const%7Bmx%2Cmy%2Cmw%2Cmh%2CtoX%2CtoY%7D%3Daxes2D(p.title%2Cp.xlabel%2Clx%2Cly%2Clw%2Clh%2Cxd%2Cyd)%3B%0A%20%20%20%20ctx.save()%3Bctx.beginPath()%3Bctx.rect(mx%2Cmy%2Cmw%2Cmh)%3Bctx.clip()%3B%0A%20%20%20%20const%20colors%3Dp.colors%7C%7C%5B%22%237dd3fc%22%2C%22%23f97316%22%2C%22%23a594f2%22%2C%22%23e0635a%22%5D%3B%0A%20%20%20%20const%20upTo%3DMath.min(f%2B1%2Cn)%3B%0A%20%20%20%20for(let%20si%3D0%3Bsi%3Cys_arr.length%3Bsi%2B%2B)%7B%0A%20%20%20%20%20%20const%20ya%3Dys_arr%5Bsi%5D%3Bctx.strokeStyle%3Dcolors%5Bsi%25colors.length%5D%3Bctx.lineWidth%3D2*DPR%3Bctx.beginPath()%3B%0A%20%20%20%20%20%20for(let%20k%3D0%3Bk%3CupTo%3Bk%2B%2B)%7Bconst%20x%3DtoX(ts%5Bk%5D)%2Cy%3DtoY(ya%5Bk%5D)%3Bk%3D%3D%3D0%3Fctx.moveTo(x%2Cy)%3Actx.lineTo(x%2Cy)%3B%7D%0A%20%20%20%20%20%20ctx.stroke()%3B%0A%20%20%20%20%7D%0A%20%20%20%20ctx.restore()%3B%0A%20%20%20%20if(p.labels)legend(lx%2Cly%2Clw%2Clh%2Cp.labels.map((l%2Cj)%3D%3E%5Bl%2C(p.colors%7C%7Ccolors)%5Bj%25colors.length%5D%5D))%3B%0A%20%20%20%20%20%20%7D%0A%20%20%20%20%20%20ctx.restore()%3B%0A%20%20%20%20%7D%0A%20%20%20%20function%20draw()%7B%0A%20%20%20%20%20%20const%20W%3Dcv.width%2CH%3Dcv.height%3B%0A%20%20%20%20%20%20ctx.fillStyle%3D%22%230f1117%22%3Bctx.fillRect(0%2C0%2CW%2CH)%3B%0A%20%20%20%20%20%20for(let%20i%3D0%3Bi%3CS.panels.length%3Bi%2B%2B)%7B%0A%20%20%20%20const%20lay%3DS.layout%5Bi%5D%3B%0A%20%20%20%20drawPanel(S.panels%5Bi%5D%2Ci%2Clay%5B0%5D*W%2Clay%5B1%5D*H%2Clay%5B2%5D*W%2Clay%5B3%5D*H)%3B%0A%20%20%20%20%20%20%7D%0A%20%20%20%20%20%20ctx.fillStyle%3D%22%232a3448%22%3Bctx.font%3D(9*DPR)%2B%22px%20monospace%22%3Bctx.textAlign%3D%22left%22%3B%0A%20%20%20%20%20%20ctx.fillText(%22frame%20%22%2Bframe%2B%22%2F%22%2BS.frames%2B%22%20%20click%20to%20pause%22%2C6*DPR%2CH-5*DPR)%3B%0A%20%20%20%20%7D%0A%20%20%20%20function%20loop(t)%7Bif(playing%26%26t-last%3ES.dt)%7Bframe%3D(frame%2B1)%25S.frames%3Blast%3Dt%3B%7Dtry%7Bdraw()%3B%7Dcatch(e)%7B%7DrequestAnimationFrame(loop)%3B%7D%0A%20%20%20%20cv.addEventListener(%22click%22%2C()%3D%3E%7Bplaying%3D!playing%3B%7D)%3B%0A%20%20%20%20let%20_s%3Dfalse%3B%0A%20%20%20%20const%20_ro%3Dnew%20ResizeObserver(()%3D%3E%7Bfit()%3Bif(!_s%26%26cv.width%3E0)%7B_s%3Dtrue%3BrequestAnimationFrame(loop)%3B%7D%7D)%3B%0A%20%20%20%20_ro.observe(cv)%3BsetTimeout(()%3D%3E%7Bfit()%3Bif(!_s%26%26cv.width%3E0)%7B_s%3Dtrue%3BrequestAnimationFrame(loop)%3B%7D%7D%2C100)%3B%0A%20%20%20%20%3C%2Fscript%3E%3C%2Fbody%3E%3C%2Fhtml%3E%22%22%22%0A%0A%20%20%20%20_th%20%3D%20sim%5B%22th%22%5D%0A%20%20%20%20_spec%20%3D%20%7B%0A%20%20%20%20%20%20%20%20%22frames%22%3A%20sim%5B%22n%22%5D%2C%20%22dt%22%3A%2020%2C%0A%20%20%20%20%20%20%20%20%22panels%22%3A%20%5B%0A%20%20%20%20%20%20%20%20%20%20%20%20%7B%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22kind%22%3A%20%22pendulum_chain%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22L%22%3A%20sim%5B%22L%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22th%22%3A%20_th%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22colors%22%3A%20%5B%22%237dd3fc%22%2C%20%22%23a594f2%22%2C%20%22%23f97316%22%2C%20%22%23e0635a%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%7D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%7B%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22kind%22%3A%20%22timeseries%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22title%22%3A%20%22%CE%B8%E2%82%81%20and%20%CE%B8%E2%82%82%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22xlabel%22%3A%20%22frame%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22ys%22%3A%20%5B%5Brow%5B0%5D%20for%20row%20in%20_th%5D%2C%20%5Brow%5B1%5D%20for%20row%20in%20_th%5D%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22colors%22%3A%20%5B%22%237dd3fc%22%2C%20%22%23a594f2%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22labels%22%3A%20%5B%22%CE%B8%E2%82%81%22%2C%20%22%CE%B8%E2%82%82%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%7D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%7B%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22kind%22%3A%20%22timeseries%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22title%22%3A%20%22%CE%B8%E2%82%83%20and%20%CE%B8%E2%82%84%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22xlabel%22%3A%20%22frame%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22ys%22%3A%20%5B%5Brow%5B2%5D%20for%20row%20in%20_th%5D%2C%20%5Brow%5B3%5D%20for%20row%20in%20_th%5D%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22colors%22%3A%20%5B%22%23f97316%22%2C%20%22%23e0635a%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22labels%22%3A%20%5B%22%CE%B8%E2%82%83%22%2C%20%22%CE%B8%E2%82%84%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%7D%2C%0A%20%20%20%20%20%20%20%20%5D%2C%0A%20%20%20%20%20%20%20%20%22layout%22%3A%20%5B%0A%20%20%20%20%20%20%20%20%20%20%20%20%5B0.00%2C%200.00%2C%200.40%2C%201.00%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%5B0.40%2C%200.00%2C%200.60%2C%200.50%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%5B0.40%2C%200.50%2C%200.60%2C%200.50%5D%2C%0A%20%20%20%20%20%20%20%20%5D%2C%0A%20%20%20%20%7D%0A%20%20%20%20mo.vstack(%5B%0A%20%20%20%20%20%20%20%20mo.md(%22**Quadruple%20pendulum%2C%20live.**%20Click%20to%20pause.%22)%2C%0A%20%20%20%20%20%20%20%20mo.iframe(_JS.replace(%22__SPEC__%22%2C%20json.dumps(_spec))%2C%20height%3D%22520px%22)%2C%0A%20%20%20%20%5D)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(b_sl%2C%20mo%2C%20omega_sl%2C%20sim%2C%20tau0_sl)%3A%0A%20%20%20%20_last%20%3D%20sim%5B%22th%22%5D%5B-1%5D%20if%20sim%5B%22th%22%5D%20else%20%5B0%2C%200%2C%200%2C%200%5D%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Drive%3A**%20%CF%84%E2%82%80%20%3D%20%7Btau0_sl.value%3A.2f%7D%20N%C2%B7m%2C%20%CF%89_d%20%3D%20%7Bomega_sl.value%3A.2f%7D%20rad%2Fs%2C%20b%20%3D%20%7Bb_sl.value%3A.2f%7D.%0A%20%20%20%20Final%20angles%3A%20%CE%B8%E2%82%81%20%3D%20%7B_last%5B0%5D%3A.2f%7D%2C%20%CE%B8%E2%82%82%20%3D%20%7B_last%5B1%5D%3A.2f%7D%2C%20%CE%B8%E2%82%83%20%3D%20%7B_last%5B2%5D%3A.2f%7D%2C%20%CE%B8%E2%82%84%20%3D%20%7B_last%5B3%5D%3A.2f%7D%20rad.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Equations%20of%20motion.**%20The%20Lagrangian%20gives%20a%20%244%5Ctimes4%24%20mass-matrix%20system%3A%0A%0A%20%20%20%20%24%24M(%5Ctheta)%5C%2C%5Cddot%5Ctheta%20%3D%20F(%5Ctheta%2C%5Cdot%5Ctheta%2Ct)%24%24%0A%0A%20%20%20%20where%0A%0A%20%20%20%20%24%24M_%7Bij%7D%20%3D%20%5CBigl(%5Csum_%7Bk%5Cgeq%5Cmax(i%2Cj)%7D%20m_k%5CBigr)%20L_i%20L_j%20%5Ccos(%5Ctheta_i-%5Ctheta_j)%24%24%0A%0A%20%20%20%20and%20%24F_i%24%20collects%20gravity%20(%24-M_%7Bii%7Dg%2FL_i%5Ccdot%5Csin%5Ctheta_i%24)%2C%20centrifugal%20coupling%0A%20%20%20%20(off-diagonal%20angular-velocity%20terms)%2C%20the%20drive%20torque%20%24%5Ctau_0%5Ccos(%5Comega_d%20t)%24%20at%20%24i%3D1%24%0A%20%20%20%20only%2C%20and%20viscous%20damping%20%24-b%5Cdot%5Ctheta_i%24.%0A%20%20%20%20RK4%20at%20%24%5CDelta%20t%20%3D%2010%5E%7B-2%7D%24%20s%2C%2012%20s%20trajectory.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A