实验2 线性方程组的直接解法
数值计算实验报告
![数值计算实验报告](https://img.taocdn.com/s3/m/53c141c0cf84b9d529ea7a60.png)
本科实验报告课程名称:计算机数值方法实验项目:方程求根、线性方程组的直接解法、线性方程组的迭代解法、代数插值实验地点:专业班级:学生姓名:指导教师:实验一方程求根}五、实验结果与分析二分法实验结果迭代法实验结果结果分析:本题目求根区间为[1,2],精度满足|x*-x n|<0.5×10-5,故二分法用公式|x*-x n|<(b-a)/ 2n,可求得二分次数并输出每次结果。
对迭代法首先要求建立迭代格式。
迭代格式经计算已输入程序之中,故直接给初值便可利用迭代法求出精度下的解。
六、讨论、心得每次的实验都是对已学过的理论知识的一种实战。
通过本次实验,我将二分法与迭代法的思路清晰化并且将其变成计算机设计语言编写出来,运用到了实际解决问题上感觉很好。
我自认为本次跟其他同学比较的优点在于我在二分法实现的时候首先利用换底公式将需要的二分次输输出,如此便很清晰明了的知道接下来每一步的意思。
迭代法给我的感觉便是高度的便捷简化,仅用几行代码便可以同样解决问题。
相比较二分法来说,我更喜欢迭代的思路。
实验二线性方程组的直接解法for(k=n-2;k>=0;k--){sum=0;for(j=k+1;j<n;j++)sum=sum+a[k][j]*x[j];x[k]=(b[k]-sum)/a[k][k];}for(i=0;i<n;i++)printf("x[%d]=%f ",i,x[i]); printf("\n"); //输出解向量x}五、实验结果与分析结果结果分析:如上图所示,输入线性方程组元数n=3,则会要求输入3*3的系数矩阵A与向量b构成的增广矩阵。
根据算法需要将系数矩阵A消元成上三角矩阵。
随后根据矩阵乘法公式变形做对应的回代。
六、讨论、心得本次实验在编写时候感觉还好,感觉将思路变成了程序设计语言,得以实现题目的要求。
但是在运行以及结果分析的时候,感觉到了本实验的一些不足之处:就是我的实验虽然可以实现不同的元数的线性方程组求解,但是缺少了分析初始条件——主元素不能为零。
求解线性方程组的直接解法
![求解线性方程组的直接解法](https://img.taocdn.com/s3/m/e3b5563ab52acfc789ebc9a8.png)
求解线性方程组的直接解法5.2LU分解① Gauss消去法实现了LU分解顺序消元结束时的上三角矩阵U和所用的乘数,严格下三角矩阵。
将下三角矩阵的对角元改成1,记为L,则有A=LU,这事实是一般的,我们不难从消去的第k个元素时的矩阵k行及k列元素的历史得到这一点.因为从消元的历史有u kj=a kj-m k1u1j- m k2u2j -…- m k,k-1u k-1,j, j=k,k+1,…,nm ik=(a ik-m i1u1k- m i2u2k -…-m i,k-1u k-1,k>/u kk i=k+1,k+2,…,n于是a kj=m k1u1j+m k2u2j+…+m k,k-1u k-1,j+u kj, j=k,k+1,…,na ik=m i1u1k+m i2u2k+…+m i,k-1u k-1,k+m ik u kk i=k+1,k+2,…,n从前面两个式子我们可以直接计算L和U(见下段>.将矩阵分解为单位下三角矩阵和上三角矩阵之积称为矩阵的LU分解.顺序消元实现了LU分解,同时还求出了g, Lg=b的解.②直接LU分解上段我们得到(l ij=m ij>u kj=a kj-l k1u1j-l k2u2j -…- l k,k-1u k-1,j, j=k,k+1,…,nl ik=(a ik-l i1u1k-l i2u2k -…-l i,k-1u k-1,k>/u kk i=k+1,k+2,…,n2诸元素对应乘积,只不过算L的元素时还要除以同列对角元.这一规律很容易记住.可写成算法(L和U可存放于A>:for k=1:n-1for j=k:nu kj=a kj-l k1u1j-l k2u2j -…- l k,k-1u k-1,jendfor i=k+1:nl ik=(a ik-l i1u1k-l i2u2k -…-l i,k-1u k-1,k>/u kkendend这一算法也叫Gauss消去法的紧凑格式,可一次算得L,U的元素,不需逐步计算存储.考察上面的表格会发现还可安排其它计算次序,只要在这一次序下每个元素左边的L的元素与上方的U的元素已计算在先。
数值线性代数_MathCAD实验
![数值线性代数_MathCAD实验](https://img.taocdn.com/s3/m/9ea9b36f011ca300a6c390a6.png)
,其中 Lk
=
⎜ ⎜ ⎜ ⎜⎜⎝
%
1 −lk +1,k ... −ln,k
%
⎟
⎟ ⎟
,
li,k
=
ai,k akk
,
。
% 1 ⎟⎟⎠ i = k +1,k + 2,...,n
则
A
=
L1−1Βιβλιοθήκη "L−1 n−2
L−n1−1U
=
LU
求解 Ax = b 等价于求解
{ LUx = b ⇔
先解 Ly = b, 再解 Ux = y。
""矩阵的三角分解
定理 若 A∈ Rnn 的顺序主子阵 Ai ,i = 1, 2,..., n 均非奇异,则存在唯一单位下三角阵 L 和上三角阵U 使得 A = LU 。
注: li,k
=
ai,k akk
中分母 akk
称为主元。
用高斯变换计算 A 三角分解的代码如下:
Gauss2( A) :=
"没有优化的LU分解"
for j ∈ 1 .. n
U1, j ← A1, j ⊕ "先计算U的第一行"
for i∈ 2 .. n
Li , 1
←
Ai , 1 U1,1
⊕
"再计算L的第一列"
"下一行中i不能为n,否则j要出界"
"因此Unn要单独计算"
for i∈ 2 .. n − 1
for j ∈ i.. n
i−1
∑ Ui , j ← Ai , j −
n
+ 1, n
,由于U
是n
求解线性方程组的直接解法范文
![求解线性方程组的直接解法范文](https://img.taocdn.com/s3/m/32647ba38762caaedc33d417.png)
求解线性方程组的直接解法5.2 LU 分解① Gauss 消去法实现了LU 分解顺序消元结束时的上三角矩阵U 和所用的乘数,严格下三角矩阵。
将下三角矩阵的对角元改成1,记为L ,则有A =LU ,⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-613322121121542774322这事实是一般的,我们不难从消去的第k 个元素时的矩阵k 行及k 列元素的历史得到这一点.因为从消元的历史有 u kj =a kj -m k 1u 1j - m k 2u 2j -…- m k ,k-1u k-1,j , j=k ,k+1,…,n m ik =(a ik -m i 1u 1k - m i 2u 2k -…-m i ,k-1u k-1,k )/u kk i=k+1,k+2,…,n 于是 a kj =m k 1u 1j +m k 2u 2j +…+m k ,k-1u k-1,j +u kj , j=k ,k+1,…,n a ik =m i 1u 1k +m i 2u 2k +…+m i ,k-1u k-1,k +m ik u kk i=k+1,k+2,…,n 从前面两个式子我们可以直接计算L 和U (见下段).将矩阵分解为单位下三角矩阵和上三角矩阵之积称为矩阵的LU 分解.顺序消元实现了LU 分解,同时还求出了g , Lg =b 的解.② 直接LU 分解上段我们得到(l ij =m ij ) u kj =a kj -l k 1u 1j -l k 2u 2j -…- l k ,k-1u k-1,j , j=k ,k+1,…,n l ik =(a ik -l i 1u 1k -l i 2u 2k -…-l i ,k-1u k-1,k )/u kk i=k +1,k+2,…,n2诸元素对应乘积,只不过算L 的元素时还要除以同列对角元.这一规律很容易记住.可写成算法(L 和U 可存放于A ): for k =1:n -1 for j=k :n u kj =a kj -l k 1u 1j -l k 2u 2j -…- l k ,k-1u k-1,jendfor i=k+1:nl ik =(a ik -l i 1u 1k -l i 2u 2k -…-l i ,k-1u k-1,k )/u kk end end这一算法也叫Gauss 消去法的紧凑格式,可一次算得L ,U 的元素,不需逐步计算存储.考察上面的表格会发现还可安排其它计算次序,只要在这一次序下每个元素左边的L 的元素与上方的U 的元素已计算在先。
解线性方程组的直接解法
![解线性方程组的直接解法](https://img.taocdn.com/s3/m/32aa430703d8ce2f0066232d.png)
解线性方程组的直接解法一、实验目的及要求关于线性方程组的数值解法一般分为两大类:直接法与迭代法。
直接法是在没有舍入误差的情况下,通过有限步运算来求方程组解的方法。
通过本次试验的学习,应该掌握各种直接法,如:高斯列主元消去法,LU分解法和平方根法等算法的基本思想和原理,了解它们各自的优缺点及适用范围。
二、相关理论知识求解线性方程组的直接方法有以下几种:1、利用左除运算符直接求解线性方程组为bx\=即可。
AAx=,则输入b2、列主元的高斯消元法程序流程图:输入系数矩阵A,向量b,输出线性方程组的解x。
根据矩阵的秩判断是否有解,若无解停止;否则,顺序进行;对于1p:1-=n选择第p列中最大元,并且交换行;消元计算;回代求解。
(此部分可以参看课本第150页相关算法)3、利用矩阵的分解求解线性方程组(1)LU分解调用matlab中的函数lu即可,调用格式如下:[L,U]=lu(A)注意:L往往不是一个下三角,但是可以经过行的变换化为单位下三角。
(2)平方根法调用matlab 中的函数chol 即可,调用格式如下:R=chol (A )输出的是一个上三角矩阵R ,使得R R A T =。
三、研究、解答以下问题问题1、先将矩阵A 进行楚列斯基分解,然后解方程组b Ax =(即利用平方根法求解线性方程组,直接调用函数):⎪⎪⎪⎪⎪⎭⎫ ⎝⎛--------=19631699723723312312A ,⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-=71636b 解答:程序:A=[12 -3 2 1;-3 23 -7 -3;2 -7 99 -6;1 -3 -6 19];R=chol(A)b=[6 3 -16 7]';y=inv(R')*b %y=R'\bx=inv(R)*y %x=R\y结果:R =3.4641 -0.8660 0.5774 0.28870 4.7170 -1.3780 -0.58300 0 9.8371 -0.70850 0 0 4.2514y =1.73210.9540-1.59451.3940x =0.54630.2023-0.13850.3279问题 2、先将矩阵A 进行LU 分解,然后解方程组b Ax =(直接调用函数):⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----=8162517623158765211331056897031354376231A ,⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-=715513252b解答:程序:A=[1/3 -2 76 3/4 5;3 1/sqrt(3) 0 -7 89;56 0 -1 3 13;21 65 -7 8 15;23 76 51 62 81];b=[2/sqrt(5);-2;3;51;5/sqrt(71)];[L,U]=lu(A)y=inv(L)*bx=inv(U)*y结果:L = 0.0060 -0.0263 1.0000 0 00.0536 0.0076 -0.0044 0.1747 1.00001.0000 0 0 0 00.3750 0.8553 -0.6540 1.0000 00.4107 1.0000 0 0 0U =56.0000 0 -1.0000 3.0000 13.00000 76.0000 51.4107 60.7679 75.66070 0 77.3589 2.3313 6.91370 0 0 -43.5728 -50.06310 0 0 0 96.5050y =3.0000-0.63880.859850.9836-11.0590x =0.13670.90040.0526-1.0384-0.1146问题3、利用列主元的高斯消去法,求解下列方程组:⎪⎪⎩⎪⎪⎨⎧=+--=--+=-+-=+-+01002010100511.030520001.0204321432143214321x x x x x x x x x x x x x x x x解答:程序:function [RA,RB,n,X]=liezhu(A,b)B=[A b];n=length(b);RA=rank(A);RB=rank(B);zhica=RB-RA;if zhica>0disp('Çë×¢Ò⣺RA~=RB£¬ËùÒÔ´Ë·½³Ì×éÎ޽⡣')returnendif RA==RBif RA==ndisp('Çë×¢Ò⣺ÒòΪRA=RB=n,ËùÒÔ´Ë·½³Ì×éÓÐΨһ½â¡£')X=zeros(n,1);C=zeros(1,n+1);for p=1:n-1[Y ,j]=max(abs(B(p:n,p)));C=B(p,:);for k=p+1:nm=B(k,p)/B(p,p);B(k,p:n+1)=B(k,p:n+1)-m*B(p,p:n+1)endendb=B(1:n,n+1);A=B(1:n,1:n);X(n)=b(n)/A(n,n);for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);endelsedisp('Çë×¢Ò⣺ÒòΪRA=RB¡´n£¬ËùÒÔ´Ë·½³ÌÓÐÎÞÇî¶à½â¡£') endend键入A=[1 20 -1 0.0012 -5 30 -0.15 1 -100 -102 -100 -1 1];b=[0;1;0;0];[RA,RB,n,X]=liezhu(A,b)结果:请注意:因为RA=RB=n,所以此方程组有唯一解。
线性方程组的直接解法
![线性方程组的直接解法](https://img.taocdn.com/s3/m/ec9c1252cf84b9d528ea7a98.png)
2、 实验结果 求得������1 = 2; ������2 = 2; ������3 = 3
五、实验结论
1、 在计算机上直接高斯消去法是一种比较高效率的算法, 但是要求矩阵各顺 序主子式不能为 0 以保证主元不为 0. 2、通过选主元的方式可以保证高斯消去法的顺利进行。 3、满足对角占优或对称正定的矩阵方程可以更高效地求解。
三、实验分析 Problem1:
������1 + 2������2 + 3������3 = 14, 用列主元高斯消去法解方程组:{2������1 + 5������2 + 2������3 = 18, 3������1 + ������2 + 5������3 = 20. 1、 实验步骤 1 2 3 14 利用附录程序,输入增广矩阵[2 5 2 18],运行程序得如图结果 3 1 5 20
[附录]迭代法源程序(C 语言)
① 列主元 Gauss 消去法 #include<stdio.h> #include<math.h> #define M 4 void select(int k, double x[][M+1]) { double t=x[k][k],e[M][M+1]; int i,j,m; for(i=k;i<M;i++)//第 k 次选主元 { if(fabs(x[i][k])>fabs(t)) { t=x[i][k]; m=i; } else m=k; } for(j=k;j<M+1;j++) { e[k][j]=x[k][j]; x[k][j]=x[m][j]; x[m][j]=e[k][j]; } } int main() { int i,j,k; double t=0,a[M][M+1],l[M][M+1],x[M]; printf("Input the matrix:\n"); for(i=1;i<M;i++) { for(j=1;j<M+1;j++) scanf("%lf",&a[i][j]); } for(k=1;k<M-1;k++)//消元过程 { select(k,a);//第 k 次选主元 for(i=k+1;i<M;i++)//消元 { l[i][k]=a[i][k]/a[k][k];//消元因子 for(j=k+1;j<M+1;j++) a[i][j]=a[i][j]-l[i][k]*a[k][j];
实验2 线性方程组的直接解法
![实验2 线性方程组的直接解法](https://img.taocdn.com/s3/m/813d725b3c1ec5da50e27041.png)
实验2 线性方程组解的直接解法【实验目的】1、验证高斯消去法和三角分解法2、掌握直接求解线性方程组的常用算法:列主元高斯消去法、LU 分解法等3、记录运行结果,回答问题,完成实验报告【实验要求】根据题目要求,用Matlab 完成下列实验内容。
书写实验报告。
【实验内容】1、利用Matlab 求解线性方程组考虑下面给出的线性方程组形式为Ax =B 。
其中,A 和B 均为给定矩阵1111n1,n nn n a a b A B a a b ⎛⎫⎛⎫ ⎪ ⎪== ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭ 可以由给定的A 和B 矩阵构造出解的判定矩阵C (可以用C=[A B]命令来生成C )11111n n nn n a a b C a a b ⎛⎫ ⎪= ⎪ ⎪⎝⎭ 这里先可以给出线性方程组有解的判定定理:(2)当rank(A )=rank(C )=r <n 时,方程组有无穷多解。
具体求解方法留给同学们实验课后自学。
(3)当rank(A )<rank(C )时,方程组为矛盾方程,只能解出方程的最小二乘解。
这部分内容将在第四章学习。
【问题1】:求解下面给出的线性方程组。
记录所用命令和求解结果。
12345432141324341322x ⎛⎫⎛⎫ ⎪ ⎪ ⎪ ⎪= ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭ 【思考1】:Matlab 对线性方程组求解是否要考虑主元的选取?为什么?可以用下列例子试试 12120.000112x x x x +=⎧⎨+=⎩2、用Matlab 编写高斯消去法具体实现的函数要求输入参数包含方程组的系数矩阵及常数项向量,输出方程组的解,并使用问题1的例子进行调试。
高斯顺序消去法:设线性方程组AX =b 。
1.消元过程设0)(≠k kk a ,对1,,2,1-=n k 计算()()(1)()()(1)()()/,1,2,,k k ik ik kk k k k ij ij ik kj k k k i i ik k m a a a a m a i j k k n b b m b ++⎧=⎪=-=++⎨⎪=-⎩K 其中⒉回代过程()()()()()1/1,,2,1()/n n n n nn n i i i ii ij j ii j i x b a i n x b a x a =+⎧=⎪=-⎨=-⎪⎩∑K 其中。
线性方程组的直接解法实验报告
![线性方程组的直接解法实验报告](https://img.taocdn.com/s3/m/2aa0d4e148d7c1c709a1454e.png)
本科实验报告
课程名称:数值计算方法B
实验项目:线性方程组的直接解法
最小二乘拟合多项式
实验地点:ZSA401
专业班级:学号:201000
学生姓名:
指导教师:李志
2012年4月13日
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
printf("%lf\t",A[i][j]);
printf("\n");
}
double answer[N];
Gauss_eliminate(n,answer);
/*输出解*/
for(i=1;i<=n;i++)
printf("a[%d]=%lf\t",i-1,answer[i]);
getchar();
getchar();
}
四、实验结果与讨论、心得
讨论、心得:
刚开始调试代码的时候有时候就是很小的错误导致整个程序不能运行,需要我们一步一步慢慢来,经过无数次的检查程序错误的原因,以及在老师的帮助下,完成了这次实验。
这段时间的实验课提高了我的分析问题,解决问题的能力,特别提高了对一个程序的整。
解线性方程组的直接方法实验报告
![解线性方程组的直接方法实验报告](https://img.taocdn.com/s3/m/657bd9b0700abb68a982fbc3.png)
解线性方程组的直接方法实验报告解线性方程组的直接方法实验报告1.实验目的:1、通过该课题的实验,体会模块化结构程序设计方法的优点;2、运用所学的计算方法,解决各类线性方程组的直接算法;3、提高分析和解决问题的能力,做到学以致用;4、通过三对角形线性方程组的解法,体会稀疏线性方程组解法的特点。
2.实验过程:实验代码:#include "stdio.h"#include "math.h"#includeusing namespace std;//Gauss法void lzy(double **a,double *b,int n){int i,j,k;double l,x[10],temp;for(k=0;k<n-1;k++){for(j=k,i=k;j<n;j++){if(j==k)temp=fabs(a[j][k]);else if(temp<fabs(a[j][k])) {temp=fabs(a[j][k]);i=j;}}if(temp==0){cout<<"无解" ; return;}else{for(j=k;j<n;j++){temp=a[k][j];a[k][j]=a[i][j];a[i][j]=temp;}temp=b[k];b[k]=b[i];b[i]=temp;}for(i=k+1;i<n;i++){l=a[i][k]/a[k][k];for(j=k;j<n;j++)a[i][j]=a[i][j]-l*a[k][j];b[i]=b[i]-l*b[k];}}if(a[n-1][n-1]==0){cout<<"无解" ; return;}x[n-1]=b[n-1]/a[n-1][n-1]; for(i=n-2;i>=0;i--){temp=0;for(j=i+1;j<n;j++)temp=temp+a[i][j]*x[j];x[i]=(b[i]-temp)/a[i][i];}for(i=0;i<n;i++){printf("x%d=%lf ",i+1,x[i]);printf(" ");}}//平方根法void pfg(double **a,double *b,int n) {int i,k,m;double x[8],y[8],temp;for(k=0;k<n;k++){temp=0;for(m=0;m<k;m++)temp=temp+pow(a[k][m],2);if(a[k][k]<temp)return;a[k][k]=pow((a[k][k]-temp),1.0/2.0);for(i=k+1;i<n;i++){temp=0;for(m=0;m<k;m++)temp=temp+a[i][m]*a[k][m]; a[i][k]=(a[i][k]-temp)/a[k][k]; }temp=0;for(m=0;m<k;m++)temp=temp+a[k][m]*y[m];y[k]=(b[k]-temp)/a[k][k];}x[n-1]=y[n-1]/a[n-1][n-1];for(k=n-2;k>=0;k--){temp=0;for(m=k+1;m<n;m++)temp=temp+a[m][k]*x[m];x[k]=(y[k]-temp)/a[k][k];}for(i=0;i<n;i++){printf("x%d=%lf ",i+1,x[i]);printf(" ");}}//追赶法void zgf(double **a,double *b,int n){int i;double a0[10],c[10],d[10],a1[10],b1[10],x[10],y[10]; for(i=0;i<n;i++) {a0[i]=a[i][i];if(i<n-1)c[i]=a[i][i+1];if(i>0)d[i-1]=a[i][i-1];}a1[0]=a0[0];for(i=0;i<n-1;i++){b1[i]=c[i]/a1[i];a1[i+1]=a0[i+1]-d[i+1]*b1[i];}y[0]=b[0]/a1[0];for(i=1;i<n;i++)y[i]=(b[i]-d[i]*y[i-1])/a1[i];x[n-1]=y[n-1];for(i=n-2;i>=0;i--)x[i]=y[i]-b1[i]*x[i+1];for(i=0;i<n;i++){printf("x%d=%lf ",i+1,x[i]);printf(" ");}}int main{int n,i,j;double **A,**B,**C,*B1,*B2,*B3;A=(double **)malloc(n*sizeof(double)); B=(double **)malloc(n*sizeof(double));C=(double **)malloc(n*sizeof(double));B1=(double *)malloc(n*sizeof(double));B2=(double *)malloc(n*sizeof(double));B3=(double *)malloc(n*sizeof(double));for(i=0;i<n;i++){A[i]=(double *)malloc((n)*sizeof(double)); B[i]=(double *)malloc((n)*sizeof(double)); C[i]=(double *)malloc((n)*sizeof(double)); } cout<<"第一题(Gauss列主元消去法):"<<endl<<endl; cout<<"请输入阶数n:"<<endl;cin>>n;cout<<" 请输入系数矩阵: ";for(i=0;i<n;i++)for(j=0;j<n;j++){。
数值分析实验报告
![数值分析实验报告](https://img.taocdn.com/s3/m/f442da13cc175527072208d6.png)
%消元过程
fori=k+1:n
m=A(i,k)/A(k,k);
forj=k+1:n
A(i,j)=A(i,j)-m*A(k,j);
end
b(i)=b(i)-m*b(k);
end
det=det*A(k,k);
end
det=det*A(n,n);
%回代过程
ifabs(A(n,n))<1e-10
flag='failure';return;
*x=(x0,x1….,xn),插值节点
*y=(y0,y1,…,yn);被插函数f(x)在插值节点处的函数值
*t求插值函数Pn(x)在t处的函数值
*返回值 插值函数Pn(x)在t处的函数值
*/
procedureNewton
forj=0to n
d1jyj;
endfor
forj=1to n
fori=j to n
[n,m]=size(A);nb=length(b)
%当方程组行与列的维数不相等时,停止计算,并输出出错信息
ifn~=m
error('The row and columns of matrix A must beepual!');
return;
end
%当方程组与右端项的维数不匹配时,停止计算,并输出错误信息
clear
fprintf('gauss-seidel迭代法')
x1_(1)=0;
x2_(1)=0;
x3_(1)=0;
fori=1:9
x1_(i+1)=7.2+0.1*x2_(i)+0.2*x3_(i);
数值分析线性方程组的直接解法
![数值分析线性方程组的直接解法](https://img.taocdn.com/s3/m/89b1cb11aa00b52acec7ca50.png)
数值分析课程实验报告实验名称线性方程组的直接解法_____________________实验目的①掌握高斯消去法的基本思路和迭代步骤;②了解高斯消去法可能遇到的困难。
用文字或图表记录实验过程和结果列主元高斯消去法算法描述将方程组用增广矩阵B=[A:b]=(a j \心申)表示。
步骤1:消兀过程,对k=12|j|, n—1(1)选主元,找i k亡{k,k+1,川,n}使得k卜maxi a ikai k,(2)如果a i k,k = 0 ,则矩阵A奇异,程序结束;否则执行(3)。
(3)如果ik^k,则交换第k行与第i k行对应兀素位置,aq㈠a i k j,j=k,IH, n + 1。
(4)消兀,对i = k +1」H,n,计算m k=a k / a kk,对j = k +1,川,n +1,计算a j = a ij — m ik a^.步骤2:回代过程:(1)右a nn -0,则矩阵奇异,程序结束;否则执行(2)。
厲(2)nX n =a ng/a nn;对i = n—1川,2,1,计算X j = a,n 出一》a j X j /a H< j4 丿三、练习与思考题分析解答1、解方程组0.10伙2.304X2 3.555X3 =1.183-1.347为3.712X2 4.623X3 = 2.137-2.835X, 1.072X25.643X^3.035(1)编程用顺序高斯消去法求解上述方程组,记下解向量,验证所得到的解向量是否是原方程组的解,若不是原方程组的解,试分析原因,并证实你的分析的正确性!解:采用顺序消元法求得如下结果:请输入一个3行矩阵0.101 2.304 3.555 1.183-1.347 3.712 4.623 2.137-2.835 1.072 5.643 3.0350.101 2.304 3.555 1.1830 34.4396 52.0347 17.91420 0 6.09738 2.0435最后计算得到x =(-0.3982,0.0138,0.3351) T,代入原方程验证可知解向量是原方程组的解。
实验二 线性方程组的求解
![实验二 线性方程组的求解](https://img.taocdn.com/s3/m/5164705e3b3567ec102d8a34.png)
3.线性方程组解的情况
1)当 时,方程组无解;
2)当 时,方程组有唯一解;此时可采用Cramer法则或左除法来求解;
3)当 时,方程组有无穷多解;首先要求出特解,以及对应的齐次线性方程组得基础解系。
3.求线性方程组及矩阵方程所用的Matlab命令:
计算方阵A的行列式:det(A)
增广矩阵:[A,b],b是向量
>> %c1,c2为任意常数.
【练习与思考】
1.求解下列线性方程组
1) ;2) .
2.求解下列矩阵方程
1) ;2) .
3.求解下列方程:
1). 2).
4.求下列齐次线性方程组的基础解系并用它表示出全部解
5.讨论 取什么值时下列方程组有解,并求解
1)
2)
Matlab实验题(大作业)
>> b=[1;1;-1];%常数b
>> rank(A)%系数矩阵的秩
ans =
2
>> rank([A,b])%增广矩阵的秩
ans =
2
>> %由上可知系数矩阵的秩等于增广矩阵的秩,并且小于未知量的个数。所以方程有无穷多解.
>> %以下对应几种方法:(1)化行最简形,用rref命令
>> rref([A,b])
>> C1=A;C2=A;C3=A;C4=A;b=[8;9;-5;0];%将A赋值给不同变量
>> C1(:,1)=b;D1=det(C1);x1=D1/D%求解x1
x1 =
3
>> C2(:,2)=b;D2=det(C2);x2=D2/D%求解x2
x2 =
计算方法第三章线性方程组的直接解法
![计算方法第三章线性方程组的直接解法](https://img.taocdn.com/s3/m/12811e80195f312b3169a5f1.png)
5 3
3 1
r3
r1 6
6 1 18 2
1 0
4 5 1 3
3 1
r3 r225
1 0
4 1
5 3
3 1
0 25 48 16
0 0 27 9
林龙
计算方法
6
化原方程组为三角方程组的过程为消元过程. 解三角方程组的过程为回代过程.
也可将上边的增广矩阵进一步化简.
1 4 5 3
1 0 7 1
xi
Di D
(i
1, 2,3,
),由于方程含有n 1个
行列式.如对每个行列式按展开定理来计算.
用克莱姆法则求解,所需要的乘除运算量为
n!(n2 1) n次,若n 20用每秒一千万次的
计算机要三百万年,所以并不是凡直接法都
可以用来做实际运算.
林龙
计算方法
4
设有
§3.1直接法
a11x1 a12 x2 a21x1 a22 x2
解 : 10
7
0
7
r1 r2
5 1 5 6
林龙
计算方法
16
10 3 5
7 2 1
0 6 5
7 4 6
r2
3 10
r1
r3
5 10
r1
10
0
0
7 0.1 2.5
0 7 6 6.1 5 2.5
r2 r3
r3
1 25
r2
10 7 0 7 x3 1
0
2.5
5
2.5
x2
2.5 5x
nn
a11 a12 .... a1n 1 0 0
a21
a22
第3章 线性方程组求解的直接解法
![第3章 线性方程组求解的直接解法](https://img.taocdn.com/s3/m/52f9da94ad51f01dc281f1e8.png)
线性方程组求解的直接法5.2线性方程组直接解法概述直接解法就是利用一系列公式进行有限步计算,直接得到方程组的精确解的方法.当然,实际计算结果仍有误差,譬如舍入误差,而且舍入误差的积累有时甚至会严重影响解的精度.这是一个众所周知的古老方法,但用在计算机上仍然十分有效.求解线性方程组最基本的一种直接法是消去法.消去法的基本思想是,通过将一个方程乘以或除以某个常数,以及将两个方程相加减这两种手段,逐步减少方程中的变元的数目,最终使每个方程仅含一个变元,从而得出所求的解.高斯(Gauss )消去法是其中广泛应用的方法,其求解过程分为消元过程和回代过程两个环节.消元过程将所给的方程组加工成上三角方程组,所归结的方程组再通过回代过程得出它的解.Gauss 消去法由于添加了回代的过程,算法结构稍复杂,但这种改进的算法明显减少了计算量.直接法比较适用于中小型方程组.对高阶方程组,即使系数矩阵是稀疏的,但在运算中很难保持稀疏性,因而有存储量大,程序复杂等不足.5.3直接解法5.3.1Gauss 消去法Gauss 消去法是一个古老的求解线性方程组的方法,由它改进而来的选主元法是目前计算机上常用的有效的求解低阶稠密矩阵线性方程组的方法.例5.1用Gauss 消去法解方程组1231231232221(5.3.1)1324 (5.3.2)2539(5.3.3)2x x x x x x x x x ⎧++=⎪⎪++=⎨⎪++=⎪⎩解〖JP4〗第1步,式35.3.12⨯-()()加到式(5.3.2)上,式()15.3.1()2⨯-加到式(5.3.3)上,得到等价方程组123232322211(5.4.4)282(5.4.5)x x x x x x x ⎧++=⎪⎪-+=-⎨⎪⎪+=⎩第2步,式()2⨯5.3.4加到式(5.3.5)上得等价的方程组12323322211100(5.3.6)x x x x x x ++=⎧⎪-+=-⎨⎪=⎩第3步,回代法求解方程组(5.3.6),即可求得该方程组的解为32110,1,.2x x x ===-.用矩阵描述其约化过程即为233(2)22221011100100r r r ⨯+⇒⎡⎤⎢⎥--⎢⎥⎢⎥⎣⎦→[]122133(1)3()21()222212221,3241/201111395/20282r r r r r r A b ⨯-+⇒⨯-+⇒⎡⎤⎡⎤⎢⎥⎢⎥=--⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦→.这种求解过程称为具有回代的Gauss 消去法.由此例可见,Gauss 消去法的基本思想是:用矩阵的初等行变换将系数矩阵A 化为具有简单形式的矩阵(如上三角阵、单位矩阵等),而三角形方程组是很容易回代求解的.一般地,设有n 个未知数的线性方程组为11112211211222221122n n n n n n nn n na x a x a xb a x a x a x b a x a x a x b +++=⎧⎪+++=⎪⎨⎪⎪++=⎩L L MM M L (5.3.7)1212)(,,)(,,)T T ij n n n n A a X x x x b b b b ⨯===L L (,,,则方程组(5.3.7)化为AX b =.方便起见,记()(1)det 0A AA ==≠,(1)b b =,且()1A的元素记为()()11,ij a b ,的元素记为()1i b ,则消去法的步骤如下:第1步:1110a≠(),,计算(1)11(1)11(2,3,4),i i a m i n a ==L 用()1i m -乘方程组(5.3.7)中的第1个方程加到第i个方程中()2,3,i n =L ,即进行行初等变换()112,3,i i i R m R R i n -⋅→=L ,消去第2个到第n个方程中的未知数1,x ,得等价方程组111121(2)(2)(2)22222(2)(2)(2)2inn n n nn n x a a b x a a b ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦LMM LM M L (5.3.8)记为(2)(2)A X b =,其中(2)(1)(1)(2)(1)(1)1111(,2,3),2,3,ij ij i j i i i a a m a i j n b b m b i n =-==-=L L ,,第k 步()1,2,1k n =-L:继续上述消元过程.第1步到第1k -步计算已完成,且得到与原方程组等价的方程组(1)(1)(1)(1)1112111(2)(2)(2)222223()()()()()()nn k k k kkkn k n k k k nk nn n a a a b x a a b xx aa b x a a b ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎣⎦⎣⎦L L LLOM L M MMM L(5.3.9)记为()(()K k A X b =,进行第k 步消元:设()0k kka≠,计算乘数()()(1,)k ikk ik kka m k k n a ==+L ,用ik m -乘方程组(5.3.9)中第k 个方程加到第i 1)i k n =+L (,,,个方程上消去方程组(5.3.9)中第i 1)i k n =+L (,,个方程的未知数k x ,得到与原方程组等价的方程组:(1)()()(1)()()(1)(1)()(,1,)( 1.)k k k ij ij ik kj k k k i i ik k k k k k a a m a i j k n b b m b i k n A A k b b k ++++⎧=-=+⎪=-=+⎨⎪⎩L L ()与前行元素相同,与前个元素相同 (5.3.10) 记为(1)(1)k k A X b ++=其中(1)(1,k k A b ++)中元素计算公式为(1)()()(1)()()(1)(1)()(,1,)( 1.)k k k ij ij ik kj k k k i i ik k k k k k a a m a i j k n b b m b i k n A A k b b k ++++⎧=-=+⎪=-=+⎨⎪⎩L L ()与前行元素相同,与前个元素相同 (5.3.11)重复上述过程,且设()0(1,2,1)k kk a k n ≠=-L ,共完成1n -步消元计算,得到与方程组(5.3.7)等价的三角形方程组1111211(2)(2)(2)22222()()n n n n n nn n x a a b x a b ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦LMOM M (5.3.12)再用回代法求方程组(5.3.12)的解,计算公式为()()()()1()(),(1,2,1)n n n nn n i i i ij j j i i i ii b x a b a x x i n n a =+⎧=⎪⎪⎨-⎪==--⎪⎩∑L (5.3.13)元素()k kka 称为约化的主元素.将方程组(5.3.7)化为方程组(5.3.12)的过程称为消元过程.方程组(5.3.12)的求解过程(5.3.13)称为回代过程.由消元过程和回代过程求解线性方程组的方法称为Gauss 消去法.定理5.1(Gauss 消去法)设AX b =。
线性方程组的直接解法迭代解法
![线性方程组的直接解法迭代解法](https://img.taocdn.com/s3/m/6f66dad6b52acfc788ebc97d.png)
广东金融学院实验报告课程名称:数值分析实验目的及要求实验目的:题一:通过数值实验,从中体会解线性方程组选主元的必要性和LU分解法的优点,以及方程组系数矩阵和右端向最的微小变化对解向最的影响。
比较各种直接接法在解线性方程组中的效果;题二:认识齐种迭代法收敛的含义、影响齐迭代法收敛速度的因素。
实验要求:题一:(1)在MATLAB中编写程序用列主元高斯消去法和LU分解求解上述方程组,输出曲b中矩阵A 及向量b和A二LU分解中的L及U, detA及解向量X.(2)将方程组中的2. 099999改为2. 1, 5. 900001改为5. 9,用列主元高斯消去法求解变换后的方程组,输出解向最x及detA,并与(1)中的结果比较。
(3)用MATLAB的内部函数inv求出系数矩阵的逆矩阵,再输入命令x=inv(A)*b,即可求出方程组的解。
请与列主元高斯消公法和LU分解法求出的解进行比较,体会选主元的方法具有良好的数值稳定性。
用MATLAB的内部曲数det求出系数行列式的值,并与(1)、(2)中输出的系数行列式的值进行比较。
(4)比较以上各种直接解法在解线性方程组中的效果。
题二:(1)选取不同的初始向M:X(0)及右端向最b,给泄迭代误差要求,用Jacobi迭代法和Gauss-Seidel迭代法求解,观察得到的序列是否收敛?若收敛,记录迭代次数,分析计算结果并得出你的结论。
列岀算法清单。
(2)用SOR迭代法求上述方程组的解,松弛系数血取1<69<2的不同的三个值,在< 10"5时停止迭代,记录迭代次数,分析计算结呆与松弛系数血的关系并得出你的结论。
(3)用MATLAB的内部函数inv求出系数矩阵的逆矩阵.再输入命令^inv(A)*b>即可求出上述各个方程组的解.并与上述三种方法求出的解进行比较。
请将比较结果列入卜表。
方程组的解X1 Xr■迭代次数误差楮确解Jacibi解法Gause・seidel 解法SOR 解法00= 60= 60=实验环境及相关情况(包含使用软件、实验设备、主要仪器及材料等)1. Win72. Mat lab 7.0实验内容及步骤(包含简要的实验步骤流程) 实验内容:题一:解卜列线性方程组'10 -7‘X 】、(8、-3 2.099999 62Xr5.9000015-1 5 -15、12> 0< 1 >题二研究解线性方程组 做=b 迭代法的收敛性、收敛速度以及SOR 方法中/佳松弛因子的选取问题, 用迭代法求解做二b,其中・4 -1r■7 A=4 -81 ,b =-21-2 ■1515实验结果(包括程序或图表、结论陈述.数据记录及分析等,可附页)题一:直接解法解线性方程组(1)列主兀高斯消去法与LU 分解求解列主元高斯消去法:编写matalab 程序(见附录gaosi.m ),输出矩阵10.000 -7.000 0.000= 0.000 2.5000-5.000一 0.000 0.0006.0000020.000 0.000 0.000向量8 b =1 8.300 L5.0800J解向量:X = (0 ・-1 , 1 r I )7 其中系数行列式的值det (A )=762.00009LU 分解求解:编J matalab 程序(见附录zhjLU. m 和LU ・m ),执行输出:-1.5 2.300 5.080-3.0001.000000.00000.5000 -25000001.0000 0.2000 -24000000.9600 10.0000 -7.0000 0.0000 1.0000n = 0.0000-0.0000010.0000 2.3000 —0.0000 0.000015000000 57500000.0000 0.0000 0.0000 5.0800在matlab 命令窗II 输入L*U ,可以得到A 二L*U ,即分解结果正确。
线性方程组直接解法实验
![线性方程组直接解法实验](https://img.taocdn.com/s3/m/12ae67479b89680203d825d2.png)
实验一 线性方程组直接解法实验一、实验目的1.运用matlab 软件完成线性方程组的直接实验;2.通过实验,了解Doolittle 分解方法和列主元消去法解方程组的过程,并比较两种方法的优点。
二、实验题目分别用Doolittle 分解方法和列主元消去法解方程组123410-7018-3 2.09999962 5.9000015-15-1521021⎛⎫⎛⎫⎛⎫⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭x x x x . 输出A ,b ;Doolittle 分解方法的L 和U ;解向量x,det A ;列主元方法的行交换次序,解向量x,det A ;比较两种方法所得的结果。
三、实验原理1) Doolittle 分解方法的原理算法原理:应用高斯消去法解n 阶线性方程Ax b =经过1n -步消去后得出一个等价的上三角形方程组()()n n A x b =,对上三角形方程组用逐步回代就可以求出解来。
这个过程也可通过矩阵分解来实现。
将非奇异阵分解成一个下三角阵L 和上三角阵U 的乘积称为对矩阵A 的三角分解,又称LU 分解。
根据LU 分解,将Ax b =分解为Ly bUx y =⎧⎨=⎩形式,简化了求解问题。
程序框图:变量说明:ij a 为系数矩阵元素,i b 为常数矩阵系数,,ij ij l u 分别为下、上三角矩阵元素。
2)列主元消去法解方程组的原理算法原理:列选主元是当变换到第k步时,从k列的kk a及以下的各元素中选取绝对值a的位置上,然后再进行消元过程。
交换系数矩阵中最大的元素,通过行交换将其交换到kk的两行(包括常数项),相当于两个方程的位置交换了。
程序框图:Array变量说明:k表示消元到a为消元第k步时第k步,kk主对角线元素3)四、实验过程及结果1)Doolittle分解方法的输出结果----------计算实习题----------Page64 第1题用Doolittle分解方法解方程组A =10.0000 -7.0000 0 1.0000-3.0000 2.1000 6.0000 2.00005.0000 -1.0000 5.0000 -1.00002.0000 1.0000 0 2.0000b =8.00005.90005.00001.0000L =1.0e+006 *0.0000 0 0 0-0.0000 0.0000 0 00.0000 -2.5000 0.0000 00.0000 -2.4000 0.0000 0.0000 U =1.0e+007 *0.0000 -0.0000 0 0.00000 -0.0000 0.0000 0.00000 0 1.5000 0.57500 0 0 0.0000 X =-0.0000-1.00001.00001.0000det(A)值为-762.00009000----------输出完毕----------2)列主元消去法输出结果----------计算实习题----------Page64 第1题列主元消去法解方程组A =10.0000 -7.0000 0 1.0000-3.0000 2.1000 6.0000 2.00005.0000 -1.0000 5.0000 -1.00002.0000 1.0000 0 2.0000b =8.00005.90005.00001.0000X =0.0000-1.00001.00001.0000detA值为-762.00009000----------输出完毕----------五、实验分析1.运用LU分解法可以成批地解方程组,且速度快.若c先求LU=A3,再解(LU)x=b,则要重新计算,计算量增加;如果按照上述方法计算,能够减少运算次数,加快运算速度.3. ⑴无论当n=10、n=100、n=1000时,x1与x2的值都相等,且随着n的增大,变化的只是解的中间部分数字,头、前后几位数都没有变化.⑵高斯消去法应用于三对角方程组得到的就是所谓的“追赶法”.追赶法不需要对零元素计算,只有6n-5次乘除法计算量,求解速度快.且当系数矩阵对角占优时数值稳定,是解三对角方程组的优秀解法.⑶用LU分解法解此方程组速度慢.顺序高斯消去法实际上就是将方程组的系数矩阵分解成单位下三角矩阵与上三角矩阵的乘积.顺序高斯消去法的消元过程相当于LU分解过程和Ly=b的求解,回代过程则相当于解线性方程组Ux=y,故其求解速度慢.六、附原程序1)Doolittle分解方法原程序fprintf('----------计算实习题----------\n')fprintf('Page64 第1题用Doolittle分解方法解方程组\n')A=[10 -7 0 1 ; -3 2.099999 6 2 ;5 -1 5 -1 ; 2 1 0 2];b=[8;5.900001;5;1];n=length(A);U=zeros(n,n);L=eye(n,n);U(1,:)=A(1,:);L(2:n,1)=A(2:n,1)/U(1,1);for i=2:n;U(i,i:n)=A(i,i:n)-L(i,1:i-1)*U(1:i-1,i:n);L(i+1:n,i)=(A(i+1:n,i)-L(i+1:n,1:i-1)*U(1:i-1,i))/U(i,i); endY=zeros(n);Y(1)=b(1);for i=2:nY(i)=b(i)-L(i,1:i-1)*Y(1:i-1,1);endX=zeros(n,1);if det(U)==0;X=0;elseX(n)=Y(n)/U(n,n);for i=n-1:-1:1X(i)=(Y(i)-U(i,i+1:n)*X(i+1:n,1))/U(i,i);endendAbLUXfprintf('det(A)值为%9.8f\n',det(A))fprintf('----------输出完毕 ----------\n')2)列主元消去法原程序fprintf('----------计算实习题----------\n')fprintf('Page64 第1题列主元消去法解方程组\n')A=[10 -7 0 1 ; -3 2.099999 6 2 ;5 -1 5 -1 ; 2 1 0 2];b=[8;5.900001;5;1];C=[A b];n=length(A);D=zeros(n,n+1);l=zeros(n,1);for i=1:nD=C;k=min(find(C(i:n,i)==max(C(i:n,i))));C(i,i:n+1)=D(k+i-1,i:n+1);C(k+i-1,i:n+1)=D(i,i:n+1);l(i+1:n,1)=C(i+1:n,i)/C(i,i);C(i+1:n,i:n+1)= C(i+1:n,i:n+1)- l(i+1:n,1)*C(i,i:n+1); endX=zeros(n,1);X(n)=C(n,n+1)/C(n,n);for i=n-1:-1:1X(i)=(C(i,n+1)-C(i,i+1:n)*X(i+1:n,1))/C(i,i); endAbXfprintf('detA值为%9.8f\n',det(A))fprintf('----------输出完毕----------\n')。
第3章 线性方程组的直接解法
![第3章 线性方程组的直接解法](https://img.taocdn.com/s3/m/bf46881614791711cc791751.png)
第三章 线性方程组的数值方法在自然科学和工程技术中很多问题的解决常常归结为求解线性代数方程组.⎪⎪⎩⎪⎪⎨⎧=++=++=++nnnn 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 22112222212111212111当它的系数行列式不为零时,由克莱姆法则可以给出方程组的唯一解,但是这一理论上完 善的结果,在实际计算中可以说没有什么用处。
因此如何建立在计算机上可以实现的有效而实用的解法,具有极其重要的意义。
这些方法大致可分为两类:一类是直接法,就是经过有限步算术运算,可求得方程组精确解的方法(如果每步计算都是精确进行的话);另一类是迭代法,就是用某种极限过程去逐步逼近其精确解的方法.本章将阐述这两类算法中最基本的高斯消元法及其变形、矩阵分解法、雅可比迭代法、高斯 -塞德尔迭代法等.第一节 高斯消元法一 回代过程设系数矩阵为n 阶上三角矩阵的线性方程组(1)2222211212111⎪⎪⎩⎪⎪⎨⎧==+=++nn nn n n n n b x a b x a x a b x a x a x a如果a ,...a ,a nn 2211都(1)自下而上可以逐次求出x ,...x ,x n n 11-为()21211⎪⎪⎩⎪⎪⎨⎧--=∑-==+=,...n ,n k ,a x a b x a b x kk jn k j kj k k nn n n按上述公式求方程组(1) 解的过程称为回代过程.不难看出,解方程组(1)共需()121+n n 次加法和()121-n n 次乘法.这恰好是用一个n 阶三角方阵乘n 维向量所需的运算次数。
当n 较大时,n n >>2,同时加法运算速度远快于乘法的运算速度,所以,可用n 221次乘法来近似表示回代过程的运算量.二 消元过程设有线性方程组()322112222212*********⎪⎪⎩⎪⎪⎨⎧=++=++=++n n 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 为了符号统一,记()()()n ...,j ;n ,...,i b a ,a a i in ij ij2121010====+, 则原方程组改写成()()()()()()()()()()()()()301020*******2022*********10120121011'nnnnn n n n )n nn n a x a x a x a a x a x a x a a x a x a x a ⎪⎪⎩⎪⎪⎨⎧=++=++=+++++如果()0011≠a ,那么就可以保留其中第一个方程并利用它分别与其余方程消去第一个未知量.令()()n,..,,i ,a a l i i 32011011==()4则以l i 1-乘第一个方程加到第i个方程中,就把方程组()3'化为()()()()()()()()()()()51112121121221220110120121011⎪⎪⎩⎪⎪⎨⎧=+=+=+++++a x a x a a x a x a a x a x a x a nn nnn n n nn n n n其中()()()1323201101+==-=n ,...,j ,n ,...,i ,a l a a j i ij ij ()6由方程组(3′)化为(5)的过程中,元素()a 011起着特殊的作用,特把元素()a 011称为主元素.如果方程组(5)中()0122≠a , 则以()a 122为主元素,并利用类似的方法消去第n ,...,43个方程中的第二个未知量,即令()()n ,..,,i ,a a l i i 43122122== ()7则以l i 2-乘以第二个方程加到第i 个方程中,于是得到新的方程组()()()()()()()()()()()()()()()()8212323213233233112123123212201101301320121011⎪⎪⎪⎩⎪⎪⎪⎨⎧=+=++=+++=++++++++a x a x a a x a ....x a a x a x a x a a x a x a x a x a nn n nn n n n n n n n n n n其中()()()13312212+==-=n ,...j ,n ,...i ,a l a a j i ij ij ()9重复上述过程1-n 步后,我们得到原方程组等价的系数矩阵为三角形方阵的方程组()()()()()()()()()()()()()()()10111213233233112123123212201101301320121011⎪⎪⎪⎩⎪⎪⎪⎨⎧==++=+++=++++-+-+++a x a a x a ....x a a x a x a x a a x a x a x a x a n nn n n nn n n n n n n n n n其中()()a a l k kk k ikik 11--= ()11()()()⎪⎩⎪⎨⎧+++=++=-=-=--1212112111n ,...k ,k j ;n ,...k ,k i n ,...,k a l a a k kj ik k ij k ij()12把方程组(3)逐步化为方程组(10)的过程称为消元过程.最后,由回代过程可求得原方程组的解为()()()()()⎪⎪⎩⎪⎪⎨⎧--=∑-==-+=--+--+122111111111,,...n .n k ,a x a a x a a x k kk n k j j k kjk kn k n nnn nn n ()13这种通过消元、再回代的求解方法称为高斯(Gauss)消元法(其特点是始终消去主对角线下方的元素).注意到,上标k 仅仅用来识别一次消元前后系数矩阵的变化,而()a k ij 变为()a k ij 1+后,()a k ij不再使用,所以在计算机存贮中只要用()a k ij1+冲掉()a k ij即可;另一方面,主元素所在列中主元素下面的各元素在消元过程中必然是零,而且在后面将要列出的回代过程中也不用它们,所以没有必要通过计算得到它们,从而在消元过程中j 就可从1+k 开始,这样做还可以节约计算时间.例1 用高斯消元法求解方程组⎪⎩⎪⎨⎧=++=-+=++52213614282321321321x x x x x x x x x解 用第一个方程消去后两个方程中的x 1,得 ⎪⎩⎪⎨⎧-=-=-=++9962214282232321x x x x x x 再用第二个方程消去第三个方程中的x 2,得⎪⎩⎪⎨⎧=-=-=++18962214282332321x x x x x x 最后,经过回代求得方程组的解为52281412262918321323=--==+=-=-=x x x x x x 三 高斯消元法的条件与运算量从消元过程可以看出,对于n 阶线性方程组,只要各步主元素不为零,即()01≠-a k kk ,经过1-n 步消元,就可以得到一个等价的系数矩阵为上三角形阵的方程组,然后再利用回代过程可求得原方程组的解.因此,有下面结论.定理1 如果在消元过程中A 的主元素不为零,即()01≠-a k kk (k=1,2,…,n),则可通过高斯消元法求出b Ax =的解.矩阵A 在什么条件下才能保证()01≠-a k kk ,下面的定理给出了这一条件.引理 在高斯消元过程中系数矩阵A 的主元素不为零,即()01≠-a k kk (k=1,2,…,n)的充要条件是矩阵A 的各阶顺序主子式不为零,即n...,k ,a ...a .........a ...a D a D kkk kk 32001111111=≠=≠=证明 首先利用归纳法证明引理的充分性,显然 ,当1=n 时,引理的充分性是成立的,现假设引理对1-n 时也成立,求证引理对n 也成立,由归纳法假设有()01≠-a k kk , 121-=n ...,k于是可用高斯消元法将A 化为⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛→--------)n (nn)n (n ,n )n (n n)(n )(n )()(n )(n )()(a a a a a a a a a a A 1212111211212201011012011且()()a a D 0110111==()()()()()a a a a a D 1220111220120112==()()()a ...a a a a a a a a a a D n nn)n (nn )(n )(n )()(n )(n )()(n 112201111211212201011012011----=⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=由假设,n ,...,k ,D k 210=≠所以有()01≠-a n nn .反过来,由上式可知必要性是显然的.定理2 如果n 阶矩阵A 的所有顺序主子式均不为零,即,n ,...,k ,D k 210=≠则可通过高斯消元法求出的解.下面考虑求解(3)的高斯消元法的运算量.消元过程需要除法()()()1211221-=+++-+-n n ...n n次,而需要的乘法和加法的次数都是()()()()131122112-=⋅++-⋅-+-⋅n n ...n n n n加上回代过程的运算次数,共需乘、除法的次数为()()()()133112113112122-+=++-+-n n n n n n n n n加、减法的次数为()()()5326112113122-+=-+-n n n n n n n 当n 较大时,n n 23>>,消元过程的运算量远大于回代过程,从而,高斯消元法中乘除法的次数与加减法的次数近似为n331.第二节 第二节 高斯主元素消元法一 问题的提出由高斯消元法可知,在消元过程中如果出现()01=-a k kk 的情况,这时消元法将无法进行;另一方面,即使主元素()01≠-a k kk ,但很小时,用其作除数,会导致其它元素数量级的严重增长和舍入误差的扩散,最后也使得计算结果很不可靠.例1 例1 解方程组⎩⎨⎧=+=+0000100001000010001200003000302121.x .x ..x .x .(它的精确解为323121==x ,x )解法一 用高斯消元法求解(取5位有效数字),用第一个方程消去第二个方程中的x 1得⎩⎨⎧-=-=+0666609999000120000300030221.x ..x .x . 因而再回代,得000300666700003000120666799990666612=⨯-=≈--=....x ...x 显然,这个解与精确解相差太远,不能作为方程组的近似解其原因是我们在消元过程中使用了小主元素,使得约化后的方程组中的元素量级大大增长,再经舍入使得计算中舍入误差扩散,因此经消元后得到的三角形方程组就不准确了.为了控制舍入误差,我们采用另一种消元过程. 解法二 为了避免绝对值很小的元素作为主元,先交换两个方程,得到⎩⎨⎧=+=+0001200003000300000100001000012121.x .x ..x .x .消去第二个方程中的x 1,得⎩⎨⎧==+9998199972000010000100001221.x ..x .x .再回代,解得()3333066670000010000166670999729998112....x ...x =⨯-=≈=结果与准确解非常接近.这个例子告诉我们,在采用高斯消元法解方程组时,用做除法的小主元素可能使舍入误差增加,主元素的绝对值越小,则舍入误差影响越大。
实验2_求解线性方程组直接法(完成版)
![实验2_求解线性方程组直接法(完成版)](https://img.taocdn.com/s3/m/4cd28930b7360b4c2e3f6453.png)
数值分析实验报告二求解线性方程组的直接方法(2学时)一 实验目的1.掌握求解线性方程组的高斯消元法及列主元素法; 2. 掌握求解线性方程组的克劳特法; 3. 掌握求解线性方程组的平方根法。
二 实验内容1.用高斯消元法求解方程组(精度要求为610-=ε):1231231233272212240x x x x x x x x x -+=⎧⎪-+-=-⎨⎪-+=⎩ 2.用克劳特法求解上述方程组(精度要求为610-=ε)。
3. 用平方根法求解上述方程组(精度要求为610-=ε)。
4. 用列主元素法求解方程组(精度要求为610-=ε):1231231233432222325x x x x x x x x x -+=⎧⎪-+-=⎨⎪--=-⎩ 三 实验步骤(算法)与结果1. 高斯消元法求解根据算法思想用C 语言编程(源程序见附录2.1) 编译结果如下图:2. 克劳特法求解根据算法思想用C 语言编程(源程序见附录2.2)编译结果如下图:3平方根法求解根据算法思想用C语言编程(源程序见附录2.3)编译结果如下图:4列主元素法求解根据算法思想用C语言编程(源程序见附录2.4)编译结果如下图:四实验收获与教师评语1.实验收获:对于这次实验,我可以锻炼到上机实验的能力,并且第一次感受到数学知识在现实生活中的应用,也是第一次运用计算机解决数学问题。
另外,正是因为这次上机实验,让我重温了有些遗忘的编程知识。
2.教师评语:附录:2.1至2.4源程序代码2.1 #include <stdio.h>#define N 8void main(){floatsum,a[N][N]={0},u[N][N]={0},l[N][ N]={0},z[N]={0},x[N]={0};int n,i,j,k;printf("input the number of roots:");/*方程的个数小于8*/scanf("%d",&n);printf("input xishu matrix:\n"); for(i=1;i<=n;i++)for(j=1;j<=n+1;j++)scanf("%f",&a[i][j]);for(i=1;i<=n;i++)l[i][i]=1;for(i=1;i<=n;i++)for(j=1;j<=n;j++){if(i>j){for(k=0,sum=0;k<=j-1;k++)sum+=l[i][k]*u[k][j];l[i][j]=(a[i][j]-sum)/u[j][j];}else{for(k=0,sum=0;k<=i-1;k++)sum+=l[i][k]*u[k][j];u[i][j]=a[i][j]-sum; }}for(i=1,j=n+1;i<=n;i++){for(k=0,sum=0;k<=i;k++)sum+=l[i][k]*z[k];z[i]=a[i][j]-sum;}for(i=n;i>=1;i--){sum=0;for(k=i+1;k<=n;k++)sum+=u[i][k]*x[k]; x[i]=(z[i]-sum)/u[i][i];}printf("changing fangcheng xi shu is:\n");for(i=1;i<=n;i++){for(j=1;j<=n;j++)printf("%-10g",u[i][j]); printf("%-10g\n",z[i]);}for(i=1;i<=n;i++)printf("x(%d)=%-g\n",i,x[i]); getch();}2.2#include <conio.h>#include <stdio.h>#define N 8void main(){floatsum,a[N][N]={0},u[N][N]={0},l[N][ N]={0},z[N]={0},x[N]={0};int n,i,j,k;printf("input the number of roots:");/*方程的个数小于8*/scanf("%d",&n);printf("input xishu matrix:\n"); for(i=1;i<=n;i++)for(j=1;j<=n+1;j++)scanf("%f",&a[i][j]);for(i=1;i<=n;i++)u[i][i]=1;for(i=1;i<=n;i++)for(j=1;j<=n;j++){if(i>=j){for(k=0,sum=0;k<=j-1;k++)sum+=l[i][k]*u[k][j];l[i][j]=a[i][j]-sum; }else{for(k=0,sum=0;k<=i-1;k++)sum+=l[i][k]*u[k][j];u[i][j]=(a[i][j]-sum)/l[i][i];}}for(i=1,j=n+1;i<=n;i++){for(k=0,sum=0;k<=i;k++)sum+=l[i][k]*z[k];z[i]=(a[i][j]-sum)/l[i][i];}for(i=n;i>=1;i--){sum=0;for(k=i+1;k<=n;k++)sum+=u[i][k]*x[k];x[i]=z[i]-sum;}printf("changing fangcheng xi shu is:\n");for(i=1;i<=n;i++){for(j=1;j<=n;j++)printf("%-10g",u[i][j]); printf("%-10g\n",z[i]);}for(i=1;i<=n;i++)printf("x(%d)=%-g\n",i,x[i]); getch();}2.3#include <conio.h>#include <stdio.h>#include <math.h>#define N 8void main(){floatsum,a[N][N]={0},u[N][N]={0},l[N][ N]={0},z[N]={0},x[N]={0};int n,i,j,k;printf("input the number of roots:");/*方程的个数小于8*/ scanf("%d",&n);printf("input xishu matrix:\n"); for(i=1;i<=n;i++)for(j=1;j<=n+1;j++)scanf("%f",&a[i][j]);for(i=1;i<=n;i++)for(j=1;j<=n;j++){if(i>j){for(k=0,sum=0;k<=j-1;k++)sum+=l[i][k]*l[j][k];l[i][j]=(a[i][j]-sum)/l[j][j];u[j][i]=l[i][j];}if(i==j){for(k=0,sum=0;k<=i-1;k++)sum+=l[i][k]*l[i][k];l[i][i]=sqrt(a[i][i]-sum);u[i][i]=l[i][i];}}for(i=1,j=n+1;i<=n;i++){for(k=0,sum=0;k<=i;k++)sum+=l[i][k]*z[k];z[i]=(a[i][j]-sum)/l[i][i];}for(i=n;i>=1;i--){sum=0;for(k=i+1;k<=n;k++)sum+=u[i][k]*x[k];x[i]=(z[i]-sum)/l[i][i];}printf("changing fangcheng xi shu is:\n");for(i=1;i<=n;i++){for(j=1;j<=n;j++)printf("%-10g",u[i][j]); printf("%-10g\n",z[i]);}for(i=1;i<=n;i++)printf("x(%d)=%-g\n",i,x[i]); getch();}2.4#include <stdio.h>#include <math.h>#define N 8void main(){ int i,j,k,n;floata[N][N],l[N][N],u[N][N],z[N],x1[N ],max,x;for(i=0;i<N;i++)for(j=0;j<N;j++)a[i][j]=u[i][j]=l[i][j]=0;z[0]=0; printf("input the number of roots:");scanf("%d",&n);printf("input xishu matrix:\n");for(i=1;i<=n;i++)for(j=1;j<=n+1;j++)scanf("%f",&a[i][j]);for(i=1;i<n;i++){ max=fabs(a[i][i]);k=i;for(j=i+1;j<=n;j++)if(fabs(a[j][i])>max){ max=fabs(a[j][i]);k=j;} for(j=i;j<=n+1;j++){ x=a[i][j];a[i][j]=a[k][j];a[k][j]=x;}for(j=i+1;j<=n;j++){ l[j][i]=a[j][i]/a[i][i];a[j][i]=0;for(k=i+1;k<=n+1;k++) a[j][k]=a[j][k]-l[j][i]*a[i][k]; }z[i]=a[i][n+1];}z[n]=a[n][n+1];for(i=1;i<=n;i++)for(j=1;j<=n;j++)u[i][j]=a[i][j];for(i=n;i>0;i--){ k=i+1;x=0;while(k<=n){ x+=u[i][k]*x1[k];k++;}x1[i]=(z[i]-x)/u[i][i];}for(i=1;i<=n;i++)for(j=1;j<=n;j++){ printf("%-10g",u[i][j]);if(j==n)printf("%-10g\n",z[i]);}for(i=1;i<=n;i++)printf("x%d=%-10g\n",i,x1[i]);getch();}。
计算方法实验:解线性方程组的直接法
![计算方法实验:解线性方程组的直接法](https://img.taocdn.com/s3/m/278c3a8b3b3567ec112d8a9b.png)
实验二解线性方程组的直接法一、实验目的用列主元素高斯消去法和三角分解法解线性方程组Ax=b。
式中,A为n阶非奇异方阵,x,b是n阶列向量,并分析选主元素的重要性。
二、实验方法(1)列主元素高斯消去法通过变换,将系数矩阵换成等价的上三角矩阵,在每步消元过程中,选列主元素。
对k=1,2,……n-1,逐次计算l ik=a ik(k-1)/a kk(k-1) (i=k+1,k+2,……,n)a ij(k)=a ij(k-1)-l ik a kj(k-1) (i,j=k+1,k+2,……,n)b i(k)=b i(k-1)-l ik b k(k-1) (i=k+1,k+2,……,n)逐步回代气的原方程组的解X n=b i(n-1)/a nn(n-1)X k=(b k(k-1)_a kj(k-1)x j)/a kk(k-1) (k=n-1,n-2, (1)(2)直接三角分解法由于两个矩阵相等就是它们的对应元素相等,因此通过比较A与LU的对应元素,即可得到直接计算L,U的元素的公式。
设A=L×U,其中U的第一行、L的第一列的元素分别为对(依次:U的第二行,L的第二列,U的第三行,L的第三列……),有由上述两种方法得到矩阵A的LU分解后,求解Ly=b与Ux=y的计算公式为∑+=n1kj三、实验内容解下列方程组·=四、实验程序(1)列主元素高斯消去法(2)直接三角分解法0147.06721.109998.42371.13142.17643.89217.44129.35435.15330.27875.15301.04017.31651.18326.31348.14321xxxx9237.164231.183941.65342.9五、实验结果(仅供参考)精确解为:(1,1,1,1)T六、结果分析实验的数学原理很容易理解,也容易上手。
把运算的结果带入原方程组,可以发现符合的还是比较好。
这说明列主元消去法计算这类方程的有效性。
线性方程组的直接解法实验报告
![线性方程组的直接解法实验报告](https://img.taocdn.com/s3/m/2aa0d4e148d7c1c709a1454e.png)
本科实验报告
课程名称:数值计算方法B
实验项目:线性方程组的直接解法
最小二乘拟合多项式
实验地点:ZSA401
专业班级:学号:201000
学生姓名:
指导教师:李志
2012年4月13日
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
printf("%lf\t",A[i][j]);
printf("\n");
}
double answer[N];
Gauss_eliminate(n,answer);
/*输出解*/
for(i=1;i<=n;i++)
printf("a[%d]=%lf\t",i-1,answer[i]);
getchar();
getchar();
}
四、实验结果与讨论、心得
讨论、心得:
刚开始调试代码的时候有时候就是很小的错误导致整个程序不能运行,需要我们一步一步慢慢来,经过无数次的检查程序错误的原因,以及在老师的帮助下,完成了这次实验。
这段时间的实验课提高了我的分析问题,解决问题的能力,特别提高了对一个程序的整。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
实验2 线性方程组解的直接解法
【实验目的】
1、验证高斯消去法和三角分解法
2、掌握直接求解线性方程组的常用算法:列主元高斯消去法、LU 分解法等
3、记录运行结果,回答问题,完成实验报告
【实验要求】
根据题目要求,用Matlab 完成下列实验内容。
书写实验报告。
【实验内容】
1、利用Matlab 求解线性方程组
考虑下面给出的线性方程组形式为Ax =B 。
其中,A 和B 均为给定矩阵
1111n1,n nn n a a b A B a a b ⎛⎫⎛⎫ ⎪ ⎪== ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭
可以由给定的A 和B 矩阵构造出解的判定矩阵C (可以用C=[A B]命令来生成C )
11111
n n nn n a a b C a a b ⎛⎫ ⎪= ⎪ ⎪⎝⎭
这里先可以给出线性方程组有解的判定定理:
A=[1 2 3 ;2 -1 -1;1 -2 -2;1 -1 -1];B=[2;1;-1;0];
C=[A,b];
C,
rank(A),
rank(C),
x=A\B,
x=[1;2;-1];
norm(A*x-B)
答案是:
C =
1 2 3 2
2 -1 -1 1
1 -
2 -2 -1
1 -1 -1 0
ans =
3
ans =
3
x =
1.0000
2.0000
-1.0000
ans =
(2)当rank(A )=rank(C )=r <n 时,方程组有无穷多解。
具体求解方法留给同学
们实验课后自学。
(3)当rank(A )<rank(C )时,方程组为矛盾方程,只能解出方程的最小二乘解。
这部分内容将在第四章学习。
【问题1】:求解下面给出的线性方程组。
记录所用命令和求解结果。
123454
3214132434
1322x ⎛⎫⎛⎫ ⎪ ⎪ ⎪ ⎪= ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭
输入命令:
A=[1 2 3 4;4 3 2 1;1 3 2 4;4 1 3 2];
B=[5;4;3;2];
C=[A,B];
C,
rank(A),
rank(C),
x=A\B,
x=[-1.8;1.8667;3.8667;-2.1333];
B=[5;4;3;2];
norm(A*x-B)
答案是:
C =
1 2 3 4 5
4 3 2 1 4
1 3
2 4 3
4 1 3 2 2
ans =
4
ans =
4
x =
-1.8000
1.8667
3.8667
-2.1333
ans =
5.0990e-004
【思考1】:Matlab 对线性方程组求解是否要考虑主元的选取?为什么?可以用
下列例子试试
121
20.000112x x x x +=⎧⎨+=⎩
答:
A=[0.0001 1;1 1];
B=[1;2];
C=[A,B];
C,
rank(A),
rank(C),
x=A\B,
X=[1.0001;0.9999];
norm(double(A*X-B))
输出答案是:
C =
0.0001 1.0000 1.0000
1.0000 1.0000
2.0000
ans =
2
ans =
2
x =
1.0001
0.9999
ans =
1.0000e-008
2、用Matlab 编写高斯消去法具体实现的函数
要求输入参数包含方程组的系数矩阵及常数项向量,输出方程组的解,并使用问题1的例子进行调试。
高斯顺序消去法:设线性方程组AX =b 。
1.消元过程
设0)(≠k kk a ,对1,,2,1-=n k 计算
()()(1)()()(1)()()/,1,2,,k k ik ik kk k k k ij ij ik kj k k k i i ik k m a a a a m a i j k k n b b m b ++⎧=⎪=-=++⎨⎪=-⎩K 其中
⒉回代过程
()()()()()1/1,,2,1()/n n n n nn n i i i i i ij j ii j i x b a i n x b a x a =+⎧=⎪=-⎨=-⎪⎩∑K 其中。
参考实验框架
输入:
function [x] = guass(A,b)
%开始计算,先赋初值
A=input('输入系数矩阵A的值');
b=input('输入值矩阵b的值');
n=length(A);
x = zeros(n,1);
%消元过程
for k = 1:n-1
for i=j=k+1:n
m(i:i,k:k)=A(i:i,k:k)/A(k:k,k:k);
A(i:i,j:j)=A(i:i,j:j)-m(i:i,k:k)*A(k:k,j:j);
b(i:i)=b(i:i)-m(i:i,k:k)*b(k:k);
end
end
%回代过程
for k = n:-1:1
for i=j=k+1:n
x(n:n)=b(n:n)/A(n:n,n:n);
x(i:i)=(b(i:i)-sum(A(i:i,i+1:n)*x(i+1:n))/A(i:i,i:i));
end
end
想试下用实例做的,不会做啊。
A=[1 2 3 4;4 3 2 1;1 3 2 4;4 1 3 2];
B=[5;4;3;2];
n=length(A);
x = zeros(n,1);
%消元过程
for k = 1:n-1
i=k+1;
j=j+1;
for i=k+1:n
m=A(i:i,k:k)/A(k:k,k:k);
A(i:i,k:n)=A(i:i,k:n)-m*A(k:k,k:n);
end
end
C=A
%回代过程
for k = n-1:1
i=k+1
x(n:n,n:n)=b(n:n,n:n)/A(n:n,n:n)
x(i:i)=(b(i:i)-sum(A(i:i,i+1:n)*x(i+1:n))/A(i:i,i:i)); end。