非约束优化:非线性方程组求解方法

介绍

寻找局部极小值和求解非线性方程组之间有一些密切的联系. 给定一个有 个未知变量和 个方程的方程组,寻求 解 就等于极小化平方和 ,而此时在最小值处的残差为零,因此这和 高斯-牛顿 方法有很密切的联系. 事实上,用于局部极小化的高斯-牛顿步骤和用于非线性方程组的 牛顿法 步骤是完全相同的. 此外,对于一个平滑函数,关于局部极小化的牛顿法和关于非线性方程组 的牛顿法是一样的. 毫不奇怪,这些算法的许多方面是相似的;然而,它们也有许多重要的差别.

与极小化算法具有的另一个共同点是这里我们也需要某种形式的步控制. 通常情况下,步控制是基于和极小化相同的方法的,只是它是应用于一个优值函数的,而经常这个优值函数是光滑的 2-范数平方,即 .

"Newton"使用精确雅可比或者是有限差分近似在一个局部线性模型基础上求解步骤
"Secant"通过利用过去的 个步骤建立一个雅可比的割线近似工作,而不需要导数; 需要在每个维度上都有两个初始条件
"Brent"这是在一个保持对根的包围的维度方向上的方法; 需要两个初始条件来包围一个根
"AffineCovariantNewton"优化使用最少数量的雅可比计算

FindRoot 的基本方法选项.

牛顿法

非线性方程的牛顿法基于如下线性近似

因此,通过简单地设置 ,可以找到牛顿步,

在接近方程的根的位置,牛顿法具有 -二次收敛速度,这类似于用来极小化的牛顿方法. 牛顿法被用作 FindRoot 的默认方法.

牛顿法可以通过线搜索或者信赖域这两种步控制方法来使用. 工作的时候,线搜索控制通常更快,但是信赖域法往往更强稳健.

这里加载一个包含一些功用函数的程序包:
这是作为一个 FindRoot 问题的Rosenbrock问题:
这里通过默认的线搜索方法寻找非线性系统的解.(牛顿法是 FindRoot 的默认方法.)

请注意,每个线搜索都沿着 的线开始. 这是对于此特定问题,牛顿步的特殊性质.

这里使用符号化方法计算Rosenbrock问题的雅克比和牛顿步:

当对点 施加这个步骤时,就很容易看出为什么步骤跑到线 上. 这只是这一问题的一个特点,对大部分函数来说并不典型. 因为信赖域方法不会尝试牛顿步,除非它位于该域的边界以内,所以当我们使用信赖域步控制时,这个特点不会如此强烈地显示出来.

这里使用信赖域方法寻找非线性系统的解. 该搜索与在 FindMinimum 中关于Rosenbrock目标函数的高斯-牛顿方法几乎完全相同:

当雅可比矩阵的结构是稀疏的,Wolfram 语言将使用 SparseArray 对象来计算雅可比和处理必要的数值线性代数.

当把求解非线性方程组作为一个更一般的数值程序的一部分,比如使用隐式方法求解微分方程,初始值往往是相当不错的,而完全收敛不是绝对有必要的. 往往计算牛顿步代价最高的部分是求雅可比和计算矩阵分解. 然而,当足够接近根时,有可能让雅可比在几个步骤内不进行更新(尽管这样肯定会影响收敛速度). 您可以在 Wolfram 语言中使用方法选项 "UpdateJacobian" 来这样做,这给出在更新雅可比之前所需的步骤数目. 默认值为 "UpdateJacobian"->1,因此雅可比在每一步都进行更新.

这里显示了当每隔三步更新一次雅可比时,用来寻找一个简单的平方根所需的步数目,函数计算次数,和雅可比计算次数:
这里显示了当在每一步更新一次雅可比时用来寻找一个简单的平方根所需的步数目,函数计算次数,和雅可比计算次数:

当然对于一个简单的一维根,更新雅可比在成本上是微不足道的,因此,这里使用更新仅在于演示这一理念.

选项名称
默认值
"UpdateJacobian"1在更新雅可比之前,我们所采取的步数目
"StepControl""LineSearch"步控制方法,可以是 "LineSearch""TrustRegion" 或者 None(这个是不被推荐的)

FindRootMethod->"Newton" 的方法选项.

割线法

当不能使用符号化方法计算导数时,牛顿方法将被使用,但要利用对雅可比的有限差分近似. 这会在时间和可靠性上花费成本. 和极小化中的处理一样,我们有另一种专门为不使用导数而设计的算法.

在一维空间中,割线方法的思想是使用两个连续的搜索点之间线的斜率,来计算步骤,而不是使用最新点的导数. 同样地,在 个维度上,我们使用在 个点的残差的差值来建立雅可比排序的近似. 请注意,这类似于有限差分方法,但不是为了获得尽量好的雅可比近似而试图使差异间隔变小,它实际上就像一维割线方法一样,是用一个平均导数. 最初,这 个点从在所有 个维度上不同的两个初始点构造. 在随后的步骤中,只有具有最小的优值函数值的 个点被保留下来. 有种情况,虽然罕见,但是可能,即步骤是共线的,因而雅可比的割线近似变得奇异. 在这种情况下,算法使用不同的点重新启动.

该方法需要在每个维度上有两个初始点. 事实上,如果在每个维度上给出两个初始点,割线法就是默认的方法,除了在一维空间,这时候可能会选择布伦特方法.

这里加载一个包含一些功用函数的程序包:
这里显示了使用割线法,Rosenbrock 问题的解:

请注意,与牛顿方法相比,我们需要更多的残差函数计算. 然而,该方法能够不直接使用导数信息而沿着相对狭窄的谷部移动.     

这里显示了在使用有限差分计算雅可比的牛顿方法的情况下,Rosenbrock 问题的解:

然而,当与有限差分法牛顿法比较时,残差函数计算的次数就具有可比性. 对于较大问题的稀疏雅可比矩阵,有限差分牛顿法通常会更有效率,因为割线法不以任何方式利用稀疏的特点.

布伦特方法

当寻找实值函数的一个实数单根的时候,我们有可能利用问题的特殊几何结构,及函数由负向正,或由正向负穿过轴线. 布伦特方法 [Br02] 实际上是一个有保障的割线法,即它总是能保持一个点,在该点函数是正的,以及另一个点,在该点函数是负的,以便使得根总是被包围. 在任何特定的步骤,都要在插值(正割)步骤和二分法之间做一个选择,而选择的方式是要确保最终的收敛.

如果给 FindRoot 的两个实初始条件包围一个实函数的根,那么将使用布伦特方法. 因此,如果问题是一维的,并且您可以决定包围一个根的初始条件,这往往是一个不错的想法,因为布伦特方法是可用于 FindRoot 的最稳健的算法.

虽然基本上所有求解非线性方程组和局部极小值的理论都是基于光滑函数的,布伦特方法以其足够稳健,使您甚至可以对不连续函数的零点得到一个良好的估计.

这里加载一个包含一些功用函数的程序包:
这里显示在试图寻找一个不连续函数的根的过程中,所执行的步骤和函数求值:

当根被紧密地包围,但是该方法又不能够找到一个为零的函数值的时候,它会给出了一个提示信息. 这种稳健性在非常陡峭的连续函数上也得到了很好的保持.

这里显示在寻找一个在根附近迅速变化的函数的根的过程中,所用的步骤和函数求值:

仿射协变牛顿方法

这里给出的仿射协变牛顿法是 NLEQ_ERR 算法 [Deu06] 的实现,可作为 FindRoot 方法的选项. 为了说明该算法的工作原理及其优点,本文将给出仿射协变牛顿求解器的一些简单示例. 它的主要应用领域是求解离散非线性偏微分方程. 这些方程组可以变得很大,重要的是要注意计算速度和求解这些方程组所需的内存.

首先,指定一个需要求根的方程组. 考虑这个非线性方程组:

ex-2-y=0
y2-x=0

给定变量和初值,FindRoot 会求出方程组的根.

使用 FindRoot 找打零交叉:

为了使本教程更容易理解和更统一,我们设置一个接受向量参数并计算所有方程的向量值函数 f.

设置函数:
分组变量名和初始向量:
调用 FindRoot 并用值替代变量:

为了在下面进行比较,有一个很精确的参考解很有用.

使用更高的精度计算参考解:

下面是一个帮助函数,使得对于本教程中出现的根类型调用 FindRoot 更容易.

FindRoot 的封装函数:
求根:

牛顿方法,例如 FindRoot 的默认方法会需要计算函数 f,并计算和评估函数 f 的雅可比 (Jacobian). 雅可比的评估是相当计算昂贵的. 事实是,需要评估的雅可比不是本章讨论的牛顿方法类型中可以避免的;但是,有存在的算法可以减少雅可比的计算量. 为了研究 FindRoot 所需要的评估,可以使用诸如 EvaluationMonitorStepMonitor 的监测器.

在求根过程中,监测函数评估数,获取解的步数,以及雅可比评估数是很有意义的.

研究 FindRoot 中的各种评估数目:

求解需要 7 个步骤和 10 个函数评估. 有时会使用试验步骤,增加函数评估数,但是,试验步骤随后被拒. 被拒步骤不会增加步骤数,但是增加函数评估数. 你们还看到完成求解过程需要 7 个雅可比评估.

FindRoot 找到的默认解很接近参考解:

FindRoot 有多种方法选项,其中一个是 "AffineCovariantNewton". 仿射协变牛顿方法的详细解释请参见 [Deu06].

研究 FindRoot"AffineCovariantNewton" 方法中的各种评估数:

仿射协变牛顿方法需要更多的步骤和函数计算,但是更少的雅可比计算. 在求解非线性偏微分方程组中,评估雅可比的计算非常昂贵,不经常评估雅可比对于高效求解方程组会有很大的不同.

比较参考解:

注意,仿射协变解的精度比默认的 FindRoot 方法找到的小,但仍然在 PrecisionGoal 默认的设置范围内. 在查看之前,你们检查可以是 TrueFalse"AffineCovariantNewton" 方法选项 "BroydenUpdates". Broyden 更新是在不必重构雅可比下估计雅可比更新因子的相对廉价的方法. 仿射协变牛顿求解器估计何时使用 Broyden 更新以及如何高效计算. 如果在求解过程中非线性增加,Broyden 更新不再有用,算法会切换回计算雅可比. 在 QNERR 算法 [Deu06] 中实现 Broyden 更新算法.

求没有 Broyden 更新的系统的根:

关闭 Broyden 更新,求解会使用更少的步骤和函数评估,但是更多的雅可比评估.

比较参考解:

回到前面的问题,增加求解的精度,精度的增加可以通过增加精度目标实现.

指定 PrecisionGoal 减少与牛顿解的差别:
比较参考解:

仿射协变牛顿方法是一个误差范数-控制算法,因此,没有监视 AccuracyGoal. 换句话说,算法只监视误差不是残差. 误差与 PrecisionGoal 相关联,残差与精确目标相关联.

方法选项

仿射协变牛顿方法除了 FindRoot 选项还接受其他选项. 以后会讨论这些方法选项. 首先,给出可用的选项列表. 然后用例子演示这些选项,并给出有用选项的典型场景.

以下方法选项可以出现在 FindRoot 中的 Method"AffineCovariantNewton".

选项名
默认值
"BroydenUpdates"AutomaticBroyden 更新
"InitialDampingFactor"Automatic第一步的阻尼因子
"InverseJacobian"Automatic如果已知指定雅可比的逆
"LinearSolver"Automatic指定线性求解器函数或选项
"MinimalDampingFactor"Automatic阻尼因子可以变得多小

FindRootMethod"AffineCovariantNewton" 的方法选项.

"BroydenUpdates"

Broyden 更新是重用先前计算的雅可比的一种方法. 重用先前计算的雅可比要比重新计算新的雅可比有效得多. 缺点是,Broyden 更新可能不如新计算的雅可比精确. 仿射协变牛顿法算法实现了一种机制,在安全时切换到 Broyden 更新,在 Broyden 更新不再可用时也可以切换回计算雅可比. 不过,可以禁用 Broyden 更新.

求没有 Broyden 更新的根:
"MinimalDampingFactor"

在求解过程的第 步,该方程被解:. 雅可比 被反转获取新的增量 . 该增量然后被加入到之前的解以获取一个更新的解 . 是介于 之间的阻尼因子,它限制步长的大小. 对于极端非线性,可以减少默认最小阻尼因子 . 默认值是 . 通常,这是不必要的.

"InitialDampingFactor"

当知道非线性不是太极端,可以选择更大的初始阻尼因子. 默认的初始阻尼因子是 1/100. 选择的值可以比 "MinimalDampingFactor" 更大,但是小于等于 1:.

指定一个更大的 "InitialDampingFactor"

更大的初始阻尼因子可以减少雅可比评估.

"InverseJacobian"

牛顿算法实际上需要反转雅可比来完成. 该反转通过创建 LinearSolveFunction 来计算. 在某些情况下,逆转雅可比或许已知或易于计算. 在这种情况下,反转雅可比可以通过方法选项传递给 "AffineCovariantNewton" 方法.

构建一个函数在评估点 ep 返回已评估的反转雅可比:
为问题函数 构建雅可比反转:
求带有指定雅可比逆的根:

指定雅可比的逆不会影响雅可比评估的数目. 事实是,雅可比评估数显示为 0 是雅可比逆与 FindRoot 交互的一个弱点.

"LinearSolver"

可以通过使用 LinearSolve 改变线性求解器或选项来构建雅可比的逆. 这可以使用 "LinearSolver" 方法选项实现.

LinearSolve 指定迭代的 "Krylov" 方法:

当指定像 "Krylov" 这样的迭代方法作为 LinearSolve 方法时,然后雅可比因式分解在可能的情况下不会被存储,而是在每次需要的时候被重新计算.

可以重新定义线性求解器. 这允许你使用自定义的线性求解器. 重定义线性求解器允许你提取 Jacobian.

重定义线性求解器,将 Sow 应用于 Jacobian:
Reap 应用于 Jacobians: