泛函是一种映射，它以函数为自变量，数为因变量，把函数映射为数。

例如

F[y]=∫abydx F[y] = \int_a^b y \mathrm{d}x

就是一个 C[a,b]→RC[a,b] \to \mathbb{R} 泛函，其中 C[a,b]C[a,b] 表示区间 [a,b][a,b] 上所有连续函数组成的集合。

若一个泛函依赖函数的前 nn 阶导数，则通常将其定义域写为 Cn[a,b]C^{n}[a,b]，这表示区间 [a,b][a,b] 上所有 nn 阶连续可导函数组成的集合。许多重要的 Cn[a,b]→RC^{n}[a,b] \to \mathbb{R} 泛函都可以用定积分表示为

F[y]=∫abL(x,y,y′,y′′,⋯ ,y(n))dx F[y] = \int_a^b \mathcal{L}(x,y,y',y'',\cdots,y^{(n)}) \mathrm{d}x

本文将推导 L=L(x,y,y′)\mathcal{L}=\mathcal{L}(x, y, y') 这类泛函的导数，并用得到的结果求解最速降线问题与狄多问题。

泛函的微分：变分

变分就是泛函的微分。类比于多元函数的微分，变分定义为自变量（函数）沿某一方向的微小扰动所引起因变量（数）变化的线性主部。

设

yε(x)=y(x)+εη(x) y_\varepsilon(x) = y(x) + \varepsilon \eta(x)

于是自变量的微小变化为

δy(x)=εη(x) \delta y(x) = \varepsilon \eta(x)

从而因变量变为

F[yε]=F[y+εη]=∫abL(x,y+εη,y′+εη′)dx F[y_\varepsilon] = F[y + \varepsilon\eta] = \int_{a}^{b} \mathcal{L}(x, y + \varepsilon\eta, y' + \varepsilon\eta') \mathrm{d}x

定义关于 ε\varepsilon 的函数 Φ\Phi

Φ(ε)=F[y+εη] \Phi(\varepsilon) = F[y + \varepsilon\eta]

从而因变量的变化量为

Φ(ε)−Φ(0)=Φ′(0)ε+o(ε) \Phi(\varepsilon) - \Phi(0) = \Phi'(0)\varepsilon + o(\varepsilon)

其线性主部，即变分为

δF=Φ′(0)ε \delta F = \Phi'(0)\varepsilon

于是对 Φ\Phi 求导

Φ′(ε)=ddεF[y+εη]=∫ab∂∂εL(x,y+εη,y′+εη′)dx \Phi'(\varepsilon) = \frac{\mathrm{d}}{\mathrm{d} \varepsilon} F[y + \varepsilon \eta] = \int_{a}^{b} \frac{\partial}{\partial \varepsilon} \mathcal{L}(x, y + \varepsilon\eta, y'+ \varepsilon\eta') \mathrm{d}x

由链式法则

∂L∂ε=∂L∂x∂x∂ε+∂L∂yε∂yε∂ε+∂L∂yε′∂yε′∂ε \frac{\partial \mathcal{L}}{\partial \varepsilon} =\frac{\partial \mathcal{L}}{\partial x} \frac{\partial x}{\partial \varepsilon} +\frac{\partial \mathcal{L}}{\partial y_\varepsilon} \frac{\partial y_\varepsilon}{\partial \varepsilon} +\frac{\partial \mathcal{L}}{\partial y_\varepsilon'} \frac{\partial y_\varepsilon'}{\partial \varepsilon}

代入

∂x∂ε=0∂yε∂ε=η∂yε′∂ε=η′ \begin{align*} \frac{\partial x}{\partial \varepsilon} &= 0 \\ \frac{\partial y_\varepsilon}{\partial \varepsilon} &= \eta \\ \frac{\partial y_\varepsilon'}{\partial \varepsilon} &= \eta' \end{align*}

得到

∂L∂ε=∂L∂yεη+∂L∂yε′η′ \frac{\partial \mathcal{L}}{\partial \varepsilon} = \frac{\partial \mathcal{L}}{\partial y_\varepsilon}\eta + \frac{\partial \mathcal{L}}{\partial y_\varepsilon'}\eta'

从而有

Φ′(ε)=∫ab(∂L∂yεη+∂L∂yε′η′)dx \Phi'(\varepsilon) = \int_{a}^{b} \left( \frac{\partial \mathcal{L}}{\partial y_\varepsilon}\eta + \frac{\partial \mathcal{L}}{\partial y_\varepsilon'}\eta' \right)\mathrm{d}x

令 ε→0\varepsilon \to 0 得到

Φ′(0)=∫ab(∂L∂yη+∂L∂y′η′)dxΦ′(0)ε=∫ab(∂L∂yεη+∂L∂y′εη′)dx \begin{align*} \Phi'(0) &= \int_{a}^{b} \left( \frac{\partial \mathcal{L}}{\partial y}\eta + \frac{\partial \mathcal{L}}{\partial y'}\eta' \right)\mathrm{d}x \\ \Phi'(0)\varepsilon &= \int_{a}^{b} \left( \frac{\partial \mathcal{L}}{\partial y}\varepsilon\eta + \frac{\partial \mathcal{L}}{\partial y'}\varepsilon\eta' \right)\mathrm{d}x \end{align*}

代入 δy(x)=εη(x)\delta y(x) = \varepsilon \eta(x) 得到最终的变分

δF=∫ab(∂L∂yδy+∂L∂y′δy′)dx \boxed{ \delta F = \int_a^b \left( \frac{\partial \mathcal{L}}{\partial y} \delta y + \frac{\partial \mathcal{L}}{\partial y'} \delta y' \right) \mathrm{d}x }

泛函的极值点：欧拉-拉格朗日方程

泛函取到极值的一个必要条件就是变分为零。但上述变分的积分项中同时包含了 δy\delta y 和 δy′\delta y'，不方便判断取到零的条件。因此接下来要尝试消去 δy′\delta y'。

首先将 δF\delta F 展开为两个积分项

δF=∫ab∂L∂yδydx+∫ab∂L∂y′δy′dx \delta F = \int_a^b \frac{\partial \mathcal{L}}{\partial y} \delta y\mathrm{d}x + \int_a^b \frac{\partial \mathcal{L}}{\partial y'} \delta y' \mathrm{d}x

然后对含 δy′\delta y' 的项进行分部积分

∫ab∂L∂y′δy′dx=[∂L∂y′δy]ab−∫abddx(∂L∂y′)δydx \int_a^b \frac{\partial \mathcal{L}}{\partial y'} \delta y' \mathrm{d}x = \left[ \frac{\partial \mathcal{L}}{\partial y'} \delta y \right]_a^b - \int_a^b \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) \delta y \mathrm{d}x

代回原式

δF=∫ab∂L∂yδydx+[∂L∂y′δy]ab−∫abddx(∂L∂y′)δydx=∫ab[∂L∂y−ddx(∂L∂y′)]δydx+[∂L∂y′δy]ab \begin{align*} \delta F &= \int_a^b \frac{\partial \mathcal{L}}{\partial y} \delta y\mathrm{d}x + \left[ \frac{\partial \mathcal{L}}{\partial y'} \delta y \right]_a^b - \int_a^b \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) \delta y \mathrm{d}x \\ &= \int_a^b \left[ \frac{\partial \mathcal{L}}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) \right] \delta y \mathrm{d}x + \left[ \frac{\partial \mathcal{L}}{\partial y'} \delta y \right]_a^b \end{align*}

其中

[∂L∂y′δy]ab=∂L∂y′(b)δy(b)−∂L∂y′(a)δy(a) \left[ \frac{\partial \mathcal{L}}{\partial y'} \delta y \right]_a^b = \frac{\partial \mathcal{L}}{\partial y'}(b) \delta y(b) - \frac{\partial \mathcal{L}}{\partial y'}(a) \delta y(a)

固定边界变分

对于常见的变分问题，自变量函数在区间边界上固定，即满足固定边界条件 δy(a)=0\delta y(a)=0 和 δy(b)=0\delta y(b)=0，从而

[∂L∂y′δy]ab=0 \left[ \frac{\partial \mathcal{L}}{\partial y'} \delta y \right]_a^b = 0

于是

δF=∫ab[∂L∂y−ddx(∂L∂y′)]δydx \boxed{ \delta F = \int_a^b \left[ \frac{\partial \mathcal{L}}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) \right] \delta y \mathrm{d}x }

由于区间 (a,b)(a,b) 内的 δy\delta y 可以任意变化，因此 δF=0\delta F=0 要求和 δy\delta y 相乘的部分必须为零（这就是变分法基本引理），从而得到欧拉-拉格朗日方程

∂L∂y−ddx(∂L∂y′)=0 \boxed{ \frac{\partial \mathcal{L}}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) = 0 }

类比于普通微积分中的导数定义，这里可以推广符号，将上式左边定义为泛函的导数

δFδy=∂L∂y−ddx(∂L∂y′) \frac{\delta F}{\delta y} = \frac{\partial \mathcal{L}}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right)

这样泛函的变分就可视为泛函变分密度的积分

δF=∫abδFδyδydx \delta F = \int_a^b \frac{\delta F}{\delta y} \delta y \mathrm{d}x

类似于多变量函数，其全微分为各自变量的偏微分之和

dF=∑i=1n∂F∂yidyi \mathrm{d} F = \sum_{i=1}^n \frac{\partial F}{\partial y_i} \mathrm{d}y_i

只不过泛函的自变量是无穷维的、需要用连续变量标记的向量，因此不同于有限维向量中的求和，需要使用积分。

自由边界变分

如果不假设边界固定，让右端点 x=bx=b 可以自由变化（包括横向和纵向的变化），那么此时问题会更复杂一点。

由于右端点的横坐标 bb 可以变化，因此泛函需要改写为

F[y,b]=∫abL(x,y,y′)dx F[y, b] = \int_a^b \mathcal{L}(x,y,y') \mathrm{d}x

然后类似地定义 bb 的变化量

bε=b+εζ b_\varepsilon = b + \varepsilon\zeta

从而

δb=εζ \delta b=\varepsilon\zeta

泛函变为

F[y,bε]=F[y,b+εζ]=∫ab+εζL(x,y,y′)dx F[y, b_\varepsilon] = F[y, b + \varepsilon\zeta] = \int_a^{b + \varepsilon\zeta} \mathcal{L}(x,y,y') \mathrm{d}x

同样定义一个辅助函数 Φ\Phi

Φ(ε)=F[y,b+εζ] \Phi(\varepsilon) = F[y, b + \varepsilon\zeta]

泛函的变化量为

F[y,b+εζ]−F[y,b]=∫ab+εζL(x,y,y′)dx−∫abL(x,y,y′)dxΦ(ε)−Φ(0)=∫bb+εζL(x,y,y′)dx \begin{align*} F[y, b + \varepsilon\zeta] - F[y, b] &= \int_a^{b + \varepsilon\zeta} \mathcal{L}(x,y,y') \mathrm{d}x - \int_a^b \mathcal{L}(x,y,y') \mathrm{d}x \\ \Phi(\varepsilon) - \Phi(0) &= \int_b^{b + \varepsilon\zeta} \mathcal{L}(x,y,y') \mathrm{d}x \end{align*}

展开 Φ(ε)−Φ(0)\Phi(\varepsilon) - \Phi(0) 得到线性主部

Φ(ε)−Φ(0)=Φ′(0)ε+o(ε) \Phi(\varepsilon) - \Phi(0) = \Phi'(0)\varepsilon + o(\varepsilon)

而

Φ′(ε)=ddεF[y,b+εζ]=∂∂ε∫bb+εζL(x,y,y′)dx=ζL(b+εζ,y(b+εζ),y′(b+εζ)) \begin{align*} \Phi'(\varepsilon) &= \frac{\mathrm{d}}{\mathrm{d} \varepsilon} F[y, b + \varepsilon\zeta] \\ &= \frac{\partial}{\partial \varepsilon} \int_b^{b + \varepsilon\zeta} \mathcal{L}(x,y,y') \mathrm{d}x \\ &= \zeta\mathcal{L}(b+\varepsilon\zeta, y(b+\varepsilon\zeta), y'(b+\varepsilon\zeta)) \end{align*}

得到端点变化产生的变分

Φ′(0)ε=L(b,y(b),y′(b))δb \Phi'(0)\varepsilon = \mathcal{L}(b, y(b), y'(b))\delta b

与曲线变化产生的变分加在一起，得到总的变分

δF=L(b,y(b),y′(b))δb+∫ab(∂L∂yδy+∂L∂y′δy′)dx \delta F = \mathcal{L}(b, y(b), y'(b))\delta b + \int_a^b \left( \frac{\partial \mathcal{L}}{\partial y} \delta y + \frac{\partial \mathcal{L}}{\partial y'} \delta y' \right) \mathrm{d}x

然后进行同样的分部积分，但由于右端点不固定，所以相关项不为零，得到

δF=L(b,y(b),y′(b))δb+∂L∂y′(b)δy(b)+∫ab[∂L∂y−ddx(∂L∂y′)]δydx \boxed{ \delta F = \mathcal{L}(b, y(b), y'(b))\delta b + \frac{\partial \mathcal{L}}{\partial y'}(b) \delta y(b) + \int_a^b \left[ \frac{\partial \mathcal{L}}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) \right] \delta y \mathrm{d}x }

此时让 δF=0\delta F=0 不仅要满足欧拉-拉格朗日方程，还多出了一个额外的条件

L(b,y(b),y′(b))δb+∂L∂y′(b)δy(b)=0 \mathcal{L}(b, y(b), y'(b))\delta b + \frac{\partial \mathcal{L}}{\partial y'}(b) \delta y(b) = 0

该条件还可以进一步讨论，从而分为只允许纵坐标变化、只允许横坐标变化和可以任意变化三种情况。

经过扰动（包括端点的扰动和曲线的扰动）后，最终的右端点变为

Bε(bε,yε(bε)) B_{\varepsilon}(b_\varepsilon, y_\varepsilon(b_\varepsilon))

其纵坐标的变化量为

δyB=yε(bε)−y(b) \delta y_B = y_\varepsilon(b_\varepsilon) - y(b)

而

yε(bε)=y(b+εζ)+εη(b+εζ) y_\varepsilon(b_\varepsilon) = y(b+\varepsilon\zeta) + \varepsilon\eta(b+\varepsilon\zeta)

对 y(b+εζ)y(b+\varepsilon\zeta) 和 η(b+εζ)\eta(b+\varepsilon\zeta) 分别进行一阶展开

{y(b+εζ)=y(b)+εy′(b)ζ+o(ε)η(b+εζ)=η(b)+εη′(b)ζ+o(ε) \begin{cases} y(b+\varepsilon\zeta) = y(b) + \varepsilon y'(b)\zeta + o(\varepsilon) \\ \eta(b+\varepsilon\zeta) = \eta(b) + \varepsilon\eta'(b)\zeta + o(\varepsilon) \end{cases}

忽略二阶小量 ε2η′(b)ζ\varepsilon^2\eta'(b)\zeta 和高阶小量 o(ε)o(\varepsilon)，代回原式

δyB=y(b)+εy′(b)ζ+εη(b)−y(b)=εy′(b)ζ+εη(b)=y′(b)δb+δy(b) \begin{align*} \delta y_B &= y(b) + \varepsilon y'(b)\zeta + \varepsilon\eta(b) - y(b) \\ &= \varepsilon y'(b)\zeta + \varepsilon\eta(b) \\ &= y'(b)\delta b + \delta y(b) \end{align*}

移项得到

δy(b)=δyB−y′(b)δb \delta y(b) = \delta y_B - y'(b)\delta b

再把上式代入原来的条件

L(b,y(b),y′(b))δb+∂L∂y′(b)[δyB−y′(b)δb]=0[L(b,y(b),y′(b))−y′(b)∂L∂y′(b)]δb+∂L∂y′(b)δyB=0 \begin{align*} \mathcal{L}(b, y(b), y'(b))\delta b + \frac{\partial \mathcal{L}}{\partial y'}(b) [\delta y_B - y'(b)\delta b] &= 0 \\ \left[ \mathcal{L}(b, y(b), y'(b))-y'(b)\frac{\partial \mathcal{L}}{\partial y'}(b) \right]\delta b + \frac{\partial \mathcal{L}}{\partial y'}(b) \delta y_B &= 0 \end{align*}

因此

  • 当端点横坐标可自由变化而纵坐标固定时，与 δb\delta b 相乘的部分必须为零，得到
    L(b,y(b),y′(b))−y′(b)∂L∂y′(b)=0 \boxed{ \mathcal{L}(b, y(b), y'(b))-y'(b)\frac{\partial \mathcal{L}}{\partial y'}(b) = 0 }
  • 当端点的纵坐标可自由变化而横坐标固定时，与 δyB\delta y_B 相乘的部分必须为零，得到
    ∂L∂y′(b)=0 \boxed{ \frac{\partial \mathcal{L}}{\partial y'}(b) = 0 }
  • 当端点的横纵坐标均可自由变化时，上面两项都要为零，稍微化简后便得到
    ∂L∂y′(b)=0L(b,y(b),y′(b))=0 \boxed{ \frac{\partial \mathcal{L}}{\partial y'}(b) = 0 \qquad \mathcal{L}(b, y(b), y'(b)) = 0 }

上述条件被称作自然边界条件。

另一种情况是，端点被限制在曲线 y=ψ(x)y=\psi(x) 上，此时有

δyB=ψ′(b)δb \delta y_B = \psi'(b) \delta b

代入得到

[L(b,y(b),y′(b))−y′(b)∂L∂y′(b)]δb+ψ′(b)∂L∂y′(b)δb=0[L(b,y(b),y′(b))−y′(b)∂L∂y′(b)+ψ′(b)∂L∂y′(b)]δb=0 \begin{align*} \left[ \mathcal{L}(b, y(b), y'(b))-y'(b)\frac{\partial \mathcal{L}}{\partial y'}(b) \right]\delta b + \psi'(b)\frac{\partial \mathcal{L}}{\partial y'}(b)\delta b &= 0 \\ \left[ \mathcal{L}(b, y(b), y'(b))-y'(b)\frac{\partial \mathcal{L}}{\partial y'}(b) + \psi'(b)\frac{\partial \mathcal{L}}{\partial y'}(b) \right]\delta b &= 0 \end{align*}

从而有

L(b,y(b),y′(b))−y′(b)∂L∂y′(b)+ψ′(b)∂L∂y′(b)=0 \boxed{ \mathcal{L}(b, y(b), y'(b))-y'(b)\frac{\partial \mathcal{L}}{\partial y'}(b) + \psi'(b)\frac{\partial \mathcal{L}}{\partial y'}(b)= 0 }

这叫做横截条件。

变分法的应用：最速降线问题

最速降线问题

设一质点从原点 A=(0,0)A=(0,0) 运动到 B=(w,h)B=(w, h) 处（w,h>0w, h > 0 且 yy 轴正方向向下），初速度为零，重力加速度为 gg，忽略所有摩擦力。质点沿哪条曲线运动所需的时间最短？

这是固定边界变分问题。

将解表示为时间泛函的极值点

设曲线为

y=y(x) y = y(x)

由能量守恒得到 (x,y)(x,y) 处的速度为

v=2gy v = \sqrt{2gy}

弧长微元为

ds=1+(y′)2dx \mathrm{d}s = \sqrt{1 + (y')^2}\mathrm{d}x

从而时间可表示为一个泛函

T[y]=∫ABdsv=∫0w1+(y′)22gydx T[y] = \int_A^B \frac{\mathrm{d}s}{v} = \int_0^{w} \frac{\sqrt{1 + (y')^2}}{\sqrt{2gy}} \mathrm{d}x

因此最速降线就是使 T[y]T[y] 取极小值的函数。

求时间泛函的极值点

由于常数 2g\sqrt{2g} 不影响泛函的极值点，因此定义拉格朗日函数时忽略掉它

L(y,y′)=1+(y′)2y=y−12[1+(y′)2]12 \mathcal{L}(y, y') = \frac{\sqrt{1 + (y')^2}}{\sqrt{y}} = y^{-\frac{1}{2}} [1 + (y')^2]^{\frac{1}{2}}

显然 L\mathcal{L} 中不包含自变量 xx。下面利用这一性质计算欧拉-拉格朗日方程的首次积分——这也叫贝尔特拉米恒等式。

首先得到

∂L∂x=0 \frac{\partial \mathcal{L}}{\partial x} = 0

然后计算全导数

dLdx=∂L∂x+∂L∂ydydx+∂L∂y′dy′dx=y′∂L∂y+y′′∂L∂y′ \frac{\mathrm{d}\mathcal{L}}{\mathrm{d}x} = \frac{\partial \mathcal{L}}{\partial x} + \frac{\partial \mathcal{L}}{\partial y} \frac{\mathrm{d}y}{\mathrm{d}x} + \frac{\partial \mathcal{L}}{\partial y'} \frac{\mathrm{d}y'}{\mathrm{d}x} = y'\frac{\partial \mathcal{L}}{\partial y} + y''\frac{\partial \mathcal{L}}{\partial y'}

接着把欧拉-拉格朗日方程 ∂L∂y=ddx(∂L∂y′)\frac{\partial \mathcal{L}}{\partial y} = \frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) 代入其中

dLdx=y′ddx(∂L∂y′)+y′′∂L∂y′ \frac{\mathrm{d}\mathcal{L}}{\mathrm{d}x} = y'\frac{\mathrm{d}}{\mathrm{d}x}\left( \frac{\partial \mathcal{L}}{\partial y'} \right) + y''\frac{\partial \mathcal{L}}{\partial y'}

反向使用导数的乘积法则，得到

dLdx=ddx(y′∂L∂y′) \frac{\mathrm{d}\mathcal{L}}{\mathrm{d}x} = \frac{\mathrm{d}}{\mathrm{d}x}\left( y'\frac{\partial \mathcal{L}}{\partial y'} \right)

再移项合并求导

ddx(L−y′∂L∂y′)=0 \frac{\mathrm{d}}{\mathrm{d}x} \left( \mathcal{L} - y'\frac{\partial \mathcal{L}}{\partial y'} \right) = 0

最后对两侧同时积分，得到贝尔特拉米恒等式

L−y′∂L∂y′=C \mathcal{L} - y'\frac{\partial \mathcal{L}}{\partial y'} = C

原先的欧拉-拉格朗日方程中还有二阶导数，现在只剩下一阶导数了。接下来的任务就是求解这个一阶微分方程。

首先计算偏导数

∂L∂y′=y−1212[1+(y′)2]−12(2y′)=y′y[1+(y′)2] \frac{\partial \mathcal{L}}{\partial y'} = y^{-\frac{1}{2}} \frac{1}{2} [1 + (y')^2]^{-\frac{1}{2}} (2y') = \frac{y'}{\sqrt{y[1 + (y')^2]}}

然后把 L\mathcal{L} 和 ∂L∂y′\frac{\partial \mathcal{L}}{\partial y'} 代入贝尔特拉米恒等式

1+(y′)2y−y′y′y[1+(y′)2]=C1y[1+(y′)2]=C1y[1+(y′)2]=C2y[1+(y′)2]=1C2 \begin{align*} \frac{\sqrt{1 + (y')^2}}{\sqrt{y}} - y' \frac{y'}{\sqrt{y[1 + (y')^2]}} &= C \\ \frac{1}{\sqrt{y[1+(y')^2]}} &= C \\ \frac{1}{y[1+(y')^2]} &= C^2 \\ y[1+(y')^2] &= \frac{1}{C^2} \end{align*}

这里更换一下常数符号 1C2=2r\frac{1}{C^2} = 2r ，得到

y[1+(y′)2]=2r y[1+(y')^2] = 2r

接下来是需要注意力的三角变换。令

y′=cot⁡θ2 y'=\cot\frac{\theta}{2}

代入上式得到

y(1+cot⁡2θ2)=2rycsc⁡2θ2=2ry=2rsin⁡2θ2y=r(1−cos⁡θ) \begin{align*} y\left(1+\cot^2\frac{\theta}{2}\right) &= 2r \\ y \csc^2\frac{\theta}{2} &= 2r \\ y &= 2r\sin^2\frac{\theta}{2} \\ y &= r(1-\cos\theta) \end{align*}

现在有了 yy 关于 θ\theta 的参数方程，再得到 xx 关于 θ\theta 的参数方程就可以确定曲线了。

先求 dy\mathrm{d}y

dy=rsin⁡θdθ \mathrm{d}y = r \sin\theta \mathrm{d}\theta

然后间接计算 dx\mathrm{d}x

dx=dyy′=rsin⁡θdθcot⁡θ2=r(2sin⁡θ2cos⁡θ2sin⁡θ2cos⁡θ2)dθ=r(2sin⁡2θ2)dθ=r(1−cos⁡θ)dθ \begin{align*} \mathrm{d}x &= \frac{\mathrm{d}y}{y'} \\ &= \frac{r \sin\theta \mathrm{d}\theta}{\cot\frac{\theta}{2}} \\ &= r (2\sin\frac{\theta}{2}\cos\frac{\theta}{2}\frac{\sin\frac{\theta}{2}}{\cos\frac{\theta}{2}}) \mathrm{d}\theta \\ &= r (2\sin^2\frac{\theta}{2}) \mathrm{d}\theta \\ &= r (1 - \cos\theta) \mathrm{d}\theta \end{align*}

最后两侧同时积分

x=r(θ−sin⁡θ)+C x = r(\theta - \sin\theta) + C

代入曲线必过的原点坐标 A=(0,0)A=(0,0) 即可确定常数值

C=0 C = 0

于是最速降线为

{x=r(θ−sin⁡θ)y=r(1−cos⁡θ) \begin{cases} x = r(\theta - \sin\theta) \\ y = r(1-\cos\theta) \end{cases}

这正是摆线的参数方程。

最速降线问题的解

代入曲线经过的另一个点 B=(w,h)B=(w,h) 得到方程

{w=r(θB−sin⁡θB)h=r(1−cos⁡θB) \begin{cases} w = r(\theta_B - \sin\theta_B) \\ h = r(1-\cos\theta_B) \end{cases}

消去 rr 后得到

θB−sin⁡θB1−cos⁡θB=wh \frac{\theta_B - \sin\theta_B}{1-\cos\theta_B} = \frac{w}{h}

可惜这是一个超越方程，无法得到 θB\theta_B 的初等函数表示。

假设用数值方法求出了 0<θB<2π0 < \theta_B < 2\pi，代入便可得到

r=h1−cos⁡θB=wθB−sin⁡θB r = \frac{h}{1-\cos\theta_B} = \frac{w}{\theta_B - \sin\theta_B}

代回参数方程便得到了最终的曲线

{x=wθB−sin⁡θB(θ−sin⁡θ)y=h1−cos⁡θB(1−cos⁡θ)θ∈[0,θB] \boxed{ \begin{cases} x = \frac{w}{\theta_B - \sin\theta_B}(\theta - \sin\theta) \\ y = \frac{h}{1-\cos\theta_B}(1-\cos\theta) \\ \theta \in [0, \theta_B] \end{cases} }

需要注意的是，这条曲线并不严格满足 C1[0,w]C^1[0,w]，因为在 θ=0\theta=0 处斜率无穷大，不可求导。但变分法的实际要求更为宽松，因此结论仍然有效。

最短的时间也可以计算出来。首先是速度

v=2gy=2gr(1−cos⁡θ)=2grsin⁡θ2 v = \sqrt{2gy} = \sqrt{2gr(1-\cos\theta)} = 2\sqrt{gr}\sin\frac{\theta}{2}

然后求弧长微元

ds=dx2+dy2=(dxdθ)2+(dydθ)2dθ \mathrm{d}s = \sqrt{\mathrm{d}x^2 + \mathrm{d}y^2} = \sqrt{\left(\frac{\mathrm{d}x}{\mathrm{d}\theta}\right)^2 + \left( \frac{\mathrm{d}y}{\mathrm{d}\theta} \right)^2 }\mathrm{d}\theta

分别得到

{dxdθ=r(1−cos⁡θ)dydθ=rsin⁡θ \begin{cases} \frac{\mathrm{d}x}{\mathrm{d}\theta} = r(1-\cos\theta) \\ \frac{\mathrm{d}y}{\mathrm{d}\theta} = r\sin\theta \end{cases}

于是

ds=r2(1−cos⁡θ)2+r2sin⁡2θdθ=2rsin⁡θ2dθ \mathrm{d}s = \sqrt{r^2(1-\cos\theta)^2+r^2\sin^2\theta}\mathrm{d}\theta = 2r\sin\frac{\theta}{2}\mathrm{d}\theta

从而

Tmin=∫ABdsv=∫0θBrgdθ=θBrg T_{\text{min}} = \int_{A}^{B} \frac{\mathrm{d}s}{v} = \int_{0}^{\theta_B} \sqrt{\frac{r}{g}}\mathrm{d}\theta = \theta_B\sqrt{\frac{r}{g}}

代入 r=h1−cos⁡θBr=\frac{h}{1-\cos\theta_B} 得到最短时间

Tmin=θBhg(1−cos⁡θB) \boxed{ T_{\text{min}} = \theta_B\sqrt{\frac{h}{g(1-\cos\theta_B)}} }

可以对比一下直线下降所需的时间。此时质点做匀加速直线运动，初速度为 00，终速度为 2gh\sqrt{2gh}，然后由平均速度乘时间等于路程得到

2gh2Tline=w2+h2 \frac{\sqrt{2gh}}{2}T_{\text{line}} = \sqrt{w^2 + h^2}

即

Tline=2w2+2h2gh \boxed{ T_{\text{line}} = \sqrt{ \frac{2w^2 + 2h^2}{gh} } }

另外，本例其实只证明了取极值的必要条件，并没有证明充分性。因此实际的结论是：如果存在最速降线，那么一定是前文求出的摆线，且时间为前文求出的时间。

  const board = JXG.JSXGraph.initBoard(BOARDID, {
  boundingbox: [-1, -1, 10, 10],
  axis: true,
  showCopyright: false,
  showNavigation: false,
});

const g = 9.8

function theta(w, h) {
  let lastW;
  let lastH;
  let lastTheta;

  if (w === lastW && h === lastH) {
    return lastTheta;
  }

  lastW = w;
  lastH = h;
  lastTheta = JXG.Math.Numerics.fzero(
    (t) => (t - Math.sin(t)) - w / h * (1 - Math.cos(t)),
    [1e-2, 2 * Math.PI - 1e-2]
  );

  return lastTheta;
}

function r(w, h) {
  return h / (1 - Math.cos(theta(w, h)));
}

function timeCycloid(w, h) {
  return theta(w, h) * Math.sqrt(r(w, h) / g);
}

function timeLine(w, h) {
  return Math.sqrt(2 * (w * w + h * h) / (g * h))
}

const A = board.create("point", [0, 0], {
  name: "A",
  fixed: true
});
const B = board.create("point", [2, 1], { name: "B" });

board.create('line', [A, B], {
  straightFirst: false,
  straightLast: false,
  strokeColor: 'green',
  strokeWidth: 3
})
board.create('curve', [
  (t) => r(B.X(), B.Y()) * (t - Math.sin(t)),
  (t) => r(B.X(), B.Y()) * (1 - Math.cos(t)),
  0, () => theta(B.X(), B.Y())
], {
  strokeWidth: 3
});

board.create("text", [
  () => B.X() + 0.25,
  () => B.Y() + 0.25,
  () => `(${B.X().toFixed(2)}, ${B.Y().toFixed(2)})`
])
board.create("text", [
  0.5,
  8.5,
  () => `
    摆线时间 ${timeCycloid(B.X(), B.Y())}<br>
    直线时间 ${timeLine(B.X(), B.Y())}
  `
]);

变分法的应用：狄多问题

狄多问题

设一曲线从原点 A=(0,0)A=(0, 0) 出发，经过上半平面回到 xx 轴 B=(d,0)B=(d, 0) 处（d>0d>0 且可任取），其总长度为 ll。哪条曲线与 xx 轴围成的面积最大？

这是右端点横坐标可变化的自由边界变分问题。

将解表示为面积泛函的极值点

设曲线为

y=y(x) y=y(x)

狄多问题有一个约束条件，即长度需满足

l=∫ABds=∫0d1+(y′)2dx l = \int_A^B \mathrm{d}s = \int_0^d \sqrt{1+(y')^2} \mathrm{d}x

围成的面积可表示为

S[y,d]=∫0dydx S[y, d] = \int_0^d y \mathrm{d}x

因此与 xx 轴围成面积最大的曲线就是使 S[y,d]S[y, d] 取极大值，且满足约束条件 ∫0d1+(y′)2dx−l=0\int_0^d \sqrt{1+(y')^2} \mathrm{d}x - l = 0 的函数。

求面积泛函的极值点

由于带有约束条件，因此需要用拉格朗日乘子构造一个辅助泛函

J[y,d]=S[y,d]+λ(∫0d1+(y′)2dx−l)=∫0dydx+∫0dλ1+(y′)2dx−λl=∫0d(y+λ1+(y′)2)dx−λl \begin{align*} J[y, d] &= S[y, d] + \lambda \left( \int_0^d \sqrt{1+(y')^2} \mathrm{d}x - l \right) \\ &= \int_0^d y \mathrm{d}x + \int_0^d \lambda\sqrt{1+(y')^2} \mathrm{d}x - \lambda l\\ &= \int_0^d \left( y + \lambda\sqrt{1 + (y')^2} \right)\mathrm{d}x - \lambda l \end{align*}

再对 y,dy, d 变分时 λ\lambda 视为常数，因此 λl\lambda l 不改变泛函的极值条件，定义拉格朗日函数时忽略掉它

L(y,y′)=y+λ1+(y′)2=y+λ[1+(y′)2]12 \mathcal{L}(y, y') = y + \lambda\sqrt{1 + (y')^2} = y + \lambda[1 + (y')^2]^{\frac{1}{2}}

由于问题是端点横坐标可自由变化的变分问题，因此极值点需要满足欧拉-拉格朗日方程和自然边界条件。

先考虑欧拉-拉格朗日方程。由于拉格朗日函数中不包含自变量 xx，因此复用之前求出的贝尔特拉米恒等式，得到

L−y′∂L∂y′=C \mathcal{L} - y'\frac{\partial \mathcal{L}}{\partial y'} = C

计算偏导数

∂L∂y′=λ12[1+(y′)2]−12(2y′)=λy′1+(y′)2 \frac{\partial \mathcal{L}}{\partial y'} = \lambda\frac{1}{2}[1 + (y')^2]^{-\frac{1}{2}}(2y') = \lambda \frac{y'}{\sqrt{1+(y')^2}}

代入 L\mathcal{L} 和 ∂L∂y′\frac{\partial \mathcal{L}}{\partial y'}

y+λ1+(y′)2−y′λy′1+(y′)2=Cy+λ11+(y′)2=C \begin{align*} y + \lambda\sqrt{1 + (y')^2} - y'\lambda \frac{y'}{\sqrt{1+(y')^2}} &= C \\ y + \lambda\frac{1}{\sqrt{1+(y')^2}} &= C \end{align*}

然后化简

λ1+(y′)2=C−yλ21+(y′)2=(C−y)21+(y′)2=λ2(C−y)2(y′)2=λ2−(C−y)2(C−y)2y′=±λ2−(C−y)2C−y \begin{align*} \frac{\lambda}{\sqrt{1+(y')^2}} &= C - y \\ \frac{\lambda^2}{1+(y')^2} &= (C - y)^2 \\ 1+(y')^2 &= \frac{\lambda^2}{(C-y)^2} \\ (y')^2 &= \frac{\lambda^2 - (C-y)^2}{(C-y)^2} \\ y' &= \pm \frac{\sqrt{\lambda^2 - (C-y)^2}}{C-y} \end{align*}

接着分离变量

dydx=±λ2−(C−y)2C−yC−yλ2−(C−y)2dy=±dx \begin{align*} \frac{\mathrm{d}y}{\mathrm{d}x} &= \pm \frac{\sqrt{\lambda^2 - (C-y)^2}}{C-y} \\ \frac{C-y}{\sqrt{\lambda^2 - (C-y)^2}} \mathrm{d}y &= \pm \mathrm{d}x \end{align*}

再进行积分变换。令 u=C−yu = C-y，于是

du=−dy \mathrm{d}u = -\mathrm{d}y

代入得到

−uλ2−u2du=±dx \frac{-u}{\sqrt{\lambda^2 - u^2}} \mathrm{d}u = \pm \mathrm{d}x

两侧同时积分

λ2−u2=±(x−C′)λ2−u2=(x−C′)2 \begin{align*} \sqrt{\lambda^2 - u^2} &= \pm (x - C') \\ \lambda^2 - u^2 &= (x - C')^2 \end{align*}

代回 u=C−yu = C-y 得到

λ2−(C−y)2=(x−C′)2 \lambda^2 - (C-y)^2 = (x - C')^2

令 C=y0,C′=x0,r=∣λ∣C=y_0, C'=x_0, r=|\lambda|，便得到了标准的圆方程

(x−x0)2+(y−y0)2=r2 (x - x_0)^2 + (y-y_0)^2 = r^2

因此曲线是过上半平面的一段圆弧

y(x)=y0+r2−(x−x0)2 y(x) = y_0 + \sqrt{r^2-(x-x_0)^2}

再来考虑自然边界条件。由于只允许右端点的横坐标变化，因此

(L−y′∂L∂y′)∣x=d=0 \left. \left(\mathcal{L}-y'\frac{\partial \mathcal{L}}{\partial y'} \right)\right|_{x=d} = 0

而贝尔特拉米恒等式要求 L−y′∂L∂y′=C\mathcal{L}-y'\frac{\partial \mathcal{L}}{\partial y'} = C，因此有

C=0 C = 0

由于此前令 C=y0C=y_0，于是可知圆心在 xx 轴上

y(x)=r2−(x−x0)2 y(x) = \sqrt{r^2-(x-x_0)^2}

代入曲线必过的原点 A(0,0)A(0,0)，并由 r,x0>0r, x_0>0 得到

0=r2−x02r=x0 \begin{align*} 0 &= \sqrt{r^2 - x_0^2} \\ r &= x_0 \end{align*}

从而有

y(x)=r2−(x−r)2=2rx−x2 y(x) = \sqrt{r^2-(x-r)^2} = \sqrt{2rx - x^2}

因此最终的圆弧一定是个以 ABAB 为直径的上半圆。

狄多问题的解

已知半圆的弧长为

l=πr l = \pi r

从而

r=lπ r = \frac{l}{\pi}

因此右端点的横坐标为

d=2r=2lπ d = 2r = \frac{2l}{\pi}

得到最终的曲线表达式

y(x)=2lπx−x2x∈[0,2lπ] \boxed{ y(x) = \sqrt{\frac{2l}{\pi}x-x^2} \qquad x\in\left[0,\frac{2l}{\pi}\right] }

同样地，上述曲线并不满足 C1[0,d]C^1[0,d]，因为在 x=0,x=dx=0,x=d 处斜率都为无穷大。但变分法的适用条件更为宽松，因此结论仍然有效。

还可以顺便求出最大的面积

Smax=12πr2=l22π S_{\text{max}} = \frac{1}{2}\pi r^2 = \frac{l^2}{2\pi}

即

Smax=l22π \boxed{ S_{\text{max}} = \frac{l^2}{2\pi} }

同样地，本例只证明了极值曲线的必要条件。因此，实际的结论是：如果存在与 xx 轴围成面积最大的曲线，那么这条曲线一定是前文求出的半圆，且面积为前文求出的面积。

  const board = JXG.JSXGraph.initBoard(BOARDID, {
  boundingbox: [-2, 10, 20, -1],
  axis: true,
  showCopyright: false,
  showNavigation: false,
  keepAspectRatio: true
});

function r(l) {
  return l / Math.PI
}

function d(l) {
  return 2 * l / Math.PI
}

function S(l) {
  return l * l / (2 * Math.PI)
}

const l = board.create('slider', [
  [6, 9],
  [16, 9],
  [0, 10, 20]
], {
  name: 'l'
});

const A = board.create("point", [0, 0], {
  name: "A",
  fixed: true
});
const B = board.create("point", [
  () => d(l.Value()),
  0
], {
  name: "B",
  fixed: true
});

board.create('semicircle', [A, B])
board.create("text", [
  10,
  8,
  () => `面积 ${S(l.Value())}`
]);