[Elementary Mechanics Using Python-02]Feather in tornado

Problem

9.17 Feather in tornado.

In this project you will learn to use Newton’s laws and the force model for air resistance in a wind field to address the motion of a light object in strong winds. We start from a simple model without wind and gradually add complexity to the model, until we finally address the motion in a tornado. Motion without wind:

First, we address the motion of the feather without wind.

(a) Identify the forces acting on a feather while it is falling and draw a free-body diagram for the feather.

(b) Introduce quantitative force models for the forces, and find an expression for the acceleration of the feather. You may assume a quadratic law for air resistance.

(c) If you release the feather from rest, its velocity will tend asymptotically toward the terminal velocity, \(v_T\) . Show that the terminal velocity is \(v_T = −(mg/D)^{1/2}\), where D is the constant in the air resistance model.

(d) We release the feather from a distance h above the floor and measure the time t until the feather hits the floor. You may assume that the feather falls with a constant velocity equal to the terminal velocity. Show how you can determine D/mg by measuring the time t. Estimate \(D/mg\) when you release the feather from a height of 2.4 m above the floor and it takes 4.8 s until it hits the floor.

(e) We will now develop a more precise model where we do not assume that the velocity is constant. You release the feather from the height h at the time t = 0 s. Find the equation you have to solve to find the position of the feather as a function of time. What are the initial conditions?

(f) Write a program that solves this equation to find the velocity and position as a function of time t. Use the parameters you determined above, and test the program by ensuring that it produces the correct terminal velocity.

(g) Fig. 9.19 shows the position and velocity calculated with the program using the parameters found above. Was the approximation in part (d) reasonable? Explain your answer. Model with wind: We have now found a model that can be used to find the motion of the feather. We will now find the motion of the feather in three dimensions while it is blowing. The velocity of the wind varies in space, so that the wind velocity w is a function of the position r. We write this as w = w(r).

(h) Find an expression for the acceleration of the feather. The expression may include the wind velocity w(r). Let the z-axis correspond to the vertical direction.

(i) Assume that the feather is moving in an approximately horizontal plane—that is you may assume that the vertical acceleration is negligible. How does the wind have to blow in order for the feather to move in a circular orbit of radiusr0 with a constant speed v0?

Motion in a tornado: For a tornado with a center at the origin, the wind velocity is expected to be approximately given by the model:

\[w(r) = u_0re^{−r/R}\widehat{u}_θ=u_0(−y,x,0)e^{−r/R} \]

where u0 is a characteristic velocity for the wind, R is the radius of the tornado, and \(\widehat{u}_θ\) is a tangential unit vector in the horizontal plane (\(\widehat{u}_θ\) is normal to r). Here, r = (x, y,z), and r = |r|. The velocity field is illustrated in Fig. 9.20.

(j) Is is possible to choose an appropriate set of initial conditions so that the feather moves in a circular path in the tornado? Explain your answer.

(k) Rewrite your program to find the velocity and position of the feather as a function of time. For the tornado you may use the values u0 = 100 m/s and R = 20 m.

(l) Find the trajectory for the feather if it is released from rest from a height of 2.4m, and in a position corresponding to .

\[r = −R i + hk. \]

(m) The trajectory of the feather is shown in Figs. 9.20 and 9.21. Compare the results with what happended when you dropped the feather without wind. Why does the feather now take longer to reach the ground?

Solution

(\(a\))

易得

(\(b\))

易得

\[ma=-Dv|v|-mg \]

\(\Longrightarrow\)

\[a=-g-Dv|v|/m \]

(\(c\))

当\(a=0\)时,可得

\[v_T|v_T|=-mg/D \]

\(\Longrightarrow\)

\[v_T=-(mg/D)^{1/2} \]

(\(d\))

易得

\[\begin{cases} x=-v_Tt\\ v_T=-(mg/D)^{1/2} \end{cases} \]

解得

\[D/mg=4 \]

(\(e\))

\[y(0)=h \]

\[v(0)=0 \]

(\(f\))

from pylab import *

g = 9.8
m = 0.001
D=m*g*4
t = 5
figure(figsize=[7,7])
dt = 0.0001
n = int(t / dt)
t1 = arange(0, 5, dt)
v = zeros(n)
x = zeros(n)
x[0]=2.4
for i in range(n - 1):
    f = v[i] * v[i] * D
    a = (m * g - f) / m
    v[i + 1] = v[i] - a * dt
    x[i + 1] = x[i] + v[i] * dt
subplot(211)
ylim(0,3)
xlim(0,5)
xlabel("t [s]")
ylabel("x [m]")
plot(t1, x)
ax2 = subplot(212)
xlabel("t [s]")
ylabel("v [m/s]")
plot(t1, v)
show()

(\(h\))

\[ma=-D(v-w)|v-w|-mg \]

\(\Longrightarrow\)

\[a=-D(v-w)|v-w|/m-g \]

(\(i\))

由题意易得可以忽略重力加速度的影响,当羽毛以稳定速度匀速圆周运动时,其状态为平衡状态.

故易得

\[w_T=v_T=v_0 \]

(\(j\))

No.

羽毛还受重力影响.

(\(k\))

from pylab import *

# 初始数值
k = array([0, 0, 1])
R = 20
H = 2.4
g = 9.81
m = 0.001
D = m * g * 4
dt = 0.00001
t1 = 25
u0 = 100


# 计算长度
def clen(x):
    return sqrt(x[0] ** 2 + x[1] ** 2 + x[2] ** 2)


# 计算风速(会报错但无伤大雅）
def wind(x, y, z=0):
    rlen = clen([x, y, z])
    return array([u0 * exp(-rlen / R) * (-y), u0 * exp(-rlen / R) * x, 0])


# 创建网格
x, y = meshgrid(
    linspace(-100, 100, 100),
    linspace(-100, 100, 100)
)
# 计算风场
u1, v1, w1 = wind(x, y)
M = hypot(u1, v1)
# 绘制风场
quiver(x, y, u1, v1, M, pivot='mid', width=0.001)
# 初始量
n = int(t1 / dt)
t = zeros(n, float)
r = zeros((n, 3), float)
v = zeros((n, 3), float)
a = zeros((n, 3), float)
t[0] = 0
r[0, :] = array([-1 * R, 0, 2.4])
v[0, :] = array([0, 0, 0])
for i in range(n - 1):
    if r[i, 2] == 0:
        print("羽毛下落时间为%.6fs",t[i])
        break
    else:
        w = wind(r[i, 0], r[i, 1], r[i, 2])
        a[i] = -(D / (m * g) * (v[i] - w) * clen(v[i] - w)) * g
        a[i][2] -= g
        v[i + 1] = v[i] + a[i] * dt
        r[i + 1] = r[i] + v[i + 1] * dt
        t[i + 1] = dt + t[i]
imax = n - 1
i = range(imax)
i1 = range(0, imax, 50000)
plot(r[i, 0], r[i, 1], '-', r[i1, 0], r[i1, 1], 'o')

axis('equal')
xlabel("x [m]")
ylabel("y [m]")
savefig("1.jpg", dpi=3000)
show()

易得羽毛下落时间为20.61322s,之所以下落时间变多是因为空气阻力正相关与速度的平方。

标签：Mechanics,02,mg,Python,tornado,velocity,feather,dt,wind
来源： https://www.cnblogs.com/Foreze/p/14534633.html