一、理论知识
解线性方程组的直接方法–低阶稠密矩阵三角分解法–Courant分解学习matlab实现
二、代码实现
function [L_matrix,U_matrix,y_matrix,x_matrix] = LU_separetion(A_matrix, B_matirx)
% LU系数矩阵分解
% AX=B A=LU LUX=B
% -> LY=B ->YX=B
% UX=Y
% inputs:
% A_matrix:输入的系数矩阵,尺寸为[n,n]
% B_matrix:输入的乘积矩阵,尺寸为[n,1]
% outputs:
% L_matrix:下三角阵,尺寸为[n,n]
% U_matrix:上三角阵,尺寸为[n,n]
% y_matrix:中间矩阵,尺寸为[n,1]
% x_matrix:结果矩阵,尺寸为[n,1]
%% 第一步:初始化
% 获取n值
[row_a, col_a] = size(A_matrix);
% 初始化上三角阵的第一行
for r = 1:col_a
U_matrix(1,r) = A_matrix(1,r);
end % for-1
% 初始化下三角阵的第一列
L_matrix(1,1) = 1;
for i = 2:row_a
L_matrix(i,1) = A_matrix(i,1)/U_matrix(1,1);
end
%% 第二步:求L和U
for r = 2:row_a % U当前行,L当前列
for i = r:col_a % U当前列,L当前行
%按每一行的每一列逐个求出每一个u(i,r)
temp_sum = 0; % 上面的乘积
for k = 1:r-1
temp_sum = temp_sum + L_matrix(r,k)*U_matrix(k,i);
end
U_matrix(r,i) = A_matrix(r,i) - temp_sum;
temp_sum_1 = 0; %左面的乘积
for k = 1:r-1
temp_sum_1 = temp_sum_1 + L_matrix(i,k)*U_matrix(k,r);
end
L_matrix(i,r) = (A_matrix(i,r) - temp_sum_1)/U_matrix(r,r);
end
end
%% 第三步:回代计算y LY=B
x_matrix = zeros(row_a,1);
% 计算中间矩阵Y
y_matrix(1,1) = B_matirx(1,1);
for i = 2:row_a
temp_sum_2 = 0; %上面的乘积
for r = 1:i-1
temp_sum_2 = temp_sum_2 + L_matrix(i,r)*y_matrix(r,1);
end
y_matrix(i,1) = B_matirx(i) - temp_sum_2;
end
%% 第四步:回代计算x UX=Y
% 计算结果矩阵X
x_matrix(row_a,1) = y_matrix(row_a,1)/U_matrix(row_a,col_a);
for i=row_a-1:-1:1
temp_sum_3 = 0; %下面的乘积
for r = i+1:row_a
temp_sum_3 = temp_sum_3 + U_matrix(i,r)*x_matrix(r,1);
end
x_matrix(i,1) = (y_matrix(i,1) - temp_sum_3)/U_matrix(i,i);
end
end
%% 测试案例
% A =[1 2 3;
% 2 5 2;
% 3 1 5];
% B=[14;18;20];
% 结果:X=[1,2,3]'