最速降线问题探讨

最速降线及相关问题讨论

最速降线求解

• 最速降线问题的提出：

​ 考虑如下运动：

​ 质点从一个定点运动到不在其垂直下方的任意一个点（即排除自由落体运动），整个运动过程只受重力作用，不计摩擦 力，那么问质点应该沿什么样的曲线下滑，才能使得运动时间最少？!

• 最速降线问题推导:

​ 下面我们开始推导最速绛线方程，如下图运动过程，假设物体从原点O点开始出发，沿任意曲线运动到A点：

​

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-XXec7ue9-1679124156195)(E:\MarkDown\images\屏幕截图 2023-02-26 121745.jpg)]

​ 考虑整个过程，只有重力做功：

​ 整个运动过程只有重力做功，根据动能守护可以得到最终末速度:
$$\frac{1}{2}m{v}^2=mgy\Rightarrow v=\sqrt{2gy}$$
​ 取物体运动过程中的任意的弧长微分：
$$ds=\sqrt{1+y^{{}^{\prime }2}}dx$$
​ 则物体运动时间微分：
$$dt=\frac{ds}{v}=\sqrt{\frac{1+y^{{}^{\prime }2}}{2gy}}dx$$

​ 设O点坐标为(x₁,y₁)，A点坐标为(x₂,y₂)，则物体运动的总时间为：
$$T={\int }_{x1}^{x2}dt={\int }_{x1}^{x2}\sqrt{\frac{1+y^{{}^{\prime }2}}{2gy}}dx\text{\hspace{1em}(*)}$$

至此我们得到了物体沿任意曲线下滑的时间积分函数表达式，下面就是求解这个积分函数。但是上述积分里包含y和y’两个未知数，而y和y’都是关于x的函数，也就是说函数的自变量包含了另外两个函数，为此需要用到泛函的知识来求解

泛函的基本概念

• 泛函的定义：

​ 区别于我们之前学过的函数，泛函是广义的函数，是定义了从函数域到数域的一个映射，即可表示为：
$$I=F(x,f(x),g(x),f^{{}^{\prime }}(x),g^{{}^{\prime }}(x))$$
将定义域D扩展到一个函数集合，通过映射法则F，得到一组值的数域，即完成函数空间到数域的映射关系，这种映射称为泛函

• (*)式泛函表述：

由上述定义，我们可得出时间表达式为一泛函，即：
$$\begin{gathered} T={\int }_{x1}^{x2}dt={\int }_{x1}^{x2}\sqrt{\frac{1+y^{{}^{\prime }2}}{2gy}}dx={\int }_{x1}^{x2}F(x,y,y^{{}^{\prime }}) \\ F(x,y,y^{{}^{\prime }})=\sqrt{\frac{1+y^{{}^{\prime }2}}{2gy}} \end{gathered}$$
而解决泛函的问题，我们就要用到变分法。

变分法求解最速降线问题

• 变分的思想：

  如图所示，红色曲线表示最速下降线，而其余则表示物体沿任意曲线下滑，我们的问题是：如何从这一簇曲线中找到下降时间最短的一条呢？

  [外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-XG5aXSvn-1679124156196)(E:\MarkDown\images\屏幕截图 2023-02-26 151138.jpg)]

我们用$\bar{F}(x)$表示实验曲线簇—实验泛函，用$F(x)$表示最佳解泛函（最速降线），定义D(x)为不包含最佳解的曲线簇，即有：
$$D(x)=\bar{F}(x)-F(x)\text{\hspace{1em}(1)}$$
​ 我们引入一个任意实数的控制因子$\epsilon$来控制$D(x)$的变化：
$$D(x)=ϵ\frac{D(x)}{ϵ}=ϵ\eta (x)$$

​ 则（1）式可改写为：
$$\bar{F}(x)=F(x)+ϵ\eta (x)$$

用图阐述上式关系为：

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-FXr67jmt-1679124156196)(E:\MarkDown\images\屏幕截图 2023-02-26 155106.jpg)]

（1）dx和dy的意义是：当自变量x发生左右伸缩变化时，最佳解泛函的因变量y也发生上下伸缩变化

（2）δy的意义是：自变量x并没有发生任何变化，而最佳解泛函和实验泛函在因变量y上发生上下伸缩变化

（3）δy'同δy的意义一致

dx和dy叫做泛函的微分，而δy和δy’叫做泛函的变分，其表示要求解的最佳解泛函到任意实验泛函的差异。用数学语言翻译为：
$$F(x,y+ϵ\eta (x),y^{{}^{\prime }}+ϵ{\eta }^{{}^{\prime }}(x))=F(x,y,y^{{}^{\prime }})+ϵ\eta (x)$$
尽管扰动函数$\eta (x)$是任意的，但我们只需使控制因子趋于0，即可让所有的实验泛函完全收敛到最佳解泛函上，这就是变分法的核心，用数学语言表述为：
$$\underset{ϵ\to 0}{\lim }F(x,y+ϵ\eta (x),y^{{}^{\prime }}+ϵ{\eta }^{{}^{\prime }}(x))=F(x,y,y^{{}^{\prime }})+ϵ\eta (x)$$
至此，我们得出：

我们要寻找的最佳解泛函，就是实验泛函出现极值的情况，而实验泛函出现极值的充要条件就是控制因子趋近于0

​

• 欧—拉方程推导：

  ​ 方式（1）：

​ 利用变分法的思想，我们来求解积分函数$I$
$$I(ϵ)={\int }_{x1}^{x2}F(x,y+ϵ\eta (x),y^{{}^{\prime }}+ϵ{\eta }^{{}^{\prime }}(x))dx$$
​ 利用换元法，我们令：
$$u=y+ϵ\eta (x);v=y^{{}^{\prime }}+ϵ{\eta }^{{}^{\prime }}(x)$$
​ 则积分函数为：
$$I(ϵ)={\int }_{x1}^{x2}F(x,u,v)dx$$
​ 由变分法的极值条件：
$$\frac{d}{dϵ}{\mid }_{ϵ\to 0}{\int }_{x1}^{x2}F(x,u,v)dx=0$$
​ 即为：
$$\frac{d}{dϵ}{\mid }_{ϵ\to 0}{\int }_{x1}^{x2}F(x,u,v)dx={\int }_{x1}^{x2}\frac{\partial F}{\partial ϵ}dx$$
​ 计算偏微分：
$$\begin{gathered} \frac{\partial F}{\partial ϵ}=\frac{\partial F}{\partial x}\frac{\partial x}{\partial ϵ}+\frac{\partial F}{\partial u}\frac{\partial u}{\partial ϵ}+\frac{\partial F}{\partial v}\frac{\partial v}{\partial ϵ} \\ =\frac{\partial F}{\partial u}\eta (x)+\frac{\partial F}{\partial v}\eta {}^{\prime }(x) \end{gathered}$$
​ 根据变分法的极值条件，控制因子趋于0:
$$\begin{gathered} \underset{ϵ\to 0}{\lim }u=\underset{ϵ\to 0}{\lim }y+ϵ\eta (x)=y \\ \underset{ϵ\to 0}{\lim }v=\underset{ϵ\to 0}{\lim }{y}^{{}^{\prime }}+ϵ{\eta }^{{}^{\prime }}(x)=y^{{}^{\prime }} \end{gathered}$$
​ 代入积分函数$I$:
$$\begin{gathered} \frac{d}{dϵ}{\mid }_{ϵ\to 0}{\int }_{x1}^{x2}F(x,u,v)dx=0 \\ \Rightarrow {\int }_{x1}^{x2}(\frac{\partial F}{\partial u}\eta (x)+\frac{\partial F}{\partial v}\eta {}^{\prime }(x))=0 \\ \Rightarrow {\int }_{x1}^{x2}(\frac{\partial F}{\partial y}\eta (x)+\frac{\partial F}{\partial y^{{}^{\prime }}}\eta {}^{\prime }(x))=0 \end{gathered}$$
​ 又有：
$$\eta \prime (x)=\frac{d\eta }{dx}$$
​ 代入：
$${\int }_{x1}^{x2}(\frac{\partial F}{\partial y}\eta (x)+\frac{\partial F}{\partial y^{{}^{\prime }}}\frac{d\eta }{dx})=0$$
​ 下面我们研究上式第二部分，即：
$${\int }_{x1}^{x2}\frac{\partial F}{\partial y^{{}^{\prime }}}\frac{d\eta }{dx}dx$$
​ 利用分布积分法，令$u=\frac{\partial F}{\partial y\prime }$，$v=\eta$，$dv=\frac{d\eta }{dx}dx$:
$$\begin{gathered} {\int }_{x1}^{x2}\frac{\partial F}{\partial y^{{}^{\prime }}}\frac{d\eta }{dx}dx={\int }_{x1}^{x2}udv \\ =\frac{\partial F}{\partial y\prime }[\eta (x_2)-\eta (x_1)]-{\int }_{x1}^{x2}\eta \frac{d}{dx}(\frac{\partial F}{\partial y\prime })dx \end{gathered}$$

​ 易知： $\eta (x_1)$=$\eta (x_2)$

​ 因此：
$$\begin{gathered} {\int }_{x1}^{x2}\frac{\partial F}{\partial y^{{}^{\prime }}}\frac{d\eta }{dx}dx={\int }_{x1}^{x2}udv \\ =\frac{\partial F}{\partial y\prime }[\eta (x_2)-\eta (x_1)]-{\int }_{x1}^{x2}\eta \frac{d}{dx}(\frac{\partial F}{\partial y\prime })dx \\ =-{\int }_{x1}^{x2}\eta \frac{d}{dx}(\frac{\partial F}{\partial y\prime })dx \end{gathered}$$
​ 代回积分$I$:
$$\begin{gathered} {\int }_{x1}^{x2}(\frac{\partial F}{\partial y}\eta (x)+\frac{\partial F}{\partial y^{{}^{\prime }}}\frac{d\eta }{dx})=0 \\ \Rightarrow {\int }_{x1}^{x2}(\frac{\partial F}{\partial y}-\frac{d}{dx}(\frac{\partial F}{\partial y\prime }))\eta dx=0 \end{gathered}$$
​ 要使其恒成立：
$$\frac{\partial F}{\partial y}-\frac{d}{dx}(\frac{\partial F}{\partial y\prime })=0$$
​ 以上就是欧—拉方程的第一种形式。

​ 方式（2）：

​ 欧—拉方程的第二种形式在于，在某些条件下，我们可以简化问题，其推导过程如下：

​ 计算泛函对$F$对$x$的全微分：
$$\begin{gathered} \frac{d}{dx}F(x,y,y\prime )=\frac{\partial F}{\partial x}\frac{dx}{dx}+\frac{\partial F}{\partial y}\frac{dy}{dx}+\frac{\partial F}{\partial y\prime }\frac{dy\prime }{dx} \\ =\frac{\partial F}{\partial x}+\frac{\partial F}{\partial y}y\prime +\frac{\partial F}{\partial y\prime }y\prime \prime \end{gathered}$$
​ 得到：
$$\frac{\partial F}{\partial y\prime }y\prime \prime =\frac{dF}{dx}-\frac{\partial F}{\partial x}-\frac{\partial F}{\partial y}y\prime \text{\hspace{1em}(*)}$$
​ 为化简上式，我们求微分算子$\frac{d}{dx}(y\prime \frac{\partial F}{\partial y\prime })$:
$$\frac{d}{dx}(y\prime \frac{\partial F}{\partial y\prime })=y\prime \prime \frac{\partial F}{\partial y\prime }+y\prime \frac{d}{dx}(\frac{\partial F}{\partial y\prime })\text{\hspace{1em}(**)}$$
​

​ 将(*)式与(**)式联立，得：
$$\begin{gathered} \frac{d}{dx}(y\prime \frac{\partial F}{\partial y\prime })=\frac{dF}{dx}-\frac{\partial F}{\partial x}-\frac{\partial F}{\partial y}y\prime +y\prime \frac{d}{dx}(\frac{\partial F}{\partial y\prime }) \\ =\frac{dF}{dx}-\frac{\partial F}{\partial x}+y\prime (\frac{d}{dx}(\frac{\partial F}{\partial y\prime })-\frac{\partial F}{\partial y}) \end{gathered}$$
​

​ 上式中出现了欧—拉方程的第一种形式，继续化简得：
$$\frac{d}{dx}(y\prime \frac{\partial F}{\partial y\prime })=\frac{dF}{dx}-\frac{\partial F}{\partial x}$$
​ 整理得：
$$\frac{\partial F}{\partial x}-\frac{d}{dx}(F-y\prime \frac{\partial F}{\partial y\prime })=0$$
​ 这就是欧—拉方程的第二种形式，特别的，当$\frac{\partial F}{\partial x}=0$时，有：
$$\frac{d}{dx}(F-y\prime \frac{\partial F}{\partial y\prime })=0$$
​ 也即：
$$F-y\prime \frac{\partial F}{\partial y\prime }=C$$

问题的求解

​ 由上述推导，我们得到了欧—拉方程以及时间的泛函，现在我们利用欧—拉方程求解。
$$\begin{gathered} F-y\prime \frac{\partial F}{\partial y\prime }=C \\ \Rightarrow \frac{\sqrt{1+y{\prime }^2}}{\sqrt{2gy}}-\frac{y{\prime }^2}{\sqrt{2gy}\sqrt{1+y{\prime }^2}}=C \\ \Rightarrow 2gy(1+y{\prime }^2)=\frac{1}{C^2} \\ \Rightarrow y(1+y{\prime }^2)=C(常数) \end{gathered}$$
​ 下面进行求解：

​ 如图：

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-MUHWDpRq-1679124156196)(E:\MarkDown\images\20200215110425743.png)]

​ 易知，$y\prime$=$tan\beta$=$cot\alpha$,因为质点一直向下滑(如果向上，肯定不是最优曲线),所以$0\le \alpha \le \frac{\pi }{2}$;

​ 令$2R=C$,令$y\prime =\frac{dy}{dx}=cot\frac{\theta }{2}$,则$\theta =2\alpha$,所以$0\le \theta \le \pi$.

​ 设$x$为$\theta$的函数，即$x=x(\theta )$,我们的目标是求$x(\theta )$.

​ 由：
$$1+co{t}^2\frac{\theta }{2}=\frac{2}{1-cos\theta }$$
​ 得：
$$\begin{gathered} y(1+y{\prime }^2)=y(1+co{t}^2\frac{\theta }{2})=\frac{2y}{1-cos\theta }=2R \\ \Rightarrow y(\theta )=R(1-cos\theta ) \end{gathered}$$
​ 又有：
$$\frac{dy}{dx}=cot\frac{\theta }{2}=\frac{cos\frac{\theta }{2}}{sin\frac{\theta }{2}},\frac{dy}{d\theta }=Rsin\theta =2Rsin\frac{\theta }{2}cos\frac{\theta }{2}$$
​ 得：
$$x\prime (\theta )=\frac{dx}{d\theta }=\frac{\frac{dy}{d\theta }}{\frac{dx}{d\theta }}=2Rsi{n}^2\frac{\theta }{2}=R(1-cos\theta )$$
​ 所以：
$$x(\theta )=\int x\prime (\theta )d\theta =\int R(1-cos\theta )d\theta =R(\theta -sin\theta )+C_1$$
​ 得到参数方程：
$$\{\begin{array}{l}x(\theta )=R(\theta -sin\theta )+C_1\\ y(\theta )=R(1-cos\theta )\end{array}$$
​ 解得$\theta$为：
$$\theta =arccos(1-\frac{y}{R}),0\le \theta \le \pi$$
​ 则：
$$x=Rarccos(1-\frac{y}{R})-\sqrt{y(2R-y)}+C_1$$
​ 边界条件有：
$$x(0)=0,x(H)=c$$
​ 代入方程得：
$$\begin{gathered} x(0)=Rarccos(1)+C_1=C_1=0; \\ x(H)=Rarccos(1-\frac{H}{R})-\sqrt{H(2R-H)}+C_1=c \end{gathered}$$
​ 即$R$的值由下落点坐标$H、c$共同决定。

​ 综上，我们得到了最速降线的轨迹方程：
$$\{\begin{array}{l}x(\theta )=R(\theta -sin\theta )\\ y(\theta )=R(1-cos\theta )\end{array}$$

最速降线的弧长和运动时间

​ 在得到最速降线轨迹方程后，我们可得其弧长：
$$\begin{gathered} x\prime (\theta )=\frac{dx}{d\theta }=\frac{\frac{dy}{d\theta }}{\frac{dx}{d\theta }}=2Rsi{n}^2\frac{\theta }{2}=R(1-cos\theta ) \\ y\prime (\theta )=Rsin\theta \end{gathered}$$
​ 则弧长微元为：
$$\begin{gathered} ds=\sqrt{(\frac{dx}{d\theta }{)}^2+(\frac{dy}{d\theta }{)}^2} \\ =R\sqrt{2-2cos\theta }d\theta \\ =2Rsin\frac{\theta }{2}d\theta \end{gathered}$$
​ 弧长$s$为：
$$s={\int }_{\theta =0}^{\theta =\pi }2Rsin\theta d\theta =4R$$
​ 时间$T$:
$$\begin{gathered} T={\int }_{\theta =0}^{\theta =\pi }\frac{ds}{v}d\theta \\ =\sqrt{\frac{R}{g}}{\int }_{\theta =0}^{\theta =\pi }d\theta =\pi \sqrt{\frac{R}{g}} \end{gathered}$$
​ 可见，最速降线的运动时间T是常量，与物体的初始位置无关。

基于MATLAB仿真对最速降线上点下落时间的研究

​ 我们给出几种不同轨道的曲线：

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-fBZj4fY0-1679124156197)(E:\MarkDown\images\不同轨道曲线.png)]

​ 在给定横坐标位移 $x$ 下，我们可得出不同轨道下落所需的时间曲线：

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-y8CkNIzV-1679124156197)(E:\MarkDown\images\时间.png)]

​ 可见，同等条件下，旋轮线的下降时间最短。速度方面：

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-JMJC4ZV9-1679124156197)(E:\MarkDown\images\速度.png)]

​ 可以看出，旋轮线上速度最快。

​ MATLAB仿真代码如下：

clc;
clear;
close all;
%% 参数设置
dx = 1*10^-6;
x = 0:dx:0.3;
d_theta = 1*10^-5;
theta = 0:d_theta:3.3;
gravity = 10;
%% 轨道生成
circle_y = sqrt(0.04-((x-0.2).^2));%以（0，0.2）为圆心，0.2为半径的圆弧
line_y = (sqrt(0.03)/3).*x; %直线
parabola_y = -1.9245.*(x-0.3).^2+0.17320528; %用时最短的抛物线
cycloid_x = 0.08715.*(theta-sin(theta));
cycloid_y = 0.08715.*(1-cos(theta));
%% 绘制轨迹图
figure();
plot(x,-circle_y,'-.',x,-line_y,':',x,-parabola_y,'--',cycloid_x,-cycloid_y,'LineWidth',3);
grid on;
title('不同的轨道曲线','FontSize',16);
xlabel('\itx\rm/m','FontSize',16);
ylabel('\ity\rm/m','FontSize',16);
axis equal;
legend('圆弧','直线','抛物线','旋轮线');
%% 计算弧微分
dy_circle = diff(circle_y);
dy_line = diff(line_y);
dy_parabola = diff(parabola_y);
dy_cycloid = diff(cycloid_y);
dx_cycloid = diff(cycloid_x);
ds_cirle = sqrt(dx.^2+dy_circle.^2);
ds_line = sqrt(dx.^2+dy_line.^2);
ds_parabola = sqrt(dx.^2+dy_parabola.^2);
ds_cycloid = sqrt(dx_cycloid.^2+dy_cycloid.^2);
%% 计算弧长
s_circle = sum(ds_cirle);
s_line = sum(ds_line);
s_parabola = sum(ds_parabola);
s_cycloid = sum(ds_cycloid);
%% 计算各点速度
% v=sqrt(2gh)
v_circle = sqrt(2*gravity.*circle_y);
v_line = sqrt(2*gravity.*line_y);
v_parabola = sqrt(2*gravity.*parabola_y);
v_cycloid = sqrt(2*gravity.*cycloid_y);
%% 绘制速度-位置曲线
figure();
plot(x,v_circle,'-.',x,v_line,':',x,v_parabola,'--',cycloid_x,v_cycloid,'LineWidth',3);
grid on;
title('速度与x的曲线','FontSize',16);
xlabel('\itx\rm/m','FontSize',16);
ylabel('\itv\rm/m\cdots^-^1','FontSize',16);
legend('圆弧','直线','抛物线','旋轮线');
%% 计算时间
%每段微分对应时间dt = ds/v;
dt_circle = ds_cirle./v_circle(2:length(v_circle));
dt_line = ds_line./v_line(2:length(v_line));
dt_parabola = ds_parabola./v_parabola(2:length(v_parabola));
dt_cycloid = ds_cycloid./v_cycloid(2:length(v_cycloid));
%对应各点耗时间
t_circle = cumsum(dt_circle);
t_line = cumsum(dt_line);
t_parabola = cumsum(dt_parabola);
t_cycloid = cumsum(dt_cycloid);
%总时间
t_circle(length(t_circle));
t_line(length(t_line));
t_parabola(length(t_parabola));
t_cycloid(length(t_cycloid));
%% 绘制位置与时间图
figure();
plot([0 t_circle],x,'-.',[0 t_line],x,':',[0 t_parabola],x,'--',[0 t_cycloid],cycloid_x,'LineWidth',3);
grid on;
title('x与时间t的关系','FontSize',16);
xlabel('\itx\rm/m','FontSize',16);
ylabel('\itt\rm/s','FontSize',16);
legend('圆弧','直线','抛物线','旋轮线');

总结

​ 本文利用变分法以及MATLAB仿真对最速降线问题进行了探讨，可以看初在给定位置下，旋轮线为下落时间最短的曲线，而解决这一问题的变分法也是后续解决很多问题的强有力工具。与我们所学习的拉格朗日方程不同的是，变分法以及其对应的欧-拉方程引入了变分这一概念，利用泛函的思想求解，而拉格朗日方程则基于我们所学习的虚位移原理以及达朗贝尔原理来得出动力学普遍方程，再对其推导得出两类不同形式的拉格朗日方程，进而求解动力学相关问题。后者也是当下众多仿真软件的理论依据。
t_cycloid],cycloid_x,‘LineWidth’,3);
grid on;
title(‘x与时间t的关系’,‘FontSize’,16);
xlabel(‘\itx\rm/m’,‘FontSize’,16);
ylabel(‘\itt\rm/s’,‘FontSize’,16);
legend(‘圆弧’,‘直线’,‘抛物线’,‘旋轮线’);

## 总结

​		本文利用变分法以及MATLAB仿真对最速降线问题进行了探讨，可以看初在给定位置下，旋轮线为下落时间最短的曲线，而解决这一问题的变分法也是后续解决很多问题的强有力工具。与我们所学习的拉格朗日方程不同的是，变分法以及其对应的欧-拉方程引入了变分这一概念，利用泛函的思想求解，而拉格朗日方程则基于我们所学习的虚位移原理以及达朗贝尔原理来得出动力学普遍方程，再对其推导得出两类不同形式的拉格朗日方程，进而求解动力学相关问题。后者也是当下众多仿真软件的理论依据。