考虑如下形式绝对值方程

$ \begin{equation}\label{eq1} Ax-|Bx-c|=b, \end{equation} $

(1.1)

其中 $A, B\in R^{n\times n}$, $b, c\in R^n$, 绝对值 $|x|$依分量形式选取.

绝对值方程首次由Rohn于1989年在文献[1]中提出, 此后许多学者对其理论与算法进行了广泛的研究.目前国内外研究的主要绝对值方程形式有以下两种

$ \begin{align} Ax-|x|&=b, \end{align} $

(1.2)

$\begin{align} Ax+B|x|&=b. \end{align} $

(1.3)

对于绝对值方程(1.2), 文献[2]指出方程(1.2) 等价于双线性规划、广义线性互补问题, 证明了当1不是矩阵 $A$的特征值时, 方程(1.2) 可以转化为线性互补问题, 给出了方程(1.2) 有解、无解、唯一解、非负解及个解的充分条件; 文献[3]中给出了在无任何假设条件下将方程(1.2) 转化为线性互补问题的一种方法, 研究了方程(1.2) 的解的凸性; 文献[4]给出了一种具有线性收敛速度的广义牛顿法; 文献[5]给出了求解绝对值方程(1.2) 的一个光滑牛顿法, 并证明了该方法具有二次收敛性; 中国矿业大学王海军研究小组在文献[6-8]结合区间数学方法, 给出了求解绝对值方程(1.2) 的区间算法; 文献[9-11]基于极大熵函数光滑化处理, 给出了求解绝对值方程(1.2) 的极大熵自适应微粒群混合算法、牛顿法及和声搜索算法; 当绝对值方程(1.2) 中矩阵是对称矩阵时, 文献[12, 13]结合优化技术给出了一类新的迭代算法, 并对算法收敛性及收敛速度进行了分析比较.

对于绝对值方程(1.3), 文献[14]给出了方程(1.3) 的择一定理及几个等价形式; 文献[15]证明了该形式绝对值方程的求解是NP-hard问题; 文献[16]给出了该形式方程唯一可解性的条件, 即 $A$的最小奇异值大于 $B$的最大奇异值; 文献[17]利用广义Jacobian矩阵及广义牛顿法对绝对值方程(1.3) 进行求解, 并且证明了该方程与二阶锥互补问题的等价性, 进而延伸到求解二阶锥互补问题; 文献[18]利用极大熵函数直接对 $|x|$进行光滑化处理, 建立求解绝对值方程(1.3) 的光滑牛顿算法.

此外, 文献[19, 20]对广义绝对值方程(1.1) 进行了较为详细的理论研究, 给出了相关的等价性理论.目前国内外学者主要对绝对值方程(1.2) 和(1.3) 进行了大量理论和算法方面的研究, 对广义绝对值方程(1.1) 的研究则很少, 且目前绝对值方程的现有算法多是基于半光滑或光滑化处理技术下的有效算法及逐次线性化方法, 然而很少有文章通过将绝对值方程等价转化为其他相关问题进行求解.基于此, 本文首先利用等价性转化将广义绝对值方程(1.1) 转化为互补问题, 构造光滑函数, 然后利用牛顿法进行求解.

引理1^[20]  广义绝对值方程(1.1) 等价于以下垂直线性互补问题

$ \begin{equation}\label{eq4} 0\leq(Mx+q)\perp(Nx+p)\geq 0, \end{equation} $

(2.1)

其中 $M=\frac{1}{2}(A+B), N=\frac{1}{2}(A-B), q=-\frac{1}{2}(b+c), p=-\frac{1}{2}(b-c)$.

记 $\sigma_{\min}(A)$和 $\sigma_{\max}(B)$分别为矩阵 $A$的最小奇异值和矩阵 $B$的最大奇异值, 下面给出广义绝对值方程(1.1) 的唯一可解性定理及推论.

引理2^[20]  若 $\sigma_{\max}(B)<\sigma_{\min}(A)$, 则对任意给定的向量 $b, c\in R^n$, 广义绝对值方程(1.1) 有唯一解.

推论1^[20]  若 $\sigma_{\max}(B)<\sigma_{\min}(A)$, 则由广义绝对值方程(1.1) 中系数矩阵 $A$、 $B$构成的矩阵 $(A+B)(A-B)^{-1}$为正定矩阵.

接下来, 将垂直线性互补问题(2.1) 转化为方程组进行求解.

定义1^[21]  函数 $\Psi:R^{2}\rightarrow R^1$被称为“NCP(nonlinear complementarity problem)函数”, 如果对任意 $(a, b)^T\in R^2$, $\Psi(a, b)=0$当且仅当 $a\geq 0$, $b\geq 0$, $ab=0$.

从事优化理论与算法研究的学者提出了许多不同的NCP函数, 其中较为常用的有以下两个NCP函数:

(1) Fischer-Burmeister函数:

$ \begin{equation}\label{eq5} \Psi_{FB}(a, b)=\sqrt{a^2+b^2}-(a+b), \forall(a, b)^T\in R^2. \end{equation} $

(2.2)

(2) min函数:

$ \begin{equation}\label{eq6} \Psi_{\min}(a, b)=\frac{1}{2}(a+b-\sqrt{(a-b)^2}), \forall(a, b)^T\in R^2. \end{equation} $

(2.3)

由引理1^[20], 定义1^[21]及文献[21]相关理论知识可知, 垂直线性互补问题(2.1) 可以等价地转化为以下方程组

$ \begin{equation}\label{eq7} \Phi(x)=\left( \begin{array}{c} \Psi(H_{1}, G_{1})\\ \vdots\\ \Psi(H_{n}, G_{n}) \\ \end{array}\right)=0, \end{equation} $

(2.4)

其中 $H(x)=Mx+q=\frac{1}{2}(A+B)x-\frac{1}{2}(b+c)$, $G(x)=Nx+p=\frac{1}{2}(A-B)x-\frac{1}{2}(b-c)$, $H_i$、 $G_i (i=1, 2, \cdots, n )$分别代表 $H(x)$及 $G(x)$的第 $i$个分式.

本文中利用Fischer-Burmeister函数代替(2.4) 式中 $\Psi$函数, 则(2.4) 式可转化为

$ \begin{equation}\label{eq8} \Phi(x)=\left( \begin{array}{c} \Psi_{FB}(H_{1}, G_{1})\\ \vdots\\ \Psi_{FB}(H_{n}, G_{n}) \\ \end{array}\right)=0. \end{equation} $

(2.5)

由以上等价性转化易得以下定理:

定理1  非线性方程组(2.5) 的解即为广义绝对值方程(1.1) 的解.

易知Fischer-Burmeister函数并不是处处可微的, 为构造光滑函数方程给出如下定义:

定义2^[21]  给定函数 $\Psi:R^n\rightarrow R^n$, 称光滑函数 $\Psi_\mu:R^n\rightarrow R^n(\mu>0)$为 $\Psi$的光滑逼近函数, 如果对任意的 $x\in R^n$, 存在 $\kappa>0$, 使得

$ \begin{equation}\label{eq9} \|\Psi_\mu(x)-\Psi(x)\|\leq\kappa\mu, \forall\mu>0, \end{equation} $

(2.6)

如果 $\kappa$不依赖于 $x$, 则称 $\Psi_\mu$为 $\Psi$的一致光滑逼近函数.

参考文献[21]引进Fischer-Burmeister函数的光滑函数:

$ \begin{equation}\label{eq10} \Psi_{FB}(\mu, a, b)=\sqrt{a^2+b^2+\mu^2}-(a+b), \mu>0. \end{equation} $

(2.7)

利用光滑函数(2.7) 代替方程组(2.5) 中Fischer-Burmeister函数得到以下光滑方程组

$ \begin{equation}\label{eq11} \Phi_\mu(x)=\left( \begin{array}{c} \Psi_{FB}(\mu, H_{1}, G_{1})\\ \vdots\\ \Psi_{FB}(\mu, H_{n}, G_{n}) \\ \end{array}\right)=0. \end{equation} $

(2.8)

引理3^[21]  对于任意的 $\mu_1\geq0$, $\mu_2\geq0$, $(a, b)\in R^2$, 有以下不等式存在

$ \begin{equation}\label{eq12} |\Psi_{FB}(\mu_1, a, b)-\Psi_{FB}(\mu_2, a, b)|\leq|\mu_1-\mu_2|. \end{equation} $

(2.9)

由引理3^[21]及方程组(2.5) 及(2.8) 可以得到如下结论:

定理2  光滑方程组 $\Phi_\mu(x)$为方程组 $\Phi(x)$的一致光滑逼近函数.

证  由引理3^[21], 定义2^[21]知

$ \|\Phi_\mu(x)-\Phi(x)\|\leq\sqrt{n}\max|\Phi(\mu, H_i, G_i)-\Phi(H_i, G_i)|\leq\sqrt{n}\mu~(1\leq i\leq n), $

该式的右项与 $x$无关, 由定义2^[21]知方程组 $\Phi_\mu(x)$为方程组 $\Phi(x)$的一致光滑逼近函数.

所以可以通过求光滑方程组(2.8) 的解来构造算法近似求解方程组(2.5) 的解即广义绝对值方程(1.1) 的解.

为保证广义绝对值方程(1.1) 的唯一可解性, 本文中总是假设矩阵 $A$及矩阵 $B$的奇异值满足 $\sigma_{\max}(B)<\sigma_{\min}(A)$这一前提.接下来给出算法的具体构造及其收敛性分析.

算法3.1

1) 给定一个初始向量 $x\in R^n$, 容许误差 $\varepsilon=10^{-9}$及参数 $\mu>0$, 构造光滑方程组 $\Phi_\mu(x)$, 置 $k=0$.

2) 计算 $x^{k+1}$, 通过 $ \left\{ \begin{aligned} &x^k=x^k+\Delta x^k, \\ &\Phi_\mu'(x^k)\Delta x^k+\Phi_\mu(x^k)=0, \end{aligned} \right. $ $k=0, 1, 2, \cdots$.

3) 若 $\|x^{k+1}-x^{k}\|\leq\varepsilon$, 则停止, 得到解 $x^*=x^{k}$; 否则转入步骤2.

为分析算法2.1的收敛性, 给出以下相关知识:

定义3  矩阵 $W\in R^{n\times n}$称为 $P_0$矩阵, 如果对任意向量 $x\in R^n$, $x\neq0$, 存在一个分量 $x_i\neq0$, 使得 $x_i(Wx)_i\geq0$; 矩阵 $W\in R^{n\times n}$称为 $P$矩阵, 如果对任意向量 $x\in R^n$, $x\neq0$, 存在一个分量 $x_i\neq0$, 使得 $x_i(Wx)_i>0$.

定义4^[22]  矩阵 $W\in R^{n\times n}$称为广义正定矩阵, 如果对任意向量 $x\in R^n$, $x\neq0$, 都有 $x^TWx>0$.

根据正定矩阵定义, 半正定矩阵定义, 定义3及定义4^[22]易得以下推论.

推论2  正定矩阵(广义正定矩阵)必是 $P$矩阵, 半正定矩阵必是 $P_0$矩阵.

定理3  若 $\sigma_{\max}(B)<\sigma_{\min}(A)$, 则式(2.8) 中 $\Phi_\mu(x)$的Jacobian矩阵 $\Phi_\mu'(x)$对任意向量 $x\in R^n$是非奇异的.

证  由引理2^[20]及推论1^[20]可知, 当 $\sigma_{\max}(B)<\sigma_{\min}(A)$时 $(A-B)^{-1}$存在且 $(A+B)(A-B)^{-1}$为正定矩阵, 令 $z=G(x)=Nx+p=\frac{1}{2}(A-B)x-\frac{1}{2}(b-c)$, 则

$ x=2(A-B)^{-1}z+(A-B)^{-1}(b-c), $

将其代人 $H(x)=Mx+q=\frac{1}{2}(A+B)x-\frac{1}{2}(b+c)$, 得

$ \tilde{H}(z)=(A+B)(A-B)^{-1}z+\frac{1}{2}(A+B)(A-B)^{-1}(b-c)-\frac{1}{2}(b+c), $

则 $H(x)\perp G(x)\Leftrightarrow z\perp\tilde{H}(z)$, 将 $z$及 $\tilde{H}(z)$代人式(2.8) 得

$ \Phi_\mu(z)=\left( \begin{array}{c} \Psi_{FB}(\mu, z_{1}, \tilde{H}_{1}(z))\\ \vdots\\ \Psi_{FB}(\mu, z_{n}, \tilde{H}_{n}(z)) \\ \end{array}\right)=0, $

其中 $z_i$, $\tilde{H}_{i}(z)$分别代表向量 $z$及方程组 $\tilde{H}(z)$的分量.

所以 $\Phi_\mu'(z)=D_1(z)+D_2(z)\tilde{H}'(z)$, 其中

$ D_1(z)={\hbox{diag}}\{a_1(z), a_2(z), \cdots, a_n(z)\}, ~~ D_2(z)={\hbox{diag}}\{b_1(z), b_2(z), \cdots, b_n(z)\}, $

$\tilde{H}'(z)$为 $\tilde{H}(z)$的Jacobian矩阵

$ \tilde{H}'(z)=(A+B)(A-B)^{-1}, a_i(z)=\frac{z_i}{\sqrt{z_i^2+\tilde{H}_{i}(z)^2+\mu^2}}-1, b_i(z)=\frac{\tilde{H}_{i}(z)}{\sqrt{z_i^2+\tilde{H}_{i}(z)^2+\mu^2}}-1. $

因为 $\Psi_{FB}(\mu, a, b)=\sqrt{a^2+b^2+\mu^2}-(a+b)$, 所以

$ \Psi_a'=\frac{a}{\sqrt{a^2+b^2+\mu^2}}-1, \Psi_b'=\frac{b}{\sqrt{a^2+b^2+\mu^2}}-1, $

对于任意的 $\mu>0$, $\Psi_a'$, $\Psi_b'$是连续的且 $-2 < \Psi_a' < 0$, $-2 < \Psi_b' < 0$.所以 $-2 < a_i(z) < 0$, $-2 < b_i(z) < 0$.所以, 对角矩阵 $-D_1(z)$及 $-D_2(z)$的对角元素恒为正数, 又

$ \tilde{H}'(z)=(A+B)(A-B)^{-1} $

为正定矩阵, 由推论2可知矩阵 $\tilde{H}'(z)=(A+B)(A-B)^{-1}$为 $P$矩阵, 由正定矩阵(或 $P$矩阵)的性质及文献[23]中引理3可知 $-D_2(z)\tilde{H}'(z)$是正定的, 从而 $-D_1(z)-D_2(z)\tilde{H}'(z)$也是正定的, 所以 $-D_1(z)-D_2(z)\tilde{H}'(z)$是非奇异的, 所以

$ \Phi_\mu'(z)=D_1(z)+D_2(z)\tilde{H}'(z)=-[-D_1(z)-D_2(z)\tilde{H}'(z)] $

也是非奇异的, 进而可知 $\Phi_\mu'(x)$也是非奇异的, 证毕.

引理4^[24]  (Ostrowski)设映象 $G:D\subset R^n\rightarrow R^n$有一不动点 $x^*\in int(D)$, 且在 $x^*$处为F可导, $Q'(x^*)$的谱半径 $\rho(Q'(x^*))<1$, 则存在开球 $S=S(x^*, \delta)\subset D$, 对任意初值 $x^0\subset S$, $x^*$是迭代序列 $x^{k+1}=Q(x^k)$, $k=0, 1, \cdots$的一个吸引点.

定理4  若 $\sigma_{\max}(B)<\sigma_{\min}(A)$则算法3.1是适定的, 设 $x^*$是方程组(2.8) 的解, 则算法3.1产生的序列 $\{x^k\}$超线性收敛到 $x^*$.

证  如果 $\sigma_{\max}(B)<\sigma_{\min}(A)$则由定理3可知 $\Phi_\mu'(x)$是非奇异的, 则 $\Phi_\mu'(x)$可逆, 所以算法3.1是适定的; 在引理4^[24]中取 $Q(x^k)=x^k-[\Phi_\mu'(x^k)]^{-1}\Phi_\mu(x^k)$, 则

$ Q'(x^*)=I-[\Phi_\mu'(x^*)]^{-1}\Phi_\mu(x^*)=0, $

这里 $I$为单位矩阵, 所以 $\rho(Q'(x^*))<1$, 由引理4^[24]知算法3.1产生的序列 $\{x^k\}$收敛到 $x^*$, 且当 $x^k\neq x^*$时有 $\lim\limits_{k\rightarrow\infty}\frac{\|x^{k+1}-x^*\|}{\|x^k-x^*\|}=0$, 即算法3.1产生的序列 $\{x^k\}$是超线性收敛的, 证毕.

运用MATLAB 7.0进行编程计算, 参数的设定和符号说明如下:

计算精度 $\varepsilon=1.0e-9$; $k$表示迭代次数; $T$表示算法运行的CPU时间(计算单位为秒); $x^*$表示广义绝对值方程的近似解; $n$表示维数; 初始点 $x^0$从区间[-10,10]中随机选取 $x^0={\hbox{randint}}(n, 1, [-10,10])$.

算例1  随机选取矩阵 $B=\left( \begin{array}{ccc} 5&0&0\\ -3&4&-5\\ 1&3&4\\ \end{array}\right)$, 为保证方程有唯一解(即:矩阵 $A$的最小奇异值大于矩阵 $B$的最大奇异值), 构造矩阵 $A=\left( \begin{array}{ccc} 9.1229&1.6002&-2.8123\\ -1.3812&6.1472&-5.6905\\ 0.5348&6.5201&6.0355\\ \end{array}\right).$随机选取 $c=\left( \begin{array}{ccc} -1\\ 1\\ 3\\ \end{array}\right).$令 $b=A\times\left(\begin{array}{ccc} 1\\ 1\\ 1\\ \end{array}\right)-\left|B\times\left(\begin{array}{ccc} 1\\ 1\\ 1\\ \end{array}\right)-c\right|=\left(\begin{array}{ccc} 1.9108\\ -5.9246\\ 8.0905\\ \end{array}\right), $则 $(1, 1, 1)^{T}$即为广义绝对值方程 $Ax-|Bx-c|=b$的准确解.随机选取初始点进行迭代, 计算结果如表 1.

算例2  随机选取矩阵 $B=\left( \begin{array}{ccccc} -4&-3&0&2&4\\ -3&-5&-1&4&0\\ -3&3&4&-5&2\\ 1&-1&0&2&-1\\ -3&5&-3&-1&-2\\ \end{array}\right)$, 为保证方程有唯一解(即矩阵 $A$的最小奇异值大于矩阵 $B$的最大奇异值), 构造矩阵

$ A=\left( \begin{array}{ccccc} -4.6101&0.1977&-1.3355&4.6007&10.9050\\ -7.4993&-8.1556&0.0062&1.5591&-6.2505\\ -5.9196&2.0731&9.1091&-5.1312&0.1729\\ 0.9829&4.645&7.8260&11.1556&-2.6493\\ -6.8459&7.9552&-6.1789&1.4257&-3.4887\\ \end{array}\right). $

随机选取 $c=\left( \begin{array}{ccccc} -3\\ -3\\ 2\\ -2\\ 0\\ \end{array}\right)$, 令 $b=A\times\left(\begin{array}{ccccc} 1\\ 1\\ 1\\ 1\\ 1\\ \end{array}\right)-\left|B\times\left(\begin{array}{ccccc} 1\\ 1\\ 1\\ 1\\ 1\\ \end{array}\right)-c\right|=\left(\begin{array}{ccccc} 7.7579\\ -22.3401\\ -0.6957\\ 19.0797\\ -11.1326\\ \end{array}\right)$, 则 $(1, 1, 1, 1, 1)^{T}$即为广义绝对值方程 $Ax-|Bx-c|=b$的准确解.随机选取初始点进行迭代, 计算结果如下表 2.

对于更高阶的算例运算结果见下表 3, 其中

$ \begin{eqnarray*}&&B={\hbox{diag}}(1:n), A={\hbox{diag}}(n+1:2n), c={\hbox{randint}}(n, 1, [-5,5]), \\ &&b=A*{\hbox{ones}}(n, 1)-{\hbox{abs}}(B*{\hbox{ones}}(n, 1)-c, \end{eqnarray*} $

$\mu=0.001, $初始点 $x^0=(x_1, x_2, \cdots, x_n)^T$, $x_i=0.1*i~~ (i=1, 2, \cdots, n)$, 计算结果如下表 3.

算例3  (同文献[19]中算例)令

$ \begin{eqnarray*}&&A=\left( \begin{array}{cccc} 8.1537&2.6795&-0.6134&-1.7996\\ -1.9500&13.7084&-0.0169&4.6010\\ 3.7437&4.9008&6.1396&2.2663\\ -4.8499&2.8886&1.4349&8.1195\\ \end{array}\right), \\ &&B=\left( \begin{array}{cccc} 0.0158&0.0576&0.6927&0.3533\\ 0.0164&0.3676&0.0841&0.1536\\ 0.1901&0.6315&0.4544&0.6756\\ 0.5869&0.7176&0.4418&0.6992\\ \end{array}\right), \\ &&c=\left( \begin{array}{cccc} 0.7159\\ 0.8928\\ 0.2731\\ 0.2548\\ \end{array}\right), b=\left( \begin{array}{cccc} 7.3121\\ -5.3530\\ 6.0974\\ 8.1680\\ \end{array}\right).\end{eqnarray*} $

随机选取不同初始点进行迭代, 计算结果如下表 4.

通过以上表 1、表 2、表 3可以看出该方法对于求解广义绝对值方程(1.1) 是十分有效的, 而且收敛速度也比较快, 当 $\mu=0.01$时计算结果就已经十分准确, 随着 $\mu$的不断减小计算的准确性会越来越高, 而且通过初始点的随机选取也可看出算法的全局收敛性; 表 4是文献[19]中的算例, 通过结果比较不难发现本文算法迭代次数减少了近一倍, 要明显优于文献[19]中的迭代算法.

本文基于广义绝对值方程与垂直线性互补的等价性转化, 将广义绝对值方程进一步转化为光滑方程组并用牛顿法对其进行求解, 该算法不仅收敛速度较快、计算准确而且具有全局收敛性, 大量数值实验表明在计算过程中只要选取适当小的参数 $\mu$便可得出十分精确的近似解.鉴于绝对值方程求解是比较困难的数学问题, 本文算法可以作为求解广义绝对值方程的一个有效算法.有没有更好、更快、收敛性更好的算法对广义绝对值方程及其它形式的绝对值方程进行求解将成为笔者后续的研究工作.