线性方程组的数值解法
- 格式:ppt
- 大小:756.00 KB
- 文档页数:142


线性方程组数值方法Ax ,必须满足系数矩阵的对角元素非零的条件。
1、回代:用回代法求解上三角线性方程组bfunction x=backsub(A,b)%Input -A is an nxn upper-triangular nonsingular matrix% -b is an nx1 matrix%Output-x is the solution to the linear system Ax=b%Find the dimension of b and initialize xn=length(b);x=zeros(n,1);x(n)=b(n)/A(n,n);for k=n-1:-1:1x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k);end例:A=[1 2 1 4;0 -4 2 -5;0 0 -5 -7.5;0 0 0 -9];b=[13;2;-35;-18];x=backsub(A,b)2、高斯消去法(不选主元素)function x=gauss(A,b)%Input -A is an nxn nonsingular matrix% -b is an nx1 matrix%Output-x is an nx1 matrix containing the solution to Ax=b%Initialize x[n,n]=size(A);x=zeros(n,1);%Form the augmented matrix:Aug=[A|b]Aug=[A b];for k=1:n-1%Elimination process for column kfor i=k+1:nm=Aug(i,k)/Aug(k,k);Aug(i,k:n+1)=Aug(i,k:n+1)-m*Aug(k,k:n+1);endend%Output Elimination matrixAug%Back SubstitutionA=Aug(1:n,1:n);b=Aug(1:n,n+1);x(n)=b(n)/A(n,n);for k=n-1:-1:1x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k);end例: A=[4 -9 2;2 -4 6;-1 1 -3];b=[5;3;-4]; x=gauss(A,b)3、列主元消去法function x=gauss_column(A,b)%Input -A is an nxn nonsingular matrix% -b is an nx1 matrix%Output-x is an nx1 matrix containing the solution to Ax=b%Initialize x and the temporary storage matrix C[n,n]=size(A);x=zeros(n,1);C=zeros(1,n+1);%Form the augmented matrix:Aug=[A|b]Aug=[A b];for k=1:n-1%Partial pivoting for column k[y,j]=max(abs(Aug(k:n,k)));%Interchange row k and jC=Aug(k,:);Aug(k,:)=Aug(j+k-1,:);Aug(j+k-1,:)=C;if Aug(k,k)==0'A was singular,No unique solution'breakend%Elimination process for column kfor i=k+1:nm=Aug(i,k)/Aug(k,k);Aug(i,k:n+1)=Aug(i,k:n+1)-m*Aug(k,k:n+1);endend%Back Substitution on [U|y] using backsubx=backsub(Aug(1:n,1:n),Aug(1:n,n+1));例:A=[1 2 1 4;2 0 4 3;4 2 2 1;-3 1 3 2];b=[13;28;20;6];x=gauss_column(A,b)4、LU分解法(不选主元素)function x=lufact(A,b)%Input -A is an nxn matrix% -b is an nx1 matrix%Output-x is an nx1 matrix containing the solution to Ax=b%Initilize x,y,[n,n]=size(A);x=zeros(n,1);y=zeros(n,1);for k=1:n-1%Calculate multiplier and place in subdiagonal portion of A for i=k+1:nmult=A(i,k)/A(k,k);A(i,k)=mult;A(i,k+1:n)=A(i,k+1:n)-mult*A(k,k+1:n);endend%Output LU factorization compact matrixA%solve for yy(1)=b(1);for i=2:ny(i)=b(i)-A(i,1:i-1)*y(1:i-1);end%Solve for xx(n)=y(n)/A(n,n);for i=n-1:-1:1x(i)=(y(i)-A(i,i+1:n)*x(i+1:n))/A(i,i);end例:A=[4 -9 2;2 -4 6;-1 1 -3];b=[5;3;-4];x=lufact(A,b)5、列主元LU分解法(PA=LU)A是非奇异矩阵function x=lu_column(A,b)%Input -A is an nxn matrix% -b is an nx1 matrix%Output-x is an nx1 matrix containing the solution to Ax=b%Initilize x,y,the temporary storage matrix C,and the row%permutation information matrix R[n,n]=size(A);x=zeros(n,1);y=zeros(n,1);C=zeros(1,n);R=1:n;for k=1:n-1%Find the pivot row for column k[maxl,j]=max(abs(A(k:n,k)));%Interchange row k and jC=A(k,:);A(k,:)=A(j+k-1,:);A(j+k-1,:)=C;d=R(k);R(k)=R(j+k-1);R(j+k-1)=d;if A(k,k)==0'A is singular,No unique solution'breakend%Calculate multiplier and place in subdiagonal portion of A for i=k+1:nmult=A(i,k)/A(k,k);A(i,k)=mult;A(i,k+1:n)=A(i,k+1:n)-mult*A(k,k+1:n);endend%Output LU factorization compact matrixA%solve for yy(1)=b(R(1));for i=2:ny(i)=b(R(i))-A(i,1:i-1)*y(1:i-1);end%Solve for xx(n)=y(n)/A(n,n);for i=n-1:-1:1x(i)=(y(i)-A(i,i+1:n)*x(i+1:n))/A(i,i);end例:A=[1 2 6;4 8 -1;-2 3 -5];b=[1;2;3];x=lu_column(A,b)6、Jacobi迭代:程序可用的充分条件是A具有严格对角优势。
线性代数方程组的数值解法讨论解线性方程组的方法,主要分为直接方法和迭代方法两种。
直接法是在没有舍入误差的假设下能在预定的运算次数内求得精确解。
而实际上,原始数据的误差和运算的舍入误差是不可以避免的,实际上获得的也是近似解。
迭代法是构造一定的递推格式,产生逼近精确解的序列。
对于高阶方程组,如一些偏微分方程数值求解中出现的方程组,采用直接法计算代价比较高,迭代法则简单又实用,因此比较受工程人员青睐。
小组成员本着工程应用,讨论将学习的理论知识转变为matlab 代码。
讨论的成果也以各种代码的形式在下面展现。
1 Jacobi 迭代法使用Jacobi 迭代法,首先必须给定初始值,其计算过程可以用以下步骤描述: 步骤1 输入系数矩阵A ,常熟向量b ,初值(0)x ,误差限ε,正整数N ,令1k =.步骤2 (0)11ni i ij jj ii j i x b a x a =≠⎡⎤⎢⎥=-⎢⎥⎢⎥⎣⎦∑,(0)j x 代表(0)x 的第j 个分量。
步骤3 计算11ni i ij j j ii j i y b a x a =≠⎡⎤⎢⎥=-⎢⎥⎢⎥⎣⎦∑,判断1max i i i n x y ε≤≤-<,如果是,则结束迭代,转入步骤5;否则,转入步骤4。
步骤4 判断k N =?如果是,则输出失败标志;否则,置1k k =+,i i x y ⇐,1,2,,i n =,转入步骤2。
步骤5 输出12,,n y y y 。
雅可比迭代代码function [x,k]=Fjacobi(A,b,x0,tol)% jacobi 迭代法 计算线性方程组% tol 为输入误差容限,x0为迭代初值max1= 300; %默认最多迭代300,超过要300次给出警告 D=diag(diag(A)); L=-tril(A,-1);U=-triu(A,1); B=D\(L+U); f=D\b; x=B*x0+f;k=1; %迭代次数while norm(x-x0)>=tol x0=x;x=B*x0+f; k=k+1;if(k>=max1)disp('迭代超过300次,方程组可能不收敛'); return; end%[k x'] %显示每一步迭代的结果 End2 高斯赛德尔迭代由Jacobi 迭代法中,每一次的迭代只用到前一次的迭代值,若每一次迭代充分利用当前最新的迭代值,即在计算第i 个分量(1)k i x +时,用最新分量11()k x +,12()k x +…(1)1k i x +-代替旧分量)1(k x ', )2(k x …)3(k x 就得到高斯赛德尔迭代格式,其数学表达式为:1(1)(1)()111(1,2,,)i n k k k ii ij j ij j j j i ii xb a x a x i n a -++==+⎛⎫=--= ⎪⎝⎭∑∑具体形式如下:()()()(1)()()()11221331111(1)(1)()()22112332222(1)(1)(1)(1)(1)112233,11111k k k k n n k k k k n n k k k k k n n n n n n n n nnx a x a x a x b a x a x a x a x b a x a x a x a x a x b a ++++++++--=----+=----+⋯⋯⋯⋯⋯⋯=-----+矩阵形式表示为:()(1)1(1)()(0,1,2,,),k k k k n +-+=++=x D Lx Ux b将(1)(1)()(0,1,2,,)k k k k n ++=++=Dx Lx Ux b 移项整理得: (1)1()1()()(0,1,2,,))k k x D L Ux D L b k n +--=-+-=记11(),()--=-=-M D L U g D L b ,则(1)()k k x x +=+M g高斯塞德尔迭代function [x,k]=Fgseid(A,b,x0,tol)%高斯-塞德尔迭代法 计算线性方程组 % tol 为误差容限max1= 300; %默认最高迭代300次D=diag(diag(A)); L=-tril(A,-1); U=-triu(A,1); G=(D-L)\U; f=(D-L)\b; x=G*x0+f;k=1; while norm(x-x0)>=tol x0=x;x=G*x0+f; k=k+1;if(k>=max1)disp('迭代次数太多,可能不收敛'); return; end% [k,x'] %显示每一步迭代结果 End3 超松弛迭代法在工程中最常遇到的问题便是线性代数方程组的求解,而线性代数方程组的求解一般可以分为两类,一类是直接法(精确法),包括克莱姆法则方法、LD 分解法等,另一类是迭代法(近似法),包括雅克比迭代法、高斯迭代法、超松弛迭代法等。
高斯顺序消去法的条件-回复高斯顺序消去法,也被称为高斯消元法或高斯消元算法,是一种用于解线性方程组的数值解法。
它的基本思想是通过一系列的行变换将线性方程组转化为一个上三角矩阵,然后通过回代求解出方程的解。
高斯顺序消去法适用于满足以下条件的线性方程组:1. 方程个数等于未知数个数:高斯顺序消去法只适用于方程个数等于未知数个数的情况。
假设方程个数为n,那么未知数个数也应为n。
2. 矩阵的行列式不为0:在进行高斯顺序消去法时,需要将系数矩阵进行行变换,其中包括将某一行乘以一个非零常数,或将某一行加到另一行上。
若矩阵的行列式为0,说明矩阵的行线性相关,无法进行行变换得到上三角矩阵。
3. 未知数的个数必须大于或等于方程组的秩:方程组的秩指的是矩阵的秩,即线性方程组的系数矩阵经过高斯消元法得到的上三角矩阵中非零行的个数。
未知数个数必须大于或等于方程组的秩,否则方程组存在自由变量,无法使用高斯顺序消去法求解。
4. 系数矩阵的系数不应过大或过小:在进行高斯顺序消去法时,系数矩阵中的系数大小对计算结果的精确性有影响。
若系数值过大或过小,会导致计算过程中出现大数吃小数或小数吃小数的情况,从而降低计算结果的准确性。
基于上述条件,以下是高斯顺序消去法的具体步骤:步骤1:构造增广矩阵将线性方程组的系数矩阵与常数项矩阵按列连接起来,形成一个增广矩阵。
步骤2:主元确定找到增广矩阵的第一列中绝对值最大的元素(称为主元),将其所在行与第一行进行交换。
步骤3:主元归一化将主元所在行的元素都乘以一个适当的常数,使主元变为1。
步骤4:消元过程从第一行的下一行开始,将各行的第一列元素消为零。
具体操作是,将下一行第一列的元素乘以一个适当的常数,然后将得到的结果加到第一行的对应的元素上。
重复这一步骤,直到第一列的所有元素都变为零。
步骤5:重新确定主元在步骤4的过程中可能会存在一些舍入误差,导致第一列中的某些元素不为零。
因此,在消元过程结束后,重新找到第二列中的主元,将其所在行与第二行进行交换,并将主元归一化。