沈家銘

(四川大學 數(shù)學學院， 四川 成都 610065)

一種判斷LR分解是否存在的優(yōu)化算法

沈家銘

(四川大學 數(shù)學學院， 四川 成都 610065)

使用定理直接判斷方陣A是否存在LR分解的計算難度系數(shù)為O(n3)，文中在此基礎上提出了一個改進算法,將計算難度降為O(n2)。給出了設計思路及推導方法，并用Matlab實現(xiàn)了兩種算法，通過例題驗證了計算效率的提升。

LR分解; 方陣; 算法; Matlab

0 引 言

在有應用背景的數(shù)學問題中，很多求解問題最終都可歸結為線性方程組的求解。因此，解線性方程組在科學與工程計算中有著十分重要的作用。對其系數(shù)矩陣進行三角分解，是一種解線性方程組的有效方法，LR分解是求三角分解的一個強有力的算法,受到了人們的極大關注。在不進行行列置換的前提下，為了判斷方陣A的LR分解是否存在，是進行LR分解的前提[1-4]。文中給出了一種快速判斷方陣LR分解是否存在，同時計算LR分解的算法,大大提升了計算速度。

1 問題描述

定義：LR分解：

即把矩陣A分解為一個下三角矩陣L和一個上三角矩陣R

而如何求證方陣A是否存在唯一的LR分解，有以下定理:

定理1 設A為n階方陣，A的k階主子式記為

則A的LR分解存在唯一的充分必要條件是

記

2 計算方法

為此，需要計算n-1個矩陣行列式的值。而容易證明，矩陣行列式的值就等于矩陣LR分解得到的R矩陣的對角元的乘積。

同時， d1≠0時可以得到A2的LR分解存在，由此算得d2，而若d2≠0，又已知d1≠0，所以A3的LR分解存在，以此類推，可以得到dk(k=1,2,…,n-1)的值。

由這個思路，可以得到算法1的計算步驟如下：

1)對Ak用Doolittle算法做LR分解；

2)由U計算dk；

3)若dk≠0，則k=k+1，否則返回錯誤，不能進行LR分解。按此思路即可得到算法1的Matlab程序代碼如下：

function[L,R,pp] =LR_P(A)

tic;

[n,m]=size(A);

R=zeros(n);

L=eye(n);

suM=0;

R(1,1)=A(1,1);

pan=R(1,1);

pp=0;

ifpan==0;

pp=-1;

return

end

fori=2:n;

forj=1:i-1;

fork=1:j-1

suM=suM+L(j,k)*R(k,i);

end

R(j,i)=A(j,i)-suM;

suM=0;

fork=1:j-1

suM=suM+R(k,j)*L(i,k);

end

L(i,j)=(A(i,j)-suM)/R(j,j);

end

R(i,i)=A(i,i)-L(i,1:i-1)*R(1:i-1,i);

pan=pan*R(i,i);

ifpan==0&&i～=n;

pp=-1;

return;

end

end

toc;

end

Doolittle分解算法代碼:

function[L,R]=LR_L(A)

[n,m]=size(A);

R=zeros(n);

L=eye(n);suM=0;

R(1,1)=A(1,1);

fori=2:n;

forj=1:i-1;

fork=1:j-1

suM=suM+L(j,k)*R(k,i);

end

R(j,i)=A(j,i)-suM;

suM=0;

fork=1:j-1

suM=suM+R(k,j)*L(i,k);

end

L(i,j)=(A(i,j)-suM)/R(j,j);

end

R(i,i)=A(i,i)-L(i,1:i-1)*R(1:i-1,i);

end

end

3 計算方法的改進

n次LR分解的計算量非常大，是一個十分耗費時間的過程，其時間難度為O(n3)，為此，尋找是否有優(yōu)化算法使得計算量減少[5-7]。

通過上述算法可以知道，每次LR分解的矩陣之間是有聯(lián)系的。它們之間的關系是：前一個矩陣增加一行一列就是后一個矩陣[8]，故進行嘗試，已知矩陣A的LR分解式，是否能通過LR計算在末尾增加一行一列的矩陣A*的LR分解式L*，R*。

3.1改進算法的推導

設矩陣A的LR分解式已知，不妨設增加一行一列以后的矩陣為：

其中

使用待定系數(shù)法，設

其中

則有

所以

即

可解得

然后可以得到

則通過A的LR分解式得到了在末尾添加一行一列A′的LR分解式A′=L′R′。

通過這種算法，使得在判斷A是否存在LR分解式，以及同時計算A分解式的LR的時間難度從O(n3)變?yōu)榱薕(n2)。

觀察解出來的β,γ的形式可以發(fā)現(xiàn)，每次對A加一行一列以后做LR分解，并不影響之前的分解結果，而且在算β,γ的時候推出來的遞推公式與Doolittle算法中對LR分解的算法遞推公式形式幾乎一致。通過比較得出了以下結論:改進的算法(以上)與原算法區(qū)別只是計算的次序不同。

改進的算法是從主對角元開始，一圈一圈向外計算，也就是說a11開始，計算r11,然后計算r11周邊的元素r12,r22,l21,…。

每算完一圈，利用∏rii是否為0判斷能否進行LR分解，若進行到某一步時，該式為0，則該矩陣不能進行LR分解，反之，則可以。

3.2改進算法的實現(xiàn)

我們得到如下改進的LR分解算法,稱之為算法2：

算法2的Matlab代碼如下：

function [L,R,pp] = LR_P2(A)

tic;

p=1;

[n,～]=size(A);

for i=1:n

if p==0

pp=-1;

break;

end

p=1;

A_=A(1:i,1:i);

[L,R]=LR_L(A_);//Dolittle算法，for j=1:i

p=p*R(i,i);

end

pp=0;

toc;

end

end

4 計算效率驗證

使用算法1和算法2的代碼做了一個簡單的例題來驗證計算，檢驗計算速度是否有提高。

對一個1 024×1 024的矩陣進行LR分解，矩陣來源是一張1 024×768的點陣圖片。

在平臺為intel 2640qm,內存8 GB，Matlab環(huán)境下，統(tǒng)一計算平臺下對1 024×1 024的矩陣用兩個算法進行對比。

算法1：Elapsed time is 2 331 seconds。

算法2：Elapsed time is 12 seconds。

可見效率提升了194.25倍。

5 結 語

在定理的基礎上設計實現(xiàn)了一種高效的判斷LR分解是否存在以及計算LR分解式的算法。算例測試證明了算法可行性和效率的提高。

[1] 高惠璇.統(tǒng)計計算[M].北京：北京大學出版社,1997.

[2] 徐曉飛,曹祥玉,姚旭,等.一種基于Doolittle LU分解的線性方程組并行求解方法[J].電子與信息學報,2010(8):2019-2022.

[3] G H Golub. Matrix computation 4th edition[M].北京：人民郵電出版社,2014.

[4] 王群英.矩陣分解方法的探究[J].長春工業(yè)大學學報:自然科學版，2011，32(1):95-101.

[5] 戴華.矩陣論[M].北京：科學出版社，2001：110-145.

[6] 徐泰燕，郝玉龍.非負矩陣分解及其應用現(xiàn)狀分析[J].武漢工業(yè)學院學報，2010(1):110-115.

[7] 黃麗嫦.矩陣的LU并行遞歸分解算法的設計研究[J].科學技術與工程，2012(5)：3626-3629.

[8] 周康，陳金.基于部分基變量的LP問題矩陣算法[J].運籌學學報，2012(6)：121-126.

AnoptimizationalgorithmtodeterminetheexistenceofLRdecomposition

SHENJia-ming

(SchoolofMathematics,SichuanUniversity,Chengdu610065,China)

SincethecalculatingdifficultycoefficientofLRdecompositionisO(n3),byusingthetheoremtojudgewhetherphalanxAexists,hereweputforwardanimprovedalgorithmtoreducethecalculatingdifficultycoefficienttoO(n2).Designandderivationmethodaregiven,andtwoalgorithmsareimplementedwithmatlab.Thecalculationefficiencyimprovementisverifiedwithexamples.

LRdecomposition;phalanx;algorithm;Matlab.

2014-07-18

沈家銘(1993-)，男，漢族，云南昆明人，主要從事統(tǒng)計計算方向研究,E-mail:493879714@qq.com.

O 151.21

A

1674-1374(2014)05-0593-04