抛体运动那些事儿

不难写出阻力的表达式：

\boldsymbol{f}=-mkv\boldsymbol{\hat{v}}\\

其中 \boldsymbol{\hat{v}}\equiv\frac{\boldsymbol{v}}{v} ， \boldsymbol{f} 是空气阻力，m 为物体质量， k 为与阻力系数， v 为速率。

以水平方向和竖直方向建立平面直角坐标系。

可以列出坐标分量形式微分方程：

\begin{cases} \dot{x}(0)=v_0cos\theta,\ddot{x}(t)=-k\dot{x}(t)\\[2ex] \dot{y}(0)=v_0sin\theta,\ddot{y}(t)=-k\dot{y}(t)-g\\[2ex] \end{cases}

可以看得出来关于 y 是一个二阶常系数非齐次线性微分方程，关于 x 是一个二阶常系数齐次线性微分方程，解得参数方程：

x(t)=\frac{v_0cos\theta-v_0cos\theta\cdot e^{-kt}}{k}，y(t)=\frac{(kv_0sin\theta+g)-(kv_0sin\theta+g)\cdot e^{-kt}-gkt}{k^2}

回归正题，这个参数方程我用几何画板画出来的图像长这个样子（如图）：

如果硬要解出轨迹方程 y(x) 也是可以的。我们知道因变量对参量的导数与自变量对参量的导数之比是因变量对自变量的导数。如下：

\frac{\text dy}{\text dx}=\frac{\frac{\text dy}{\text dt}}{\frac{\text dx}{\text dt}}=\frac{\dot y}{\dot x}\\

那么我们对参数方程的时间导数求得

{\begin{cases} \dot x(t)=v_0cos\theta\cdot e^{-kt}\\[2ex] \dot y(t)=v_0sin\theta\cdot e^{-kt}+\frac{g}{k}(e^{-kt}-1) \end{cases}}\

得到比值为

\frac{\text dy}{\text dx}=-\frac{g}{kv_0cos\theta}(e^{kt}-1)+tan\theta\\

不难发现，我们还是没有消去参量 t ，所以我们得找到 e^{-kt} 与 x 的关系。根据参数方程 x(t) 我们可以得到

e^{-kt}=1-\frac{k}{v_0cos\theta}x\\

代入再求积分就可以得到

{\begin{align} y(x)&=\int(-\frac{g}{kv_0cos\theta}((1-\frac{k}{v_0cos\theta}x)^{-1}-1)+tan\theta)\text dx\\[2ex] &=\int(\frac{g}{k^2x-kv_ocos\theta}+\frac{g}{kv_ocos\theta}+tan\theta)\text dx\\[2ex] &=\frac{g}{k^2}ln(1-\frac{k}{v_0cos\theta}x)+(\frac{g}{kv_0cos\theta}+tan\theta )x+C \end{align}}\\

因为 y(0)=0 ,所以代入可以求得常量 C ，最终得到解析式：

\color{red}{y=\frac{g}{k^2}ln(1-\frac{k}{v_0cos\theta}x)+(\frac{g}{kv_0cos\theta}+tan\theta )x}\\

但我们要求得是落地点得水平距离，那么我们平移坐标，讲坐标从抛射点平移到的人得脚底位置。那么在新坐标系下解析式为

y=\frac{g}{k^2}ln(1-\frac{k}{v_0cos\theta}x)+(\frac{g}{kv_0cos\theta}+tan\theta )x+h\\

落地点就是 y=0 的非原点解，得到方程：

\frac{g}{k^2}ln(1-\frac{k}{v_0cos\theta}x)+(\frac{g}{kv_0cos\theta}+tan\theta )x+h=0\\

这是一个超越方程，没有解析解，但可以通过数值分析的方法大致计算解的值。

如果我们令 k\rightarrow0 ,利用洛必达法则看看这个曲线方程是否会退化为抛物线方程。

{\begin{align} \lim_{k\to 0}{y} &=\lim_{k\to 0}\frac{g\ln(1-\frac{x}{v_0\cos\theta}k)+\frac{gx}{v_0\cos\theta}\cdot k+x\tan\theta \cdot k ^2 }{k^2}\\ &=\lim_{k\to 0}\frac{g\frac{\frac{x}{v_0\cos\theta}}{\frac{x}{v_0\cos\theta}\cdot k-1}+\frac{x}{v_0\cos\theta}g+2x\tan\theta\cdot k}{2k}\\ &=\lim_{k\to 0}\frac{-\frac{gx}{v_0\cos\theta}\frac{\frac{x}{v_0\cos\theta}}{(\frac{x}{v_0\cos\theta}\cdot k-1)^2}+2x\tan\theta}{2}\\[3ex] &=\tan\theta\cdot x-\frac{g}{2{v_0^2}\cos^2\theta}\cdot x^2 \end{align}}\\

确实退化为了斜抛运动抛物线方程。

这一部分内容呢，不是原创的解法，这是该解法的原回答

同样按照先前的思路：为了简便，我们设空气阻力大小与速度大小的平方的比例系数除以物体质量后得到的新的系数为 \mu ，其余的参数和之前的设定相同。

即阻力的表达式：

\boldsymbol{f}=-m\mu v^2\boldsymbol{\hat{v}}\\

写出坐标分量形式的微分方程：

\begin{cases} \ddot x=-\mu\sqrt{\dot x^2+\dot y ^2}\dot x,x(0)=0,\dot x(0)=v_0\cos\theta\\ \ddot y=-\mu\sqrt{\dot x^2+\dot y ^2}\dot y-g,y(0)=0,\dot y(0)=v_0\sin\theta\\ \end{cases}\

两式相除得到：

\frac{\ddot{y} }{\ddot{x} }=\frac{\mu \sqrt{\dot{x}^2+\dot{y}^2 }\dot{y} +g}{\mu \sqrt{\dot{x}^2+\dot{y}^2 }\dot{x}}=\frac{\dot{y} }{\dot{x} }+\frac{g}{\mu\sqrt{\dot{x}^2+\dot{y}^2 }\dot{x} }\\

令 \frac{\dot{y} }{\dot{x} }\equiv m= m(t)

那么又有

\frac{\ddot{y} }{\ddot{x}}=\frac{\text{d}\dot{y}/\text{d}t }{\text{d}\dot{x}/\text{d}t }=\frac{\text{d}\dot{y} }{\text{d}\dot{x}}=\frac{\text{d}(m\dot{x}) }{\text{d}\dot{x}}=\frac{\text{d}m }{\text{d}\dot{x} }\dot{x}+m\\

即：

\frac{\text{d}m }{\text{d}\dot{x} }\dot{x}+m =m+\frac{g}{\mu\dot{x}^2\sqrt{1+m^2}}\\

得到 \dot{x} 和 m 的微分关系：

\sqrt{1+m^2} \text{d}m =\frac{g\text{d}\dot{x}}{\mu\dot{x}^3}\\

进而定积分有

\fbox{$\int_{\tan\theta}^{m}\sqrt{1+m^2} \text{d}m =\int_{v_0\cos\theta}^{\dot x}\frac{g\text{d}\dot{x}}{\mu\dot{x}^3}$}\\

积分再化简得到

\ln(m+\sqrt{1+m^2})+m\sqrt{1+m^2}+\frac{g}{\mu \dot{x}^2}=\ln(\tan\theta+|\sec\theta|)+\tan\theta|\sec\theta|+\frac{g\sec^2\theta}{\mu v_0^2}\\ 可见等式右边是一个定值，由最初的参数决定。令其为 h ，即定义：

G\equiv\ln(\tan\theta+|\sec\theta|)+\tan\theta|\sec\theta|+\frac{g\sec^2\theta}{\mu v_0^2}\\ 所以上式又化简为

\ln(m+\sqrt{1+m^2})+m\sqrt{1+m^2}+\frac{g}{\mu \dot{x}^2}=G\\ 再移项得到

\dot x=\sqrt{\frac{g}{\mu(G-\ln(m+\sqrt{1+m^2})-m\sqrt{1+m^2})}}\\

令

f(m)\equiv\dot x=\sqrt{\frac{g}{\mu(G-\ln(m+\sqrt{1+m^2})-m\sqrt{1+m^2})}}\\

因为 \frac{\dot{y} }{\dot{x} }\equiv m ，所以 \dot x=f(m),\dot y=mf(m) 而想要积分得到 x(m)，y(m) ，那么就必须得到 x，y 对 m 的导数。有

\frac{\text dx}{\text dm}=\frac{\text dx/\text dt}{\text dm/\text dt}=\frac{\dot x}{\dot m},\frac{\text dy}{\text dm}=\frac{\text dy/\text dt}{\text dm/\text dt}=\frac{\dot y}{\dot m}\\

其中未知的只有 \dot m ，可将它表示为：

\dot m=\frac{\text d}{\text dt}(\frac{\dot{y} }{\dot{x} })=\frac{\ddot y\dot x-\dot y\ddot x}{\dot x^2}\\

观察原微分方程可以看出

\frac{\ddot x}{\dot x}=\frac{\ddot y+g}{\dot y}=-\mu\sqrt{\dot x^2+\dot y^2}\\

得到交叉相乘有

\ddot y\dot x-\dot y\ddot x=-g\dot x\\

所以带回有\dot m=-\frac{g}{\dot x}\\ 最后得到

\frac{\text d x}{\text dm}=-\frac{f^2(m)}{g},\frac{\text dy}{\text dm}=-\frac{mf^2(m)}{g}\\

所以分别定积分有

x(m)=-\frac1g \int_{\tan\theta}^{m}f^2(m)\text dm，y(m)=-\frac1g\int_{\tan\theta}^{m}mf^2(m)\text dm\\

即完整的参数方程形式：

\fbox{$\color{darkblue}{x(i)=\int_{\tan\theta}^{i}\frac{\text dz}{\mu(\operatorname{arsh}z+z\sqrt{1+z^2}-G)}\\ y(i)=\int_{\tan\theta}^{i}\frac{z\text dz}{\mu(\operatorname{arsh}z+z\sqrt{1+z^2}-G)}\\}$}\\

\fbox{$G\equiv\ln(\tan\theta+|\sec\theta|)+\tan\theta|\sec\theta|+\frac{g\sec^2\theta}{\mu v_0^2}$}\\

用Desmos画出图像：

这个积分是不能用初等函数表示的，还要注意 i 是参数不是时间。

至此，轨迹方程的参数方程以及得到了，但是为了一写特殊需求，我们不妨尝试得到关于时间的参数方程。由于关于 i 的参数方程已经得到了，得到关于时间的参数方程只需要得到 i(t) 即可。

之前得到\dot i=\dot m=-\frac{g}{\dot x} 即

f(i)\text di=-g\text dt\\

积分得 t 和 m 的关系：

t=-\frac{1}{g}\int_{\tan\theta}^{i}f(u)\text du\\

综上所述，我们已经得到 x(t) 和 y(t) 啦！很显然这两个函数都不是初等函数，只能得到数值解。用计算机带入积分式算还是不成问题的。