[关闭]
@zy-0815 2016-10-23T14:05:56.000000Z 字数 9011 阅读 886

张梓桐计算物理第六次作业

1.Problem

2.10

2.Abstract

In this passage,we are going to discuss the track of a cannon shell with adbiatic effect and wind resistance(headwind,or tailwind).Additionally,a calculation of minimun velocity with regard to a specific target is presented too.

3.Background Introduction

With the previous knowledge of Euler Method,we don't have to illustrate again the basic knowledge.What is of great importance to us is the additional wind resistance.Here we set the velocity of wind to be ,and the veloocity of cannon to be .
The components of the drag force can be represented as:



Noticing that the drag coefficient regarding to the velocity of cannon is presented in the follow figure.

http___i1.piimg.com_567571_5e3234d6f661c3f0.png-60.2kB
The analyatic form of can be expressed as:


Here and are constant, and .
On top of that,there will be a effect of adbiatic air effect with regard to the altitude:

In the situation of adbiatic effect:

which means there will be a additional coefficient plus the wind resistance force.


Having discussed all of the concerning factors,we finally deduce the core calculation equation:



4.Main body

4.1 The track of the cannon shell.

  1. import math
  2. import pylab as pl
  3. a=1
  4. class cannon:
  5. def __init__(self,intial_velocity,intial_angle,velocity_of_wind,angle_of_wind):#Variable that can be set#
  6. self.x=[0]
  7. self.y=[0]
  8. time_step=0.01
  9. mass=10
  10. intial_vx=intial_velocity*math.cos(int(intial_angle)*math.pi/180)
  11. intial_vy=intial_velocity*math.sin(int(intial_angle)*math.pi/180)
  12. self.vx=[intial_vx]
  13. self.vy=[intial_vy]
  14. self.v=[math.sqrt(intial_vx**2+intial_vy**2)]
  15. self.dt=time_step
  16. self.intial_angle=intial_angle
  17. self.m=mass
  18. self.angle_of_wind=angle_of_wind
  19. self.velocity_of_wind=velocity_of_wind
  20. self.B2=[self.m*(0.0039+0.0058/(1+math.exp((self.v[0]-35)/5)))]
  21. self.delta_v_mold=[math.sqrt((intial_vx-self.velocity_of_wind*math.cos(self.angle_of_wind*math.pi/180))**2+(intial_vy-self.velocity_of_wind*math.sin(self.angle_of_wind*math.pi/180))**2)]
  22. def run(self):
  23. i=0
  24. while self.y[i]>=0:
  25. self.B2.append(self.m*(0.0039+0.0058/(1+math.exp((self.v[i]-35)/5))))
  26. self.v.append(math.sqrt(abs(self.vx[i])**2+abs(self.vy[i])**2))
  27. self.delta_v_mold.append(math.sqrt((self.vx[i]-abs(self.velocity_of_wind*math.cos(self.angle_of_wind*math.pi/180)))**2+(self.vy[i]-self.velocity_of_wind*math.sin(self.angle_of_wind*math.pi/180))**2))
  28. self.x.append(self.x[i]+self.vx[i]*self.dt)
  29. self.y.append(self.y[i]+self.vy[i]*self.dt)
  30. self.vx.append(self.vx[i]-pow(1-0.0065*self.y[i]/273,2.5)*self.B2[i]*self.delta_v_mold[i]*(self.vx[i]-self.velocity_of_wind*math.cos(self.angle_of_wind*math.pi/180))/self.m*self.dt)
  31. self.vy.append(self.vy[i]-10*self.dt-pow(1-0.0065*self.y[i]/273,2.5)*self.B2[i]*self.delta_v_mold[i]*(self.vy[i]-self.velocity_of_wind*math.sin(self.angle_of_wind*math.pi/180))/self.m*self.dt)
  32. i+=1
  33. self.x[-1]=(self.x[-1] * self.y[-2] - self.x[-2] * self.y[-1])/ (self.y[-2] - self.y[-1])figure_2.png-74.6kB
  34. def show(self):
  35. pl.plot(self.x,self.y,label='$\\theta=$'+str(self.intial_angle)+'$^{\circ}$'+'.x='+str(int(self.x[-1]))+',Wind velocity '+str(self.velocity_of_wind)+'\n'+'Angle of wind is '+str(self.angle_of_wind)+'$^{\circ}$')
  36. pl.ylim(0,300)
  37. pl.xlim(0,400)
  38. pl.xlabel('x(m)')
  39. pl.ylabel('y(m)')
  40. pl.show()
  41. pl.legend(fontsize='xx-small')
  42. class show_trajectory_of_various_angle(cannon):
  43. def run2(self):
  44. theta=0
  45. while theta<90:
  46. a=cannon(intial_angle=theta,intial_velocity=110,velocity_of_wind=10,angle_of_wind=135)
  47. theta+=1
  48. a.run()
  49. a.show()
  50. class show_trajectory_of_various_wind_velocity(cannon):
  51. def run3(self):
  52. v=0
  53. while v<30:
  54. b=cannon(intial_velocity=110,intial_angle=45,velocity_of_wind=v,angle_of_wind=135)
  55. v+=1
  56. b.run()
  57. b.show()
  58. class show_trajectory_of_various_angle_of_wind(cannon):
  59. def run4(self):
  60. theta2=45
  61. while theta2 in [45,135]:
  62. c=cannon(intial_velocity=244,angle_of_wind=theta2,intial_angle=35,velocity_of_wind=22.27)
  63. theta2+=135
  64. c.run()
  65. c.show()
  66. a=show_trajectory_of_various_angle(110,45,10,135)
  67. a.run2()

4.2 Find the minimun firing velocity

  1. import math
  2. import pylab as pl
  3. import time
  4. a=1
  5. class cannon:
  6. def __init__(self,intial_velocity,intial_angle,velocity_of_wind,angle_of_wind):#Variable that can be set#
  7. self.x=[0]
  8. self.y=[0]
  9. time_step=0.01
  10. mass=10
  11. intial_vx=intial_velocity*math.cos(int(intial_angle)*math.pi/180)
  12. intial_vy=intial_velocity*math.sin(int(intial_angle)*math.pi/180)
  13. self.vx=[intial_vx]
  14. self.vy=[intial_vy]
  15. self.v=[math.sqrt(intial_vx**2+intial_vy**2)]
  16. self.dt=time_step
  17. self.intial_angle=intial_angle
  18. self.m=mass
  19. self.angle_of_wind=angle_of_wind
  20. self.velocity_of_wind=velocity_of_wind
  21. self.B2=[self.m*(0.0039+0.0058/(1+math.exp((self.v[0]-35)/5)))]
  22. self.delta_v_mold=[math.sqrt((intial_vx-self.velocity_of_wind*math.cos(self.angle_of_wind*math.pi/180))**2+(intial_vy-self.velocity_of_wind*math.sin(self.angle_of_wind*math.pi/180))**2)]
  23. def run(self):
  24. i=0
  25. while self.y[i]>=-100:
  26. self.B2.append(self.m*(0.0039+0.0058/(1+math.exp((self.v[i]-35)/5))))
  27. self.v.append(math.sqrt(abs(self.vx[i])**2+abs(self.vy[i])**2))
  28. self.delta_v_mold.append(math.sqrt((self.vx[i]-abs(self.velocity_of_wind*math.cos(self.angle_of_wind*math.pi/180)))**2+(self.vy[i]-self.velocity_of_wind*math.sin(self.angle_of_wind*math.pi/180))**2))
  29. self.x.append(self.x[i]+self.vx[i]*self.dt)
  30. self.y.append(self.y[i]+self.vy[i]*self.dt)
  31. self.vx.append(self.vx[i]-pow(1-0.0065*self.y[i]/273,2.5)*self.B2[i]*self.delta_v_mold[i]*(self.vx[i]-self.velocity_of_wind*math.cos(self.angle_of_wind*math.pi/180))/self.m*self.dt)
  32. self.vy.append(self.vy[i]-10*self.dt-pow(1-0.0065*self.y[i]/273,2.5)*self.B2[i]*self.delta_v_mold[i]*(self.vy[i])/self.m*self.dt)
  33. i+=1
  34. self.x[-1]=((-100-self.y[-1])*(self.x[-1]-self.x[-2])/(self.y[-1]-self.y[-2]))+self.x[-1]
  35. self.y[-1]=-100
  36. global a
  37. a=self.x[-1]
  38. def show(self):
  39. pl.plot(self.x,self.y,'r',label='$\\theta=$'+str(self.intial_angle)+'$^{\circ}$'+'.x='+str(int(self.x[-1]))+',Wind velocity '+str(self.velocity_of_wind)+'\n'+'Angle of wind is '+str(self.angle_of_wind)+'$^{\circ}$')
  40. pl.ylim(-100,300)
  41. pl.xlabel('x(m)')
  42. pl.ylabel('y(m)')
  43. pl.show()
  44. pl.legend(fontsize='xx-small')
  45. class find_spot(cannon):
  46. object_x=220
  47. object_y=-100
  48. def run4(self):
  49. v=67.7230
  50. while 67.7230<=v<68:
  51. theta2=14.99
  52. while 14.99<=theta2<20:
  53. c=cannon(intial_velocity=v,angle_of_wind=135,intial_angle=theta2,velocity_of_wind=10)
  54. c.run()
  55. if 220>a:
  56. theta2+=0.0001
  57. else:
  58. print theta2,'degree',a,'m',v,'m/s'
  59. break
  60. v+=0.00001
  61. if 220<a:
  62. break
  63. a=find_spot(110,45,10,135)
  64. a.run4()

5.Track Result

5.1 Result with different intial angles

Here we set the intial velocity all to be .The angle of wind is (headwind),and the velocity of wind is

Now we change the intial firing angles from with

figure_1.png-312kB

As this figure is not too compact and difficult to distinguish,we now set
figure_2.png-74.6kB

Here we can see,around ,the land point reachs max.
Now we take a closer look at around .
figure_3.png-64.9kB
Here we can see, at around ,the land point reaches the maxium.And .
figure_4.png-55.6kB

5.2 Results with different wind angles

Here we set the intial velocity all to be .The angle of wind is (from tailwind to headwind),,and the velocity of wind is .Intial angles is .
figure_5.png-160.2kB
Here we can see,as the angles of wind is changed from tailwind to headwind, the land point is shifting left,which means reach a lesser distance.

To illustrate the difference with or withoud wind,I am going to attach another figure.The angle of wind is and
figure_6.png-42.5kB
This figure corresponds well to our feeling in our life!

5.3 Results with different velocity of wind

Here we set the intial velocity all to be .The angle of wind is (headwind),and the velocity of wind is ,. Intial angles is .
figure_7.png-163kB
As the wind velocity is becoming larger and larger,there may be situation where the cannon shell is blown back.

6.Find the target

In this section,we are going to discuss about how to fire a cannon shell with the minimun velocity to reach the target.At first,we set the objec to be set in ,which means the target is lower than the firing places.
However,when we are about to set the scanning interval of and velocity with and .The time taken by the calculation is too long.For compromising,we set the and to identify the mighty intervals,and we find around and velocity about is the ideal result.
Step by step,we narrow the section and finally finds out that at and ,the velocity reaches the minmum.
figure_1.png-22.8kB

7.Acknowledgement

1.Yu kang has helped with the find-spot programming.
2.Matplotlib.

添加新批注
在作者公开此批注前,只有你和作者可见。
回复批注