为何我参数化的最速降线（Brachistochrone）摆线路径倒置？

问题根源分析 #

你的代码存在4个核心错误，导致路径凹性不符合预期：

1. 参数方程平移方向错误：将终点设为摆线的theta=0基准点，而非起点，导致路径方向完全颠倒
2. y轴参数符号错误：摆线的y坐标计算方向与重力势能下落方向相反，导致路径向上抬升而非下落
3. 速度物理意义错误：直接使用当前y坐标计算速度，而非初始高度与当前高度的差值，违背势能转动能的物理规律
4. 斜率匹配函数错误：目标函数的斜率逻辑与摆线的实际斜率公式不匹配，导致求解的theta值不符合要求

修正后的完整代码 #

import numpy as np
from scipy.optimize import newton
from scipy.integrate import quad
import matplotlib.pyplot as plt


class BrachistochroneCycloidPath():

    def __init__(self, initial_position, final_position, g):
        self.initial_position = initial_position
        self.final_position = final_position
        self.xi = initial_position[0]
        self.yi = initial_position[1]
        self.xf = final_position[0]
        self.yf = final_position[1]
        self.g = g
        self._initial_theta = 0
        self._final_theta = None
        self._radius = None
        self._arc_length = None
        self._elapsed_duration = None

    @property
    def initial_theta(self):
        return self._initial_theta

    @property
    def final_theta(self):
        return self._final_theta

    @property
    def radius(self):
        return self._radius

    @property
    def arc_length(self):
        return self._arc_length

    @property
    def elapsed_duration(self):
        return self._elapsed_duration

    def initialize_theta_and_radius(self, x0):
        # 修正斜率匹配函数：匹配摆线实际斜率公式
        f = lambda theta : (self.yi - self.yf)/(self.xf - self.xi) - (1 - np.cos(theta))/(theta - np.sin(theta))
        final_theta = newton(f, x0=x0)
        # 修正半径计算逻辑，基于起点到终点的高度差
        radius = (self.yi - self.yf) / (1 - np.cos(final_theta))
        self._final_theta = final_theta
        self._radius = radius

    def initialize_arc_length(self):
        integral = quad(self.ds, self.initial_theta, self.final_theta)
        self._arc_length = abs(integral[0])       

    def initialize_elapsed_duration(self):
        f = lambda theta : self.ds(theta) / self.v(theta)
        integral = quad(f, self.initial_theta, self.final_theta)
        self._elapsed_duration = abs(integral[0])       

    # 修正参数方程：以起点为theta=0基准点
    def x(self, theta):
        return self.xi + self.radius * (theta - np.sin(theta))

    def y(self, theta):
        return self.yi - self.radius * (1 - np.cos(theta))

    def dx(self, theta):
        return self.radius * (1 - np.cos(theta))

    def dy(self, theta):
        return self.radius * np.sin(theta)

    def ds(self, theta):
        dx = self.dx(theta)
        dy = self.dy(theta)
        return np.sqrt(np.square(dx) + np.square(dy))

    # 修正速度计算：基于下落高度（初始高度-当前高度）
    def v(self, theta):
        return np.sqrt(2 * self.g * (self.yi - self.y(theta)))


if __name__ == "__main__":
    ## 初始化输入参数
    initial_position = (1, 10)
    final_position = (10, 1)
    x0 = np.pi  # 使用正数初始猜测，匹配摆线的theta范围
    g = 9.8

    ## 初始化路径
    brachistochrone_path = BrachistochroneCycloidPath(
        initial_position=initial_position,
        final_position=final_position,
        g=g)
    brachistochrone_path.initialize_theta_and_radius(x0=x0)
    brachistochrone_path.initialize_arc_length()
    brachistochrone_path.initialize_elapsed_duration()

    ## 打印验证
    print(f"\n\tArc-Length [m]:\n{brachistochrone_path.arc_length:.4f}\n")
    print(f"\n\tElapsed Duration [s]:\n{brachistochrone_path.elapsed_duration:.4f}\n")

    ## 绘制路径
    theta = np.linspace(
        brachistochrone_path.initial_theta,
        brachistochrone_path.final_theta,
        100)

    path_x = brachistochrone_path.x(theta=theta)
    path_y = brachistochrone_path.y(theta=theta)
    path_label = f"最速降线摆线 (半径 = {brachistochrone_path.radius:.4f})"

    positions_x = [initial_position[0], final_position[0]]
    positions_y = [initial_position[1], final_position[1]]
    positions_label = "起点/终点"

    fig, ax = plt.subplots()
    ax.plot(path_x, path_y, color="darkorange", label=path_label)
    ax.scatter(positions_x, positions_y, color="black", label=positions_label)
    fig.subplots_adjust(bottom=0.2)
    fig.legend(mode="expand", loc="lower center", ncol=2)
    ax.grid(True)
    plt.show()
    plt.close()

关键修正说明 #

1. 参数方程调整：将起点设为theta=0的基准点，x坐标随theta增大向右延伸，y坐标随theta增大先快速下降再平缓过渡，符合最速降线的凹性要求
2. 斜率函数修正：使用(yi-yf)/(xf-xi)作为目标斜率，匹配摆线的(1-cosθ)/(θ-sinθ)斜率公式，确保求解的theta值正确
3. 速度计算修正：基于下落高度yi - y(theta)计算速度，符合重力势能转化为动能的物理规律
4. 初始猜测值调整：使用正数np.pi作为theta的初始猜测，确保求解的theta在合理范围内（正数，对应摆线的正常滚动方向）

内容的提问来源于stack exchange，提问作者user23358153