Skip to content

Adding two dimensions #1

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
Rickisterr wants to merge 1 commit into
base: testbranch
Choose a base branch
from
Open

Adding two dimensions #1

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
89 changes: 53 additions & 36 deletions Rocket_Simulator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@
#EVERYTHING WORKS EXCEPT FOR THE MAX ACCELERATION. MAX ACC IS ALMOST 10M/S^2 LESS THAN OPEN ROCKET

from this import d
from turtle import color
# from turtle import color
import matplotlib.pyplot as plt
import csv
import time as timee
Expand Down Expand Up @@ -40,22 +40,22 @@
G = 6.67408 * pow(10,-11) # Gravitational constant
ME = 5.972 * pow(10,24) # Earth mass in Kg
Re = 6371 * pow(10,3) # Radius of earth in m
cg = 0.528 #in cm
cp = 0.728 #in cm
cg = 0.528 #in cm (distance from rocket's nose to center of pressure)
cp = 0.728 #in cm (distance from rocket's nose to center of gravity)
mr = 0.618 #Mass without casing or propellant in kg
me = 0.043
mp = 0.060
im = 49.61 #in N-s
th = 14.5 #Avg thrust in N(use this)
mth = 25.3 #Max Thrust in N
tb = 3.42 #in s (im/th)
rho = 1.22 #in Kg/m^3, Drops by 10% for each 1000m of altitude for the first 14km or we can use the formula = 1.22 * 0.9^(Altitude/1000)
cd = 0.75 #Avg Cd, but we can calculate it based on our rocket dimensions.
dr = 0.045 #in m
me = 0.043 #Engine case mass
mp = 0.012 #Propellent mass defaullt=0.060
im = 8 #in N-s default=49.61
th = 6.7 #Avg thrust in N(use this) default=14.5
mth = 12 #Max Thrust in N default=25.3
tb = 3.42 #in s (im/th)
rho = 1.22 #(density of air) in Kg/m^3, Drops by 10% for each 1000m of altitude for the first 14km or we can use the formula = 1.22 * 0.9^(Altitude/1000)
cd = 0.75 #Avg Cd, but we can calculate it based on our rocket dimensions (drag coeff)
dr = 0.045 #in m (rocket diameter)
a = 3.14 * pow(dr/2,2) #Area of the rocket
dt = 0.01 #in s
dm = (mp/tb)*dt#mass flow rate
g = 9.8 #m/s
g = 9.8 #m/s2
tmass = mr + me + mp #Total Mass in kg
alt = []
vel = []
Expand All @@ -72,7 +72,7 @@
prevHeight = 0
t2apogee = 0

theta = 20 #in degrees from x axis.
theta = 20 #in degrees from x axis. Given 20
thetaV = 0 # Velocity directions

# Currently, we are considering that the angle remains the same throughout the entire flight.
Expand Down Expand Up @@ -121,11 +121,15 @@ def getThrust(t):

def getDrag(v):
d = 0.5 * rho * cd * a * pow(v,2)
return d
d_x = d * math.cos(math.radians(theta))
d_y = d * math.sin(math.radians(theta))
return [d, d_x, d_y]

def getParachuteDrag(v):
d = (cd * rho * 3.14 * pow(diameter,2) * pow(v,2))/8 #(diameter /2) ^ 2 ==> (4) and then the entire formula is 1/2 * .... ==> 1/8
return d
d_x = d * math.cos(math.radians(theta))
d_y = d * math.sin(math.radians(theta))
return [d, d_x, d_y]

def getAcc(tmass,th,v,velY,altY,t,h,parachuteFlag,ascentFlag):
# In the beginning, acceleration is -9.8(Pointing downwards), and then it gradually increases. Until, it reaches a +ve value,
Expand All @@ -134,41 +138,54 @@ def getAcc(tmass,th,v,velY,altY,t,h,parachuteFlag,ascentFlag):
if parachuteFlag == 0 and len(velY) > 0 and velY[-1]<0 and len(altY) > 0 and altY[-1]>0:
parachuteFlag = 1
if parachuteFlag == 1: # After parachute is deployed
d = getParachuteDrag(v) + getDrag(v)
d = getParachuteDrag(v)[0] # + getDrag(v)
dx = getParachuteDrag(v)[1]
dy = getParachuteDrag(v)[2]
thh = getThrust(t)
acceleration = (d-(tmass*g))/tmass # Sin theta not considered because, once the parachute is deployed, the rocket is going to fall straight down. So the angle is not going to be considered.
thh_x = thh*math.cos(math.radians(theta))
thh_y = thh*math.sin(math.radians(theta))
acceleration_x = 0 # wind velocity requires adding here
acceleration_y = ((d - (tmass * g)) / tmass) # Sin theta not considered because, once the parachute is deployed, the rocket is going to fall straight down. So the angle is not going to be considered.
acceleration = acceleration_y # for now, since acc along x is 0 without wind
print("Status : Parachute Deployed(Descent)","Time : ",t,"Thrust : ",thh,"Drag : ", d,"Weight : ",tmass*g,"Force : ",(d - (tmass*g) ),"Mass : ",tmass,"Acceleration",acceleration)
else: # Before Parachute is deployed
# If parachute hasn't deployed, then either the rocket is still on the launch pad, or it is in the ascent phase.
thh = getThrust(t)
d = getDrag(v)
acceleration = (thh - ((tmass*g)*math.sin(math.radians(theta))) - d)/tmass # math.sin considers radians, hence theta is first converted to radians, then passed into math.sin
thh_x = thh*math.cos(math.radians(theta))
thh_y = thh*math.sin(math.radians(theta))
d = getDrag(v)[0]
dx = getDrag(v)[1]
dy = getDrag(v)[2]

acceleration_x = ((thh_x - dx)/tmass) # -+ Wind can also be added
acceleration_y = ((thh_y - dy - (tmass * g))/tmass) # ERROR HERE : ACC_Y coming as negative of Acc for some reason
acceleration = ((thh - ((tmass * g) * math.sin(math.radians(theta))) - d) / tmass)# math.sin considers radians, hence theta is first converted to radians, then passed into math.sin
# The reason, why we are considering mg.sin(theta) here is because, we are considering that the rocket is going to be launched at a launch angle theta. So, when calculating the net force, at that angle, mg.sin(theta) is used.
print("Status : Ascent, Theta : ",theta,"Time : ",t,"Thrust : ",thh,"Drag : ", d,"Weight : ",((tmass*g)*math.sin(math.radians(theta))),"Force : ",(thh - ((tmass*g)*math.sin(math.radians(theta))) - d),"Mass : ",tmass,"Acceleration",acceleration)

if ascentFlag == 0 and acceleration >= 0: # If the rocket hasn't begun its ascent phase and the acceleration is +ve, meaning ascent phase is beginning
ascentFlag = 1 # Ascent has begun
elif ascentFlag == 1 and acceleration<0: #If the rocket has finished its ascent phase and the acceleration becomes -ve, meaning pointing downwards, means that the ascent phase has ended and the descent phase has begun
elif (ascentFlag == 1) and (acceleration < 0): #If the rocket has finished its ascent phase and the acceleration becomes -ve, meaning pointing downwards, means that the ascent phase has ended and the descent phase has begun
ascentFlag = 2 # Ascent has ended and Descent has begun
return [acceleration,ascentFlag]
return [acceleration, ascentFlag, acceleration_x, acceleration_y]

def isCondition(altY,t,parachuteFlag): # Condition to check whether the rocket has landed or not
if len(altY) > 0 and altY[-1]<=0 and t > 0 and parachuteFlag == 1:
def isCondition(altY, t, parachuteFlag): # Condition to check whether the rocket has landed or not
if (len(altY) > 0) and (altY[-1] <= 0) and (t > 0) and (parachuteFlag == 1):
return 0
return 1

while isCondition(altY,t,parachuteFlag1):
while isCondition(altY, t, parachuteFlag1):
print("*(***************************",len(velY),len(altY))
if len(velY) != 0 and len(altY) != 0:
print("VelY : ", velY[-1], " AltY : ", altY[-1],"Parachute Flag : ", parachuteFlag1)
if parachuteFlag1 == 0 and velY[-1]<0 and altY[-1]>0: # Apogee condition, Checking if the last vertical velocity < 0 , and last vertical altitude is > 0
if parachuteFlag1 == 0 and velY[-1] < 0 and altY[-1] > 0: # Apogee condition, Checking if the last vertical velocity < 0 , and last vertical altitude is > 0
parachuteFlag1 = 1 #Parachute is deployed
t2apogee = t # Time to apogee set at apogee

# Gravity modelling(Change in acceleration due to gravity)
# g = (G*ME)/(h + Re)
try:
theta = math.degrees(math.atan(velY[-1]/velX[-1])) # Theta is the angle at which the rocket is going to be launched.
theta = math.degrees(math.atan(velY[-1] / velX[-1])) # Theta is the angle at which the rocket is going to be launched.
print("theta:",theta)
except:
print("Exception in theta calculations")
Expand All @@ -177,11 +194,12 @@ def isCondition(altY,t,parachuteFlag): # Condition to check whether the rocket h
# Accel, Vel, Height calculations
inp = getAcc(tmass,th,v,velY,altY,t,h,parachuteFlag1,ascentFlag)
acc = float(inp[0])
acc_x = float(inp[2])
acc_y = float(inp[3])
ascentFlag = int(inp[1])
if ascentFlag == 0:
v = 0
h = 0

velY.append(v)
accelY.append(acc)
else:
Expand All @@ -190,8 +208,6 @@ def isCondition(altY,t,parachuteFlag): # Condition to check whether the rocket h
# Breaking acc into its x and y vectors
theta = 90 # Considering that when the parachute is deployed, the rocket will fall straight down. So theta is 90 degrees from x axis and 0 from y
# acc_x = acc * math.cos(math.radians(theta))
acc_x = 0
acc_y = acc * math.sin(math.radians(theta))

print("COS 90 : " , math.cos(math.radians(90)),"SIN 90 : ", math.sin(math.radians(90)),"ACC_X : ",acc_x,"ACC_Y : ",acc_y)

Expand All @@ -210,14 +226,14 @@ def isCondition(altY,t,parachuteFlag): # Condition to check whether the rocket h
# vy = v_y + acc_y * dt

v += acc_y * dt
v_y = v
v_y = v

print("++++++++++++++++++v_y : ",vy,"acc_y : ",acc_y ,"dt : ",dt,"acc_y * dt : ",acc_y * dt,"v_y + acc_y * dt : ",v_y + acc_y * dt)

# v = math.sqrt(pow(vx,2) + pow(vy,2))
# print(acc,dt,pow(dt,2),acc*pow(dt,2))

h += v_y*dt + (0.5*acc_y*pow(dt,2)) # Because altitude is basically distance in the y axis
h += v_y * dt + (0.5 * acc_y * pow(dt,2)) # Because altitude is basically distance in the y axis
hx += v_x * dt

altX.append(hx)
Expand All @@ -238,8 +254,8 @@ def isCondition(altY,t,parachuteFlag): # Condition to check whether the rocket h

else: # Parachute hasn't been deployed yet
# Breaking acc into its x and y vectors
acc_x = acc * math.cos(math.radians(theta))
acc_y = acc * math.sin(math.radians(theta))
# acc_x = acc * math.cos(math.radians(theta))
# acc_y = acc * math.sin(math.radians(theta))

# Breaking previous velocity and altitude into its x and y vectors
v_x = v * math.cos(math.radians(theta)) # Previous velocity x components
Expand All @@ -249,8 +265,9 @@ def isCondition(altY,t,parachuteFlag): # Condition to check whether the rocket h
vx = v_x + acc_x * dt
vy = v_y + acc_y * dt
v = math.sqrt(pow(vx,2) + pow(vy,2))
h += v_y*dt + (0.5*acc_y*pow(dt,2)) # Because altitude is basically distance in the y axis
hx += v_x * dt
h += v_y * dt + (0.5 * acc_y * pow(dt, 2)) # Because altitude is basically distance in the y axis
hx += v_x * dt + (0.5 * acc_x * pow(dt, 2))


altX.append(hx)
altY.append(h)
Expand Down
1 change: 1 addition & 0 deletions tempCodeRunnerFile.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
,