matlab求线性方程组的解

matlab求线性方程组的解
求解线性方程分为两种方法–直接法和迭代法常见的方法一共有8种直接法Gauss消去法Cholesky分解法迭代法Jacobi迭代法Gauss-Seidel迭代法超松弛迭代法共轭梯度法Bicg迭代法Bicgstab迭代法
这里我从计算代码的角度来解释一下，代码按以下顺序给出。

把方程组直接带入已知条件，就可以得到答案。

适用条件Gauss消去法：求解中小规模线性方程（阶数不过1000），一般用于求系数矩阵稠密而且没有任何特殊结构的线性方程组
Cholesky分解法：对称正定方程优先使用，系数矩阵A是n 阶对称正定矩阵
Jacobi迭代法非奇异线性方程组，分量的计算顺序没有关系
Gauss-Seidel迭代法与Jacobi迭代法相似，但计算的分量不能改变
超松弛迭代法Jacobi迭代法和Gauss-Seidel迭代法的加速版，由Gauss-Seidel迭代法改进而来，速度较快
共轭梯度法需要确定松弛参数w,只有系数矩阵具有较好的性质时才可以找到最佳松弛因子。

但好处是不用确定任何参数，他是对称正定线性方程组的方法也是求解大型稀疏线性方程组最热门的方法
Bicg迭代法本质是用双共轭梯度求解线性方程组的方法，对求解的方程没有正定性要求
Bicgstab迭代法本质是用稳定双共轭梯度求解线性方程组的方法，对求解的方程没有正定性要求
Gauss消去法第一、二个函数ltri、utri是一定要掌握的，后面的几乎每个函数都要用到ltri简单来说，当Ly=bb,L(非奇异下三角矩阵）已知求y
function y =ltri(L,b)
n=size(b,1);
y=zeros(n,1);for j =1:n-
1y(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-y(j)*L(j+1:n,j); end
y(n)=b(n)/L(n,n);
utri简单来说，当Ux=yy,U(非奇异上三角矩阵）已知求x
function x =utri(U,y)
n=size(y,1);
x=zeros(n,1);for j = n:-1:2x(j)=y(j)/U(j,j);y(1:j-
1)=y(1:j-1)-x(j)*U(1:j-1,j);
end
x(1)=y(1)/U(1,1);
gauss算法，计算时粘贴过去就好
function[L,U]=gauss(A)
n=size(A,1);for k =1:n-
1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k +1:n)-A(k+1:n,k)*A(k,k+1:n);
end
L=tril(A,-1)+eye(n);U=triu(A);
使用例子已经知道一个线性方程组，这里我就不写出数学形式了，A是系数矩阵，直接把上面写好的函数复制过来在运算就
可以。

主要是把这两个矩阵写进去。

A = [1 -1 1 -3;0 -1 -1 1;2 -2 -4 6;1 -2 -4 1];b = [1;0;-1;-1];
A=[1-11-3;0-1-11;2-2-46;1-2-41];
b =[1;0;-1;-1];[L,U]=gauss(A);
y =ltri(L,b);
x =utri(U,y)
Cholesky分解法这里是要用到ltri和utri，直接把上面的粘贴过来就好这里的算法比较简单，用L = chol(A,‘lower’);这里A是系数矩阵
使用例子用Cholesky分解法求接对称正定方程Ax=b，b随机选取，系数矩阵A是100阶的矩阵n = 100;A = 10 * eye(n) + diag(ones(n-1,1),-1) + diag(ones(n-1,1),1);
n =100;A=10*eye(n)+diag(ones(n-1,1),-1)+diag(ones(n-1,1),1);
b =rand(n,1);L=chol(A,'lower');
y =ltri(L,b);
x =utri(L',y)
Jacobi迭代法这个就是Jacobi迭代法的算法，用时复制。

这里一共要输入四个条件，系数矩阵A，矩阵b，x0=zeros(列数,1),tol=1.0e-6
function x =jacobi(A,b,x0,tol)D=diag(diag(A));L=-
tril(A,-1);U=-triu(A,1);B=D \ (L+U);
g =D \ b;
x =B* x0 + g;
n =1;whilenorm(x - x0)> tol
x0 = x;
x=B* x0 + g;
n = n +1;
end
x,n
使用例子用Jacobi迭代法求解下列方程A = [10 -1 0;-1 10 -2;0 -2 10];b = [9;7;6];
A=[10-10;-110-2;0-210];
b =[9;7;6];
x =zeros(3,1);
tol =1.0e-6;
x =jacobi(A,b,x,tol);
Gauss-Seidel迭代法算法函数，输入四个条件，和上面相同
function x =seidel(A,b,x0,tol)D=diag(diag(A));L=-
tril(A,-1);U=-triu(A,1);B=(D-L) \ U;
g =(D-L) \ b;
x =B* x0 +g;
n=1;whilenorm(x - x0)> tol
x0 = x;
x=B* x0 + g;
n = n +1;
end
x,n
用这个方法解上面Jacobi迭代法的例子
A=[10-10;-110-2;0-210];
b =[9;7;6];
x =zeros(3,1);
tol =1.0e-6;
x =seidel(A,b,x,tol);
超松弛迭代法算法函数，这里需要输入五个条件：系数矩阵A，矩阵b，x0=zeros(列数,1),w(松弛因子一到二之间）
=1.025，tol=1.0e-6
function x =sor(A,b,x0,w,tol)D=diag(diag(A));L=-
tril(A,-1);U=-triu(A,1);B=(D- w *L) \ ((1- w)*D+ w *U);
g =(D- w *L) \ b * w;
x =B* x0 +g;
n=1;whilenorm(x - x0)>= tol
x0 = x;
x=B* x0 + g;
n = n +1;
end
x,n
使用例子求这个方程组A = [10 -1 0;-1 10 -2;0 -2 10];b = [9;7;6];
A=[10-10;-110-2;0-210];
b =[9;7;6];
x =zeros(3,1);
tol =1.0e-6;
w =1.025;
x =sor(A,b,x,w,tol);
共轭梯度法
A=[16484;41084;881210;441012];
b =[32;26;38;30];
tol =1.0e-6;
maxit =1000;
x =pcg(A,b,tol,maxit)
Bicg迭代法这个比较简单，不需要额外写函数，直接用x = bicg(A,b,tol,maxit)maxit取1000
使用例子
n =100;A=10*eye(n)+diag(ones(n-1,1),-1)+diag(ones(n-1,1),1);
b =rand(n,1);
tol =1.0e-6;
x =bicg(A,b,tol,1000)
Bicgstab迭代法和上面一样的输入条件x =
bicgstab(A,b,tol,1000)
n =100;A=10*eye(n)+diag(ones(n-1,1),-1)+diag(ones(n-1,1),1);
b =rand(n,1);
tol =1.0e-6;
x =bicgstab(A,b,tol,1000)。