胡晶晶 ,柯藝芬 ,馬昌鳳

(1.福建師范大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,福建 福州 350117;2.福建師范大學(xué)分析數(shù)學(xué)及應(yīng)用教育部重點(diǎn)實(shí)驗(yàn)室,福建 福州 350117;3.福建省應(yīng)用數(shù)學(xué)中心,福建 福州 350117)

1.引言

四元數(shù)經(jīng)常出現(xiàn)在許多領(lǐng)域,如量子物理、信號(hào)處理、彩色圖像處理和計(jì)算機(jī)科學(xué)等[1-4].由于四元數(shù)域比復(fù)數(shù)域提供了更多的自由度,因此在處理多維問(wèn)題時(shí),四元數(shù)比實(shí)數(shù)和復(fù)數(shù)具有更多的通用性和靈活性.因此,越來(lái)越多的學(xué)者對(duì)四元數(shù)的理論性質(zhì)和計(jì)算問(wèn)題感興趣,并取得了許多有價(jià)值的成果.例如,文[5-6]通過(guò)復(fù)表示方法獲得了四元數(shù)矩陣方程AXB+CY D=E和AHXA+BHY B=C的極小范數(shù)η-(反)-Hermitian最小二乘解的顯示形式.文[7]通過(guò)實(shí)表示方法推導(dǎo)出了四元數(shù)矩陣方程AXB+CY D=E的η-Hermitian和η-反-Hermitian解.文[8]建立了一種有效的迭代方法去獲得四元數(shù)矩陣方程AXB+CY D=E的極小范數(shù)η-Hermitian和η-反-Hermitian最小二乘解.文[9]考慮了四元數(shù)Sylvester張量方程的最小二乘問(wèn)題,通過(guò)求解其等價(jià)形式從而得到了四元數(shù)Sylvester張量方程的最小二乘解.

本文采用以下符號(hào): R表示實(shí)數(shù)集,Rn表示n維實(shí)向量空間,Rn×n表示n×n階實(shí)矩陣集合,表示N階I1×I2×···×IN維實(shí)張量的集合,Q表示四元數(shù)體,Qn×n表示n×n階四元數(shù)矩陣集合,表示N階I1×I2×···×IN維四元數(shù)張量集合;AT表示矩陣A的轉(zhuǎn)置,R(A)表示矩陣A的列空間.

2.預(yù)備知識(shí)

定義2.1四元數(shù)a可以表示為

其中as ∈R,s=1,2,3,4,且i,j,k滿足

對(duì)于兩個(gè)四元數(shù)a=a1+a2i+a3j+a4k∈Q和b=b1+b2i+b3j+b4k∈Q,它們的乘積定義為

定義2.2四元數(shù)張量A可以表示為

其中As ∈,s=1,2,3,4.四元數(shù)張量A的共軛可被定義為

類似地,對(duì)于四元數(shù)矩陣A,有著相似的表示.

例2.1[10]張量A ∈R3×4×2的兩個(gè)正面切片分別為

張量A沿模-1展開(kāi)成矩陣A(1),其中

定義2.4[11]對(duì)于張量A ∈,算子‘vec’將張量的列疊加成一個(gè)向量,定義vec(A)=vec(A(1)),其中A(1)是張量A沿模-1的展開(kāi)矩陣.

例2.2若張量A ∈R3×4×2如例2.1所示,則vec(A)=[1,2,···,24]T.

證四元數(shù)Sylvester張量方程(2.2)可寫(xiě)作

將式(2.7)和式(2.8)代入式(2.6),根據(jù)定義2.7,四元數(shù)張量相等意味著實(shí)部和虛部分別對(duì)應(yīng)相等,從而可得方程組(2.3).

進(jìn)一步,根據(jù)Kronecker積和向量化算子,從而可得線性方程組(2.4).

注2.1值得注意的是當(dāng)線性方程組(2.4)的系數(shù)矩陣規(guī)模非常大時(shí),尋找其數(shù)值解將面臨很大的計(jì)算挑戰(zhàn).文[12]表明,基于張量格式的算法總體上比經(jīng)典形式的算法效率更高.下面我們將構(gòu)造一類基于張量形式的高效迭代算法求解方程組(2.3),并建立迭代算法的收斂結(jié)果.

注2.2記

3.全局收斂性

本節(jié)主要考慮四元數(shù)張量方程(2.2)的求解.利用定理2.1,可以將四元數(shù)張量方程轉(zhuǎn)化為實(shí)張量方程組(2.3),考慮用張量形式的修正共軛梯度(Modify Conjugate Gradient,MCG)算法求解(2.3),從而得到方程(2.2)的解.

為了方便起見(jiàn),我們引入下列線性算子:

從而四元數(shù)Sylvester張量方程(2.2)可寫(xiě)成

步4 置k:=k+1,返回步2.

定理3.1假設(shè)序列{Rs(k)}和{Qt(k)}(s=1,2,3,4,t=1,2,3,4)由算法3.1生成,則下列等式成立

證首先,我們證明式(3.3)和式(3.4)對(duì)于0≤u

從而,當(dāng)k=1時(shí),式(3.3)和式(3.4)成立.

現(xiàn)假設(shè)式(3.3)和式(3.4)對(duì)于k=l成立.對(duì)于k=l+1,根據(jù)算法3.1和式(3.2),我們有

最后,我們證明當(dāng)u

從而,當(dāng)k=l+1時(shí),式(3.3)和式(3.4)也成立,則對(duì)所有的0≤uv時(shí),式(3.3)和式(3.4)也成立.這樣對(duì)所有的uv,結(jié)論得證.

定理3.2記(t=1,2,3,4)為張量方程組(3.1)的解,則對(duì)于任給的初始張量Xt(0)(t=1,2,3,4),由算法3.1生成的序列{Rs(k)}和{Qt(k)}(s=1,2,3,4,t=1,2,3,4)滿足

證用數(shù)學(xué)歸納法.當(dāng)k=0時(shí),有

假設(shè)當(dāng)k=l(l ≥1)時(shí),式(3.5)成立.對(duì)于k=l+1,我們可得

從而當(dāng)k=l+1時(shí),式(3.5)成立.由數(shù)學(xué)歸納法可得,式(3.5)對(duì)所有的k=0,1,2,···都成立.□

定理3.3若實(shí)張量方程組(3.1)是相容的,對(duì)任意的初始張量[X1(0),X2(0),X3(0),X4(0)],在不計(jì)舍入誤差的情況下,其精確解可通過(guò)算法3.1在有限迭代步內(nèi)獲得.

通過(guò)算法3.1可得到Xt(4I)和Rs(4I),令R(k)=diag{R1(k)(1),R2(k)(1),R3(k)(1),R4(k)(1)},k=0,1,2,···,4I-1,其中Rs(k)(1)表示張量Rs(k)(s=1,2,3,4)沿模-1的展開(kāi)矩陣,根據(jù)定理3.1,有

則R(0),R(1),···,R(4I-1)可構(gòu)成如下線性空間的一組正交基,

根據(jù)定理3.1,有

故Rs(4I)=O.則[X1(4I),X2(4I),X3(4I),X4(4I)]是實(shí)張量方程組(3.1)的解.

4.極小范數(shù)解

引理4.1[13]若A ∈Rm×n,b ∈Rm,線性方程組Ax=b有解x?且x?∈R(AT),則x?是Ax=b的唯一極小范數(shù)解.

定理4.1若實(shí)張量方程組(3.1)是相容的,如果選擇初始張量為

由于線性方程組(2.4)和實(shí)張量方程組(3.1)等價(jià),根據(jù)引理4.1可得x?為(2.4)的唯一極小范數(shù)解,從而由算法3.1生成的解(t=1,2,3,4)是實(shí)張量方程組(3.1)的唯一極小Frobenius范數(shù)解.

5.數(shù)值實(shí)驗(yàn)

為了驗(yàn)證本文算法的可行性,本節(jié)通過(guò)三個(gè)例子來(lái)說(shuō)明所提出算法的效果.所有的程序都是基于Bader和Kolda[14]開(kāi)發(fā)的MATLAB張量工具箱,實(shí)驗(yàn)中,選取N=3,當(dāng)滿足

迭代終止.所有算法均在帶有MATLAB R2018a的個(gè)人PC機(jī)上運(yùn)行.此外,定義殘差范數(shù)

考慮初始張量Xt(0)=O(t=1,2,3,4),當(dāng)滿足式(5.1)時(shí)迭代終止.根據(jù)算法3.1,我們可獲得式(5.2)的極小Frobenius范數(shù)解,數(shù)值結(jié)果如圖1所示.

圖1 例5.1的數(shù)值結(jié)果

例5.2考慮四元數(shù)Sylvester張量方程(5.2),設(shè)

令精確解

從而可得式(5.2)的右端項(xiàng)張量C.在這個(gè)例子中,張量方程(5.2)是相容的,通過(guò)選擇初始張量Xt(0)=O(t=1,2,3,4),對(duì)于n=10時(shí),執(zhí)行算法3.1得到的收斂曲線如圖2所示.

圖2 例5.2的數(shù)值結(jié)果

下面的例子通過(guò)求解稀疏四元數(shù)Sylvester張量方程來(lái)檢測(cè)算法3.1的有效性.

例5.3考慮四元數(shù)Sylvester張量方程

其中系數(shù)矩陣滿足

與Ahmadi-Asl和Beik[8]工作稍有不同,我們考慮以下四種情況:

情況1A11=S1,A12=pivtol,A21=S2,A13=S3,A14=S4;

情況2A11=S1,A12=gre-343,A21=S2,A13=S3,A14=S4;

情況3A11=S1,A12=mesh3em5,A21=S2,A13=S3,A14=S4;

情況4A11=S1,A12=sphere3,A21=S2,A13=S3,A14=S4.

測(cè)試矩陣來(lái)自于Davis收集的矩陣[15].數(shù)據(jù)來(lái)源于HB(Harwell-Boeing) group.在這個(gè)例子中出現(xiàn)的這四種測(cè)試矩陣(pivtol,gre-343,mesh3em5和sphere3)經(jīng)常出現(xiàn)在統(tǒng)計(jì)、有向加權(quán)圖和結(jié)構(gòu)問(wèn)題中.這些矩陣的稀疏結(jié)構(gòu)如圖5.1-5.4所示,測(cè)試矩陣的性質(zhì)見(jiàn)表5.1,其中‘nz’表示非零項(xiàng)的個(gè)數(shù).需要注意的是,這里的稀疏矩陣的階數(shù)較大,因此四元數(shù)Sylvester張量方程(5.3)的解的維數(shù)非常大.令

表5.1 例5.3測(cè)試矩陣的性質(zhì)

圖5.1 第一個(gè)測(cè)試矩陣的稀疏結(jié)構(gòu)

圖5.2 第二個(gè)測(cè)試矩陣的稀疏結(jié)構(gòu)

圖5.3 第三個(gè)測(cè)試矩陣的稀疏結(jié)構(gòu)

圖5.4 第四個(gè)測(cè)試矩陣的稀疏結(jié)構(gòu)

其中對(duì)于情況1-4,n分別取為102,343,289和258作為其精確解.

通過(guò)選擇初始張量Xt(0)=O,應(yīng)用算法3.1去求解張量方程(5.3),當(dāng)滿足式(5.1)時(shí)迭代終止,相應(yīng)的收斂結(jié)果如圖5.5所示,結(jié)果證實(shí)了算法3.1的收斂性,可以看出殘差范數(shù)是顯著下降的.因此,算法3.1對(duì)于求解稀疏張量方程是有效的.

圖5.5 例5.3的收斂結(jié)果

6.結(jié)論

本文通過(guò)將四元數(shù)Sylvester張量方程等價(jià)地轉(zhuǎn)化為實(shí)數(shù)域上的張量方程組并引入線性算子,構(gòu)造了基于張量形式的修正共軛梯度算法求其等價(jià)形式.理論分析表明對(duì)任給的初始張量,由算法3.1生成的迭代序列在不計(jì)舍入誤差的情況下,可在有限迭代步內(nèi)收斂到方程組的精確解.進(jìn)一步,通過(guò)選擇特殊類型的初始張量,可獲得唯一極小Frobenius范數(shù)解.最后,數(shù)值例子證明了該算法的有效性.