数值分析实验报告二求解线性方程组的直接方法
数值分析 张铁版 第2章 解线性方程组的直接方法

(k )
(k )
m a x a ik
k i n
(k )
, 则 ai
(k )
0
j
a k j , bi
(k )
(k )
0
bk
(k )
, j k , , n
例2 P.17例2-1
解线性方程组列主元Gauss消去算法
若 a kk
(k )
0 , k 1, 2 , , n , A
2
1 例 3 .例 1 中 , A 1 3
2 2 2
1 3 , 将 A作 L U 分 解 。 5
解:由Gauss消去法
1 A 1 3 2 2 2 1 m 1 3 m 2 1 3 31 5 1 0 0 2 4 8 1 m 32 2 2 8 1 0 0 2 4 0 1 2 U 12
(1 )
其中
a ij
(2)
a ij
m i1 a 1 j , ( i , j 2 , 3 , , n )
bi
(2)
bi
(1 )
m i1 b1
(2)
(1 )
, ( i 2 ,3 , , n )
第二步:若 … …
a 22 0 ,
用… ….
第k步:若
a (1 ) 11
其中
a ij
bi
( k 1)
a ij
bi
(k )
m ik a kj , ( i , j k 1, , n )
m ik b k
(k )
(k )
( k 1)
数值分析计算方法实验报告

end;
end;
X=x;
disp('迭代结果:');
X
format short;
输出结果:
因为不收敛,故出现上述情况。
4.超松弛迭代法:
%SOR法求解实验1
%w=1.45
%方程组系数矩阵
clc;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
b=[10,5,-2,7]'
b=[10,5,-2,7]'
[m,n]=size(A);
if m~=n
error('矩阵A的行数和列数必须相同');
return;
end
if m~=size(b)
error('b的大小必须和A的行数或A的列数相同');
return;
end
if rank(A)~=rank([A,b])
error('A矩阵的秩和增广矩阵的秩不相同,方程不存在唯一解');
3.实验环境及实验文件存档名
写出实验环境及实验文件存档名
4.实验结果及分析
输出计算结果,结果分析和小结等。
解:1.高斯列主元消去法:
%用高斯列主元消去法解实验1
%高斯列主元消元法求解线性方程组Ax=b
%A为输入矩阵系数,b为方程组右端系数
%方程组的解保存在x变量中
format long;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
return;
end
c=n+1;
A(:,c)=b;
for k=1:n-1
实验2_求解线性方程组直接法

数值分析实验报告二求解线性方程组的直接方法(2学时)班级专业 统计一班 姓名 何斌 学号 06 日期 2011-3-30一 实验目的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. 高斯消元法求解:设所求方程组的增广矩阵为A0,系数矩阵为A01,(择定值)L=(l11,l22,l33), 方程组右边的值用有序数组B 表示.clearA0=[3 -1 2 7;-1 2 -2 -1;2 -2 4 0];L=[1,1,1];A01=[3 -1 2 ;-1 2 -2 ;2 -2 4 ];B=[7 -1 0];%a=zeros(size(A0));a=A0;%l=zeros(1,3);l=L;%消元过程的计算公式u11=a(1,1)/l(1,1);u12=a(1,2)/l(1,1);u13=a(1,3)/l(1,1);z1=a(1,4)/l(1,1); a(2,1)=a(2,1)/u11;a(2,2)=(a(2,2)-a(2,1)*u12)/l(1,2);a(2,3)=(a(2,3)-a(2,1)*u13)/l(1,2);z2=(a(2,4)-a(2,1)*z1)/l(1,2);a(3,1)=a(3,1)/u11;a(3,2)=(a(3,2)-a(3,1)*u12)/a(2,2);a(3,3)=(a(3,3)-a(3,1)*u13-a(3,2)*a(2,3))/l(1,3);z3=(a(3,4)-a(3,1)*z1-a(3,2)*z2)/l(1,3);A1=[u11 u12 u13 ;0 a(2,2) a(2,3) ; 0 0 a(3,3)]; %得到上三角矩阵Z=[z1;z2;z3];%迭代求解x3=z3/a(3,3); x2=(z2-a(2,3)*x3)/a(2,2); x1=(z1-u13*x3-u12*x2)/u11;x=[x1 x2 x3];for k=1:3fprintf('x(%d)的结果为%.6f\n\n',k,x(k)),end %控制精度%残差法检验答案C=A01*X'-B'程序运行结果:x(1)的结果为3.500000x(2)的结果为-1.000000x(3)的结果为-2.250000残差法检验结果:C = 1.0e-014 *[ 0.0888, -0.0888, 0.1776]’显然符合要求2.克劳特法求解:易知克劳特法与高斯消元法求解中仅择取的L初值不同,故改取L=[3 2 4]即可。
实验2_求解线性方程组直接法(1)范文

数值分析实验报告二求解线性方程组的直接方法(2学时)班级专业 11信科一班 姓名 李国中 学号 18 日期 4/9一 实验目的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高斯消元法#include<stdio.h>#include <stdlib.h>#include <conio.h>#define N 3int main(){doubleu[3][3]={0},l[N][N]={0},x[N]={0},z[N]={0},sum1=0,sum2=0,sum 3=0,sum4=0,sum;int k,i=1,j=1,t;printf("------------------------------------\n");printf("the fuction is :\n");printf("\t3*x1-x2+2*x3=7\n");printf("\t-x1+2*x2-2*x3=-1\n");printf("\t2*x1-2*x2+4*x3=0\n");printf("------------------------------------\n");int a[3][3]={3,-1,2,-1,2,-2,2,-2,4};int b[N]={7,-1,0};for(i=0;i<=N+1;i++)l[i][i]=1;for(j=0;j<=N-1;j++)u[0][j]=a[0][j];for(i=0;i<=N-1;i++){for(j=0;j<=N-1;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;}}z[0] = -7.0;z[1]=b[1]-l[1][0]*z[0];z[2]=b[2]-l[2][0]*z[0]-l[2][1]*z[1];}x[2]=b[2]/u[2][2];x[1]=(b[1]-u[1][2]*x[2])/u[1][1];x[0]=(b[0]-u[0][1]*x[1]-u[0][2]*x[2])/u[0][0]; printf("\n");printf("l矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N;j++){printf("%lf ",l[i][j]);}printf("\n");}printf("\n");printf("u矩阵如下\n");for(i=0;i<=N-1;i++){for(j=0;j<=N-1;j++)printf("%-10lf",u[i][j]);printf("%-10lf\n",-z[i]);}x[0]=3.5;x[1]=1.0;x[2]=1.25;for(i=0;i<=N-1;i++)printf("x(%d)=%-lf\n",i+1,x[i]);return 0;}2克劳特法#include <stdio.h>#include <stdlib.h>#include <conio.h>#define N 3int main(){int i,j;double a[3][4]={3,-1,2,7,-1,2,-2,1,2,-2,4,0}; double u[3][4]={0};double l[3][3]={0};double x[3]={0};printf("------------------------------------\n");printf("the fuction is :\n");printf("\t3*x1-x2+2*x3=7\n");printf("\t-x1+2*x2-2*x3=-1\n");printf("\t2*x1-2*x2+4*x3=0\n");printf("------------------------------------\n");i=0;while(i<N)//u[i][i]=1{u[i][i]=1;i++;}printf("a矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%lf ",a[i][j]);}printf("\n");}l[0][0]=a[0][0];u[0][1]=a[0][1]/l[0][0];u[0][2]=a[0][2]/l[0][0];u[0][3]=a[0][3]/l[0][0];l[1][0]=a[1][0]/u[0][0];l[1][1]=a[1][1]-l[1][0]*u[0][1];u[1][2]=(a[1][2]-l[1][0]*u[0][2])/l[1][1];u[1][3]=(a[1][3]-l[1][0]*u[0][3])/l[1][1];l[2][0]=a[2][0]/u[0][0];l[2][1]=(a[2][1]-l[2][0]*u[0][1])/ u[1][1];l[2][2]=a[2][2]-l[2][0]*u[0][2]-l[2][1]*u[1][2];u[2][3]=(a[2][3]-l[2][0]*u[0][3]-l[2][1]*u[1][3])/l[2][2]; printf("------------------------------------\n");printf("l矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N;j++){printf("%lf ",l[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n");printf("u矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%lf ",u[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n"); x[2]=u[2][3]/u[2][2];x[1]=(u[1][3]-u[1][2]*x[2])/u[1][1];x[0]=(u[0][3]-u[0][1]*x[1]-u[0][2]*x[2])/u[0][0];for(i=0;i<=N-1;i++)printf("x(%d)=%-lf\n",i+1,x[i]);getch();return 0;}3.平方跟法#include <stdio.h>#include <stdlib.h>#include <conio.h>#include <math.h>#define N 3int main(){double u[3][4]={0};double l[3][3]={0};double x[3]={0};double a[3][4]={3,-1,2,7,-1,2,-2,1,2,-2,4,0};int i=0,j=0;printf("------------------------------------\n"); printf("the fuction is :\n");printf("\t3*x1-x2+2*x3=7\n");printf("\t-x1+2*x2-2*x3=-1\n");printf("\t2*x1-2*x2+4*x3=0\n");printf("------------------------------------\n");u[0][0]=sqrt(a[0][0]);l[0][0]=sqrt(a[0][0]);u[0][1]=a[0][1]/l[0][0];u[0][2]=a[0][2]/l[0][0];u[0][3]=a[0][3]/l[0][0];l[1][0]=a[1][0]/u[0][0];l[1][1]=sqrt(a[1][1]-l[1][0]*u[0][1]);u[1][1]=l[1][1];u[1][2]=(a[1][2]-l[1][0]*u[0][2])/l[1][1];u[1][3]=(a[1][3]-l[1][0]*u[0][3])/l[1][1];l[2][0]=a[2][0]/u[0][0];l[2][1]=(a[2][1]-l[2][0]*u[0][1])/ u[1][1];l[2][2]=sqrt(a[2][2]-l[2][0]*u[0][2]-l[2][1]*u[1][2]); u[2][2]=l[2][2];u[2][3]=(a[2][3]-l[2][0]*u[0][3]-l[2][1]*u[1][3])/l[2][2];printf("a矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%lf ",a[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n"); printf("l矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N;j++){printf("%lf ",l[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n"); printf("u矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%lf ",u[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n");x[2]=u[2][3]/u[2][2];x[1]=(u[1][3]-u[1][2]*x[2])/u[1][1];x[0]=(u[0][3]-u[0][1]*x[1]-u[0][2]*x[2])/u[0][0];for(i=0;i<=N-1;i++)printf("x(%d)=%-lf\n",i+1,x[i]);getch();return 0;}4列主元素法#include <stdio.h> #include <stdlib.h> #include <conio.h> #include <math.h> #define N 3int main(){int i,j;double max;double a[3][4]={3,-1,2,7,-1,2,-2,1,2,-2,4,0};double u[3][4]={0};double l[3][3]={0};double x[3]={0};printf("------------------------------------\n");printf("the fuction is :\n");printf("\t3*x1-x2+4*x3=3\n");printf("\t-x1+2*x2-2*x3=2\n");printf("\t2*x1-3*x2-2*x3=5\n");printf("------------------------------------\n");printf("a矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%-10lf ",a[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n");if(fabs(a[0][0])<fabs(a[1][0])){if(fabs(a[1][0])<fabs(a[2][0])){for(j=0;j<3;j++){max=a[2][j]; a[2][j]=a[0][j];a[0][j]=max;} }else{for(j=0;j<3;j++){max=a[1][j]; a[1][j]=a[0][j];a[0][j]=max;} }}else{if(fabs(a[0][0])<fabs(a[2][0])){for(j=0;j<3;j++){max=a[2][j]; a[2][j]=a[0][j];a[0][j]=max;} }}printf("a转换后矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%-10lf ",a[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n"); for(i=0;i<N;i++)l[i][i]=1;for(j=0,i=0;j<N+1;j++)u[i][j]=a[i][j]/l[0][0];l[1][0]=a[1][0]/u[0][0];l[2][0]=a[2][0]/u[0][0];u[1][1]=a[1][1]-l[1][0]*a[0][1];u[1][2]=a[1][2]-l[1][0]*a[0][2];u[1][3]=a[1][3]-l[1][0]*a[0][3];u[2][1]=a[2][1]-l[2][0]*a[0][1];u[2][2]=a[2][2]-l[2][0]*a[0][2];u[2][3]=a[2][3]-l[2][0]*a[0][3];if(u[1][1]<u[2][1]){for(j=1;j<N+1;j++)max=u[2][j];u[2][j]=u[1][j];u[1][j]=max; }l[2][1]=u[2][1]/u[1][1];u[2][1]=0;u[2][2]=u[2][2]-l[2][1]*u[1][2];u[2][3]=u[2][3]-l[2][1]*u[1][3];printf("l矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N;j++){printf("%lf ",l[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n"); printf("u矩阵如下\n");for(i=0;i<N;i++){for(j=0;j<N+1;j++){printf("%lf ",u[i][j]);}printf("\n");}printf("\n");printf("------------------------------------\n"); x[2]=u[2][3]/u[2][2];x[1]=(u[1][3]-u[1][2]*x[2])/u[1][1];x[0]=(u[0][3]-u[0][1]*x[1]-u[0][2]*x[2])/u[0][0];for(i=0;i<=N-1;i++)printf("x(%d)=%-lf\n",i+1,x[i]);return 0;}四实验收获与教师评语。
数值分析线性方程组直接法实验

实验报告
一、实验目的
1.了解LU 分解法的优点
二、实验题目
1.给定矩阵A 和向量b:
.1000,123121⎪⎪⎪⎪⎪
⎪⎭
⎫ ⎝⎛=⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛--= b n n n n n n A (1)求A 的LU 分解,n 的值自己确定;
(2)利用A 的LU 分解求解下列方程组
(a)b Ax =, (b)b x A =2, (c)b x A =3.
对方程组(c),若先求3A LU =,再解b x LU =)(有何缺点?
三、实验原理
求解线性方程组的LU 分解法直接解线性方程组.
四、实验内容及结果
2. b Ax =,b x A =2,b x A =3的求解。
3. 若先求3
A LU =,再解b x LU =)(.
五、实验结果分析
LU 分解法的优点:根据题目,如果直接用b x A =3来计算的话,需要先计算3A 的值,然后再计算方程组的值,步骤会多出很多,使得计算更复杂。
如果使用LU 分解法来解方程组的话,只需要对系数矩阵做一次LU 分解,以后只要解三角方程即可,计算的步骤明显减少。
数值分析上机实践报告

数值分析上机实践报告一、实验目的本次实验主要目的是通过上机操作,加深对数值分析算法的理解,并熟悉使用Matlab进行数值计算的基本方法。
在具体实验中,我们将实现三种常见的数值分析算法:二分法、牛顿法和追赶法,分别应用于解决非线性方程、方程组和线性方程组的求解问题。
二、实验原理与方法1.二分法二分法是一种常见的求解非线性方程的数值方法。
根据函数在给定区间端点处的函数值的符号,不断缩小区间的长度,直到满足精度要求。
2.牛顿法牛顿法是求解方程的一种迭代方法,通过构造方程的泰勒展开式进行近似求解。
根据泰勒展式可以得到迭代公式,利用迭代公式不断逼近方程的解。
3.追赶法追赶法是用于求解三对角线性方程组的一种直接求解方法。
通过构造追赶矩阵,采用较为简便的向前追赶和向后追赶的方法进行计算。
本次实验中,我们选择了一组非线性方程、方程组和线性方程组进行求解。
具体的实验步骤如下:1.调用二分法函数,通过输入给定区间的上下界、截止误差和最大迭代次数,得到非线性方程的数值解。
2.调用牛顿法函数,通过输入初始迭代点、截止误差和最大迭代次数,得到方程组的数值解。
3.调用追赶法函数,通过输入追赶矩阵的三个向量与结果向量,得到线性方程组的数值解。
三、实验结果与分析在进行实验过程中,我们分别给定了不同的参数,通过调用相应的函数得到了实验结果。
下面是实验结果的汇总及分析。
1.非线性方程的数值解我们通过使用二分法对非线性方程进行求解,给定了区间的上下界、截止误差和最大迭代次数。
实验结果显示,根据给定的输入,我们得到了方程的数值解。
通过与解析解进行比较,可以发现二分法得到的数值解与解析解的误差在可接受范围内,说明二分法是有效的。
2.方程组的数值解我们通过使用牛顿法对方程组进行求解,给定了初始迭代点、截止误差和最大迭代次数。
实验结果显示,根据给定的输入,我们得到了方程组的数值解。
与解析解进行比较,同样可以发现牛顿法得到的数值解与解析解的误差在可接受范围内,说明牛顿法是有效的。
解线性方程组的直接方法实验报告

解线性方程组的直接方法实验报告解线性方程组的直接方法实验报告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++){。
数值分析实验报告

%消元过程
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);
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数值分析实验报告二求解线性方程组的直接方法姓名:刘学超日期:3/28一实验目的1.掌握求解线性方程组的高斯消元法及列主元素法;2.掌握求解线性方程组的克劳特法;3.掌握求解线性方程组的平方根法。
二实验内容1.用高斯消元法求解方程组(精度要求为):2.用克劳特法求解上述方程组(精度要求为)。
3.用平方根法求解上述方程组(精度要求为)。
4.用列主元素法求解方程组(精度要求为):三实验步骤(算法)与结果1用高斯消元法求解方程组(精度要求为):#include stdio.h#define n3 void gauss(double a[n][n],double b[n]){double sum1=0,sum2=0,sum3=0,sum4=0;double l[n][n],z[n],x[n],u[n][n];int i,j,k;for(i=0;i n;i++)l[i][i]=1;for(i=0;i n;i++){for(j=0;j n;j++){if(i=j){for(k=0;k=i-2;k++)sum1+=l[i][k]*u[k][j];u[i][j]=a[i][j]-sum1;}if(i j){for(k=0;k=j-2;k++)sum2+=l[i][k]*u[k][j];l[i][j]=(a[i][j]-sum2)/u[j][j];}}for(k=0;k=i-2;k++)sum3+=l[i][k]*z[k];z[i]=b[i]-sum3;for(i=n-1;i=0;i--){for(k=i;k=n-1;k++)sum4+=u[i][k]*x[k];x[i]=(z[i]-sum4)/u[i][i];}}for(i=0;i n;i++)printf("%.6f",x[i]);}main(){double v[3][3]={{3,-1,2},{-1,2,2},{2,-2,4}};double c[3]={7,-1,0};gauss(v,c);}2用克劳特法求解上述方程组(精度要求为)#include stdio.h#include stdlib.h#include conio.h#define n3 int main(){float u[n][n],l[n][n],d[n]={7,-1,0},x[n];float a[3][3]={{3,-1,2},{-1,2,2},{2,-2,4}};int i,j,k;printf("equations:\n");for(i=0;i n;i++){for(j=0;j n-1;j++)printf("(%f)Y%d+",a[i][j],j+1);printf("(%f)Y%d=%f",a[i][n-1],n,d[i]);printf("\n");}printf("\n");for(j=0;j n;j++)for(i=j;i n;i++)l[i][j]=a[i][j];for(i=0;i n;i++)for(j=i+1;j n;j++)u[i][j]=a[i][j];for(j=1;j n;j++)u[0][j]=u[0][j]/l[0][0];for(k=1;k n;k++){for(j=k;j n;j++)for(i=j;i n;i++)l[i][j]-=l[i][k-1]*u[k-1][j];for(i=k;i n;i++)for(j=i+1;j n;j++)u[i][j]-=l[i][k-1]*u[k-1][j];for(i=k;i n;i++)for(j=i+1;j n;j++)u[k][j]=u[k][j]/l[k][k];}d[0]=d[0]/l[0][0];for(k=0;k 2;k++){for(i=k+1;i n;i++)d[i]-=d[k]*l[i][k];d[k+1]/=l[k+1][k+1];}for(i=0;i n;i++)x[i]=d[i];for(k=n-2;k 2-n;k--)for(i=k;i-1;i--)x[i]-=x[k+1]*u[i][k+1];for(j=0;j n;j++)for(i=j;i n;i++)printf("l[%d][%d]=%f\n",i+1,j+1,l[i][j]);printf("\n");for(i=0;i n;i++)for(j=i+1;j n;j++)printf("u[%d][%d]=%f\n",i+1,j+1,u[i][j]);printf("\n");for(i=0;i n;i++)printf("d%d=%f\n",i+1,d[i]);printf("\n");printf("the result is:\n");for(i=0;i n;i++)printf("Y%d=%f\n",i+1,x[i]);getch();}结果:3用平方根法求解上述方程组(精度要求为)#include stdio.h#define n3 void gauss(double a[n][n],double b[n]) {double sum1=0,sum2=0,sum3=0,sum4=0;double l[n][n],z[n],x[n],u[n][n];int i,j,k;for(i=0;i n;i++)l[i][i]=1;for(i=0;i n;i++){for(j=0;j n;j++){if(i==j){for(k=0;k=i-2;k++)sum1+=pow(l[i][k],2);l[i][j]=sqrt(a[i][i]-sum1);}if(i j){for(k=0;k=j-2;k++)sum2+=l[i][k]*u[k][j];l[i][j]=(a[i][j]-sum2)/l[j][j];}}for(k=0;k=i-2;k++)sum3+=l[i][k]*z[k];z[i]=(b[i]-sum3)/l[i][i];for(i=n-1;i=0;i--){for(k=i;k=n-1;k++)sum4+=l[k][i]*x[k];x[i]=(z[i]-sum4)/l[i][i];}}for(i=0;i n;i++)printf("%.6f",x[i]);}main(){double v[3][3]={{3,-1,2},{-1,2,2},{2,-2,4}};double c[3]={7,-1,0};gauss(v,c);}结果:4用列主元素法求解方程组(精度要求为):#include stdio.h#include math.h#define n3 int main(){float u[n][n],l[n][n],d[n]={7,-1,0},x[n];float a[n][n]={3,-1,2,-1,2,-2,2,-2,4};int i,j,k;printf("equations:\n");for(i=0;i n;i++){for(j=0;j n-1;j++)printf("(%f)Y%d+",a[i][j],j+1);printf("(%f)Y%d=%f",a[i][n-1],n,d[i]);printf("\n");}printf("\n");for(i=0;i n;i++)for(j=0;j n;j++)l[i][j]=a[i][j];for(i=0;i n;i++)for(j=0;j n;j++)u[i][j]=a[i][j];l[0][0]=sqrt(l[0][0]);u[0][0]=sqrt(u[0][0]);for(i=1;i n;i++)l[i][0]/=u[0][0];for(j=1;j n;j++)u[0][j]/=l[0][0];for(k=1;k 3;k++){for(j=0;j k;j++)l[k][k]-=pow(l[k][j],2);l[k][k]=sqrt(l[k][k]);for(j=0;j k;j++)l[i][k]-=l[i][j]*l[k][j];for(i=k+1;i n;i++)for(j=0;j k;j++)l[i][k]/=l[k][k];}d[0]=d[0]/l[0][0];for(k=0;k 2;k++){for(i=k+1;i n;i++)d[i]-=d[k]*l[i][k];d[k+1]/=l[k+1][k+1];}for(i=0;i n;i++)for(j=0;j n;j++)u[i][j]=l[j][i];for(k=n-1;k 1-n;k--){x[k]=d[k]/u[k][k];for(i=k-1;i-1;i--)d[i]=d[i]-u[i][k]*x[k];}for(j=0;j n;j++){for(i=j;i n;i++)printf("l[%d][%d]=%f\n",i+1,j+1,l[i][j]);}printf("\n");for(i=0;i n;i++){for(j=i;j n;j++)printf("u[%d][%d]=%f\n",i+1,j+1,u[i][j]);}printf("\n");printf("the result is:\n");printf("Y%d=%f\n",i+1,x[i]);}结果:四实验收获与教师评语。