如何在重力,浮力和空气阻力的作用下绘制弹丸的运动? [英] How to plot the motion of a projectile under the effect of gravity, buoyancy and air resistance?

查看:64
本文介绍了如何在重力,浮力和空气阻力的作用下绘制弹丸的运动?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我试图绘制在重力,浮力和阻力作用下的质点的弹丸运动图.基本上,我想显示浮力和阻力对飞行距离,飞行时间和速度变化对绘图的影响.

I am trying to make a plot of a projectile motion of a mass which is under the effect of gravitational, buoyancy and drag force. Basically, I want to show effects of the buoyancy and drag force on flight distance, flight time and velocity change on plotting.

import matplotlib.pyplot as plt
import numpy as np

V_initial = 30 # m/s
theta = np.pi/6 # 30
g = 3.711
m =1
C = 0.47
r = 0.5
S = np.pi*pow(r, 2)
ro_mars = 0.0175
t_flight = 2*(V_initial*np.sin(theta)/g)
t = np.linspace(0, t_flight, 200)

# Drag force
Ft = 0.5*C*S*ro_mars*pow(V_initial, 2)

# Buoyant Force
Fb = ro_mars*g*(4/3*np.pi*pow(r, 3))

x_loc = []
y_loc = []

for time in t:
    x = V_initial*time*np.cos(theta)
        y = V_initial*time*np.sin(theta) - (1/2)*g*pow(time, 2)
    x_loc.append(x)
    y_loc.append(y)

x_vel = []
y_vel = []
for time in t:
    vx = V_initial*np.cos(theta)
    vy = V_initial*np.sin(theta) - g*time
    x_vel.append(vx)
    y_vel.append(vy)


v_ch = [pow(i**2+ii**2, 0.5) for i in x_vel for ii in y_vel]

tau = []
for velocity in v_ch:
    Ft = 0.5*C*S*ro_mars*pow(velocity, 2)
        tau.append(Ft)

buoy = []
for velocity in v_ch:
    Fb = ro_mars*g*(4/3*np.pi*pow(r, 3))
    buoy.append(Fb)

在此之后,我无法弄清楚如何在这种力下绘制弹丸运动.换句话说,我正在尝试比较三种情况下质量的弹丸运动

after this point, I couldn't figure out how to plot to a projectile motion under this forces. In other words, I'm trying to compare the projectile motion of the mass under three circumstances

  1. 仅在重力作用下的质量
  2. 在重力和空气阻力作用下的质量
  3. 在重力,空气阻力和浮力作用下的质量

推荐答案

您必须基于给定时间的力之和来计算每个位置.为此,最好从随时计算净力开始,然后使用该净力计算加速度,速度,然后计算位置.对于以下计算,假定浮力和重力是恒定的(实际上并不正确,但在这种情况下其可变性的影响可忽略不计),还假定初始位置为(0,0),尽管可以轻松地更改为任何初始位置.

You must calculate each location based on the sum of forces at the given time. For this it is better to start from calculating the net force at any time and using this to calculate the acceleration, velocity and then position. For the following calculations, it is assumed that buoyancy and gravity are constant (which is not true in reality but the effect of their variability is negligible in this case), it is also assumed that the initial position is (0,0) though this can be trivially changed to any initial position.

F_x = tau_x 
F_y = tau_y + bouyancy + gravity

其中tau_xtau_y分别是在xy方向上的拖曳力.速度v_xv_y然后由

Where tau_x and tau_y are the drag forces in the x and y directions, respectively. The velocities, v_x and v_y, are then given by

v_x = v_x + (F_x / (2 * m)) * dt
v_y = v_y + (F_y / (2 * m)) * dt

因此,在任何时候txy位置r_xr_y

So the x and y positions, r_x and r_y, at any time t are given by the summation of

r_x = r_x + v_x * dt
r_y = r_y + v_y * dt

在两种情况下,对于某些dt,都必须从0t进行评估,如果n是求和的步数,则必须对dt * n = t进行评估.

In both cases this must be evaluated from 0 to t for some dt where dt * n = t if n is the number of steps in summation.

r_x = r_x + V_i * np.cos(theta) * dt + (F_x / (2 * m)) * dt**2
r_y = r_y + V_i * np.sin(theta) * dt + (F_y / (2 * m)) * dt**2

整个计算实际上可以分两行完成

The entire calculation can actually be done in two lines,

r_x = r_x + V_i * np.cos(theta) * dt + (tau_x / (2 * m)) * dt**2
r_y = r_y + V_i * np.sin(theta) * dt + ((tau_y + bouyancy + gravity) / (2 * m)) * dt**2

除了v_xv_y需要在每个时间步进行更新.要循环遍历并计算一段时间内的xy位置,您可以简单地遵循以下示例(已编辑).

Except that v_x and v_y require updating at every time step. To loop over this and calculate the x and y positions across a range of times you can simply follow the below (edited) example.

以下代码包含一些纠正措施,可防止出现y负位置,因为g的给定值适用于表面或火星,所以我认为这是适当的-当您击中零y并尝试继续前进时,您可能会最终就像我们物理学家所说的那样,发生了计划外的快速拆卸.

The following code includes corrections to prevent negative y positions, since the given value of g is for the surface or Mars I assume this is appropriate - when you hit zero y and try to keep going you may end up with a rapid unscheduled disassembly, as we physicists call it.

为响应已编辑的问题,对以下示例进行了修改,以绘制请求的所有三种情况-重力,重力加阻力以及重力加阻力和浮力.情节设置代码也已添加

In response to the edited question, the following example has been modified to plot all three cases requested - gravity, gravity plus drag, and gravity plus drag and buoyancy. Plot setup code has also been added

import numpy as np
import matplotlib.pyplot as plt

def projectile(V_initial, theta, bouyancy=True, drag=True):
    g = 9.81
    m = 1
    C = 0.47
    r = 0.5
    S = np.pi*pow(r, 2)
    ro_mars = 0.0175

    time = np.linspace(0, 100, 10000)
    tof = 0.0
    dt = time[1] - time[0]
    bouy = ro_mars*g*(4/3*np.pi*pow(r, 3))
    gravity = -g * m
    V_ix = V_initial * np.cos(theta)
    V_iy = V_initial * np.sin(theta)
    v_x = V_ix
    v_y = V_iy
    r_x = 0.0
    r_y = 0.0
    r_xs = list()
    r_ys = list()
    r_xs.append(r_x)
    r_ys.append(r_y)
    # This gets a bit 'hand-wavy' but as dt -> 0 it approaches the analytical solution.
    # Just make sure you use sufficiently small dt (dt is change in time between steps)
    for t in time:
        F_x = 0.0
        F_y = 0.0
        if (bouyancy == True):
            F_y = F_y + bouy
        if (drag == True):
            F_y = F_y - 0.5*C*S*ro_mars*pow(v_y, 2)
            F_x = F_x - 0.5*C*S*ro_mars*pow(v_x, 2) * np.sign(v_y)
        F_y = F_y + gravity

        r_x = r_x + v_x * dt + (F_x / (2 * m)) * dt**2
        r_y = r_y + v_y * dt + (F_y / (2 * m)) * dt**2
        v_x = v_x + (F_x / m) * dt
        v_y = v_y + (F_y / m) * dt
        if (r_y >= 0.0):
            r_xs.append(r_x)
            r_ys.append(r_y)
        else:
            tof = t
            r_xs.append(r_x)
            r_ys.append(r_y)
            break

    return r_xs, r_ys, tof

v = 30
theta = np.pi/4

fig = plt.figure(figsize=(8,4), dpi=300)
r_xs, r_ys, tof = projectile(v, theta, True, True)
plt.plot(r_xs, r_ys, 'g:', label="Gravity, Buoyancy, and Drag")
r_xs, r_ys, tof = projectile(v, theta, False, True)
plt.plot(r_xs, r_ys, 'b:', label="Gravity and Drag")
r_xs, r_ys, tof = projectile(v, theta, False, False)
plt.plot(r_xs, r_ys, 'k:', label="Gravity")
plt.title("Trajectory", fontsize=14)
plt.xlabel("Displacement in x-direction (m)")
plt.ylabel("Displacement in y-direction (m)")
plt.ylim(bottom=0.0)
plt.legend()
plt.show()

请注意,这将保留并返回变量tof中的飞行时间.

Note that this preserves and returns the time-of-flight in the variable tof.

这篇关于如何在重力,浮力和空气阻力的作用下绘制弹丸的运动?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆