撰文
| 丁玖(美国南密西西比大学数学系教授)
记得有个大科学家曾经说过,这个世界是非线性的。然而非线性方程一般不能直接求解,即解析解虽可证存在却无具体表达式,因而迭代法几乎是唯一可行的办法。比如说,从介值定理可知,方程x = cos x在区间(0, 1)内定有一解,但没有一步到位的法子找到它,人们只能用基于介值定理的二分法或基于切线逼近的牛顿法,来求得此方程的迭代近似解。这样,从最古老的巴比伦平方根迭代法,到今日非线性方程组数值解的最重要方法——牛顿迭代法,人们一直热衷于迭代法的理论探索和创新实践。
既然有非线性方程,就必然有线性方程,而且后者在科学技术和日常生活中同样随处可见。电话公司的客户安排、航空公司的航班分配、运输公司的运货路线……都可以用线性规划来优化成本。四十年前,当我在南京大学数学系时,读到一篇关于线性规划的综述性文章,其中写道:全世界的计算机大约有百分之六十的时间,都在求解线性规划问题!
另一方面,微分概念告诉我们,非线性函数可以局部地用线性函数来近似,构造牛顿迭代法本质上就基于这个简单的观察。只要非线性函数f在一个点a存在导数,那么在a点的附近,函数值f(x)就约等于线性函数值f(a) + f’(a)(x-a),其几何直观则是函数图像的切线在切点附近与曲线相差无几。
此外,在物理科学中,许多描述自然规律的常微分方程或偏微分方程对于未知函数以及它们的导数或偏导数都呈现出线性关系,如热传导方程或波动方程。当用计算数学中的差商逼近导数技术或有限元思想来离散化这些连续方程时,获得的差分方程也是线性的。它们可以写成一般形式Ax = b,其中A是一个n行、n列的系数矩阵,b是一个n维的已知列向量,x是一个n维的未知列向量。这个线性方程看上去像一元一次方程ax = b一样简单,但如果按照矩阵乘法的法则将方程左边每个分量的代数表达式全部写出来,结果就是一组含有n个未知数x
1
, x
2
, …, x
n
的n个n元一次方程。如果将方矩阵A中第i行、第j列的元素记为a
ij
,将列向量b的第i个分量记为b
i
,那么线性方程组Ax = b展开后的第i个方程为
a
i1
x
1
+ a
i2
x
2
+ … + a
in
x
n
= b
i
, i = 1, 2, …, n。
计算数学中的一个子学科“数值线性代数”有个基本的职责:怎样有效求解上面的线性方程组?当n = 2或3,最多不超过4时,我们早已在中学课堂上学会了怎样求解二元一次方程组、三元一次方程组,等等。采用的方法无非是两种,一种叫做代入法,另一种称为加减法。它们都是基于同一个思想:通过消去方程中的一个或几个未知数,设法获得只含有一个未知数的方程,比方说-2y = 3,从而解出这个未知数。
中学生所学的上述方法,实际上是大名鼎鼎的高斯消去法的简单实现,发明者正是德国天才高斯(Carl Friedrich Gauss,1777-1855)
。高斯消去法的基本思想还是如上所说:将含有多个未知数的方程化简到只剩下一个未知数,因为只含一个变元的线性方程如-2y = 3总是可以一步到位解出y = -3/2的。基于这个非常简单的思路,高斯想出了一个点子,将线性方程组的系数方阵程式化系统性地化约成主对角线下方的元素全为零的上三角矩阵,这样新的等价的线性方程组的最后一个方程只含一个未知数,倒数第二个方程含有两个未知数,如此等等。因而,通过所谓的“回代法”,从最后一个方程解出唯一的未知数,然后将它的值代入到倒数第二个方程,结果该方程也只含一个未知数,马上可以解出,再将这两个未知数的值代入到倒数第三个方程,同法解出下一个的未知数。如此这般,就可以依次解出所有的未知数。
高斯消去法在经过有限次代数运算后就可以完全找出给定线性方程组的精确解(如果每步运算的误差为零的话)。所以它属于求解线性方程组的第一类数值方法——直接法,这是解线性方程组优越于解非线性方程组的一个地方。然而,对有巨量(如成千上万个或更大数量级)变元个数或者特别多的系数为零(比如系数矩阵为三对角或其他稀疏类型)
的那种方程组,直接法常常没有第二类方法——迭代法有效。对线性偏微分方程直接用差商代替偏导数的差分方法,或将微分方程问题变为数学上等价的变分问题而进行数值优化的有限元方法,所得到的大型线性方程组一般都带有特定结构的稀疏系数矩阵,这时,迭代法几乎是最有效的算法了。
我们下面只讲解线性方程组的迭代法。让我们回忆求解单变量非线性方程的迭代法,一般形式是x
n
= f(x
n-1
), n = 1, 2, …,其中f是将定义域区间映到自身的一个函数,x
0
是迭代所取的初始点。然而,对于线性迭代法,迭代函数不再是一个自变量的线性函数,而是有n个自变量的线性向量函数。由于字母n现在另有他用,我们将用字母k代表迭代次数的下标,而将多变量线性向量函数用y = L(x)表示,其中L(x)的表达式是Mx + c,M是一个有n行和n列的给定矩阵
(也称为n阶正方矩阵或n阶方阵)
,c是一个给定的n维列向量,其分量是c
1
, c
2
, …, c
n
,x = (x
1
, x
2
, …, x
n
)
T
是n维变元列向量,其中的上标字母T表示转置运算,即一般的有m行和n列的矩阵A的转置矩阵A
T
有n行和m列,其第j行的元素是原矩阵第j列的对应元素。
本章将要面对的线性迭代程序可以写成
x
k
= L(x
k-1
) = Mx
k-1
+ c, k = 1, 2, …; x
0
是给定的。
我在以前的科普文章里写过,单变量非线性方程的迭代法收敛的一个充分条件,是迭代函数在不动点的导数绝对值小于1。那么,如果被迭代的是有多个变量的线性向量函数,什么样的条件可以充分到保证迭代法收敛呢?所谓序列的收敛,本质上是用距离这个概念来定义的,这个概念我们在以前提到巴拿赫不动点定理或压缩映射定理时已经见到过。我们再来回忆有关内容。在一个抽象的距离空间(X, d)里,一个由X中的元素所组成的无穷序列{x
k
}被说成是收敛到X中的一个元素x,指的是由距离函数d所确定的非负数列{d(x
k
, x)}当k趋向于无穷大时趋向于0。对于单变量函数迭代法的收敛性问题,我们实际上选取了一个简单而具体的距离空间,它就是所有实数全体构成的一维欧几里得空间R,其任意两个实数之间的距离就是实数轴上两点x和y之间的通常距离|x – y|。
上面两数之间的距离可用线段长度的说法来等价地定义。任意一个实数a的绝对值|a|在几何上的意义是以数轴上代表a的那个点和原点0为两端点的线段的长度。不难看出,两个数x和y之间的距离就是数x – y所对应的那条线段的长度。
对于现在要考察的高维欧几里得空间R
n
,两个n维向量之间的距离可如法炮制地定义。对于二维平面R
2
,设第一个行向量为x = (a, b),第二个行向量是y = (c, d),则它们之间的距离就是平面上这两点之间的通常距离
, 或等价地,它是平面向量x – y的长度。在通常的三维空间R
3
里,设向量x 的三个分量是a, b, c,向量y的三个分量是u, v, w,则x和y之间的距离也为这两个空间点之间的通常距离
。它也就是空间向量x – y的实际长度。
三维以上的欧几里得空间,我们的眼力无法感知,但我们的想象力可以穿越界限,所以将上面直观的距离公式推广开来,就得到n维空间R
n
的欧几里得距离定义:设向量x的分量为a
1
, a
2
, …, a
n
,向量y的分量为b
1
, b
2
, …, b
n
,则它们之间的距离是
d(x, y) = {(a
1
- b
1
)
2
+ … + (a
n
- b
n
)
2
}
1/2
。
我们将三维以上空间中的向量“长度”改称为“欧几里得范数”,简称范数或2-范数,则n维向量x的范数||x||
2
被定义为x的所有分量的平方之和再开平方根。这样,x和y之间的距离就等于向量x – y的欧几里得范数||x – y||
2
。
作为长度概念的推广,范数继承了长度的三个基本性质: