求解张量绝对值方程的两步非精确Levenberg-Marquardt算法

2023-06-05马昌凤谢亚君

山西大学学报（自然科学版） 2023年3期

关键词：收敛性张量步长

马昌凤，谢亚君

（福州外语外贸学院大数据学院数字技术与智能计算重点实验室，福建福州 350202）

0 引言

考虑如下的张量绝对方程（TAVE）：寻找向量x∈Rn满足

其中A∈T(m，n)且m为偶数，b∈Rn为已知向量。这里T(m，n)表示m阶n维实张量的集合，向量|x|[m-1]定义为

当m=2 时，方程（1）退化为下面的（向量）绝对值方程（AVE）：

绝对值方程在物理和经济均衡问题等应用科学技术中有着广泛的应用［1］。关于方程（2）的数值算法，已经存在大量研究工作，可参见文献［2-8］。

本文将证明张量绝对值方程（1）等价于广义张量互补问题，然后通过引入互补函数将张量互补问题转化为等价的非线性方程组。我们将在有关研究工作的基础上提出两步非精确Levenberg-Marquardt 算法（简记为LM 算法）来求解该非线性方程组。在适当的假设条件下证明该方法的全局收敛性，并讨论它的收敛率，利用数值实验验证所做的理论分析。

1 两步非精确LM算法

本节，我们首先证明TAVE（1）等价于广义张量互补问题。首先介绍如下定义。

定义1一个m阶n维张量A称为对称张量，若成立

式中Πm表示m个指标i1i2…im的置换群。

定义2设A是m阶n维实张量，x，b∈Rn。定义

这里，I是m阶n维单位张量，广义张量互补问题就是寻找向量x∈Rn满足

定理1设A∈T(m，n)，b∈Rn，其中m为偶数，那么TAVE（1）等价于广义张量互补问题（3）。

证明事实上，由于

因此，由（1）可得：

这就意味着x是（3）的可行解。因为

这就完成了定理证明。

利用Fischer-Burmeister 函数φFB，可以将（3）再生为如下方程：

这里，x是（3）的解当且仅当H(x)=0。进一步，由于H(x)各个分量函数是强半光滑的，所以H(x)是强半光滑的（参见文献［9］）。

定理2设A是m阶n维对称张量，则函数H(x) 是强半光滑的。进一步对任意的x∈Rn，有

其中Da(x)=diag(ai(x))，Db(x)=diag(bi(x))是n×n对角矩阵，对应元素

这里∂φFB(Fi(x)，Gi(x))是∂φFB(a，b)中的(a，b)由(Fi(x)，Gi(x))替换，并且JF(x)和JG(x)由下式给出：

为了给出求解H(x)=0 的算法，定义如下的价值函数：

接下来给出价值函数的一个性质：

性质1设A是m阶n维对称张量，则价值函数Ψ(x) 是连续可微的并且对任意的Q∈∂H(x)，有

最近，许多专家学者在LM 算法的试探步和LM 参数的选择上进行研究。例如文献［10］中提出了一种修正LM 方法，使用了如下的LM近似步：

其中Fk=F(xk)，。在局部误差界条件下，可证明其具有三次收敛性。

Argyros 和Chen 提出了Chebyshev 方法（文献［11］），在每步迭代计算：

其中p∈(0，1]是一个常数，那么可得到

结合（4）可导出如下修正Chebyshev 方法：

注意到，当Jacobi 矩阵是奇异或接近奇异，（4）、（6）中的dk和k不能很好地定义。为了克服这个困难，Levenberg-Marquardt 算法计算如下的步长：

其中LM 参数λk是非负的，且当接近奇异时能阻止步长过大。目前已经有很多学者研究参数λk，并分析算法的收敛性。例如，Yamashita和Fukushima 选择λk=||Fk||2，在局部误差有界的条件下证明了其二次收敛性［13］。Fan 和Yuan选择λk=||Fk||δ，δ∈[1，2]，证明了算法仍保留二次收敛性，并且在λk=||Fk|| 时实验效果最好［14］。但是，当迭代序列距离方程的解较远时，该算法的迭代结果并不令人满意。为避免该缺陷，文献［15］进一步提出λk=μk||Fk||，μk＞0，利用信赖域技术证明了它的收敛性，数值结果表明其具有较好的收敛性质。但当迭代序列远离方程的解集时，λk=μk||Fk||的值可能会比较大，这会导致LM 尝试步长dk非常小，以致算法的效率可能很低。基于上述考虑，文献［16］考虑了新的参数：

这种选择方式避免了上面提到的不足。目前，LM 方法及其变形修正在求解非线性方程、非线性互补问题及无约束优化和约束优化问题等方面均有广泛应用。

受Chebyshev 方法的启发，本文提出两步LM 方法用于求解张量绝对值方程，在每步计算

来获得步长dk，μk＞0 随着步长的更新而更新。不难发现，步长dk也是如下凸极小化问题的极小解：

则易证明dk是如下信赖域子问题的解

考虑如下两步迭代

其中yk=xk+dk。为了提高数值性能，用Qk代替Q(yk)来避免计算新的Jacobi 矩阵Q(yk)，并提出Hk和H(zk)的凸组合代替H(yk)，即

其中θ=1/α且α≥1。因此可得到

预估下降量要求非负，根据文献［17］可得

结合（10）和（11），定义如下两步预下降

基于以上分析，提出以下两步非精确LM算法。

算法1（两步非精确LM 算法）

步0 给定初始点x0，令μ0＞m＞0，0 ＜γ0≤γ1＜γ2＜1，α＞1，ε＞0，置k≔0。

步1 若||H(xk)||≤ε，停算。否则，计算Qk∈∂H(xk)。

步2 计算LM 参数

步3 非精确求解方程组

得到dk，其中k为非精确求解的控制阈值。置zk=xk+αdk。

步4 计算

步5 非精确求解方程组

步6 计算下降比

其中Aredk和Predk分别如（9）和（12）所定义。

步7 选择μk+1为

步8 置k≔k+1，转步1。

下面分析算法1 的全局收敛性，假设算法1生成无限序列｛xk｝。类似于文献［18］定理2.1和定理2.2 的证明，可以得到如下的两个定理。

定理3设{xk}是由算法1 生成的序列，如果残差向量k，k满足如下条件：

其中η∈(0，1)，νk=o(dist(xk，X*))，其中X*是信赖域问题的解集，δ∈[ 1，2 )，那么对任意{xk}的聚点是价值函数Ψ 的稳定点。进一步，如果{xk}的聚点x*是H(x)=0 的解，则{dist(xk，X*)}超线性收敛到0。

定理4设{xk}是由算法1 生成的序列，δ=2。若残差向量k，满足如下条件：

其中η∈(0，1)，νk=O(dist(xk，X*)2)，则{xk}的任意聚点都是价值函数Ψ 的稳定点。进一步，如果{xk} 的聚点x*是H(x)=0 的解，则{dist(xk，X*)}二次收敛到0。

下面再对算法1 做一些说明。在算法1 的执行过程中，主要计算（13）、（14）当k=k=0 时的近似值。注意到方程总是可解的。事实上，由于λk＞0，那么矩阵是对称正定矩阵，因此（13）和（14）一定可解。现在考虑如何计算Qk∈∂H(xk)，我们选择在第k步迭代，通过定理2，矩阵Qk∈∂H(xk)的元素可以通过如下公式计算。令

是一个退化指数集合，定义z∈Rn是一个向量，若i∈Λ，其分量zi=1；否则zi=0，那么矩阵Qk可以如下定义

其中A(xk)和B(xk)是n×n对角矩阵，其第i个对角元素可分别如下给出

后面实验部分计算Qk用这个公式。

2 数值实验

本节，将用一个数值算例来说明用算法1求解张量绝对值方程（TAVE）的可行性及有效性。数值实验在MATLAB R2018b 软件上运行，PC 机的运行环境为2.40 GHz 中央处理器，8.0 GB 内存以及Windows 10 操作系统。

在整个计算实验中，算法1 使用的参数，取ε=10-4，γ0=10-5，γ1=0.25，γ2=0.75，m=10-5，μ0=10.0。终止条件为||H(xk)||≤ε。

例1首先随机生成一个四阶四维的张量，其元素属于[0，1]，然后将其对称化得到非负张量B。令

其中e=(1，1，1，1)T，a∈(0，+∞)。因为，这样c的选择确保c＞ρ(B)+1，那么令A=cI-B满足A-I是强M 张量。张量B和张量A见表1 和表2。

表1 随机对称非负张量B=(bi1 i2 i3 i4)∈S(4，4)Table 1 Random symmetric non-negative tensors B=(bi1 i2 i3 i4)∈S(4,4)

表2 基于张量B的对称张量A=(ai1 i2 i3 i4)∈S(4，4)Table 2 Symmetric tensors A=(ai1 i2 i3 i4)∈S(4，4) based on tensor B

这个例子主要用于观察算法1 的迭代过程。随机生成对称张量B∈S[4，4] 和随机向量x*∈R4且张量B和向量x*元素取值均在[0，1]。为了确保方程只有唯一的解，计算b=A(x*)3-|x*|[3]，这里取a=3，从而计算出c=4.899 8。算法1 的迭代过程参见表3。从表3 可以看到||H(xk)||随着迭代次数k的增加会快速趋于0。此外，||∇Ψ(xk)||随着迭代次数k的增加亦会快速趋于0，这说明了算法具有良好的收敛性。