Matlab广义追赶法（Thomas法）

矩阵分解方法（广义追赶法）

对于一般形式的线性方程组

Ax=bA x=bAx=b

将系数矩阵A分解成下三角阵L和上三角阵U

A=LU1A=LU_1A=LU1

化归为两个三角方程组

Ly=bLy=bLy=b，Ux=yUx=yUx=y

对于3阶矩阵：

A=[a11a12a13a11a12a13a11a12a13]=[l1100l21l220l31l32l33][1u12u1301u23001]=LU1A= \begin{bmatrix} a_{11}& a_{12} & a_{13} \\ a_{11} & a_{12} & a_{13} \\ a_{11} & a_{12} & a_{13} \end{bmatrix}\\= \begin{bmatrix} l_{11} &0 & 0 \\ l_{21} & l_{22} & 0 \\ l_{31} & l_{32} & l_{33} \end{bmatrix}\begin{bmatrix} 1 &u_{12} & u_{13} \\ 0 & 1 & u_{23} \\ 0 & 0 & 1 \end{bmatrix}=LU_1A=⎣⎡a11a11a11a12a12a12a13a13a13⎦⎤=⎣⎡l11l21l310l22l3200l33⎦⎤⎣⎡100u1210u13u231⎦⎤=LU1

展开

l11=a11,u12=a12/l11,u13=a13/l11l21=a21,l22=a22−l21u12,u23=(a23−l21u13)/l22l31=a31,l32=a32−(l31u12),l33=a33−l31u13−l32u23\begin{aligned} l_{11}=a_{11}, &\quad u_{12}=a_{12}/l_{11},\qquad\qquad u_{13}=a_{13}/l_{11}\\ l_{21}=a_{21}, &\quad l_{22}=a_{22}-l_{21}u_{12},\qquad u_{23}=(a_{23}-l_{21}u_{13})/l_{22}\\ l_{31}=a_{31},& \quad l_{32}=a_{32}-(l_{31}u_{12}),\quad l_{33}=a_{33}-l_{31}u_{13}-l_{32}u_{23} \end{aligned}l11=a11,l21=a21,l31=a31,u12=a12/l11,u13=a13/l11l22=a22−l21u12,u23=(a23−l21u13)/l22l32=a32−(l31u12),l33=a33−l31u13−l32u23

求解环节

1.预处理

实现矩阵分解A=LU1A=LU_1A=LU1:对i=1,2,...,ni=1,2,...,ni=1,2,...,n计算

li,j=ai,j−∑k=1j−1li,kuk,j,j=1,2,...,il_{i,j}=a_{i,j}-\sum_{k=1}^{j-1}l_{i,k}u_{k,j},\qquad j=1,2,...,ili,j=ai,j−∑k=1j−1li,kuk,j,j=1,2,...,i

ui,j=(ai,j−∑k=1j−2li,kuk,j)/li,i,j=i+1,i+2,...,nu_{i,j}=(a_{i,j}-\sum_{k=1}^{j-2}l_{i,k}u^{k,j})/l_{i,i},\qquad j=i+1,i+2,...,nui,j=(ai,j−∑k=1j−2li,kuk,j)/li,i,j=i+1,i+2,...,n

2.追的过程

解下三角方程组Ly=bLy=bLy=b即

∑j=1ili,jyj=bi,i=1,2,...,n\sum_{j=1}^il_{i,j}y_{j}=b_{i},\qquad i=1,2,...,n∑j=1ili,jyj=bi,i=1,2,...,n

回代公式

yi=（bi−∑j=1i−1li,jyj）/li,i,i=1,2,...,ny_i=（b_{i}-\sum_{j=1}^{i-1}l_{i,j}y_{j}）/l_{i,i},\qquad i=1,2,...,nyi=（bi−∑j=1i−1li,jyj）/li,i,i=1,2,...,n

3.赶的过程

解单位上三角方程U1x=yU_{1}x=yU1x=y即

xi+∑j=i+1nui,jxj=yj,i=1,2,...,nx_{i}+\sum_{j=i+1}^{n}u_{i,j}x_{j}=y_{j},\qquad i=1,2,...,nxi+∑j=i+1nui,jxj=yj,i=1,2,...,n

回代公式为

xi=yi−∑j=i+1nui,jxj,i=n,n−1,...,1x_{i}=y_{i}-\sum_{j=i+1}^{n}u_{i,j}x_{j},\qquad i=n,n-1,...,1xi=yi−∑j=i+1nui,jxj,i=n,n−1,...,1

Matlab实现：

%% 追赶法（Thomas法）
A=[1  1  1;2 -1 -2;3  1 -4];
b=[17;1;3];
m=length(A);%矩阵大小
L=zeros(m);%下三角方程
U=zeros(m)+eye(m);%上三角方程
y=zeros(m,1);%下三角方程解初始化
x=zeros(m,1);%解向量初始化
%% 1.预处理
for i=1:mfor j=1:i%计算矩阵L的各个元素，其中列坐标j=1,2,..,iK=0;if j<=1L(i,j)=A(i,j)-K;elsefor k=1:j-1K=L(i,k)*U(k,j)+K;endL(i,j)=A(i,j)-K;end%计算矩阵U的各元素，其中列坐标jj=i+1,i+2,...,nKK=0;if i<m;for jj=i+1:mif i==1;U(i,jj)=(A(i,jj)-KK)./L(i,i);elsefor k=1:i-2KK=L(i,k)*U(k,i+1)+KK;endU(i,jj)=(A(i,i+1)-KK)./L(i,i);endendelseendend
end
%% 追的过程 解Ly=b
for i=1:mu=0;j=1;if j>i-1y(i,1)=(b(i,1)-u)./L(i,i);elsefor j=1:iu=L(i,j)*y(j,1)+u;endy(i,1)=(b(i,1)-u)./L(i,i);end
end
%% 赶的过程 解Ux=y
for k=1:mi=m-k+1;jj=i+1;u=0;if jj>mx(i,1)=y(i,1)-u;elsefor j=1:mu=U(i,j)*x(j,1)+u;endx(i,1)=y(i,1)-u;end
end
%% 显示结果
disp('解向量为')
x
disp('下三角阵为')
L
disp('上三角阵为')
U

结果：
解向量为

x =

4.4706
7.9412
4.5882

下三角阵为

L =

1.0000         0         0
2.0000   -3.0000         0
3.0000   -2.0000   -5.6667

上三角阵为

U =

1.0000    1.0000    1.00000    1.0000    0.66670         0    1.0000

可以写成函数形式:function (xxx,LLL,U1U_{1}U1)=Thomas_Mathod(AAA,bbb)

输入为AAA,bbb

输出xxx,LLL,U1U_{1}U1

Matlab广义追赶法（Thomas法）相关推荐

广义最小二乘法的基本思想是什么_解决异方差问题的方法可行广义最小二乘法fgls法.ppt...
解决异方差问题的方法可行广义最小二乘法fgls法第五章模型的建立与估计中的问题及对策;本章内容第一节误设定第二节多重共线性第三节异方差性第四节自相关; OLS估计量令人满意的性质,是根据一 ...
matlab 坐标不用科学计数法,matlab不用科学计数法
『壹』 matlab中怎么才能不是科学计数法表示结果.比如1.0e+003 * 2.7581,怎么使它显示为2758.1谢谢了,很急啊 format long (小数位14) 或 format sho ...
MATLAB之牛顿下山法
MATLAB之牛顿下山法算法原理 matlab程序算法原理上一篇博客,我介绍了牛顿法迭代法,接下来我就们接着讲解一下什么是牛顿下山法. 一.迭代公式在牛顿迭代过程中,若满足单调性|f(x(k+ ...
matlab曲面拟合的算法,用Matlab 实现移动曲面拟合法生成DEM
用Matlab实现移动曲面拟合法生成DEM 杜玉军 (武汉大学测绘工程0408班 200431610007 武汉 430079) 摘要:移动曲面拟合法是DEM格网点内插常用的一种方法,利用Matl ...
【雷达通信】基于matlab距离角度解耦法MIMO-OFDM雷达波束形成【含Matlab源码 2208期】
⛄一.获取代码方式获取代码方式1: 完整代码已上传我的资源:[雷达通信]基于matlab距离角度解耦法MIMO-OFDM雷达波束形成[含Matlab源码 2208期] 点击上面蓝色字体,直接付费下载 ...
基于Matlab的结点电压法及相关定理验证的人机交互界面！
基于Matlab的结点电压法及相关定理验证的人机交互界面! 前言本篇主要针于无储能元件的电路分析由结点电压法求得相关结点的电压进而验证戴维宁定理.叠加定理等提示:以下是本篇文章正文内容,下面案 ...
中心差分法matlab实现,动力学系统时域响应计算的六种方法Matlab源程序(Newmark,Houbolt法,中心差分法)...
动力学系统时域响应计算的六种方法Matlab源程序(Newmark,Houbolt法,中心差分法).直接使用,无需再修改. Newmark法Matlab源程序 function [acc,vel,ds ...
基于matlab的频率特性测试仪,基于Matlab控制系统频率特性分析法
基于Matlab控制系统频率特性分析法基于Matlab控制系统频率特性分析法本文主要介绍了基于Matlab控制系统的频率特性分析方法.频域稳定性判据以及开环频域性能分析,并获得频率响应曲线等.通过 ...
matlab双线性z变换法设计数字低通滤波器
matlab双线性z变换法设计数字低通滤波器双线性z变换法利用了正切函数的非线性特点,将整个jΩ轴压缩到了单位圆的一周上. 低通: clear; close all; clc; fp=100;fs= ...

Matlab广义追赶法（Thomas法）

Matlab广义追赶法（Thomas法）相关推荐

最新文章

热门文章