<jupyter_start><jupyter_text>Set up rotation matrices representing a 3-1-3 $(\psi,\theta,\phi)$ Euler angle set.<jupyter_code>aCi = rotMat(3,psi)
cCa = rotMat(1,th)
bCc = rotMat(3,ph)
aCi,cCa,bCc
bCi = bCc*cCa*aCi; bCi #3-1-3 rotation
bCi_dot = difftotalmat(bCi,t,{th:thd,psi:psid,ph:phd});
bCi_dot<jupyter_output><empty_output><jupyter_text>$\tilde{\omega} = {}^\mathcal{B}C^{\mathcal{I}}  {}^\mathcal{B}{\dot{C}}^{\mathcal{I}}$<jupyter_code>omega_tilde = bCi*bCi_dot.T; omega_tilde<jupyter_output><empty_output><jupyter_text>$\left[{}^\mathcal{I}\boldsymbol{\omega}^{\mathcal{B}}\right]_\mathcal{B} = \left[ {}^\mathcal{B}C^{\mathcal{I}}_{32} \quad  {}^\mathcal{B}C^{\mathcal{I}}_{13} \quad {}^\mathcal{B}C^{\mathcal{I}}_{21} \right]$<jupyter_code>omega = simplify(Matrix([omega_tilde[2,1],omega_tilde[0,2],omega_tilde[1,0]]))
omega
w1,w2,w3 = symbols('omega_1,omega_2,omega_3')
s0 = solve(omega - Matrix([w1,w2,w3]),[psid,thd,phd]); s0<jupyter_output><empty_output><jupyter_text>Find EOM (second derivatives of Euler Angles)<jupyter_code>I1,I2,I3 = symbols("I_1,I_2,I_3",real=True,positive=True)
iWb_B = omega
I_G_B = diag(I1,I2,I3)
I_G_B
diffmap = {th:thd,psi:psid,ph:phd,thd:thdd,psid:psidd,phd:phdd}
diffmap
t1 = I_G_B*difftotalmat(iWb_B,t,diffmap) 
t2 = skew(iWb_B)*I_G_B*iWb_B
t1,t2
dh_G_B = t1+t2
dh_G_B
t3 = expand(dh_G_B[0]*cos(ph)*I2 - dh_G_B[1]*sin(ph)*I1)
sol_thdd = simplify(solve(t3,thdd)) 
sol_thdd
t4= expand(dh_G_B[0]*sin(ph)*I2 + dh_G_B[1]*cos(ph)*I1)
t4
sol_psidd = simplify(solve(t4,psidd)) 
sol_psidd
sol_phdd = solve(dh_G_B[2],phdd)
sol_phdd<jupyter_output><empty_output><jupyter_text>Find initial orientation such that $\mathbf h$ is down-pointing<jupyter_code>h = sqrt(((I_G_B*Matrix([w1,w2,w3])).transpose()*(I_G_B*Matrix([w1,w2,w3])))[0]);h
eqs1 = simplify(bCi.transpose()*I_G_B*Matrix([w1,w2,w3]) - Matrix([0,0,-h])); eqs1 #equal 0
simplify(solve(simplify(eqs1[0]*cos(psi) + eqs1[1]*sin(psi)),ph)) #phi solution
solve(simplify(expand(simplify(-eqs1[0]*sin(psi) + eqs1[1]*cos(psi)).subs(ph,atan(I1*w1/I2/w2)))),th) #th solution
simplify(eqs1[2].subs(ph,atan(I1*w1/I2/w2)))<jupyter_output><empty_output><jupyter_text>Generate MATLAB Code<jupyter_code>out = codegen(("eom1",sol_psidd[0]), 'Octave', argument_sequence=[th,thd,psi,psid,ph,phd,I1,I2,I3]);out
codegen(("eom1",sol_thdd[0]), 'Octave', argument_sequence=[th,thd,psi,psid,ph,phd,I1,I2,I3])
codegen(("eom1",sol_phdd[0]), 'Octave', argument_sequence=[th,thd,psi,psid,ph,phd,I1,I2,I3,psidd])
codegen(("eom1",[s0[psid],s0[thd],s0[phd]]), 'Octave', argument_sequence=[w1,w2,w3,th,thd,psi,psid,ph,phd,I1,I2,I3,psidd])
codegen(("eom1",bCi), 'Octave', argument_sequence=[th,thd,psi,psid,ph,phd,I1,I2,I3,psidd])
codegen(("eom1",omega), 'Octave', argument_sequence=[w1,w2,w3,th,thd,psi,psid,ph,phd,I1,I2,I3,psidd])<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.22454771195459383, 'fraction_numerical': 0.04611564384533522, 'mean_word_length': 12.364928909952607, 'pattern_counts': {'":': 0, '<': 25, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '513238', 'license': 'MIT', 'n_tokens_mistral': 1442, 'n_tokens_neox': 1400, 'n_words': 127, 'path': 'Notebooks/Torque Free 3-1-3 Body Dynamics.ipynb', 'repo_name': 'dsavransky/MAE4060', 'size': '2819'}