快速开始¶

本教程通过一个具有解析解的小型最优控制问题介绍 pockit 的基本工作流程。请先按照安装指南安装 pockit；如需进一步了解系统和阶段对象，请参阅多阶段问题指南。

如需集中查阅解析解推导和数值验证准则，请参阅双积分器数学参考。

背景¶

当加速度作为控制指令时，双积分器可以描述一维运动。虽然它的动力学很简单，但最短时间解包含控制不连续点，因此适合用来初步了解问题定义、自由时间优化和结果验证。

Pockit 使用配点上的数值代替连续的状态和控制函数，再将动力学、边界和目标转录为有限维非线性规划问题，交由数值优化器求解。

问题定义¶

质量为一千克的物体从距离目标一米处静止出发。求有界加速度，使物体在最短时间内到达目标并恢复静止：

\[ \begin{aligned} \min_{x,v,u,t_f}\quad & t_f = \int_0^{t_f} 1\,\mathrm{d}t \\ \text{约束}\quad & \dot{x}=v, \qquad \dot{v}=u, \\ & x(0)=1, \quad v(0)=0, \\ & x(t_f)=0, \quad v(t_f)=0, \\ & -1 \le u(t) \le 1. \end{aligned} \]

解析最优解先以 \(u=-1\) 朝目标加速，在 \(t=1\) 处切换，再以 \(u=1\) 制动，因此 \(t_f=2\)。

变量与单位¶

符号	含义	单位
\(t\)	时间	s
\(x\)	相对目标的位置	m
\(v\)	速度	m/s
\(u\)	加速度指令	m/s\(^2\)
\(t_f\)	自由终止时间	s

建模选择¶

实现采用包含两个状态和一个控制的 Radau 阶段。将终止时间设为 None，即可把 \(t_f\) 作为优化变量。加速度上下界通过路径约束给出，最后一个 True 标记说明该约束具有 bang-bang 特征。归一化时间 \(0.5\) 处的网格边界允许控制在两个多项式区段之间跳变。配点边界向内收紧 \(10^{-7}\ \mathrm{m/s^2}\)；这个可忽略的转录保护量用于吸收求解器容差，使插值后的数值控制仍严格位于物理边界 \(\lvert u\rvert\le1\) 内。

运行示例¶

在仓库根目录运行：

python -m examples.double_integrator

如需非交互运行并保存图片：

python -m examples.double_integrator --save double-integrator.png --no-show

关键实现¶

首先定义阶段、动力学、边界、端点和目标：

from pockit.optimizer import ipopt
from pockit.radau import System, linear_guess

MAX_ACCELERATION = 1.0
COLLOCATION_ACCELERATION_LIMIT = MAX_ACCELERATION - 1.0e-7

system = System(0)
phase = system.new_phase(["position", "velocity"], ["acceleration"])
_, velocity = phase.x
(acceleration,) = phase.u

phase.set_dynamics([velocity, acceleration])
phase.set_integral([1.0])
phase.set_phase_constraint(
    [acceleration],
    [-COLLOCATION_ACCELERATION_LIMIT],
    [COLLOCATION_ACCELERATION_LIMIT],
    True,
)
phase.set_boundary_condition([1.0, 0.0], [0.0, 0.0], 0.0, None)
phase.set_discretization([0.0, 0.5, 1.0], [6, 6])
system.set_phase([phase])
system.set_objective(phase.I[0])

平滑的状态猜测和两段式控制猜测为优化器提供了具有物理意义的起点：

import numpy as np

guess = linear_guess(phase, 0.0)
guess.t_f = 2.2
tau = guess.t_x / guess.t_f
guess.x[0] = 1.0 - 3.0 * tau**2 + 2.0 * tau**3
guess.x[1] = (-6.0 * tau + 6.0 * tau**2) / guess.t_f
guess.u[0] = np.where(
    guess.t_u < guess.t_f / 2.0,
    -COLLOCATION_ACCELERATION_LIMIT,
    COLLOCATION_ACCELERATION_LIMIT,
)

solution, info = ipopt.solve(system, guess)

完整示例还会检查终端条件、终止时间、两段控制的符号，并在 10,001 个物理时刻检查插值后的加速度。对于切换位置未知的问题，可以在求解后使用 System.check() 和 System.refine() 更新网格。

已验证结果¶

当前示例得到的最短时间为 \(2.000000100109\ \mathrm{s}\)，估计切换时间为 \(1.000000050055\ \mathrm{s}\)。二者均在数值容差内符合解析 bang-bang 解。在 10,001 个密集点上，加速度范围为 \([-0.999999900154,\ 0.999999902382]\ \mathrm{m/s^2}\)，严格位于物理边界内。

双积分器最短时间问题的位置、速度和 bang-bang 加速度

源码¶

完整可运行示例见 examples/double_integrator.py。