【发布时间】：2020-04-20 23:23:57
【问题描述】：
等式是：
d^2 r/dt^2 = -c/m (dr/dt)+g
其中r是弹丸的位置，c是阻力系数，m是弹丸的质量，g是重力加速度。
假设组件形式只有二维，当然，读取，
d^2 X/dt^2 = -c(dX/dt)= U
d^2 Y/dt^2 = -c/m(dY/dt)+g
如果我们采用上述方法并将 X 和 Y 中的速度明确定义为，
U = dX/dt
和
V = dX/dt
那么整个耦合方程组是，
dU/dt= -c/m(U)
dV/dt= - c/m(V)+g
dX/dt= U
dY/dt = V
此 ODE 系统的参数为 c = 0.5 kgs^-1、m = 2kg 和 g = -9.81 ms^-2。
将变量初始化为(U0, V0, X0, Y0) = (173, 100, 0, 0)，它从原点以与水平方向成 ∼ 30 度角发射弹丸。
我如何使用 rk4（我想知道如何编写代码）在 python 中编写一个新函数，该函数实现了上述四个 ODE 系统，解决了 2D 抛射体运动问题......？请帮助我对 ODE 和编码非常陌生。谢谢
到目前为止，我已经得到了以下内容......它不起作用，我真的不知道如何解决这个特定问题，我也打算获得一个弹丸图......有人可以改进我的代码请感谢
import numpy as np
import matplotlib.pyplot as plt
def projectileMotion_V(t, M, g, c):
return -c/M * V0 + g
def projectileMotion_U(t, c, M):
return -c/M * U0
V0 = 100
U0 = 173
ang = 30.0
c = 0.5
dt = 0.1
M = 2.0
g = -9.81
h = 0.1
t = [0]
x = [0]
y = [0]
vx = [V0*np.cos((ang*np.pi)/180)]
vy = [U0*np.sin((ang*np.pi)/180)]
ax = [-(c*V0*np.cos((ang*np.pi)/180))/M]
ay = [g-(c*U0*np.sin((ang*np.pi)/180))/M]
def solveODEsWithR4Method(t, x, y, vx, vy, ax, ay):
t.append(t[0]+dt)
vx.append(vx[0]+dt*ax[0])
vy.append(vy[0]+dt*ay[0])
x.append(x[0]+dt*vx[0])
y.append(y[0]+dt*vy[0])
vel = np.sqrt(vx[0+1]**2 + vy[0+1]**2)
drag = c*vel
ax.append(-(drag*np.cos(ang/180*np.pi))/M)
ay.append(-g-(drag*np.sin(ang/180*np.pi)/M))
return -c/M * V0 + g
fig,ax = plt.subplots()
ax.plot(t, M, g, c)
plt.show()
【问题讨论】：
-
您好，欢迎您-您的“网络搜索”出现了什么？例如，我看到：codeproject.com/tips/792927/… 和 youtube.com/watch?v=IOkwWYaZbck 这里有两件事，编码和求解 ODE。您可能想从编写一些简单的代码开始，例如二次方程的解，然后是一个简单的 1D ODE，以获得您的“python 腿”，然后解决这个问题。
-
空气阻力是速度的二次方，您使用的是水阻力（但系数错误），请参阅stackoverflow.com/a/35009733/3088138 了解计算结果。
-
这是一个具有常系数的非齐次线性微分方程组。可以完全求解，得到精确解，不需要数值积分。
标签： python math physics numerical-methods differential-equations