线性方程组AX=B的数值计算方法实验

实验5 线性代数方程组的数值解法化工系 毕啸天 2010011811【实验目的】1. 学会用MATLAB 软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析；2. 通过实例学习用线性代数方程组解决简化的实际问题。

【实验内容】题目3已知方程组Ax=b ，其中A ，定义为⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡--------------=32/14/12/132/14/14/12/132/14/14/12/132/14/12/13A试通过迭代法求解此方程组，认识迭代法收敛的含义以及迭代初值和方程组系数矩阵性质对收敛速度的影响。

实验要求：（1）选取不同的初始向量x (0)和不同的方程组右端项向量b ，给定迭代误差要求，用雅可比迭代法和高斯-赛德尔迭代法计算，观测得到的迭代向量序列是否均收敛？若收敛，记录迭代次数，分析计算结果并得出你的结论；（2）取定右端向量b 和初始向量x (0)，将A 的主对角线元素成倍增长若干次，非主对角线元素不变，每次用雅可比迭代法计算，要求迭代误差满足 ，比较收敛速度，分析现象并得出你的结论。

3.1 模型分析选取初始向量x(0) =(1,1,…,1)T ，b=(1,1,…,1)T ，迭代要求为误差满足 ，编写雅各比、高斯-赛德尔迭代法的函数，迭代求解。

3.2 程序代码function x = Jacobi( x0,A,b,m ) D=diag(diag(A)); U=-triu(A,1); L=-tril(A,-1); B1=D\(L+U); f1=D\b; x(:,1)=x0;x(:,2)=B1*x(:,1)+f1; k=1;while norm((x(:,k+1)-x(:,k)),inf)>m x(:,k+2)=B1*x(:,k+1)+f1; k=k+1; endendfunction x = Gauss( x0,A,b,m )D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1);B2=(D-L)\U;f2=(D-L)\b;x(:,1)=x0;x(:,2)=B2*x(:,1)+f2;k=1;while norm((x(:,k+1)-x(:,k)),inf)>mx(:,k+2)=B2*x(:,k+1)+f2;k=k+1;endendA1=3.*eye(20,20);A2=sparse(1:19,2:20,-1/2,20,20);A3=sparse(1:18,3:20,-1/4,20,20);AA=A1+A2+A3+A2'+A3';A=full(AA);b=ones(20,1); %输入自选右端项向量b x0=ones(20,1); %输入自选初始向量x0 m=1e-5;x1=Jacob(x0,A,b,m);x2=Gauss(x0,A,b,m);结果输出数据：3.3.1将x0各分量初值置为0【分析】从数据中可以看出，当迭代的初值变化了，达到相同精度所需要的迭代次数也变化了。

线性方程组AX=B 的数值计算方法实验一、 实验描述：随着科学技术的发展，线性代数作为高等数学的一个重要组成部分，在科学实践中得到广泛的应用。

本实验的通过C 语言的算法设计以及编程，来实现高斯消元法、三角分解法和解线性方程组的迭代法（雅可比迭代法和高斯-赛德尔迭代法），对指定方程组进行求解。

二、 实验原理：1、高斯消去法：运用高斯消去法解方程组，通常会用到初等变换，以此来得到与原系数矩阵等价的系数矩阵，达到消元的目的。

初等变换有三种：(a)、(交换变换)对调方程组两行；(b)、用非零常数乘以方程组的某一行；(c)、将方程组的某一行乘以一个非零常数，再加到另一行。

通常利用(c)，即用一个方程乘以一个常数，再减去另一个方程来置换另一个方程。

在方程组的增广矩阵中用类似的变换，可以化简系数矩阵，求出其中一个解，然后利用回代法，就可以解出所有的解。

2、选主元：若在解方程组过程中，系数矩阵上的对角元素为零的话，会导致解出的结果不正确。

所以在解方程组过程中要避免此种情况的出现，这就需要选择行的判定条件。

经过行变换，使矩阵对角元素均不为零。

这个过程称为选主元。

选主元分平凡选主元和偏序选主元两种。

平凡选主元：如果()0p pp a ≠，不交换行；如果()0p pp a =，寻找第p 行下满足()0p pp a ≠的第一行，设行数为k ，然后交换第k 行和第p 行。

这样新主元就是非零主元。

偏序选主元：为了减小误差的传播，偏序选主元策略首先检查位于主对角线或主对角线下方第p 列的所有元素，确定行k ，它的元素绝对值最大。

然后如果k p >，则交换第k 行和第p 行。

通常用偏序选主元，可以减小计算误差。

3、三角分解法：由于求解上三角或下三角线性方程组很容易所以在解线性方程组时，可将系数矩阵分解为下三角矩阵和上三角矩阵。

其中下三角矩阵的主对角线为1，上三角矩阵的对角线元素非零。

有如下定理：如果非奇异矩阵A 可表示为下三角矩阵L 和上三角矩阵U 的乘积： A LU = (1) 则A 存在一个三角分解。

第12 章 实验十一线性方程组的数值解法实验目的：理解线性方程组基本定理，掌握求解线性方程组的各种方法，逆矩阵求法、高斯消去法、带选主元的分解法、迭代法。

12.1 线性方程组（基本定理） 已知m 个方程n 个未知量的方程组：⎪⎪⎩⎪⎪⎨⎧=+++=+++=+++m n mn m m n n n n b x a x a x a b x a x a x a b x a x a x a 22112222212********* （12.1）其中ij a 为已知系数，j x 未知，j b 为已知项，当),,2,1(n j b j =全为零时，称（12.1）方程组为线性齐次方程组，当),,2,1(n j b j =不全为零时，称（12.1）方程组为线性非齐次方程组。

上述方程组也可表示如下矩阵形式：b Ax = 其中A 、x 和y 分别定义为：⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=mn m n a a a a A LL L L L 1111⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=n x x x 1 ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=m b b b 1。

方程组（12.1）中的m 与n 有三种情况：m=n 或m<n 或m>n ，m=n 是最常见的。

当m=n 时，有系数矩阵的行列式0≠A ，方程组（12.1）有唯一解。

用矩阵的秩可以得到更满意的结果：（1）当n B R A R ==)()(时，方程组（12.1）有唯一解， （2）当n B R A R <=)()(时，方程组（12.1）有无穷多组解， （3）当)()(B R A R ≠时，方程组（12.1）无解。

12.2举例例12.1 用MATLAB 求解线性方程组：⎪⎪⎩⎪⎪⎨⎧=+++-=+++=+++=+++62332022428340213424321432143214321x x x x x x x x x x x x x x x x解法一：逆矩阵求法设 A=[1 2 1 4 ; 2 0 4 3 ;4 2 2 1 ;-3 1 3 2 ]; det(A)=-180, B=[13 28 20 6]’; 则：X=A\B % A\表示矩阵A 的逆，或用1-A 表示得到X =3.0000 -1.00004.0000 2.0000也可以用如下命令： X=B ’/A ’ 得到X =3.0000 -1.00004.0000 2.0000解法二：在MA TLAB 上，用高斯消去法求解方程组：⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-62028132313122434024214321x x x x （12.2） 首先写出方程组对应的增广矩阵：a=[ 1 2 1 4 13; ... 2 0 4 3 28 ; ... 4 2 2 1 20 ; ...-3 1 3 2 6 ]其中前四列是系数矩阵，最后一列是方程组的常数列。

线性方程组AX=B的数值计算方法实验线性方程组AX=B 的数值计算方法实验【摘要】在自然科学与工程技术中很多问题的解决常常归结为解线性代数方程组。

例如电学中的网络问题，船体数学放样中建立三次样条函数问题，用最小二乘法求实验数据的曲线拟合问题，解非线性方程组的问题，用差分法或者有限元法解常微分方程，偏微分方程边值问题等都导致求解线性方程组。

线性代数方面的计算方法就是研究求解线性方程组的一些数值解法与研究计算矩阵的特征值及特征向量的数值方法。

关于线性方程组的数值解法一般有两类：直接法和迭代法。

关 键 字 高斯消元法、三角分解法、高斯-赛德尔迭代、稀疏矩阵一、实验目的1.掌握高斯消元法、三角分解法、高斯—赛德尔迭代发的编程技巧。

2.掌握线性方程组AX=B 的数值计算方法。

3.掌握矩阵的基本编程技巧。

二、实验原理1.高斯消元法数学上，高斯消元法是线性代数规划中的一个算法，可用来为线性方程组求解。

高斯（Gauss ）夏鸥按法其实是将一般的线性方程组变换为三角形（上三角）方程组求解问题（消元法），只是步骤规范，便于编写计算机程序。

一般高斯消元法包括两过程：先把方程组化为同解的上三角形方程组，再按相反顺序求解上三角方程组。

前者称为消去或消元过程，后者称回代过程。

消去过程实际上是对增广矩阵作行初等变换。

对一般的n 阶方程组，消去过程分n-1步：第一步消去11a 下方元素。

第二步消去22a 下方元素，......，第n-1步消去1-n 1-n a ，下方元素。

即第k 步将第k 行的适当倍数加于其后各行，或可说是从k+1~n 行减去第k 行的适当倍数，使它们第k 列元素变为零，而其余列元素减去第k 行对应列元素的倍数。

2.三角分解法三角分解法是将原正方 (square)矩阵分解成一个上三角形矩阵或是排列(permuted) 的上三角形矩阵和一个 下三角形矩阵，这样的分解法又称为LU 分解法。

它的用途主要在简化一个大矩阵的行列式值的计算过程，求 反矩阵，和求解联立方程组。

线性代数方程组的直接解法实验目的：线性方程组求解的直接法编程实现，实验内容：线性方程组求解的高斯消去法算法实现线性方程组求解的主元素消去法算法实现线性方程组求解的LU 分解得方法算法实现线性方程组求解追赶法算法实现实验比较：高斯消去、主元素消去、LU 分解都用实例 ⎪⎩⎪⎨⎧-=-+-=+-=++15.3181533126321321321x x x x x x x x x 这个进行比较。

知识理论解线性方程组的方法大致分为直接法和迭代法。

直接法是指假设计算过程中不产生舍入误差，经过有限次的运算可求的方程组精确解的方法。

方程（2-1）⎪⎪⎩⎪⎪⎨⎧=+++=+++=+++nn nn n n n n n n b x a x a x a b x a x a x a b x a x a x a ...... (22112)222212111212111 记为：AX=b ；一、高斯（Gauss ）消元法(1).Gauss 消元法是最基本的一种方法。

先逐次消去变量，将方程组化成同解的上三角方程组。

消元过程：先逐次消去变量，将方程组化成同解的上三角方程组 基本思想回代过程：按方程的相反顺序求解三角形方程组，得到原方程组的解。

(2) Gauss 消去法的求解思路为：若先让第一个方程组保持不变，利用它消去其余方程组中的,使之变成一个关于变元的n-1阶方程组。

按照（1）中的思路继续运算得到更为低阶的方程组。

经过n-1步的消元后，得到一个三角方程。

利用求解公式回代得到线性方程组的解。

①消元过程：第一次消元，,0)1(11≠a设 )1(11)1(11a a l i i =记（i=1,2,3……，n ）.将（2-1）中第i 个方程减去第一个方程乘以1i l （i=1,2,3……，n ），完成第一次消元，⎪⎩⎪⎨⎧=++=++=+++)2()2(2)2(2)2(2)2(22)2(22)1(1)1(12)1(121)1(11...... ... n n nn n n n n n b x a x a b x a x a b x a x a x a 得等价方程组 （2-2）)2()2(:b x A =简记为其中：)1(11)1()2(j i ij ij a l a a -=；)1(11)1()2(b l b b i i i -=次消元后经过1-k ：⎪⎪⎪⎪⎪⎪⎪⎩⎪⎪⎪⎪⎪⎪⎪⎨⎧=+++=+++=+++=+++++=++++++++++++++++++++... ...... ......... ......)()(1)(1,)()(1)(,11)(1,1)(,1)()(1)(1,)()2(2)2(21)2(1,2)2(22)2(22)1(1)1(11)1(1,1)1(12)1(121)1(11k n n k nn k k k n k k nk k k n k n k k k k k k k k k k kn k kn k k k k k k kk n n k k k k n n k k k k b x a x a x a b x a x a x a b x a x a x a b x a x a x a x a b x a x a x a x a x a 简记为：)1()1(++=k k b x A其中：)()()1(k kjik k ij k ija l a a -=+ )()()1(k k ik k i k ib l b b -=+),,1,(n k j i +=按上述方法完成n-1次消元后得到。

线性⽅程组的⽅法-数值分析-王兵团-北京交通⼤学注解：1.线性代数中线性⽅程组的⽅法：克拉默法则。

线性⽅程组：Ax=b解：x i=D i/D如果A可逆，还可以写成：x=A-1/b⽅程组的解是：系数⾏列式某⼀项换成等式右端常数项/系数⾏列式。

既然可以有这么好的公式，那为何还要学习其它解法呢？答：好多数学的公式⼀旦⽤到计算机⾥⾯，就不⾏了。

有⼈实验过，100万/s的计算量，解算40阶的线性⽅程组的解，要算⼀年。

天⽓预报的话有⼏百万⼏千万的⽅程组，怎样快速解出来？⼤规模集成电路也需要解⼤规模⽅程组的。

⼈们需要快速求解⼤规模的线性⽅程组，这样，理论解就不⾏了。

数值分析讲的是怎样⽤计算机快速求出数学问题的解。

注解：1.如果⽅程组的数量>未知数个数，没有解，或者有最⼩⼆乘解。

2.如果如果⽅程组的数量<未知数个数，没有解，有⽆穷多解。

3.如果如果⽅程组的数量=未知数个数，有唯⼀解。

计算机做的最多的是这种情况：即⽅(阵)的情况。

4.线性⽅程组怎么得来的？答：每次实验得来的。

5.⼀个系统，有n个元器件，x i代表第i个元器件，每个元器件相当于⼀个变量x i,它们之间的变化有⼀定的关系。

每实验⼀次，得到⼀个它们之间相互关系的⽅程。

实验n次，得到n个⽅程。

通过⽅程组，求出n个元器件情况。

注解：1.计算机求的解都是近似解。

2.学⼀个东西，怎样学好？答：通过类⽐去看。

注解：1.简单迭代法是怎样做的？答：注解：1.初值可以给成：[0,0,0,...].2.x是⼀个向量。

注解：1.构造迭代格式所⽤的等价形式⼀定是有的。

2.未必都收敛的意思：⽐如，如果c给的合适，就收敛，如果不合适，就不收敛。

3.前⽂的引例就是example9.注解：1.x k代表第k步的迭代值。

x k是⼀个向量，所以ε(k)也是⼀个向量，是指第k步的迭代误差。

2.红⾊部分的等式是⼀个递推式⼦。

3.。

完全取决于迭代矩阵B，跟初值怎样选择是没有关系的。

线性方程组AX=B的数值计算方法实验学号：姓名：梁哲豪一、实验描述在自然科学和工程技术中很多问题的解决常常归结为解线性代数方程组。

例如电学中的网络问题，船体数学放样中建立三次样条函数问题，用最小二乘法求实验数据的曲线拟合问题，解非线性方程组问题，用差分法或者有限元法解常微分方程，偏微分方程边值问题等都导致求解线性方程组，而且后面几种情况常常归结为求解大型线性方程组。

线性代数方面的计算方法就是研究求解线性方程组的一些数值解法与研究计算矩阵的特征值及特征向量的数值方法。

关于线性方程组的数值解法一般有两类：直接法：若在计算过程中没有舍入误差，经过有限步算术运算，可求得方程组的精确解的方法。

迭代法：用某种极限过程去逐步逼近线性方程组精确解的方法。

迭代法具有占存储单元少，程序设计简单，原始系数矩阵在迭代过程中不变等优点，但存在收敛性及收敛速度等问题。

上三角线性方程组的求解：基本算法：高斯消元法：将原方程组化为三角形方阵的方程组：（k=1,2,…,n-1; i=k+1,k+2, …,n ;j=k+1,k+2, …,n+1）由回代过程求得原方程组的解：LU分解法：将系数矩阵A转化为A=L*U，L为单位下三角矩阵，U为普通上三角矩阵，然后通过解方程组l*y=b，u*x=y,来求解x。

二、实验内容1、许多科学应用包含的矩阵带有很多零。

在实际情况中很重要的三角形线性方程组有如下形式：……构造一个程序求解三角形线性方程组。

可假定不需要变换。

而且可用第k 行消去第k+1行的x。

k核心代码：#include<iostream.h>#include<math.h>#include<iomanip.h>#define N 4//矩阵阶数void ColPivot(double c[N][N+1],double[]);//函数声明void main(){int i,j;double x[N];double c[N][N+1]={1,3,5,7,1,2,-1,3,5,2,0,0,2,5,3,-2,-6,-3,1,4};cout<<"----------------------------------------"<<endl;cout<<"系数矩阵为: \n";for(i=0;i<N;i++){for(j=0;j<N;j++)cout<<setw(10)<<c[i][j];cout<<endl;}cout<<"右侧矩阵 y 为: \n";for(i=0;i<N;i++)cout<<setw(10)<<c[i][N];cout<<endl;cout<<"----------------------------------------"<<endl;ColPivot(c,x);//调用函数，进行高斯消去法变换cout<<"变换后得到的三角矩阵: \n";for(i=0;i<N;i++){for(j=0;j<N;j++)cout<<setw(10)<<c[i][j];cout<<endl;}cout<<"变换后的右侧矩阵 y 为: \n";for(i=0;i<N;i++)cout<<setw(10)<<c[i][N];cout<<endl;cout<<"----------------------------------------"<<endl; cout<<"方程的解为: \n";for(i=0;i<N;i++)cout<<" x["<<i<<"]= "<<x[i]<<endl;cout<<"----------------------------------------"<<endl; }void ColPivot(double c[N][N+1],double x[]){int i,j,k;double p,max;double t[N];for(i=0;i<=N-2;i++){max=0;k=i;for(j=i+1;j<N;j++)if(fabs(c[j][i])>max){k=j;max=fabs(c[j][i]);//选主元}if(k!=i)for(j=i;j<=N;j++){p=c[i][j];c[i][j]=c[k][j];//选出主元后进行交换c[k][j]=p;}for(j=i+1;j<N;j++){p=c[j][i]/c[i][i];for(k=i;k<=N;k++)c[j][k]-=p*c[i][k];//高斯消去，进行计算}}for(i=0;i<N;i++)t[i]=c[i][N];for(i=N-1;i>=0;i--)//利用回代法求最终解{for(j=N-1;j>i;j--)t[i]-=c[i][j]*x[j];x[i]=t[i]/c[i][i];}}运行结果：（以具体方程组为例）2、（PA=LU：带选主元的分解法）求解线性方程组AX=B,其中：A=B=核心代码：#include <stdio.h>#include <math.h>#define L 30double a[L][L],b[L],l[L][L],u[L][L],x[L],y[L];int main(){int n,i,j,k,r;printf("请输入矩阵元次：\n");scanf("%d",&n);printf("请输入矩阵各项：\n");for(i=1;i<=n;++i){for(j=1;j<=n;++j){scanf("%lf",&a[i][j]);}}printf("请输入方程组的常数项：\n");for(i=1;i<=n;++i){scanf("%lf",&b[i]);}for(i=1;i<=n;++i){for(j=1;j<=n;++j){l[i][j]=0;u[i][j]=0.0;}}for(k=1;k<=n;++k){for(j=k;j<=n;++j){u[k][j]=a[k][j];for(r=1;r<k;++r){u[k][j]-=l[k][r]*u[r][j];}}for(i=k+1;i<=n;++i){l[i][k]=a[i][k];for(r=1;r<k;++r){l[i][k]-=l[i][r]*u[r][k];}l[i][k]/= u[k][k];}l[k][k]=1.0;}for(i=1;i<=n;++i){y[i]= b[i];for(j=1;j<i;++j){y[i]-=l[i][j]*y[j];}}for(i=n;i>0;--i){x[i]= y[i];for(j=i+1;j<=n;++j){x[i]-=u[i][j]*x[j];}x[i]/= u[i][i];}for(i=1;i<=n;++i){printf("%0.2lf\n",x[i]);}return0;}运行结果：3、使用程序3.3求解线性方程组AX=B，其中，A= [a ij] N×N= i j-1,而且B=[b ij] N×1, b11=N,当i≥2时，b i1=(i N-1)/(i-1),对N=3,7,11的情况分别求解。

实验5 线性代数方程组的数值解法分1 黄浩 43一、实验目的1.学会用MATLAB软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析；2.通过实例学习用线性代数方程组解决简化的实际问题。

二、For personal use only in study and research; not forcommercial use三、四、实验内容1.《数学实验》第二版（问题1）问题叙述：通过求解线性方程组，理解条件数的意义和方程组性态对解的影响，其中是n阶范德蒙矩阵，即是n阶希尔伯特矩阵，b1，b2分别是的行和。

（1）编程构造（可直接用命令产生）和b1，b2；你能预先知道方程组和的解吗？令n=5，用左除命令求解（用预先知道的解可验证程序）。

（2）令n=5,7,9，…，计算和的条件数。

为观察他们是否病态，做以下试验：b1，b2不变，和的元素，分别加扰动后求解；和不变，b1，b2的分量b1(n),b2(n)分别加扰动后求解。

分析A与b的微小扰动对解的影响。

取10^-10,10^-8,10^-6。

（3）经扰动得到的解记做，计算误差，与用条件数估计的误差相比较。

模型转换及实验过程：（1）小题.由b1,b2为,的行和，可知方程组和的精确解均为n 行全1的列向量。

在n=5的情况下，用matlab编程（程序见四.1），构造,和b1，b2，使用高斯消去法得到的解x1,x2及其相对误差e1,e2（使用excel计算而得）为：由上表可见，当n=5时，所得的解都接近真值，误差在10^-12的量级左右。

（2）小题分别取n=5,7,9,11,13,15，计算和的条件数c1和c2，（程序见四.2），结果如下：由上表可见，二者的条件数都比较大，可能是病态的。

为证实和是否为病态，先保持b不变，对做扰动，得到该情况下的高斯消元解，（程序见四.3），结果如下：(为使结果清晰简洁，在此仅列出n=5,9,13的情况，n=7,11,15略去)=10^-10时：=10^-8时：=10^-6时：由上表可见：a)对于希尔伯特阵，随着阶数的增加，微小扰动对解带来的影响越来越大，到了n=9时，已经有了6倍误差的解，到了n=13时，甚至出现了22倍误差的解元素；而随着的增加，解的偏差似乎也有增加的趋势，但仅凭上述表格无法具体判断（在下一小题中具体叙述）。