弹丸运动：结果显示有空气阻力的线和没有空气阻力的线之间存在很大差异

Question

我在 Spyder 上用 python 模拟弹丸运动，一个有空气阻力，一个没有空气阻力。我开始关注this exercise。 它引导我在空气阻力下进行抛射运动。为了进行比较，我尝试用我自己的相同参数在同一张图中绘制没有空气阻力的弹丸运动。然而，结果图让我感到惊讶，因为它看起来有点不对劲。就范围而言，这两条线之间存在很大差异，因为我在网上遇到的其他类似图表仅显示出细微的差异。

我的图表：

（小的是有拖拽的，大的是没有拖拽的）

我的问题是：

1. 这是什么原因造成的？是不是因为我放的阻力只包括阻力系数和速度，所以图比还包括横截面积和rho的图要小？

2. 我做的方式是不是很复杂，你会如何修改它以使其更高效和整洁？

如果你也愿意指出我犯的任何错误，那就太好了

import numpy as np
import matplotlib.pyplot as plt

# Model parameters
M = 1.0          # Mass of projectile in kg
g = 9.8          # Acceleration due to gravity (m/s^2)
V = 80           # Initial velocity in m/s
ang = 60.0       # Angle of initial velocity in degree
Cd = 0.005       # Drag coefficient
dt = 0.5         # time step in s

# Set up the lists to store variables
# Start by putting the initial velocities at t=0
t = [0]                         # list to keep track of time
vx = [V*np.cos(ang/180*np.pi)]  # list for velocity x and y components
vy = [V*np.sin(ang/180*np.pi)]

# parameters for the projectile motion without drag force
t1=0
vx_nodrag=V*np.cos(ang/180*np.pi)
vy_nodrag=V*np.sin(ang/180*np.pi)

while (t1 < 100):
    x_nodrag=vx_nodrag*t1
    y_nodrag=vy_nodrag*t1+(0.5*-9.8*t1**2)
    plt.ylim([0,500])
    plt.xlim([0,570])
    plt.scatter(x_nodrag, y_nodrag)
    print(x_nodrag,y_nodrag)
    t1=t1+dt

# Drag force
drag = Cd*V**2                      # drag force 

# Create the lists for acceleration components
ax = [-(drag*np.cos(ang/180*np.pi))/M]        
ay = [-g-(drag*np.sin(ang/180*np.pi)/M)]

# Use Euler method to update variables
counter = 0
while (counter < 100):
    t.append(t[counter]+dt)            # increment by dt and add to the list of time 
    vx.append(vx[counter]+dt*ax[counter])  
    vy.append(vy[counter]+dt*ay[counter])  

    # With the new velocity calculate the drag force
    vel = np.sqrt(vx[counter+1]**2 + vy[counter+1]**2)   
    drag = Cd*vel**2                                  
    ax.append(-(drag*np.cos(ang/180*np.pi))/M)    
    ay.append(-g-(drag*np.sin(ang/180*np.pi)/M))    
    
    # Increment the counter by 1
    counter = counter +1

x=[0]#creating a list for x
y=[0]#creating a list for y

counter1=0
while (counter1<50):

    #t.append(t[counter1]+dt),t already has a list.
    x.append(x[counter1]+dt*vx[counter1])    
    y.append(y[counter1]+dt*vy[counter1])  
    plt.ylim([0,500])
    plt.xlim([0,570])
    plt.plot(x,y)
    #print(x,y)
    counter1=1+counter1

# Let's plot the trajectory
plt.plot(x,y,'ro')
plt.ylabel("height")
plt.xlabel("range")
print("Range of projectile is {:3.1f} m".format(x[counter]))

Answer 1

1. 是的，您肯定必须包括射弹的所有参数，因为总行进距离对此非常敏感。
2. 通常尝试使用可以调用 for 循环的函数。这使您的代码更易于阅读和更可重用。例如，实现一个计算新时间步长的函数。还可以将 numpy 数组用作向量，而不是分别使用 vx 和 vy。

此外，您的计算中存在错误。在计算新速度时，您还必须计算更新的角度。这是因为在每一步你都有新的角度（即在飞行结束时，y 方向的阻力与重力的符号相反）

例如：

    vel = np.sqrt(vx[counter+1]**2 + vy[counter+1]**2)   
    drag = Cd*vel**2


    # calculate new angle for force separation
    a, v = np.array([1,0]), np.array([vx[counter+1], vy[counter+1]])
    if vy[counter+1] > 0:
        ang = np.arccos( np.dot(a,v) / (1*np.linalg.norm(v)) )
    else:
        ang = 2*np.pi - np.arccos( np.dot(a,v) / (1*np.linalg.norm(v)) )   

                           
    ax.append(-(drag*np.cos(ang))/M)    
    ay.append(-g-(drag*np.sin(ang)/M))