解线性方程组的基本思想文档格式.docx
- 文档编号:22735605
- 上传时间:2023-02-05
- 格式:DOCX
- 页数:20
- 大小:95.60KB
解线性方程组的基本思想文档格式.docx
《解线性方程组的基本思想文档格式.docx》由会员分享,可在线阅读,更多相关《解线性方程组的基本思想文档格式.docx(20页珍藏版)》请在冰豆网上搜索。
这样可以大大提高运算速度。
命令
[L,U]=lu(A)
在matlab中可以编如下通用m文件:
在Matlab中建立M文件如下
%exp1.m
A;
b;
[L,U]=lu(A);
X=U\(L\b)
II)Cholesky分解
若A为对称正定矩阵,则Cholesky分解可将矩阵A分解成上三角矩阵和其转置的乘积,即:
其中R为上三角阵。
方程
A*X=b
变成
所以
%exp2.m
[R’,R]=chol(A);
X=R\(R’\b)
III)QR分解
对于任何长方矩阵A,都可以进行QR分解,其中Q为正交矩阵,R为上三角矩阵的初等变换形
式,即:
A=QR
变形成
QRX=b
所以
X=R\(Q\b)
上例中
[Q,R]=qr(A)
X=R\(Q\B)
%exp3.m
[Q,R]=qr(A);
X=R\(Q\b)
2.求线性齐次方程组的通解(A*X=0)
在Matlab中,函数null用来求解零空间,即满足A&
#8226;
X=0的解空间,实际上是求出解空
间的一组基(基础解系)。
%exp4.m
formatrat
%指定有理式格式输出
b=0;
r=rank(A);
bs=null(A,‘r’);
%一组基含(n-r)个列向量
%k,k,……,k
%X=k*bs(:
1)+k*bs(:
2)+……+k*bs(:
n-r)方程组的通解
pretty(X)
%让通解表达式更加精美
3
求非齐次线性方程组的通解(A*X=b)
非齐次线性方程组需要先判断方程组是否有解,若有解,再去求通解。
因此,步骤为:
第一步:
判断AX=b是否有解,(利用基本思路的第一条)
若有解则进行第二步
第二步:
求AX=b的一个特解
第三步:
求AX=0的通解
第四步:
AX=b的通解为:
AX=0的通解加上AX=b的一个特解。
clearall
%输入矩阵A,b
[m,n]=size(A);
R=rank(A);
B=[Ab];
Rr=rank(B);
formatrat
ifR==Rr&
R==n
%n为未知数的个数,判断是否有唯一解
x=A\b;
elseifR==Rr&
R<
n
%判断是否有无穷解
x=A\b
%求特解
C=null(A,r)
%求AX=0的基础解系,所得C为n-R列矩阵,这n-R列即为对
%应的基础解系
%这种情形方程组通解xx=k(p)*C(:
P)(p=1…n-R)
elseX=Nosolution!
%判断是否无解
end
第3章线性方程组的迭代解法
3.1实验目的
理解线性方程组计算机解法中的迭代解法的求解过程和特点,学习科学计算的方法和简单的编程技术。
3.2
概念与结论
1.
n阶线性方程组
如果未知量的个数为n,而且关于这些未知量x1,x2,…,xn的幂次都是一次的(线性的)那末,n个方程
a11x1+a12x2+…+a1nxn=b1
┆┆┆
(1)
an1x1+an2x2+…+annxn=bn
构成一个含n个未知量的线性方程组,称为n阶线性方程组。
其中,系数a11,…,a1n,a21,…,a2n,…,an1,…,ann和b1,…,bn都是给定的常数。
方程组
(1)也常用矩阵的形式表示,写为
Ax=b
其中,A是由系数按次序排列构成的一个n阶矩阵,称为方程组的系数矩阵,x和b都是n维向量,b称为方程组的右端向量。
2.n阶线性方程组的解
使方程组
(1)中每一个方程都成立的一组数x1*,x2*,…,xn*称为式
(1)的解,把它记为向量的形式,称为解向量.
3.向量范数的三种常用范数
4.矩阵的四种常用范数
5.谱半径
设nn阶矩阵A的特征值为i(i=1,2,3……n),则称
(A)=MAX|i|
1in
为矩阵A的谱半径.
矩阵范数与谱半径之间的关系为:
(A)||A||
6.严格(行)对角占优阵A
如果矩阵A=(aij)满足
n
|aii|>
|aij|i=1,2,……n,
j=1,ji
则称方阵A是严格(行)对角占优的.
7.收敛定理
对任意初始向量x(0)及任意右端向量g,由迭代x(k+1)=Bx(k)+g产生的迭代向量序列{x(k)}收敛的充要条件是谱半径(B)<
1
8.收敛判别条件
判别条件1:
若||B||<
1,则迭代x(k+1)=Bx(k)+g对任何初始向量x(0)都收敛.
判别条件2:
如果A为严格对角占优阵,则其Jacobi迭代和Seidel迭代对任何初始向量x(0)都收敛。
判别条件3:
如果A为对称正定阵,则其Seidel迭代对任何初始向量x(0)都收敛。
9.迭代法的误差估计
1,则对迭代格式x(k+1)=Bx(k)+g 有
3.3程序中Mathematica语句解释
a*matrix数a与矩阵matrix相乘
matrix1+matrix2矩阵matrix1和矩阵matrix2相加(注意矩阵的大小相同)
matrix1.matrix2矩阵matrix1和矩阵matrix2相乘(注意矩阵乘法的规则)
Transpose[matrix]求矩阵matrix转置
Inverse[matrix]求矩阵(方阵)matrix的逆
DiagonalMatrix[list]使用列表list中的元素生成一个对角矩阵.
IdentityMatrix[n]生成n阶单位矩阵
Max[x]求向量x中元素的最大值
3.4方法、程序、实验
解线性方程组的迭代法是将线性方程组Ax=b化为等价线性方程组
x=Bx+f
再由矩阵迭代格式
x(k+1)=Bx(k)+f
构造向量序列{x(k)}来求线性方程组解的。
如果得出的向量序列{x(k)}收敛至某个向量x*,则可得该向量x*就是所求方程组Ax=b的准确解.
线性方程组的迭代法主要有Jocobi迭代法、Seidel迭代法和超松弛(Sor)迭代法。
1.Jocobi迭代法
1)Jocobi迭代法的构造过程
假设aii0,依次在第i个方程解出xi,i=1,2,,n并令
cij=-aij/aii(ij),gi=bi/aii
就得到如下Jocobi迭代格式:
x1(k+1)=c12x2(k)+c13x3(k)++c1nxn(k)+g1
x2(k+1)=c21x1(k)+c23x3(k)++c2nxn(k)+g2
。
。
xn(k+1)=cn1x1(k)+cn2x2(k)++cn(n-1)xn-1(k)+gn
若令
则有Jocobi迭代的矩阵格式:
x(k+1)=BJx(k)+gJ
BJ称为Jocobi迭代矩阵。
Jocobi迭代可以写成如下紧凑格式:
在给定初始迭代向量x(0)后就可以进行Jocobi迭代求解了。
2)Jacobi迭代算法
1.输入变量个数n、初值向量x(0)、迭代精度eps、系数矩阵A、常数项b和迭代最大次数nmax
2Fori=1,2,…,n
2.1如果|aii|<
eps1,则输出“迭代失败”提示并终止
3.BjE-D-1A
4.gjD-1b
5.Fork=1,2,…,nmax
5.1xBj.x0+gj
5.2如果||x-x0||<
eps,输出解向量x,终止;
否则x(0)x
6.如果||x-x0||>
eps,输出迭代失败,终止。
3)Jacobi迭代法程序
Clear[a,b,x];
nmax=500;
n=Input[“线性方程组阶数n=”];
a=Input["
系数矩阵A="
];
b=Input["
常数项b="
x0=Input["
输入迭代初值向量x0"
eps1=0.000001;
eps=Input["
输入精度控制eps="
Do[If[Abs[a[[i,i]]]<
eps1,t1=1;
Return[],t1=0],{i,1,n}];
If[t1==1,
Print["
Jacobi迭代法失效"
],
d=DiagonalMatrix[Table[a[[i,i]],{i,1,n}]];
d1=Inverse[d];
bj=IdentityMatrix[n]-d1.a;
gj=d1.b;
Do[x=bj.x0+gj;
err=Max[Abs[x-x0]];
x="
x//N,"
i="
i,"
err="
err//N];
If[N[err]<
eps,Break[],x0=x],
{i,1,nmax}];
If[err>
=eps,Print["
迭代失败"
]]
]
说明本程序用于求线性方程组Ax=b的解。
程序执行后,先通过键盘输入线性方程组阶数n、系数矩阵A、常数项b、迭代初值向量x0和输入精度控制eps,程序即可给出每次迭代的次数和对应的迭代向量序列x(k),其中最后输出的结果即为所求的根。
如果迭代超出500次还没有求出满足精度的根则输出迭代失败提示,如果出现主对角线元素aii=0给出Jacobi迭代法失效提示。
程序中变量说明
x0:
存放初始向量和迭代过程中的向量x(k)
x:
存放迭代过程中的向量x(k+1)
nmax:
存放迭代允许的最大次数
err:
存放误差||x-x0||
t1:
临时变量
注:
迭代最大次数可以修改为其他数字。
4)例题与实验
例1.用Jacobi迭代法解如下线性方程组
5x1+2x2+x3=-12
-x1+4x2+2x3=20
2x1-3x2+10x3=3
要求误差||x(k+1)-x(k)||<
10-4,并用取不同初值的方法实验观察迭代收敛的情况。
解:
执行Jacobi迭代法程序后在输入的四个窗口中按提示分别输入:
3、{{5,2,1},{-1,4,2},{2,-3,10}}、{-12,20,3}、{0,0,0}、0.0001
每次输入后用鼠标点击窗口的“OK”按扭,得如下输出结果:
x={-2.4,5.,0.3}i=1err=5.
x={-4.46,4.25,2.28}i=2err=2.06
x={-4.556,2.745,2.467}i=3err=1.505
x={-3.9914,2.6275,2.0347}i=4err=0.5646
x={-3.85794,2.9848,1.88653}i=5err=0.3573
x={-3.97123,3.09225,1.96703}i=6err=0.113286
x={-4.03031,3.02368,2.02192}i=7err=0.0685705
x={-4.01386,2.98146,2.01316}i=8err=0.042216
x={-3.99522,2.98995,1.99721}i=9err=0.0186374
x={-3.99542,3.00259,1.99603}i=10err=0.0126367
x={-4.00024,3.00313,1.99986}i=11err=0.0048186
x={-4.00122,3.00001,2.00099}i=12err=0.00312067
x={-4.0002,2.9992,2.00025}i=13err=0.00102319
x={-3.99973,2.99983,1.9998}i=14err=0.000625697
x={-3.99989,3.00017,1.99989}i=15err=0.00034136
x={-4.00005,3.00008,2.00003}i=16err=0.000155236
x={-4.00004,2.99997,2.00003}i=17err=0.000106099
x={-4.,2.99997,2.}i=18err=0.0000414468
此结果说明迭代18次,求得误差为err=0.0000414468的近似解,最后显示的近似解向量为x={-4.,2.99997,2.},它表示所求解为
x1=-4,x2=2.99997,x3=2。
本题的准确解为
x1=-4,x2=3,x3=2。
如果将如上输入的初值改为{21,-18,30},执行Jacobi迭代法程序后得如下输出结果:
x={-1.2,-4.75,-9.3}i=1err=39.3
x={1.36,9.35,-0.885}i=2err=14.1
x={-5.963,5.7825,2.833}i=3err=7.323
x={-5.2796,2.09275,3.22735}i=4err=3.68975
x={-3.88257,2.06642,1.98375}i=5err=1.39703
x={-3.62332,3.03749,1.69644}i=6err=0.97106
x={-3.95428,3.24595,1.93591}i=7err=0.330963
x={-4.08556,3.04347,2.06464}i=8err=0.202475
x={-4.03032,2.94629,2.03015}i=9err=0.0971858
x={-3.98455,2.97734,1.98995}i=10err=0.0457716
x={-3.98893,3.00889,1.99011}i=11err=0.0315451
x={-4.00158,3.00771,2.00045}i=12err=0.0126504
x={-4.00318,2.99938,2.00263}i=13err=0.00833246
x={-4.00028,2.99789,2.00045}i=14err=0.00289753
x={-3.99925,2.99971,1.99942}i=15err=0.0018145
x={-3.99977,3.00048,1.99976}i=16err=0.000770763
x={-4.00014,3.00018,2.0001}i=17err=0.000375926
x={-4.00009,2.99992,2.00008}i=18err=0.000261658
x={-3.99998,2.99994,1.99999}i=19err=0.000107579
x={-3.99997,3.00001,1.99998}i=20err=0.0000714045
从计算结果可以看到虽然所取的初值不同,但迭代总是收敛相同的结果。
考察本题系数矩阵可以发现它是严格对角占优矩阵,由收敛判别法,Jacobi迭代对任意初值都是收敛的,我们通过实际计算验证了这个结果。
例2.通过Jacobi迭代序列观察用Jacobi迭代法解如下线性方程组
x1+2x2+2x3=-12
-x1+4x2+x3=20
2x1-3x2+x3=3
的收敛性。
3、{{1,2,2},{-1,4,1},{2,-3,1}}、{-12,20,3}、{0,0,0}、0.0001
x={-12.,5.,3.}i=1err=12.
x={-28.,1.25,42.}i=2err=39.
x={-98.5,-12.5,62.75}i=3err=70.5
x={-112.5,-35.3125,162.5}i=4err=99.75
…………………………………………….
{246.719,-120.176,696.719}i=8err=763.5
x={-1165.09,-107.5,-850.965}i=9err=1547.68
x={1904.93,-73.5303,2010.67}i=10err=3070.02
x={-3886.28,-21.4355,-4027.45}i=11err=6038.12
x={8085.77,40.2917,7711.26}i=12err=11972.1
x={-15515.1,98.6279,-16047.7}i=13err=23758.9
x={31886.1,138.141,31329.1}i=14err=47401.2
x={-62946.5,144.247,-63354.7}i=15err=94832.5
x={126409.,107.068,126329.}i=16err=189683.
x={-252883.,25.0775,-252494.}i=17err=379292.
x={504925.,-92.4305,505845.}i=18err=758339.
…………………………………….
通过观察迭代过程中的误差是不断变大的特点可以知道本题的Jacobi迭代序列是不收敛的,因此,本题线性方程组不能用Jacobi迭代法求解。
这个实验说明并不是每个线性方程组都能用Jacobi迭代法求解。
2.
Seidel迭代
1)Seidel迭代的构造过程
为了加快收敛速度,同时节省计算机的内存,对Jocobi迭代作如下的改进:
每算出一个分量的近似值,立即用到下一个分量的计算中去,即用迭代格式:
x1(k+1)=c12x2(k)+c13x3(k)++c1nxn(k)+g1
x2(k+1)=c21x1(k+1)+c23x3(k)++c2nxn(k)+g2
xn(k+1)=cn1x1(k+1)+cn2x2(k+1)++cn(n-1)xn-1(k+1)+gn
如上迭代可以写成如下紧凑格式:
这样所得的迭代法就称为Seidel迭代法。
利用Ax=b及A=L+D+U,其中D为对角矩阵,L,U分别为严格下,上三角矩阵.则有Seidel迭代法的矩阵形式为
x(k+1)=Bsx(k)+gs其中:
Bs=-(L+D)-1U,gs=D-1b
Seidel迭代矩阵为Bs=-(L+D)-1U。
在给定初始迭代向量x(0)后就可以进行Seidel迭代求解了。
Jacobi迭代和Seidel迭代格式可表述为统一形式:
2)Seidel迭代算法
2.Fori=1,2,…,n
3.Fori=1,2,…,nmax
3.1Fori=1,2,…,n
3.2如果||x(k+1)-x(k)||<
否则x(k)x(k+1)
4.如果||x-x0||<
3)Seidel迭代法程序
x=x0;
Do[If[Abs[a[[i,i]]]<
Seidel迭代法失效"
Do[Do[u1=Sum[a[[i,j]]*x[[j]],{j,1,i-1}];
u2=Sum[a[[i,j]]*x0[[j]],{j,i+1,n}];
x[[i]]=(b[[i]]-u1-u2)/a[[i,i]],
{i,1,n}];
k="
k,"
{k,1,50}];
]]
程序执行后,先通过键盘输入线性方程组阶数n、系数矩阵A、常数项b、迭代初值向量x0和输入精度控制eps,程序即可给出每次迭代的次数和对应的迭代向量序列x(K),其中最后输出的结果即为所求的根。
存放初始向量
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 线性方程组 基本 思想