Thursday, February 2, 2017

Function to Simulate Parabolic Shot with Drag.

Script to calculate the most important quantitative information of the drag parabolic shot in International System of Units.
If you want to use another system of units you can do it by making some simple changes.

This is a first approach to the problem and is totally perfectible.

To calculate the maximum range, a transcendental equation is solved regularly, but here we did a first run with a number of time values left over to approximate the time index for which the vertical position becomes zero.

Inputs of the function are:
a) initial velocity vo (scalar) [m/s]
b) shot angle alpha [degrees]
c) (drag coefficient / mass) = b [m^-1]

Parameter
g) gravity acceleration g = 9.81 [m/s^2]

The outputs of the function are:
T_1: ascending time [s]
H: maximum height [m]
L: maximum horizontal range [m]

We need to calculate:
x: horizontal position at time t [m]
y: vertical position at time t [m]
vox: horizontal initial velocity [m/s]
voy: vertical initial velocity [m/s]

drag_parabolic <- function(vo, alpha, b){
 
  g <- 9.81
  an <- (2*pi*alpha)/360
  vox <- vo*cos(an)
  voy <- vo*sin(an)
  T_1 <- (1/b)*log(1+(b*voy/g))
  H <- (voy/b)-((g/b^2)*log(1+(b*voy/g)))
 
 
# first run

  t <- seq(0, 25, 1/10)

  x <- vector()
  y <- vector()
 
  x <- (vox/b)*(1-exp(-b*t))
  y <- ((1/b)*((g/b)+voy)*(1-exp(-b*t)))-((g/b)*t)
 
 
# time index of max. range

  a_T <- which(y < 0)[1]-1
 
 
# second run

  t_t <- head(t, a_T)
 
  xx <- vector()
  yy <- vector()
 
  xx <- (vox/b)*(1-exp(-b*t_t))
  yy <- ((1/b)*((g/b)+voy)*(1-exp(-b*t_t)))-((g/b)*t_t)
 
  R <- round(xx[length(xx)],2)
  H <- round(H,2)
  

  plot(xx, yy, xlab="X", ylab="Y", type = "o", col = "blue", axes=F)
  

  axis(1, at = seq(0,R,R/10),labels = seq(0,R,R/10), cex.axis = 0.7)
  axis(2, at = seq(0,H,H/10),labels = seq(0,H,H/10), cex.axis = 0.7)
 
  print("Initial Velocity");print(paste(vo,"m/s"))
  print("Angle of Shot");print(paste(alpha,"degrees"))
  print("Ascending Time");print(paste(T_1,"s"))
  print("Maximum Height");print(paste(H,"m"))
  print("Aprox. Max. Range");print(paste(R,"m", "+-2%"))
 
  legend(R/3, H/2, legend = c(paste("vo =", vo, "m/s"),
                   paste("alpha =", alpha,"degrees"),
                   paste("Ascending time", paste(T_1,"s")),
                   paste("Maximum height", paste(H,"m")),
                   paste("Aprox. Max. Range", paste(R,"m","+-2%"))),
                   cex=0.7, bg = par("bg"))
  

title(main = "Drag Parabolic Shot", sub = "")
}


Let's try the function:

drag_parabolic(50, 45, 0.5)


 We compare with a simple parabolic shot with the same parameters:


 


I will be glad to receive your comments and suggestions to improve the script.
 
Get the script in:
https://github.com/pakinja/Data-R-Value 

Tuesday, January 31, 2017

Function to Simulate Simple Parabolic Shooting.

This function simulates a parabolic shot from the origin and does not take into account the friction of the air.

# Script to calculate the most important info about the parabolic
# shot without friction and in International System of Units

# if you want to use another system of units only change the value
# of g and the units in the prints and legends

# Inputs of the function are:
# a) initial velocity vo
# b) shot angle alfa

# Parameter
# c) gravity acceleration g = 9.81 m/s^2

# The outputs of the function are:
# T_1: ascending time
# T_2: descending time
#   H: maximum height
#  Tt: total time
#   L: maximum horizontal range


parabolic <- function(vo, alfa){
 
  g <- 9.81
 an <- round((2*pi*alfa)/360,2)
T_1 <- round(vo*sin(an)/g,2)
  H <- round((vo^2)*(sin(an)^2)/(2*g),2)
T_2 <- round(vo*sin(an)/g,2)
T_t <- round(2*vo*sin(an)/g,2)
  L <- round(vo^2*sin(2*an)/g,2)
 
  print("Ascending time");print(paste(T_1,"s"))
  print("Maximum height");print(paste(H,"m"))
  print("Descending time");print(paste(T_2, "s"))
  print("Total time");print(paste(T_t,"s"))
  print("Max. horizontal range");print(paste(L,"m"))
 
  y <- vector()
  p <- seq(0.0, round(L), round(L)/100)
 
  y <- (tan(an)*p)-((g/((2*vo^2)*cos(an)^2))*(p^2))
   
  plot(p, y, xlab="X", ylab="Y", type = "o", col = "red", axes=F)
  axis(1, at = seq(0,L,L/10),labels=seq(0,L,L/10),
       cex.axis=0.7)
  axis(2, at = seq(0,H,H/100),labels=seq(0,H,H/100),
       cex.axis=0.7)
  legend(L/4, H/3, legend = c(paste("vo =", vo, "m/s"),
                paste("alfa =", alfa,"degrees"),
                paste("Ascending time", paste(T_1,"s")),
                paste("Maximum height", paste(H,"m")),
                paste("Descending time", paste(T_2, "s")),
                paste("Total time", paste(T_t,"s")),
                paste("Max. horizontal range", paste(L,"m"))),
         cex=0.7, bg = par("bg"))
  title(main = "Parabolic Shot", sub = "From Origin & Frictionless")

}



We test the function with a shot of 100 m / s of initial velocity and 45 degrees:



parabolic(100,45)


In a following post we will do the simulation but with the effect of the air friction.

If you want to use the function in another system of units, you only have to change the value of the acceleration of gravity g.

We eliminate some for loop to generate the trajectory -y- for been useless and for the kind recommendation of Andrej-Nikolai Spiess.

Get the script in:
https://github.com/pakinja/Data-R-Value