The trajectory of a cannon shell

Calculate the trajectory of our cannon shell including both air drag and the reduced air density at high altitudes so that you can reproduce the results in the figure.Perform your calculation for different firing angles and determine the value of the angle that gives the maximum range.

Projectile motion

  • Newton's second law in two spatial dimensions



(1)

Euler method for second-order ODE

Resistance

  • the magnitude of the drag force is given by



    which
  • the componets of the drag force



(2)

How to decide the point of fall?

  • the last point above the ground(n) and the first point below the ground(n+1) is given by



Code(python)

  • without air resistance

import numpy as np
import math
import matplotlib.pyplot as pl
class cannon_shell:
    def __init__(self,time_step=0.1,total_time=30,gravity=9.8,initial_dis=0):
        print("enter the angle of the cannon shell ->")
        self.rad=float(input())
        self.angle=(self.rad/180)*math.pi
        print("enter the speed of the cannon shell ->")
        self.velocity=input()
        self.g=gravity
        self.dt=time_step
        self.v_x=[self.velocity*(math.cos(self.angle))]
        self.v_y=[self.velocity*(math.sin(self.angle))]
        self.x=[initial_dis]
        self.y=[initial_dis]
    def run(self):
        while(self.y[-1]>=0):
            self.x.append(self.x[-1]+self.v_x[-1]*self.dt)
            self.y.append(self.y[-1]+self.v_y[-1]*self.dt)
            self.v_x.append(self.v_x[-1])
            self.v_y.append(self.v_y[-1]-self.g*self.dt)
    def show_result(self):
        pl.plot(self.x,self.y)
        pl.title('connon shell without air resistance')
        pl.xlabel('x/R')
        pl.ylabel('y/H')
        pl.ylim(0.0)
        pl.show()
a=cannon_shell()
a.run()
a.show_result()
  • consider the resistance of air

import numpy as np
import math
import matplotlib.pyplot as pl
class cannon_shell:
    def __init__(self,time_step=0.1,resis_coeff=4e-5,total_time=30,gravity=9.8,initial_dis=0):
        print("enter the angle of the cannon shell ->")
        self.rad=float(input())
        self.angle=(self.rad/180)*math.pi
        print("enter the speed of the cannon shell ->")
        self.v=input()
        self.B2_m=resis_coeff
        self.g=gravity
        self.dt=time_step
        self.v_x=[self.v*(math.cos(self.angle))]
        self.v_y=[self.v*(math.sin(self.angle))]
        self.x=[initial_dis]
        self.y=[initial_dis]
    def run(self):
        while(self.y[-1]>=0):
            self.x.append(self.x[-1]+self.v_x[-1]*self.dt)
            self.y.append(self.y[-1]+self.v_y[-1]*self.dt)
            self.v_x.append(self.v_x[-1]-self.B2_m*self.v*self.v_x[-1]*self.dt)
            self.v_y.append(self.v_y[-1]-self.g*self.dt-\
            self.B2_m*self.v*self.v_y[-1]*self.dt)
    def show_result(self):
        pl.plot(self.x,self.y)
        pl.title('connon shell without air resistance')
        pl.xlabel('x/R')
        pl.ylabel('y/H')
        pl.ylim(0.0)
        pl.show()
a=cannon_shell()
a.run()
a.show_result()
(3)

Effects of altitude

  • isothermal (constant temperature) ideal gas

=p(0)e^{-mgy/k_BT})

which
![](http://latex.codecogs.com/png.latex?y_0=k_BT/mg \approx 1.0 \times 10^4m)

  • adiabatic approximation

    which
    ![](http://latex.codecogs.com/png.latex?a\approx6.5 \times 10^{-3}K/m,\alpha \approx2.5)
  • the drag force due to air resistance is proportional to the density


consider the reduce of air density

import numpy as np
import math
import matplotlib.pyplot as pl
class cannon_shell:
    def __init__(self,time_step=0.1,resis_coeff=4e-5,total_time=30,gravity=9.8,initial_dis=0):
        print("enter the angle of the cannon shell ->")
        self.rad=float(input())
        self.angle=(self.rad/180)*math.pi
        self.v=700
        self.B2_m=resis_coeff
        self.g=gravity
        self.dt=time_step
        self.v_x=[self.v*(math.cos(self.angle))]
        self.v_y=[self.v*(math.sin(self.angle))]
        self.x=[initial_dis]
        self.y=[initial_dis]
    def run(self):
        while(self.y[-1]>=0):
            self.rho=(1-(2.257e-5)*self.y[-1])**2.5
            self.x.append(self.x[-1]+self.v_x[-1]*self.dt)
            self.y.append(self.y[-1]+self.v_y[-1]*self.dt)
            self.v_x.append(self.v_x[-1]-self.rho*self.B2_m*self.v*self.v_x[-1]*self.dt)
            self.v_y.append(self.v_y[-1]-self.g*self.dt-\
            self.rho*self.B2_m*self.v*self.v_y[-1]*self.dt)
        labeltext=str(self.rad)
        pl.plot(self.x,self.y,label=labeltext)
def show_result():
    pl.title('connon shell without air resistance')
    pl.xlabel('x/m')
    pl.ylabel('y/m')
    pl.ylim(0.0)
    pl.legend(loc='upper left')
    pl.show()
a=cannon_shell()
b=cannon_shell()
c=cannon_shell()
d=cannon_shell()
e=cannon_shell()
f=cannon_shell()
a.run()
b.run()
c.run()
d.run()
e.run()
f.run()
show_result()

show the result

(4)
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
平台声明:文章内容(如有图片或视频亦包括在内)由作者上传并发布,文章内容仅代表作者本人观点,简书系信息发布平台,仅提供信息存储服务。

推荐阅读更多精彩内容

  • 早上起床,就看到一则泪目的新闻《10岁男童临死前,请求爸爸签遗体自愿捐让书》。 懂事,对于一个孩子来说真是一个残忍...
    钧吾戏研阅读 450评论 2 1
  • 大学放寒假的时候, 我跟几个同样喜欢写作的朋友相约去深圳看海。我们到达深圳的时候恰逢下雨,入住的酒店刚好是繁...
    寒时未央阅读 432评论 2 4
  • 刚刚经历丢手机风波。 一起丢的还有之前很多的回忆,不过这样也好,从今天开始,做个没有故事的女同学,我希望你简单,快...
    枫叶子飘飘阅读 183评论 1 1
  • 看了林更新的新作品,《枪手快手快枪手》,要说剧情也一般,有点侠义江湖情节,故事里穿插了小庄和若儿的爱情,为这部戏增...
    64b79db7c71b阅读 278评论 0 1
  • 写于16年10月25日 现象: 目前规划了10km跑步计划,在5km跑步时遇到了吃力的现象 腹部上身有明显赘肉 大...
    会咬人的春虫虫阅读 337评论 0 0