You will have to set values for the initial mass of the engine, the mass of the fuel, and the exhaust velocity.
In [1]:
%matplotlib inline
In [2]:
from __future__ import division, print_function
from ivisual import *
from math import *
import matplotlib.pyplot as plt
In [50]:
#set these values
m0engine= #kg
mfuel= #kg
vex= #m/s
######
g=9.8
mdisk=2.94 #kg
Rdisk=0.6048/2 #m
mf=m0engine-mfuel
m=m0engine
m1=m0engine
r=0.3 #m
mdot=-mfuel/2 #2 s burn time
t=0
dt=0.001
theta=0
omega=0
Idisk=1/2*mdisk*Rdisk*Rdisk
Ieng1=m1*r*r
Ieng2=m*r*r
I=Idisk+Ieng1+Ieng2
r1=r*vector(cos(theta+pi),sin(theta+pi),0)
r2=r*vector(cos(theta),sin(theta),0)
tlist=[]
thetalist=[]
omegalist=[]
mlist=[]
taulist=[]
tau=vector(0,0,0)
while t<3:
if(m>mf):
m=m+mdot*dt
Fgrav1=m1*vector(0,-g,0)
Fgrav2=m*vector(0,-g,0)
thrust=-mdot*vex
taugrav1=cross(r1,Fgrav1)
taugrav2=cross(r2,Fgrav2)
if(m>mf):
tauthrust=vector(0,0,thrust*r)
else:
tauthrust=vector(0,0,0)
tau=taugrav1+taugrav2+tauthrust
Ieng2=m*r*r
I=Idisk+Ieng1+Ieng2
omega=omega+tau.z/I*dt
theta=theta+omega*dt
r1=r*vector(cos(theta+pi),sin(theta+pi),0)
r2=r*vector(cos(theta),sin(theta),0)
t=t+dt
tlist.append(t)
thetalist.append(theta)
omegalist.append(omega)
mlist.append(m)
taulist.append(tau.z)
In [51]:
plt.title('omega vs. t')
plt.xlabel('t (s)')
plt.ylabel('omega (rad/s)')
plt.plot(tlist,omegalist,'m.')
plt.show()
In [17]:
plt.title('tau vs. t')
plt.xlabel('t (s)')
plt.ylabel('tau (N m)')
plt.plot(tlist,taulist,'y.')
plt.show()
In [ ]: