在數值分析求解線性方程組當中常常會用到高斯分解,對線性方程組的係數矩陣a
aa進行 l
ll、u
uu分解, 即:a=l
∗ua = l*u
a=l∗
u其中 l
ll 為主對角線元素均為 1的單位下三角矩陣, u
uu 為 上三角矩陣,可以將原線性方程組 ax=
bax=b
ax=b
,轉換成兩個三角矩陣方程組,即:令y=u
xy = ux
y=ux
下三角線性方程組 ly=
bly=b
ly=b
和上三角線性方程組 ux=
yux=y
ux=y
,對 ly=
bly=b
ly=b
和 ux=y
ux=y
ux=y
分別採用前代法、回代法實現快速求解的目的。但問題在於如何得到 l
ll、u
uu 矩陣,採用最簡單的待定係數法,利用矩陣結構進行運算:
問題在於:求解出來的矩陣不太對,即演算法出現了問題,**如下:
% 待定係數法求解高斯變換矩陣l 使得上三角矩陣u、下三角矩陣l ,a = l*u
clear
a =randi([
1,10]
,6,6
);n =length
(a);
l =zeros
(n);u =
zeros
(n);
i =1
;j =2;
u(1,i:n)=a
(1,i:n)
;% 首先求解 u的第一行元素 並採用向量化運算
l(j:n,1)
=a(j:n,1)
/u(1
,1);
% 求解 l 的第一列元素
l = l+
eye(n)
;for k =2:n
for j1 = k:n
sum =0;
%暫時暫存器
for i1 =
1:k-
1 sum = sum +
l(k,i1)*u
(i1,j1)
; end
u(k,j1)=a
(k,j1)
-sum;
endfor i2= k+1:n
sum =0;
%暫時暫存器
for j2 =
1:k-
1 sum = sum +
l(i2,j2)*u
(j2,k)
; end
l(i2,k)=(
a(k,j2)
-sum)/u
(k,k)
; end
end
參考文獻: 東南大學《數值分析》孫志忠 70頁
也有發現正確的解法,**如下:
l =
zeros
(n,n)
;% initial the n*n matrix m zeros
for i =1:n
l(i,i)=1
;end % let the diagonal elements be 1
for j =
1: n-1if
abs(
a(j,j)
)error
('zero pivot encountered');
endfor i = j+
1: n
mult =
a(i,j)/a
(j,j);l
(i,j)
= mult;
for k = j:n
a(i,k)=a
(i,k)
- mult*
a(j,k)
; end
endenddisp
(' l=');
disp
(l);
disp
(' u=');
disp
(a);
disp
(' lu=');
disp
(l*a)
;
靜待有緣人(隔壁老王),答疑解惑 ! 數值分析 python LU分解法
lu分解在本質上是高斯消元法的一種表達形式。實質上是將a通過初等行變換變成乙個上三角矩陣,其變換矩陣就是乙個單位下三角矩陣。這正是所謂的杜爾里特演算法 doolittle algorithm 從下至上地對矩陣a做初等行變換,將對角線左下方的元素變成零,然後再證明這些行變換的效果等同於左乘一系列單位下...
M 數值分解
description 對乙個自然數n 1 n 50 n可以分解成若干個數字 數字可以是1,2,3,9 之和,問題是如何分解能使這些數字的乘積最大。input 輸入資料有多組,每組佔一行,每行包含乙個自然數n 1 n 50 輸入檔案直到eof為止!output 對每組輸入,輸出有2行。第一行是n分解...
M 數值分解
description 對乙個自然數n 1 n 50 n可以分解成若干個數字 數字可以是1,2,3,9 之和,問題是如何分解能使這些數字的乘積最大。input 輸入資料有多組,每組佔一行,每行包含乙個自然數n 1 n 50 輸入檔案直到eof為止!output 對每組輸入,輸出有2行。第一行是n分解...