何國(guó)良， 張文星， 趙熙樂(lè)

(電子科技大學(xué) 數(shù)學(xué)科學(xué)學(xué)院，成都611731)

1 引 言

數(shù)學(xué)是一門(mén)實(shí)踐性很強(qiáng)的學(xué)科，目的是從數(shù)量關(guān)系和空間結(jié)構(gòu)來(lái)研究客觀對(duì)象的一些特殊屬性(恩格斯、布爾巴基學(xué)派、李大潛等語(yǔ)). 數(shù)學(xué)成果的發(fā)現(xiàn)通常不是簡(jiǎn)單的靈機(jī)一動(dòng)，而是需要反復(fù)實(shí)踐，是建立在直觀和多次失敗的計(jì)算分析上，尤其是現(xiàn)代數(shù)學(xué)更是如此. 但和物理、化學(xué)、生物學(xué)等需要做各種各樣的實(shí)驗(yàn)來(lái)發(fā)現(xiàn)新的問(wèn)題不一樣，數(shù)學(xué)研究主要實(shí)驗(yàn)方式是“思考”. 這種方式和其他學(xué)科的“硬核實(shí)驗(yàn)”比起來(lái)，顯得比較抽象，難以理解. 實(shí)踐表明如果讓“思考”變得直觀一些，則無(wú)論對(duì)概念、方法的理解，還是把數(shù)學(xué)作為一種工具來(lái)處理實(shí)際問(wèn)題就會(huì)更加流暢. 著名數(shù)學(xué)家F. Enriques、G. Castelnuovo、J.A. Todd、小平邦彥(菲爾茲獲得者)等人均認(rèn)為數(shù)學(xué)的發(fā)現(xiàn)需要直觀，需要通過(guò)“觀察”數(shù)學(xué)現(xiàn)象來(lái)總結(jié)規(guī)律[1].

共軛梯度法(Conjugate Gradient Method)，簡(jiǎn)稱(chēng)CG法，是應(yīng)用數(shù)學(xué)和工程領(lǐng)域中的一種重要方法，自該方法被研發(fā)出來(lái)以后，無(wú)論是在工程界還是在理論界都引起了大家的廣泛興趣和應(yīng)用[2-5]. 目前仍然是一個(gè)非?；钴S的研究話(huà)題. 人們對(duì)其不斷改進(jìn)，不斷推廣其應(yīng)用場(chǎng)景，獲得了良好的效果[6-7]，尤其是在如日中天的智能算法方面，更是起著舉足輕重的作用[8-9]. 共軛梯度法已經(jīng)成為理工類(lèi)大學(xué)本科或者是研究生課程數(shù)值分析教學(xué)里面的一個(gè)重要內(nèi)容[10-12]. 然而，由于這個(gè)方法涉及比較多、比較深的數(shù)學(xué)理論，在教學(xué)過(guò)程中我們發(fā)現(xiàn)，如果直接按照純粹代數(shù)推導(dǎo)的方法去講授，同學(xué)們普遍掌握不了抽象、繁瑣的推導(dǎo)方法，或?qū)τ谠摲椒ǖ幕舅枷肜斫獠皇呛芡笍兀陀绊懥嗽摲椒ㄔ趯?shí)踐中的應(yīng)用.

本文結(jié)合多年的教學(xué)經(jīng)驗(yàn)，從幾何的角度來(lái)討論共軛梯度法基本思想. 在此基礎(chǔ)上，通過(guò)把代數(shù)的嚴(yán)密和幾何的直觀結(jié)合起來(lái)，比較流暢地展示了該方法“自然而然”的思維發(fā)展和遞進(jìn)過(guò)程，為學(xué)生進(jìn)一步理解和掌握共軛梯度法提供了一種比較簡(jiǎn)單的途徑.

2 基本概念

極小化方法是求函數(shù)極小值的典型方法，在工程設(shè)計(jì)、人工智能、運(yùn)籌與優(yōu)化、線(xiàn)性方程組和非線(xiàn)性方程組求解等方面都有廣泛的用途[5，8，9]. 為避免引入過(guò)多的背景知識(shí)，這里以求解線(xiàn)性方程組所導(dǎo)出的非線(xiàn)性模函數(shù)(多元二次函數(shù)，參見(jiàn)下面定理1中的公式(2))的極小值為例來(lái)討論共軛梯度法，這與在其他方面導(dǎo)出的求函數(shù)極小化方法是一樣的. 下面簡(jiǎn)要介紹一下將方程組求解轉(zhuǎn)換為函數(shù)極小化問(wèn)題的簡(jiǎn)要過(guò)程.

設(shè)A是n階對(duì)稱(chēng)正定矩陣，x，y∈n，待求解的線(xiàn)性方程組為

Ax=b.

(1)

若x*為方程組(1)的唯一解，則有如下重要的結(jié)果[10]：

定理1設(shè)A= (aij)n×n為n階實(shí)對(duì)稱(chēng)正定矩陣，b，x∈n，則使下面二次函數(shù)

(2)

取極小值的x*是線(xiàn)性方程組Ax=b的解，即Ax*=b；反之亦然.

上述定理1是線(xiàn)性方程組與函數(shù)極小值的等價(jià)定理，通常把函數(shù)(2)叫模函數(shù)[10].根據(jù)該定理，求解線(xiàn)性方程組的唯一解就轉(zhuǎn)換成為求模函數(shù)極小值的問(wèn)題：

任取一個(gè)初始向量x(0)；

構(gòu)造迭代格式:

x(k+1)=x(k)+αkx(k)， (k=0，1，2，…)，

(3)

其中p(k)是搜索方向，αk是搜索步長(zhǎng)；

選擇p(k)和αk，使得

φ(x(k+1))=φ(x(k)+αkp(k)) <φ(x(k+1))，

在上述的計(jì)算過(guò)程里面，關(guān)鍵是搜索方向p(k)的選取.實(shí)際計(jì)算中k一般取不太大的有限值，只要x(k)滿(mǎn)足一定要求即可停止迭代.

3 共軛梯度法的一般特點(diǎn)

3.1 共軛梯度法的出現(xiàn)緣由

無(wú)論是出于工程方面的需求還是算法研究本身，都要求數(shù)值方法能迅速找到模函數(shù)(2)的極小值.由于模函數(shù)(2)是一個(gè)多元連續(xù)可微的函數(shù)，并且在矩陣是正定、對(duì)稱(chēng)情況下，模函數(shù)是一個(gè)凸函數(shù)，所以為迅速求模函數(shù)的最小值點(diǎn)，人們自然的想法就是把模函數(shù)的負(fù)梯度方向作為搜索方向，從而得到經(jīng)典的最速下降法[10，12]：

算法1 最速下降法

(i) 選取x(0)∈n；

(iii) 當(dāng)‖x(k+1)-x(k)‖<ε時(shí)，終止迭代.

從上面算法1可以看出，該方法計(jì)算格式簡(jiǎn)單，其搜索方向p(k)就是殘差方向r(k)，整個(gè)迭代過(guò)程只需要一些矩陣和向量作代數(shù)運(yùn)算即可.但實(shí)踐表明，該方法整體效率可能不高，甚至對(duì)于一些對(duì)稱(chēng)正定矩陣，也可能存在嚴(yán)重的鋸齒現(xiàn)象，如圖1所示.這樣一來(lái)，為了提高求解效率，人們研究發(fā)展了共軛梯度法.從上面的思路發(fā)展過(guò)程可以看出最速下降法、共軛梯度法的發(fā)生、發(fā)展情況如下：

求模函數(shù)極小值→最速下降法→需要改進(jìn)→共軛梯度法

一般來(lái)說(shuō)，教學(xué)活動(dòng)如果按照如上過(guò)程設(shè)計(jì)，到此具有此邏輯清晰、容易理解等特點(diǎn)，思維也是比較流暢的.

圖1 最速下降法的鋸齒現(xiàn)象

3.2 共軛梯度法的傳統(tǒng)教學(xué)設(shè)計(jì)

從前面關(guān)于極小值的通用計(jì)算過(guò)程可以看出，極小值計(jì)算過(guò)程中，關(guān)鍵是公式(3)中的方向向量p(k)的選取.縱觀國(guó)內(nèi)外教材，由于共軛梯度法“不直觀”，所以通常采用純代數(shù)分析的角度來(lái)展開(kāi)討論選取方案——共軛(一個(gè)抽象概念)，下面簡(jiǎn)要介紹一下傳統(tǒng)的處理方法，比如如下形式[10]：

設(shè)p(0)=r(0)，當(dāng)k≥1時(shí)，對(duì)于搜索方向p(k)的確定，不但希望使

而且希望{p(k)}的選擇使

由于x∈span{p(0)，…，p(k)}，則可將x記成

x=y+αp(k)，y∈span{p(0)，…，p(k-1)}，α∈.

所以有

(4)

為比較簡(jiǎn)單地計(jì)算(4)，可令

(Ay，p(k))=0，y∈span{p(0)，…，p(k-1)}，

即

(Ap(j)，p(k))=0，j=0，1，…，k-1.

(5)

公式(5)在計(jì)算過(guò)程(4)中起著重要作用，為此需要引出下面共軛的定義：

定義1設(shè)A是對(duì)稱(chēng)正定矩陣，若n中的向量組{p(0)，p(1)，…，p(k)}滿(mǎn)足

(Ap(i)，p(j))=0，i≠j，

這樣一來(lái)，如果將方向向量取為共軛向量組，則

(6)

對(duì)上式(6)中的第一個(gè)最小問(wèn)題的解為x(k)；對(duì)于第二問(wèn)題，若取p(k)，使得A(x(k)，p(k))=0(共軛)，則通過(guò)一些和最速下降法一樣的簡(jiǎn)單推導(dǎo)，就可以到如下的共軛梯度法[10，12]：

算法2共軛梯度法

(i) 選取x(0)∈n，r(0)=b-Ax(0)；

(iii) 當(dāng)‖x(k+1)-x(k)‖<ε時(shí)，終止迭代.

從前面關(guān)于極小值的通用計(jì)算過(guò)程可以看出，極小值計(jì)算過(guò)程中，關(guān)鍵是求公式(3)中的方向向量{p(k)}.由于模函數(shù)(公式(2))的負(fù)梯度方向在求解模函數(shù)極小值方面的意義明確清晰，并且對(duì)于對(duì)稱(chēng)正定矩陣來(lái)講，其對(duì)應(yīng)模函數(shù)負(fù)梯度方向剛好就是殘差方向，所以無(wú)論從代數(shù)角度還是從幾何角度來(lái)講，最速下降法都容易理解.但是對(duì)于從“基”及“向量空間”的構(gòu)造來(lái)介紹的共軛梯度法，雖然組織方式嚴(yán)謹(jǐn)，但是不具備良好的幾何直觀性，也無(wú)法向?qū)W生直觀展示共軛是什么，為什么要這樣選，在迭代過(guò)程中解有什么樣的變化特點(diǎn)等.根據(jù)筆者的教學(xué)經(jīng)驗(yàn)來(lái)看，以這種方式來(lái)介紹共軛梯度法方法，事倍功半，講解時(shí)間花不少，但是學(xué)反應(yīng)不佳.對(duì)于空間抽象思維比較好的學(xué)生，基本能夠順利理解該方法，但是對(duì)于數(shù)學(xué)基礎(chǔ)稍弱的學(xué)生，上面的介紹過(guò)程就很具有挑戰(zhàn)性了.

4 圖解共軛梯度法

下面簡(jiǎn)要介紹從幾何圖形的角度來(lái)的討論共軛梯度法，為簡(jiǎn)單起見(jiàn)，這里首先以二維情況來(lái)說(shuō)明.

4.1 二維情況

對(duì)求解兩個(gè)變量的對(duì)稱(chēng)正定線(xiàn)性方程組來(lái)講，其對(duì)應(yīng)的二元模函數(shù)，不妨記為φ(x1，x2)(注意這里的x1，x2表示二維向量的兩個(gè)分量)，該函數(shù)的圖形為三維空間一曲面，如下圖：

圖2 二維方程組模函數(shù)曲面及在XOY坐標(biāo)面的投影 圖3 二次搜索計(jì)算示意圖

在給定初值的情況下，為求該曲面的最小值，也選取初值的殘差作為首次搜索方向.在此基礎(chǔ)上，下面研究第二次搜索方向該如何選取——希望直接選取類(lèi)似圖3中的p(1)向量(直接指向最小值點(diǎn))，而不是負(fù)梯度向量-▽?duì)?x(1))，這樣總共只需要2次迭代就可以得到極小值點(diǎn)，從而實(shí)現(xiàn)高效計(jì)算.下面來(lái)看看是否能夠找到這樣的向量.

如果p(1)存在，希望在x(1)和p(1)的基礎(chǔ)上，通過(guò)選擇恰當(dāng)?shù)膖，就可以得到

x*=x(1)+tp(1)，t∈，

其中x*表示原來(lái)方程組的解.

下面具體來(lái)看這樣的p(1)是否存在.將上式帶入原始的線(xiàn)性方程組有

Ax*-b=Ax(1)-b+tAp(1).

當(dāng)A為對(duì)稱(chēng)正定矩陣時(shí)，由于▽?duì)?x)=Ax-b，所以上式可寫(xiě)為

▽?duì)?x*)=▽?duì)?x(1))+tAp(1).

由于▽?duì)?x*)=0，有

0=▽?duì)?x(1))+tAp(1).

在上式兩端同時(shí)左乘p(0)T，有

0=p(0)T▽?duì)?x(1))+tp(0)TAp(1).

由于p(0)和第1次迭代的殘差向量垂直，即(p(0)，▽?duì)?x(1)))=0，則上式變?yōu)?/p>

0=p(0)TAp(1).

(7)

上式說(shuō)明什么呢？上式說(shuō)明：要想兩次迭代就得到線(xiàn)性方程組的準(zhǔn)確解(忽略舍入誤差的意義下)，那么第2次迭代中理想的搜索方向應(yīng)該滿(mǎn)足(7)式，也就是前面提到的共軛(定義1).那具體如何計(jì)算這個(gè)共軛向量呢？

對(duì)于二維問(wèn)題，由于-▽?duì)?x(1))⊥p(0)，且三向量p(0)， -▽?duì)?x(1))，p(1)共面，則可將向量p(1)表示為

p(1)=-▽?duì)?x(1))+β1p(0).

上式兩邊左乘(p(0))TA，得到

p(0)TAp(1)=-p(0)T▽?duì)?x(1))+β1p(0)Tp(0).

由于p(0)和p(1)共軛，得

0=-p(0)T▽?duì)?x(1))+β1p(0)Tp(0)，

即得

這樣一來(lái)，就唯一確定出第2次的方向向量.通過(guò)這樣的過(guò)程，不僅從幾何知道理想的方向向量是存在的，而且是可計(jì)算的.

4.2 對(duì)于三維情況

三個(gè)變量的方程組對(duì)應(yīng)的模函數(shù)，不妨記為φ(x1，x2，x3)，其圖形在笛卡爾坐標(biāo)系下為四維超曲面.為直觀起見(jiàn)，這里用如下的3維空間(OXYZ)來(lái)展示.圖4中的左圖為超曲面在Z方向取不同的值得到不同三維曲面，右圖為在左圖的基礎(chǔ)上將圖形沿著XOY平面的第一象限對(duì)角線(xiàn)條逐漸平移的結(jié)果：包括曲面和在XOY平面上的等值線(xiàn)變化情況.

圖4 三維方程組模函數(shù)曲面及在坐標(biāo)面的投影

圖5為三元模函數(shù)φ(x1，x2，x3)在x3=0的投影曲面和在x2=0的投影曲面，藍(lán)色曲線(xiàn)為x2=0的投影曲面在XOZ平面的等值線(xiàn)，彩色曲線(xiàn)為x3=0的投影曲面在XOY平面的等值線(xiàn).通過(guò)圖4和圖5可以從不同的角度理解模函數(shù)的基本特點(diǎn).下面討論三元函數(shù)的共軛梯度的幾何特點(diǎn).

圖5 三元模函數(shù)投影曲面及在坐標(biāo)面的等值線(xiàn)投影

如圖6所示，三維問(wèn)題的共軛梯度法首先可以看成在超曲面的在3維空間的一個(gè)投影曲面上找一個(gè)初值x(0)，然后選取該點(diǎn)的負(fù)梯度方向作為搜索方向，沿著這個(gè)方向到達(dá)一個(gè)平面x(1)，這樣一來(lái)，三維的共軛梯度法就被降維為變成二維的共軛梯度法(投影)，這樣再運(yùn)用前面二維的討論方法就可以得出完整的三維共軛梯度法，如圖6的左圖和右圖所示.為進(jìn)一步直觀，下面通過(guò)幾個(gè)數(shù)值例子來(lái)進(jìn)一步了解共軛梯度法的一些其他特點(diǎn).

圖6 三元模函數(shù)共軛梯度法搜索方向示意圖

4.3 數(shù)值例子

例1利用共軛梯度法求AX=b的解.其中

該矩陣對(duì)應(yīng)的模函數(shù)為

共軛梯度法計(jì)算結(jié)果如下表：

表1 2階方程迭代情況

這里選的初值為x=(0，0)T，該問(wèn)題的真解為x=(-0.5，1)T.所以在四舍五入范圍內(nèi)，共軛梯度法只需要2步就可以得到非常準(zhǔn)確的結(jié)果.從計(jì)算結(jié)果來(lái)看，確實(shí)符合前面的幾何解釋.

例2利用共軛梯度法求3階線(xiàn)性方程組AX=b的解.其中

迭代計(jì)算結(jié)果為

表2 3階方程迭代情況

這里選的初值為x=(0，0，0)T，該問(wèn)題的真解為x=(1，1，1)T.所以在四舍五入范圍內(nèi)，共軛梯度法也只需要3步得到非常準(zhǔn)確的結(jié)果.表3是不同初值對(duì)求解過(guò)程的影響.

表3 不同初值對(duì)應(yīng)的共軛向量

從上表3可以看出，選取不同的初值，只需經(jīng)過(guò)3次迭代后都可以得到理想的結(jié)果，且共軛梯度法確實(shí)保持的前后向量的共軛性(在矩陣乘法意義下的正交).

表4 不同初值對(duì)應(yīng)的求解精度

不同的初值，經(jīng)過(guò)3次迭代后的殘差向量和解的絕對(duì)誤差都很小，所以都能得到精度很高的解而且計(jì)算穩(wěn)定.

例3利用共軛梯度法求AX=b的解， 其中

圖7 三元模函數(shù)共軛梯度法搜索路徑

圖7分別從三維和二維角度展現(xiàn)對(duì)于3階線(xiàn)性方程，當(dāng)初值選為x(0)=(1，1，1)T時(shí)，共軛梯度法的搜索向量是如何變化的.它能從幾何上地展現(xiàn)正交性以及近似解到達(dá)最小值點(diǎn)的變化過(guò)程——這樣就便于學(xué)生理解共軛梯度法.

圖8 不同初值對(duì)三元模函數(shù)共軛梯度法搜索路徑的影響

圖8分別從三維和二維的幾個(gè)不同角度展現(xiàn)3個(gè)初值：x(0)=(1，1，1)T，x(0)=(-1，1，-1)T，x(0)=(-1，-1，-1)T，對(duì)于3階線(xiàn)性方程的搜索向量和在近似解到達(dá)最小值點(diǎn)的變化過(guò)程.從這里可以看出，不同的初值搜索方向系不一樣.當(dāng)然，這是容易理解的，不同的初值殘差自然不同，這導(dǎo)致緊隨其后的搜索方向也不同.此外，這個(gè)例子也展示了該算法的穩(wěn)定性——對(duì)初值不敏感.

4.4 高維推廣

圖9是高維空間堆積示意圖.雖然下面這種構(gòu)成很簡(jiǎn)單，但它向我們展現(xiàn)了高維空間的一種直觀構(gòu)成形式，這種形式類(lèi)似于直角坐標(biāo)系(標(biāo)準(zhǔn)正交)的疊加構(gòu)成.從這個(gè)示意圖，可以大致理解抽象的高維空間的構(gòu)成形式及體會(huì)共軛梯度法在其中的變化過(guò)程.

圖9 1至6維空間構(gòu)成示意圖

這樣一來(lái)，高維空間的共軛梯度法，就可以比較容易理解了.比如對(duì)含5個(gè)自變量線(xiàn)性方程組，其對(duì)應(yīng)的的5元模函數(shù)為6維空間(由三位空間疊加而成)的一個(gè)超曲面.那么為求其在該空間的極小值(方程組的解)，可以直觀地理解為：首先將6維空間的超曲面投影到5維空間中——沿著6維空間的豎直方向(Z′方向)進(jìn)行投影，并在該5維超曲面上取一個(gè)初值，x(0)∈5；然后再沿著斜向方向(Y′方向)求一次最小值，得x(1)∈5，再沿著橫向方向(X′方向)求最小值x(2)∈5.這樣兩次迭代以后，最小值將落入普通的三維空間曲面(值還是5維的)，這樣就可以利用繼續(xù)前面三維空間的圖示來(lái)理解了.

利用類(lèi)似的思想可以繼續(xù)構(gòu)造或理解7維、8維、9維的空間的共軛梯度法，它們可以3的倍數(shù)進(jìn)行躍遷.這似乎和我們傳統(tǒng)文化的“三生萬(wàn)物”有異曲同工之效，所以先賢們確實(shí)具有非凡智慧和想象力.由于這種共軛梯度采取空間投影來(lái)實(shí)現(xiàn)降維的目的，共軛梯度法有時(shí)候也被歸入投影法，或子空間方法的范疇.

值得說(shuō)明的是，當(dāng)然可以采取仿射變換、非線(xiàn)性坐標(biāo)變換的形式去構(gòu)造非線(xiàn)性、非正交的空間疊加形式，只不過(guò)這種構(gòu)成比較復(fù)雜，這里就不展示了.

5 結(jié) 論

和其他學(xué)科自然學(xué)科一樣，如果能把數(shù)學(xué)方法和結(jié)果作為一種可以直接觀察的現(xiàn)象時(shí)，所產(chǎn)生的進(jìn)步是非常顯著的，這也能讓抽象的數(shù)學(xué)概念或方法成為一個(gè)有意思，容易被理解、接受的的成果.懷這樣的愿望，從二維問(wèn)題出發(fā)，討論了共軛梯度法的幾何意義和實(shí)現(xiàn)：首先通過(guò)圖形進(jìn)行直觀分析，然后再進(jìn)行數(shù)學(xué)推導(dǎo)，從而比較直觀地導(dǎo)出了共軛梯度法.通過(guò)比較豐富的數(shù)值例子和在不同空間中的推廣，本文比較好地解釋了共軛梯度法的基本思想和計(jì)算過(guò)程.近兩年的教學(xué)實(shí)踐也表明，這樣的教學(xué)設(shè)計(jì)，比起傳統(tǒng)純代數(shù)分析的教學(xué)過(guò)程，學(xué)生更容易理解，也更能體會(huì)到該方法的好處.尤其是在理解諸如“投影”“正則化”這些在數(shù)學(xué)和工程應(yīng)用中非常重要，但又很抽象概念和方法的基本思想時(shí)至關(guān)重要，這也能讓學(xué)生將這些方法廣泛地應(yīng)用科學(xué)研究和工程實(shí)踐中.

致謝作者非常感謝相關(guān)文獻(xiàn)對(duì)本文的啟發(fā)以及審稿專(zhuān)家提出的寶貴意見(jiàn)；感謝學(xué)校數(shù)學(xué)分析課程教學(xué)團(tuán)隊(duì)在教學(xué)交流過(guò)程中予以的富有啟發(fā)意義的建議；感謝《大學(xué)數(shù)學(xué)》編輯予以的耐心指導(dǎo)和細(xì)致編撰工作.