非线性规划中的原始-对偶内点法(Primal-Dual Interior-Point Method in Nonlinear Programming)
字数 5456 2025-12-18 22:08:05

非线性规划中的原始-对偶内点法(Primal-Dual Interior-Point Method in Nonlinear Programming)

非线性规划中的原始-对偶内点法是一种用于求解具有不等式约束的非线性优化问题的强大算法。它通过结合原始变量、对偶变量和障碍参数,在可行域内部生成一系列迭代点,并最终收敛到满足一阶最优性条件(KKT条件)的边界点。我会从基础概念开始,循序渐进地解释其核心思想、数学构造、算法步骤以及关键特性。

第一步:问题形式与核心思想

  1. 标准问题形式:考虑如下形式的非线性规划问题:

\[ \begin{aligned} \min_{x} \quad & f(x) \\ \text{s.t.} \quad & c_i(x) = 0, \quad i \in \mathcal{E}, \\ & c_i(x) \ge 0, \quad i \in \mathcal{I}, \end{aligned} \]

其中 \(f: \mathbb{R}^n \to \mathbb{R}\) 是目标函数,\(c_i(x)\) 是约束函数(等式约束集合 \(\mathcal{E}\) 和不等式约束集合 \(\mathcal{I}\) 通常假设是光滑的)。我们的目标是找到满足所有约束并使 \(f\) 最小的点 \(x^*\)

  1. 基本思路:原始-对偶内点法的核心是处理不等式约束 \(c_i(x) \ge 0\)。传统方法(如积极集法)在约束边界上工作。而内点法(Interior-Point Method, IPM)的关键是让迭代点始终严格满足不等式约束(即 \(c_i(x) > 0\)),保持在可行域“内部”,同时驱动迭代点趋向边界上的最优解。这需要一个机制来“阻挡”迭代点碰到边界,同时又能引导其逼近边界。

第二步:障碍函数与扰动KKT系统

  1. 障碍函数引入:为了处理不等式约束 \(c_i(x) \ge 0\),我们将其合并到目标函数中。定义一个对数障碍函数(Logarithmic Barrier Function):

\[ B(x; \mu) = f(x) - \mu \sum_{i \in \mathcal{I}} \ln(c_i(x)), \]

其中 \(\mu > 0\) 称为障碍参数(Barrier Parameter)。当 \(c_i(x) > 0\) 时,\(\ln(c_i(x))\) 是良定义的;当 \(c_i(x) \to 0^+\) 时,\(-\ln(c_i(x)) \to +\infty\),这就在边界上形成了一道“障碍”,阻止迭代点越界。因此,我们转而求解一系列依赖于 \(\mu\) 的无约束/等式约束子问题。

  1. 扰动KKT条件:原问题的KKT条件包括原始可行性对偶可行性互补松弛条件。对于对数障碍法,最优性条件对应的是障碍问题 \(\min_x B(x; \mu)\) 的驻点条件(在只考虑不等式约束时)。更一般地,对于包含等式约束的标准形式,引入松弛变量 \(s_i = c_i(x) \ge 0\) 和对偶变量(拉格朗日乘子) \(y_i\)(对应等式约束)和 \(z_i \ge 0\)(对应不等式约束),我们可以推导出扰动(或中心路径)KKT条件

\[ \begin{aligned} \nabla f(x) - A_E(x)^T y - A_I(x)^T z &= 0 \quad &\text{(梯度条件)} \\ c_i(x) &= 0, \quad i \in \mathcal{E} &\text{(等式约束)} \\ c_i(x) - s_i &= 0, \quad i \in \mathcal{I} &\text{(松弛变量定义)} \\ S Z e &= \mu e \quad &\text{(扰动互补松弛条件)} \\ s, z &\ge 0 \end{aligned} \]

其中:
  • \(A_E(x)\)\(A_I(x)\) 分别是等式和不等式约束的雅可比矩阵。
  • \(S = \text{diag}(s_i)\), \(Z = \text{diag}(z_i)\) 是以松弛变量和对偶变量为对角元素的对角矩阵。
  • \(e\) 是适当维度的全1向量。
  • 关键方程 \(S Z e = \mu e\) 是标准互补松弛条件 \(s_i z_i = 0\) 的“扰动”版本。当 \(\mu > 0\) 时,它强制 \(s_i > 0\)\(z_i > 0\),从而保证迭代点严格在内部(\(c_i(x) = s_i > 0\))且对偶变量为正,这正是“内点”的体现。当 \(\mu \to 0\) 时,这个条件恢复到精确的互补松弛。

第三步:中心路径与算法框架

  1. 中心路径:对于每一个 \(\mu > 0\),扰动KKT系统通常有唯一解(在一定的正则性条件下),记为 \((x(\mu), y(\mu), z(\mu), s(\mu))\)。当 \(\mu\) 从正数趋近于0时,这些解的轨迹称为中心路径(Central Path)。中心路径上的点有两个重要性质:
  • 内点性:所有点严格可行(\(s(\mu) > 0, z(\mu) > 0\))。
  • 导向最优:当 \(\mu \to 0\) 时,中心路径指向原问题的一个KKT点(通常是局部最优解)。
  1. 基本算法框架:原始-对偶内点法的策略不是精确求解每一个 \(\mu\) 对应的扰动KKT系统,而是:
    a. 外层迭代:生成一个单调递减趋于0的障碍参数序列 \(\{\mu_k\}\)
    b. 内层迭代:对于当前的 \(\mu_k\) 和当前迭代点 \((x_k, y_k, z_k, s_k)\),通过牛顿法(或拟牛顿法)求解扰动KKT系统得到一个搜索方向(牛顿方向),然后沿该方向进行线搜索,获得新的迭代点,使其更接近中心路径上对应于 \(\mu_k\) 的点。
    c. 参数更新:在迭代点足够接近当前中心路径后,减小 \(\mu_k\)\(\mu_{k+1}\),开始下一轮“追随”新的中心路径点。

第四步:牛顿方程与搜索方向

  1. 线性化扰动KKT系统:在迭代点 \((x, y, z, s)\) 处,将扰动KKT条件线性化(即应用牛顿法),得到关于搜索方向 \((\Delta x, \Delta y, \Delta z, \Delta s)\) 的线性方程组:

\[ \begin{bmatrix} H & -A_E^T & -A_I^T & 0 \\ A_E & 0 & 0 & 0 \\ A_I & 0 & 0 & -I \\ 0 & 0 & S & Z \end{bmatrix} \begin{bmatrix} \Delta x \\ \Delta y \\ \Delta z \\ \Delta s \end{bmatrix} = - \begin{bmatrix} \nabla f - A_E^T y - A_I^T z \\ c_E(x) \\ c_I(x) - s \\ S Z e - \sigma \mu e \end{bmatrix} \]

其中 \(H = \nabla_{xx}^2 L(x, y, z)\) 是拉格朗日函数 \(L = f(x) - y^T c_E(x) - z^T c_I(x)\) 的海森矩阵,\(\sigma \in (0, 1]\)中心参数(Centering Parameter)。右侧最后一项 \(S Z e - \sigma \mu e\) 是关键:

  • \(\sigma = 1\) 时,我们寻求直接指向当前 \(\mu\) 对应的中心路径点(称为中心化方向)。
  • \(\sigma = 0\) 时,我们寻求直接减小互补间隙 \(s^T z\)(即朝着 \(\mu=0\) 的方向,称为仿射尺度方向)。
  • 实际中常取 \(0 < \sigma < 1\),以获得一个平衡了“趋近中心路径”和“减小互补间隙”的搜索方向。
  1. 求解与简化:这个线性系统可以通过代数消元(利用第4个方程解出 \(\Delta s\),然后代入其他方程)简化成一个规模更小的、关于 \((\Delta x, \Delta y, \Delta z)\) 的对称不定系统(通常称为增广系统),再进一步简化成关于 \(\Delta x\) 的正定系统(通过舒尔补),从而可以利用高效的线性代数技术求解。

第五步:步长选择与算法步骤

  1. 步长选择:得到搜索方向后,需要选择步长 \(\alpha\) 来更新迭代点:\((x^+, y^+, z^+, s^+) = (x, y, z, s) + \alpha (\Delta x, \Delta y, \Delta z, \Delta s)\)。步长选择必须保证内点性(即 \(s^+ > 0, z^+ > 0\))不被破坏,通常采用分数到边界规则(Fraction to the Boundary Rule):

\[ \alpha_{\max} = \max \{ \alpha \in (0, 1] : s + \alpha \Delta s \ge (1-\tau)s, \; z + \alpha \Delta z \ge (1-\tau)z \}, \]

其中 \(\tau \in (0,1)\) 是一个安全因子(如0.995),确保新点不会太靠近边界。然后,实际步长 \(\alpha\) 可能在 \(\alpha_{\max}\) 基础上再进行线搜索,以充分减少某个价值函数(如原始-对偶残差的范数或障碍函数值)。

  1. 算法基本步骤
    a. 初始化:选择初始点 \((x_0, y_0, z_0, s_0)\) 满足 \(s_0 > 0, z_0 > 0\),初始障碍参数 \(\mu_0 > 0\),中心参数 \(\sigma_0\),容忍度 \(\epsilon > 0\)
    b. 主循环(对于 \(k=0,1,2,...\)):
  2. 终止检验:如果原始可行性、对偶可行性和互补间隙的某种度量(如 \(\|c_E(x_k)\|\), \(\|\nabla_x L\|\), \(s_k^T z_k\))都小于 \(\epsilon\),则停止并返回当前点作为近似解。
  3. 障碍参数更新:计算当前互补间隙 \(s_k^T z_k\),据此更新 \(\mu_k = \sigma_k \cdot (s_k^T z_k / m_I)\),其中 \(m_I\) 是不等式约束个数。\(\sigma_k\) 可以动态调整。
  4. 计算搜索方向:在当前点构建并求解上述牛顿方程,得到方向 \((\Delta x_k, \Delta y_k, \Delta z_k, \Delta s_k)\)
  5. 计算步长:使用分数到边界规则计算最大步长 \(\alpha_{\max}\),并可能通过线搜索确定最终步长 \(\alpha_k\)
  6. 更新迭代点\(x_{k+1}=x_k+\alpha_k \Delta x_k\),同理更新 \(y, z, s\)

第六步:算法特性与扩展

  1. 收敛性与复杂度:在适当的假设下(如函数二阶连续可微,约束满足某种约束规格,海森矩阵估计正定等),原始-对偶内点法具有全局收敛性,并且在局部具有超线性或二次收敛速率。对于凸规划问题,存在多项式时间复杂度的理论结果。

  2. 与线性/二次规划内点法的关系:非线性规划的原始-对偶内点法是线性规划和二次规划内点法的自然推广。在线性规划中,目标函数和约束都是线性的,海森矩阵 \(H\) 为零,算法大大简化,并催生了著名的多项式时间算法。非线性情况更复杂,但框架一致。

  3. 关键优势

    • 规模不敏感性:迭代次数对问题规模(变量和约束个数)相对不敏感,尤其适合大规模问题。
    • 边界处理能力强:能高效处理大量不等式约束,无需预先知道哪些约束是积极的。
    • 框架统一:提供了一个处理线性、二次、非线性凸甚至某些非凸问题的统一算法框架。
  4. 实际实现考虑

    • 海森矩阵计算:需要计算或近似拉格朗日函数的二阶导数。对于大规模问题,常使用拟牛顿法(如BFGS)来近似海森矩阵或其作用。
    • 线性系统求解:牛顿方程的求解是计算瓶颈。通常利用其特殊结构,采用稀疏矩阵技术、预条件共轭梯度法等。
  • 初始化:找到一个严格的初始内点 \((x_0, s_0>0, z_0>0)\) 本身可能是一个挑战,常采用同伦方法两阶段法

总结来说,非线性规划中的原始-对偶内点法通过引入对数障碍函数和扰动互补松弛条件,将不等式约束优化问题转化为一系列近似问题,并利用牛顿法在中心路径的邻域内迭代,最终收敛到最优解。它结合了原始变量和对偶变量的信息,是现代大规模非线性优化求解器的核心算法之一。

非线性规划中的原始-对偶内点法(Primal-Dual Interior-Point Method in Nonlinear Programming) 非线性规划中的原始-对偶内点法是一种用于求解具有不等式约束的非线性优化问题的强大算法。它通过结合原始变量、对偶变量和障碍参数,在可行域内部生成一系列迭代点,并最终收敛到满足一阶最优性条件(KKT条件)的边界点。我会从基础概念开始,循序渐进地解释其核心思想、数学构造、算法步骤以及关键特性。 第一步:问题形式与核心思想 标准问题形式 :考虑如下形式的非线性规划问题: \[ \begin{aligned} \min_ {x} \quad & f(x) \\ \text{s.t.} \quad & c_ i(x) = 0, \quad i \in \mathcal{E}, \\ & c_ i(x) \ge 0, \quad i \in \mathcal{I}, \end{aligned} \] 其中 \(f: \mathbb{R}^n \to \mathbb{R}\) 是目标函数,\(c_ i(x)\) 是约束函数(等式约束集合 \(\mathcal{E}\) 和不等式约束集合 \(\mathcal{I}\) 通常假设是光滑的)。我们的目标是找到满足所有约束并使 \(f\) 最小的点 \(x^* \)。 基本思路 :原始-对偶内点法的核心是处理不等式约束 \(c_ i(x) \ge 0\)。传统方法(如积极集法)在约束边界上工作。而内点法(Interior-Point Method, IPM)的关键是让迭代点始终 严格满足 不等式约束(即 \(c_ i(x) > 0\)),保持在可行域“内部”,同时驱动迭代点趋向边界上的最优解。这需要一个机制来“阻挡”迭代点碰到边界,同时又能引导其逼近边界。 第二步:障碍函数与扰动KKT系统 障碍函数引入 :为了处理不等式约束 \(c_ i(x) \ge 0\),我们将其合并到目标函数中。定义一个 对数障碍函数 (Logarithmic Barrier Function): \[ B(x; \mu) = f(x) - \mu \sum_ {i \in \mathcal{I}} \ln(c_ i(x)), \] 其中 \(\mu > 0\) 称为 障碍参数 (Barrier Parameter)。当 \(c_ i(x) > 0\) 时,\(\ln(c_ i(x))\) 是良定义的;当 \(c_ i(x) \to 0^+\) 时,\(-\ln(c_ i(x)) \to +\infty\),这就在边界上形成了一道“障碍”,阻止迭代点越界。因此,我们转而求解一系列依赖于 \(\mu\) 的无约束/等式约束子问题。 扰动KKT条件 :原问题的KKT条件包括 原始可行性 、 对偶可行性 和 互补松弛条件 。对于对数障碍法,最优性条件对应的是障碍问题 \(\min_ x B(x; \mu)\) 的驻点条件(在只考虑不等式约束时)。更一般地,对于包含等式约束的标准形式,引入松弛变量 \(s_ i = c_ i(x) \ge 0\) 和对偶变量(拉格朗日乘子) \(y_ i\)(对应等式约束)和 \(z_ i \ge 0\)(对应不等式约束),我们可以推导出 扰动(或中心路径)KKT条件 : \[ \begin{aligned} \nabla f(x) - A_ E(x)^T y - A_ I(x)^T z &= 0 \quad &\text{(梯度条件)} \\ c_ i(x) &= 0, \quad i \in \mathcal{E} &\text{(等式约束)} \\ c_ i(x) - s_ i &= 0, \quad i \in \mathcal{I} &\text{(松弛变量定义)} \\ S Z e &= \mu e \quad &\text{(扰动互补松弛条件)} \\ s, z &\ge 0 \end{aligned} \] 其中: \(A_ E(x)\) 和 \(A_ I(x)\) 分别是等式和不等式约束的雅可比矩阵。 \(S = \text{diag}(s_ i)\), \(Z = \text{diag}(z_ i)\) 是以松弛变量和对偶变量为对角元素的对角矩阵。 \(e\) 是适当维度的全1向量。 关键方程 \(S Z e = \mu e\) 是标准互补松弛条件 \(s_ i z_ i = 0\) 的“扰动”版本。当 \(\mu > 0\) 时,它强制 \(s_ i > 0\) 且 \(z_ i > 0\),从而保证迭代点严格在内部(\(c_ i(x) = s_ i > 0\))且对偶变量为正,这正是“内点”的体现。当 \(\mu \to 0\) 时,这个条件恢复到精确的互补松弛。 第三步:中心路径与算法框架 中心路径 :对于每一个 \(\mu > 0\),扰动KKT系统通常有唯一解(在一定的正则性条件下),记为 \((x(\mu), y(\mu), z(\mu), s(\mu))\)。当 \(\mu\) 从正数趋近于0时,这些解的轨迹称为 中心路径 (Central Path)。中心路径上的点有两个重要性质: 内点性 :所有点严格可行(\(s(\mu) > 0, z(\mu) > 0\))。 导向最优 :当 \(\mu \to 0\) 时,中心路径指向原问题的一个KKT点(通常是局部最优解)。 基本算法框架 :原始-对偶内点法的策略不是精确求解每一个 \(\mu\) 对应的扰动KKT系统,而是: a. 外层迭代 :生成一个单调递减趋于0的障碍参数序列 \(\{\mu_ k\}\)。 b. 内层迭代 :对于当前的 \(\mu_ k\) 和当前迭代点 \((x_ k, y_ k, z_ k, s_ k)\),通过牛顿法(或拟牛顿法)求解扰动KKT系统得到一个搜索方向(牛顿方向),然后沿该方向进行线搜索,获得新的迭代点,使其更接近中心路径上对应于 \(\mu_ k\) 的点。 c. 参数更新 :在迭代点足够接近当前中心路径后,减小 \(\mu_ k\) 到 \(\mu_ {k+1}\),开始下一轮“追随”新的中心路径点。 第四步:牛顿方程与搜索方向 线性化扰动KKT系统 :在迭代点 \((x, y, z, s)\) 处,将扰动KKT条件线性化(即应用牛顿法),得到关于搜索方向 \((\Delta x, \Delta y, \Delta z, \Delta s)\) 的线性方程组: \[ \begin{bmatrix} H & -A_ E^T & -A_ I^T & 0 \\ A_ E & 0 & 0 & 0 \\ A_ I & 0 & 0 & -I \\ 0 & 0 & S & Z \end{bmatrix} \begin{bmatrix} \Delta x \\ \Delta y \\ \Delta z \\ \Delta s \end{bmatrix} = - \begin{bmatrix} \nabla f - A_ E^T y - A_ I^T z \\ c_ E(x) \\ c_ I(x) - s \\ S Z e - \sigma \mu e \end{bmatrix} \] 其中 \(H = \nabla_ {xx}^2 L(x, y, z)\) 是拉格朗日函数 \(L = f(x) - y^T c_ E(x) - z^T c_ I(x)\) 的海森矩阵,\(\sigma \in (0, 1]\) 是 中心参数 (Centering Parameter)。右侧最后一项 \(S Z e - \sigma \mu e\) 是关键: 当 \(\sigma = 1\) 时,我们寻求直接指向当前 \(\mu\) 对应的中心路径点(称为 中心化方向 )。 当 \(\sigma = 0\) 时,我们寻求直接减小互补间隙 \(s^T z\)(即朝着 \(\mu=0\) 的方向,称为 仿射尺度方向 )。 实际中常取 \(0 < \sigma < 1\),以获得一个平衡了“趋近中心路径”和“减小互补间隙”的搜索方向。 求解与简化 :这个线性系统可以通过代数消元(利用第4个方程解出 \(\Delta s\),然后代入其他方程)简化成一个规模更小的、关于 \((\Delta x, \Delta y, \Delta z)\) 的对称不定系统(通常称为 增广系统 ),再进一步简化成关于 \(\Delta x\) 的正定系统(通过舒尔补),从而可以利用高效的线性代数技术求解。 第五步:步长选择与算法步骤 步长选择 :得到搜索方向后,需要选择步长 \(\alpha\) 来更新迭代点:\((x^+, y^+, z^+, s^+) = (x, y, z, s) + \alpha (\Delta x, \Delta y, \Delta z, \Delta s)\)。步长选择必须保证 内点性 (即 \(s^+ > 0, z^+ > 0\))不被破坏,通常采用 分数到边界规则 (Fraction to the Boundary Rule): \[ \alpha_ {\max} = \max \{ \alpha \in (0, 1 ] : s + \alpha \Delta s \ge (1-\tau)s, \; z + \alpha \Delta z \ge (1-\tau)z \}, \] 其中 \(\tau \in (0,1)\) 是一个安全因子(如0.995),确保新点不会太靠近边界。然后,实际步长 \(\alpha\) 可能在 \(\alpha_ {\max}\) 基础上再进行线搜索,以充分减少某个价值函数(如原始-对偶残差的范数或障碍函数值)。 算法基本步骤 : a. 初始化 :选择初始点 \((x_ 0, y_ 0, z_ 0, s_ 0)\) 满足 \(s_ 0 > 0, z_ 0 > 0\),初始障碍参数 \(\mu_ 0 > 0\),中心参数 \(\sigma_ 0\),容忍度 \(\epsilon > 0\)。 b. 主循环 (对于 \(k=0,1,2,...\)): 1. 终止检验 :如果原始可行性、对偶可行性和互补间隙的某种度量(如 \(\|c_ E(x_ k)\|\), \(\|\nabla_ x L\|\), \(s_ k^T z_ k\))都小于 \(\epsilon\),则停止并返回当前点作为近似解。 2. 障碍参数更新 :计算当前互补间隙 \(s_ k^T z_ k\),据此更新 \(\mu_ k = \sigma_ k \cdot (s_ k^T z_ k / m_ I)\),其中 \(m_ I\) 是不等式约束个数。\(\sigma_ k\) 可以动态调整。 3. 计算搜索方向 :在当前点构建并求解上述牛顿方程,得到方向 \((\Delta x_ k, \Delta y_ k, \Delta z_ k, \Delta s_ k)\)。 4. 计算步长 :使用分数到边界规则计算最大步长 \(\alpha_ {\max}\),并可能通过线搜索确定最终步长 \(\alpha_ k\)。 5. 更新迭代点 :\(x_ {k+1}=x_ k+\alpha_ k \Delta x_ k\),同理更新 \(y, z, s\)。 第六步:算法特性与扩展 收敛性与复杂度 :在适当的假设下(如函数二阶连续可微,约束满足某种约束规格,海森矩阵估计正定等),原始-对偶内点法具有全局收敛性,并且在局部具有超线性或二次收敛速率。对于凸规划问题,存在多项式时间复杂度的理论结果。 与线性/二次规划内点法的关系 :非线性规划的原始-对偶内点法是线性规划和二次规划内点法的自然推广。在线性规划中,目标函数和约束都是线性的,海森矩阵 \(H\) 为零,算法大大简化,并催生了著名的多项式时间算法。非线性情况更复杂,但框架一致。 关键优势 : 规模不敏感性 :迭代次数对问题规模(变量和约束个数)相对不敏感,尤其适合大规模问题。 边界处理能力强 :能高效处理大量不等式约束,无需预先知道哪些约束是积极的。 框架统一 :提供了一个处理线性、二次、非线性凸甚至某些非凸问题的统一算法框架。 实际实现考虑 : 海森矩阵计算 :需要计算或近似拉格朗日函数的二阶导数。对于大规模问题,常使用拟牛顿法(如BFGS)来近似海森矩阵或其作用。 线性系统求解 :牛顿方程的求解是计算瓶颈。通常利用其特殊结构,采用稀疏矩阵技术、预条件共轭梯度法等。 初始化 :找到一个严格的初始内点 \((x_ 0, s_ 0>0, z_ 0>0)\) 本身可能是一个挑战,常采用 同伦方法 或 两阶段法 。 总结来说, 非线性规划中的原始-对偶内点法 通过引入对数障碍函数和扰动互补松弛条件,将不等式约束优化问题转化为一系列近似问题,并利用牛顿法在中心路径的邻域内迭代,最终收敛到最优解。它结合了原始变量和对偶变量的信息,是现代大规模非线性优化求解器的核心算法之一。