数值分析报告线性方程组

合集下载

迭代法解线性方程组数值分析实验报告

迭代法解线性方程组数值分析实验报告

迭代法解线性方程组数值分析实验报告一、实验目的本次实验旨在深入研究和掌握迭代法求解线性方程组的基本原理和方法,并通过数值实验分析其性能和特点。

具体目标包括:1、理解迭代法的基本思想和迭代公式的推导过程。

2、掌握雅克比(Jacobi)迭代法、高斯赛德尔(GaussSeidel)迭代法和超松弛(SOR)迭代法的算法实现。

3、通过实验比较不同迭代法在求解不同类型线性方程组时的收敛速度和精度。

4、分析迭代法的收敛性条件和影响收敛速度的因素。

二、实验原理1、线性方程组的一般形式对于线性方程组$Ax = b$,其中$A$ 是$n×n$ 的系数矩阵,$x$ 是$n$ 维未知向量,$b$ 是$n$ 维常向量。

2、迭代法的基本思想迭代法是从一个初始向量$x^{(0)}$出发,按照某种迭代公式逐步生成近似解序列$\{x^{(k)}\}$,当迭代次数$k$ 足够大时,$x^{(k)}$逼近方程组的精确解。

3、雅克比迭代法将系数矩阵$A$ 分解为$A = D L U$,其中$D$ 是对角矩阵,$L$ 和$U$ 分别是下三角矩阵和上三角矩阵。

雅克比迭代公式为:$x^{(k+1)}= D^{-1}(b +(L + U)x^{(k)})$。

4、高斯赛德尔迭代法在雅克比迭代法的基础上,每次计算新的分量时立即使用刚得到的最新值,迭代公式为:$x_i^{(k+1)}=(b_i \sum_{j=1}^{i-1}a_{ij}x_j^{(k+1)}\sum_{j=i+1}^{n}a_{ij}x_j^{(k)})/a_{ii}$。

5、超松弛迭代法在高斯赛德尔迭代法的基础上引入松弛因子$\omega$,迭代公式为:$x_i^{(k+1)}= x_i^{(k)}+\omega((b_i \sum_{j=1}^{i-1}a_{ij}x_j^{(k+1)}\sum_{j=i}^{n}a_{ij}x_j^{(k)})/ a_{ii} x_i^{(k)})$。

数值分析计算方法实验报告

数值分析计算方法实验报告
break;
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

线性方程组求解的数值实验报告

线性方程组求解的数值实验报告

线性方程组求解的数值实验一、课题名称:双对角线线性方程组的数值实验二、班级和姓名:09101900 张争(学号 0910190161) 三、实验问题摘要:考虑一种特殊的对角线元素不为零的双对角线线性方程组(以n=7为例)⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡7665544332211d a d a d a d a d a da d ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡7654321x x x x x x x =⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡7654321b b b b b b b 写出解一般的n (奇数)阶方程组的程序(不要用消元过程,因为不用它可以十分方便的解出这个方程组)。

四、过程中的问题及处理方法:解方程组时应当先解出x1和xn,r 然后根据公式求出其他的解,由于矩阵的特殊性,应该从头尾开始解,最后解出21+n x 。

编程序时使用循环语句要分情况讨论。

五、计算公式: 1x =11d b kk k k k d x a b x 11---=(k<n/2)n x =nnd b kk k k k d x a b x 1+-=(k>n/2)21212123212121+--++++--=n n n n n n n d x a x a b x六、结构程序设计输入矩阵阶数→输入各系数→解出前(n-1)/2个解→解出后(n-1)/2个解→解出第(n+1)/2个解七、程序运行结果:这里令所有的系数等于1,实验结果如下:八、主要程序段介绍:解方程函数:由于矩阵具有好的对称性,解的形式也一样,用循环语句。

但是第一个和最后一个解以及中间一个解的形式不同,需要单独解出,程序如下:void jiefangcheng(){int i;x[1]=b[1]/d[1];for(i=2;i<(i+1)/2;i++)x[i]=(b[i]-a[i-1]*x[i-1])/d[i];x[n]=b[n]/d[n];for(i=n-1;i>(i+1)/2;i--)x[i]=(b[i]-a[i]*x[i+1])/d[i];x[(n+1)/2]=(a[(n-1)/2]*x[(n-1)/2]-a[(n+1)/2]*x[(n+3)/2])/d[(n+1)/2];}九、源程序如下://*****************张争******学号0910190161**********#include<iostream.h>#define N 10000int n;int a[N],b[N],d[N],x[N];void input(){int i;cout<<"请输入方程组的阶数:";cin>>n;cout<<"请输入方程组的系数:"<<endl;for(i=1;i<=n;i++){cout<<"d"<<i<<"=";cin>>d[i];cout<<"b"<<i<<"=";cin>>b[i];}for(i=1;i<n;i++){cout<<"a"<<i<<"=";cin>>a[i];}}void jiefangcheng(){int i;x[1]=b[1]/d[1];for(i=2;i<(i+1)/2;i++)x[i]=(b[i]-a[i-1]*x[i-1])/d[i];x[n]=b[n]/d[n];for(i=n-1;i>(i+1)/2;i--)x[i]=(b[i]-a[i]*x[i+1])/d[i];x[(n+1)/2]=(a[(n-1)/2]*x[(n-1)/2]-a[(n+1)/2]*x[(n+3)/2])/d[(n+1)/2];}void output(){ int i;cout<<"方程组的解为:"<<endl;for(i=1;i<=n;i++)cout<<"x"<<i<<"="<<x[i]<<endl;}void main(){input();jiefangcheng();output();}。

数值分析实验报告二求解线性方程组的直接方法

数值分析实验报告二求解线性方程组的直接方法

数值分析实验报告二求解线性方程组的直接方法姓名:刘学超日期: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]);}结果:四实验收获与教师评语。

线性方程组的数值解法实验报告

线性方程组的数值解法实验报告

实验报告——线性方程组的数值解法姓名:班级:学号:日期:一 实践目的1. 熟悉求解线性方程组的有关理论和方法。

2. 会编列主元消去法,全主元消去法,雅克比迭代法及高斯-赛德尔迭代法的程序。

3.通过实际计算,进一步了解各种方法的优缺点,选择合适的数值方法。

4. 进一步应用数学知识,扩展数学思维,提高编程能力。

二 问题定义及题目分析1.求解线性方程是实际中常遇到的问题。

而一般求解方法有两种,一个是直接求解法,一个是迭代法。

2.直接法是就是通过有限步四则运算求的方程准确解的方法。

但实际计算中必然存在舍入误差,因此这种方法只能得到近似解。

3.迭代法是先给一个解的初始近似值,然后按一定的法则求出更准确的解,即是用某种极限过程逐步逼近准确解的方法。

4.这次实践,将采用直接法:列主元高斯消去法,全主元高斯消去法;迭代法:雅克比迭代法和高斯-赛德尔迭代法。

三 详细设计 设有线性方程组11112211n n a x a x a x b +++=21122222n n a x a x a x b +++= 1122n n nn n n a x a x a x b +++=令A= 111212122212n n n n nn a a a aa a a a a ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦ x=12n x x x ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦ b= 12n b b b ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦一. 上三角方程组的求解很简单。

所以若能将方程组化为上三角方程组,那就很容易得到方程的解,高斯消去法则是一种简洁高效的方法,可以将方程组化为上三角方程组,但是这种方法必须满足每一步时0kk a =才能求解,还有当kk a 绝对值很小时,作分母会引起较大的舍入误差。

因此在消元过程中应尽量选择绝对值较大的系数作为主元素。

1.列主元消去法很好的解决这一问题。

将1x 的系数1(1)k a k n ≤≤中绝对值最大者作为主元素,交换第一行和此元素所在的行,然后主元消去非主元项1x 的系数,。

数值分析实验报告--解线性方程组的迭代法及其并行算法

数值分析实验报告--解线性方程组的迭代法及其并行算法

disp('请注意:高斯-塞德尔迭代的结果没有达 到给定的精度,并且迭代次数已经超过最大迭 代次数max1,方程组的精确解jX和迭代向量X 如下: ') X=X';jX=jX' end end X=X';D,U,L,jX=jX'
高斯-塞德尔的输入为:
A=[10 2 3;2 10 1;3 1 10]; b=[1;1;2]; X0=[0 0 0]'; X=gsdddy(A,b,X0,inf, 0.001,100) A=[10 2 3;2 10 1;3 1 10]; 请注意:因为对角矩阵 D 非奇异,所以此方程组有解.
0.0301 0.0758 0.1834
8.心得体会:
这已经是第三次实验了, 或多或少我已经对 MATLAB 有了更多的了 解与深入的学习。通过这次实验我了解了雅可比迭代法和高斯- 塞德尔迭代法的基本思想,虽然我们不能熟练编出程序,但还是 能看明白的。运行起来也比较容易,让我跟好的了解迭代法的多 样性,使平常手算的题能得到很好的验证。通过这次实验让我对 MATLAB 又有了更深一层的认识,使我对这门课兴趣也更加浓厚。
运行雅可比迭代程序输入: A=[10
b=[1;1;2];X0=[0 0 0]'; X=jacdd(A,b,X0,inf,0.001,100)
2 3;2 10 1;3 1 10];
结果为:
k= 1 X=
0.1000 k= 2 X= 0.0200 k= 3 X= 0.0400 k= 4 X= 0.0276 k= 5 X= 0.0314 k= 6 X= 0.0294 k= 7 X= 0.0301 k= 8 X= 0.0297
6、 设计思想:先化简,把对角线的项提到左边,其它项

数值分析报告讲义——线性方程组解法

数值分析报告讲义——线性方程组解法

数值分析讲义第三章线性方程组的解法§3.0 引言§3.1 雅可比(Jacobi)迭代法§3.2 高斯-塞德尔(Gauss-Seidel)迭代法§3.3 超松驰迭代法§3.7 三角分解法§3.4 迭代法的收敛性§3.8 追赶法§3.5 高斯消去法§3.9 其它应用§3.6 高斯主元素消去法§3.10 误差分析§3 作业讲评3 §3.11 总结§3.0 引言重要性:解线性代数方程组的有效方法在计算数学和科学计算中具有特殊的地位和作用.如弹性力学、电路分析、热传导和振动、以及社会科学及定量分析商业经济中的各种问题.分类:线性方程组的解法可分为直接法和迭代法两种方法.(a) 直接法:对于给定的方程组,在没有舍入误差的假设下,能在预定的运算次数内求得精确解.最基本的直接法是Gauss消去法,重要的直接法全都受到Gauss消去法的启发.计算代价高.(b) 迭代法:基于一定的递推格式,产生逼近方程组精确解的近似序列.收敛性是其为迭代法的前提,此外,存在收敛速度与误差估计问题.简单实用,诱人.§3.1 雅可比Jacobi 迭代法 (AX =b )1基本思想:与解f (x )=0 的不动点迭代相类似,将AX =b 改写为X =BX +f 的形式,建立雅可比方法的迭代格式:X k +1=BX (k )+f ,其中,B 称为迭代矩阵.其计算精度可控,特别适用于求解系数为大型稀疏矩阵(sparse matrices)的方程组. 2问题:(a) 如何建立迭代格式?(b) 向量序列{X k }是否收敛以及收敛条件? 3 例题分析:考虑解方程组⎪⎩⎪⎨⎧=+--=-+-=--2.453.82102.7210321321321x x x x x x x x x (1)其准确解为X *={1, 1.2, 1.3}. 建立与式(1)相等价的形式:⎪⎩⎪⎨⎧++=++=++=84.02.01.083.02.01.072.02.01.0213312321x x x x x x x x x (2) 据此建立迭代公式:⎪⎩⎪⎨⎧++=++=++=+++84.02.01.083.02.01.072.02.01.0)(2)(1)1(3)(3)(1)1(23)(2)1(1k k k k k k kk k x x x x x x x x x (3) 取迭代初值0)0(3)0(2)0(1===xxx,迭代结果如下表.JocabiMethodP31.cpp迭代次数 x1 x2 x30 0 0 01 0.72 0.83 0.842 0.971 1.07 1.153 1.057 1.1571 1.24824 1.08535 1.18534 1.282825 1.095098 1.195099 1.2941386 1.098338 1.198337 1.2980397 1.099442 1.199442 1.2993358 1.099811 1.199811 1.2997779 1.099936 1.199936 1.29992410 1.099979 1.199979 1.29997511 1.099993 1.199993 1.29999112 1.099998 1.199998 1.29999713 1.099999 1.199999 1.29999914 1.1 1.2 1.315 1.1 1.2 1.34Jocobi迭代公式:设方程组AX=b, 通过分离变量的过程建立Jocobi迭代公式,即),,2,1()(1),,2,1(0,11n i x a b a x n i a b x a n ij j j ij i iii ii ni i j ij =∑-==≠∑=≠== 由此我们可以得到Jacobi 迭代公式:),,2,1()(11)1(n i x a b a xn ij j k i ij i iik i=∑-=≠=+[Jacobi 迭代公式的算法] 1: 初始化. n , (a ij ), (b j ), (x 1) , M . 2: 执行k =1直到M 为止. ① 执行i =1直到n 为止.ii nij j j ij i i a x a b u /)(1∑-←≠= ;② 执行i =1直到n 为止.i i u x ← ;③输出k , (x i ).另外,我们也可以建立Jacobi 迭代公式的矩阵形式. 设方程组AX =b ,其中,A =(a ij )n 为非奇异阵,X =(x 1,x 2,…,x n )T , b =(b 1,b 2,…,b n )T将系数阵A 分解为: A =U +D +L ,U 为上三角矩阵,D 为对角矩阵,L 为下三角矩阵.于是AX =b 可改写为 (U +D +L )X =b⇔ X =D -1b -D -1(U +L )X由此可得矩阵形式的Jocobi 迭代公式: X k +1=BX (k )+f □§3.2 高斯-塞德尔Gauss-Seidel 迭代法注意到利用Jocobi 迭代公式计算)1(+k ix 时,已经计算好)(1)(2)(1,,,k i k k x x x - 的值,而Jocobi 迭代公式并不利用这些最新的近似值计算,仍用)(1)(2)(1,,,k i k k x x x - .这启发我们可以对其加以改进,即在每个分量的计算中尽量利用最新的迭代值,得到),,2,1()(1111)1()1(n i x a x a b a xn i j k jij i j k j ij i iik i=∑-∑-=+=-=++上式称为Gauss-Seidel 迭代法. 其矩阵形式是X =-(D +L )-1UX +(D +L )-1b , X k +1=BX (k )+f .迭代次数 x1 x2 x3 0 0 0 0 1 0.72 0.902 1.1644 2 1.04308 1.167188 1.282054 3 1.09313 1.195724 1.2977714 1.099126 1.199467 1.2997195 1.09989 1.199933 1.2999656 1.099986 1.199992 1.2999967 1.099998 1.199999 1.2999998 1.1 1.2 1.3§3.3 超松驰迭代法SOR 方法1基本思想:逐次超松弛迭代法(Successive Over Relaxation Method,简写为SOR)可以看作带参数ω的高斯-塞德尔迭代法,是G-S 方法的一种修正或加速.是求解大型稀疏矩阵方程组的有效方法之一. 2 SOR 算法的构造:设方程组AX =b , 其中,A =(a ij )n 为非奇异阵,X =(x 1,x 2,…,x n )T , b =(b 1,b 2,…,b n )T . 假设已算出x (k ),),,2,1()(1111)1()1(n i x a x a b a xn i j k j ij i j k j ij i iik i=∑-∑-=+=-=++ (1)相当于用高斯-塞德尔方法计算一个分量的公式. 若对某个参数ω,作)1(+k ix与)(k i x 加权的平均,即)()1()()1()()1()(1k i k ik i k ik ik ix xx xxx-+=+-=+++ωωω (2)其中,ω称为松弛因子.用(1)式代入(2)式,就得到解方程组AX =b 的逐次超松弛迭代公式:⎪⎩⎪⎨⎧=∑-∑-=∆∆+==-=++),,2,1()()(11)1()()1(n i x a x a b a x x x x n ij k j ij i j k j ij i iii i k i k i ω (3) 显然,当取ω=1时,式(3)就是高斯-塞德尔迭代公式. 3 例题分析:利用SOR 方法解方程组⎪⎩⎪⎨⎧=+---=-+-=--3322242024321321321x x x x x x x x x (1) 其准确解为X *={1, 1, 2}. 建立与式(1)相等价的形式:⎪⎪⎩⎪⎪⎨⎧++=-+=+=132315.05.05.025.05.021*******x x x x x x x x x (2) 据此建立迭代公式:⎪⎪⎩⎪⎪⎨⎧++=-+=+=+++132315.05.05.025.05.0)(2)(1)1(3)(3)(1)1(23)(2)1(1k k k k k k kk k x x x x x x x x x (3)利用SOR 算法,取迭代初值1)0(3)0(2)0(1===x x x ,ω=1.5,迭代结果如下表.逐次超松弛迭代法次数 x1 x2 x3 1 0.625000 0.062500 1.750000 2 0.390625 0.882813 1.468750 3 1.017578 0.516602 1.8085944 0.556885 0.880981 1.7104495 1.023712 0.743423 1.8681036 0.746250 0.908419 1.8387377 0.997715 0.860264 1.9138948 0.864050 0.936742 1.9086059 0.986259 0.922225 1.94552310 0.928110 0.958649 1.94749311 0.985242 0.955944 1.96619812 0.961661 0.973818 1.96952113 0.988103 0.974699 1.97928914 0.979206 0.983746 1.98217215 0.991521 0.985318 1.98741616 0.988509 0.990038 1.98951317 0.994341 0.991414 1.99239718 0.993538 0.993946 1.99380619 0.996367 0.994950 1.99542420 0.996313 0.996342 1.99633121 0.997724 0.997018 1.99725422 0.997871 0.997798 1.99782223 0.998596 0.998234 1.998355GS迭代法须迭代85次得到准确值X*={1, 1, 2};而SOR方法只须55次即得准确值.由此可见,适当地选择松弛因子ω,SOR法具有明显的加速收敛效果. □§3.4 迭代法的收敛性1. 向量和矩阵范数 (a) 向量范数R n 空间的向量范数 || · || ,对任意n R y x ∈,, 满足下列条件:00||||;0||||)1(=⇔=≥x x x (正定性)||||||||||)2(x x⋅=αα (齐次性)||||||||||||)3(y x y x+≤+ (三角不等式)常见的向量范数有: (1) 列范数:(2) 谱范数:(欧几里德范数或向量的长度,模)(3) 行范数:(4) p 范数:上述范数的几何意义是:∞||||x =max(|x 2-x 1|,|y 2-y 1|) ; 1||||x =|x 2-x 1|+|y 2-y 1| ;2122122)()(||||y y x x x -+-=.向量序列}{)(k x依坐标收敛于向量x * 的充要条件是向量序列}{)(k x 依范数收敛于向量x *,即0||||lim *)(=-∞→x x k k .(b) 矩阵范数n m R ⨯空间的向量范数 || ·|| ,对任意 n m R B A ⨯∈,, 满足下列条件:|||||||| || AB || (4)||||||||||||)3(||||||||||)2(00||||;0||||)1(B A B A B A A A A A A ≤+≤+⋅==⇔=≥αα常见的矩阵范数有:∑==∞≤≤nj ij a A ni 1||max ||||1 (行和范数)∑==≤≤ni ij a A nj 11||max ||||1 (列和范数))(||||m ax 2A A A T λ= (谱范数)若A 对称,则有)()(2m ax m ax A A A T λλ=.矩阵A 的谱半径记为)(||||2A A ρ=,ρ(A ) =||m ax 1i ni λ≤≤,其中λi 为A 的特征根。

数值分析线性方程组的迭代解法

数值分析线性方程组的迭代解法

数值分析课程实验报告实验名称 线性方程组的迭代解法Ax b =的系数矩阵对角线元素容许误差。

雅可比(Jacobi )迭代法解方程组的算法描述如下:任取初始向量(0)(0)1(xx =1+,并且 1,2,...,n ,计算 11(ni j ii j ib a a =≠-∑()k x ,结束;否则执行④,则不收敛,终止程序;否则转② 迭代法的算法描述)迭代法中,如果当新的分量求出后,马上用它来代替旧的分量,则可能会更快地接近方程组的准确解。

基于这种设想构造的迭代公式,n ,k = (2)算法可相应地从雅可比(Jacobi )迭代法改造得到(Gauss-Seidel)迭代得到的值进一()()()1((1k i ii k k i i x b a x x ωω==+-1,2,,n ,2,k =(3)为松弛因子(显然当1ω=塞德尔迭代公式) ()k ix 通常优于旧值(1)k ix -,在将两者加工成松弛值时,自然要求松弛因子1ω>,以尽量发挥新值的优势,这类迭代就称为逐次超松弛迭代法。

SOR 迭代的关键在于选取合适的松弛因子,松弛因子的取值对收敛速度影响很大,但如何选取最佳松弛因子的问题,至今仍未有效解决,在实际计算时,通常依据系数矩阵的特点,并结合以往的经验选取合适的松弛因子。

练习与思考题分析解答(0)(1,1,1,1)x =[ -0.999976, -0.999976, -0.999976, -0.999976]x =[ -0.99999, -0.999991, -0.999992, -0.999993]x =塞德尔迭代算法的收敛速度要比雅可比迭代算法的收敛速度快SOR 迭代实质上是高斯原理和基本方法相同。

如果选择合适的松弛因子,它能够加快收敛速度。

SOR 迭代算法更加普通,当选取一个合适的松弛因子后收敛速度明显加快。

迭代算法将前一步的结果[ -0.99999, -0.999991, -0.999992, -0.999993]x =[ -0.999992, -0.999993, -0.999994, -0.999995]x =[ -0.999993, -0.999994, -0.999995, -0.999995]x =[ -0.999992, -0.999993, -0.999994, -0.999995]x =[ -0.999999, -1.0, -1.0, -1.0]x =[ -0.999999, -1.0, -1.0, -1.0]x =因为为了保证迭代过程收敛,松弛因子1.3左右。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
0 0 0 0 0 0 1 3 2 0
0 0 0 0 0 0 0 1 3 2
0 0 0 0 0 0 0 0 1 3];
[n,n]=size(A);
L=zeros(n,n);
U=zeros(n,n);
fori=1:n
L(i,i)=1;
end
fork=1:n
forj=k:n
U(k,j)=A(k,j)-sum(L(k,1:k-1).*U(1:k-1,j)');
0 0 0 1 1 0 0 0 0 0
0 0 0 0 1 1 0 0 0 0
0 0 0 0 0 1 1 0 0 0
0 0 0 0 0 0 1 1 0 0
0 0 0 0 0 0 0 1 1 0
0 0 0 0 0 0 0 0 1 1
2、求解下面代数方程组的解
x =
1
2
3
4
5
6
7
8
9
10
3、jacobi迭代求下面方程组的解
数值分析实验报告
班级
姓名
学号
课程名称
数值分析
时间
地点
一、实验名称:
实验四、线性代数方程组直接解法、迭代法
二、实验目的:
1.高斯消去法的原理和计算步骤,理解迭代法基本原理及步骤;
2.了解顺序消去法能够实现的条件;
3.掌握矩阵的三角分解(即LU分解),平方分解的原理与步骤;
4. 利用LU分解、平方分解求解线性代数方程组.
5.利用Jacobi,Gauss-seidel迭代求解方程组
三、实验内容及要求:
1.编写LU分解matlab通用程序,并用该程序求下面矩阵的LU分解
试写出 , 阵。
2.编写平方分解matlab通用程序,并用该程序求解下面代数方程组的解:
其中
3.编写Jacobi、Gauss-seidel迭代法求解线性代数方程组通用程序,并利用Jacobi、Gauss-seidel迭代求下面方程组的解,精度控制为 (误差最大范数).
End
以Fjacobi.m保存
a=[28 1 3;2 8 1; 2 -3 15];
b=[32 11 14]';
x0=[0 0 0]';
[x,k]=Fjacobi(a,b,x0,1e-3)
以jacobi保存
再运行程序jacobi
max1=300;%默认迭代300次
D=diag(diag(A));%求a的对角矩阵
end
fori=k+1:n
L(i,k)=(A(i,k)-sum(L(i,1:k-1).*U(1:k-1,k)'))/U(k,k);
end
end
2、a=[1 2 0 0 0 0 0 0 0 0
2 5 2 0 0 0 0 0 0 0
0 2 5 2 0 0 0 0 0 0
0 0 2 5 2 0 0 0 0 0
x =
1.0e+005 *
-0.3211
-0.3203
1.1929
k =
300
Gauss-seidel迭代求下面方程组的解
x =
14.7758
-19.5433
-5.8788
k =
300
六、讨论
直接法可以求得精确解是指就计算公式而言保证得到精确解,但计算机计算过程中的舍入误差是不可避免的,这种误差对解的精度影响会不会太大,也就是计算的稳定性,是要考虑的问题。直接法普遍适用,但要求计算机有较大的存储量,迭代法要求的存储量较小,但必须在收敛性得以保证的情况下才能使用。对于迭代法,其收敛性则是要考虑的问题。
U=-tril(A,1);%求a的上角
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>=maxl)
disp('迭代超过300次,方程可能不收敛');
return;
End%[k x']如果想显示每一步的迭代结果
return;
end
%[k x']如果想显示每一步的迭代结果
End
以Fgseid.m保存
a=[28 1 3;2 8 1; 2 -3 15];
b=[32 11 14]';
x0=[0 0 0]';
[x,k]=Fgseid(a,b,x0,1e-3)
以gseid.m保存
再运行程序gseid
八、成绩评定
教师签字:年月日
L=-tril(A,-1);%求a的下三角
U=-tril(A,1);%求a的上三角
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('迭代超过300次,方程可能不收敛');
七、程序附录
1、function[L,U]=myLU(A)
A=[1 2 0 0 0 0 0 0 0 0
1 3 2 0 0 0 0 0 0 0
0 1 3 2 0 0 0 0 0 0
0 0 1 3 2 0 0 0 0 0
0 0 0 1 3 2 0 0 0 0
0 0 0 0 1 3 2 0 0 0
0 0 0 0 0 1 3 2 0 0
[l,u]=lu(a);
x=u\(l\b)
以LUfj.m保存文件
在执行LUfj
3、function [x,k]=Fjacobi(A,b,x0,tol)
%tol是输入误差容限,x0是初值
maxl=300;%默认迭代300次
D=diag(diag(A));%求a的对角矩阵
L=-tril(A,-1);%求a的下三角
四、直接法的实验步骤:
步骤1:打开matlab,存储行列式值,并选主元
步骤2:对矩阵进行消元,再回代
步骤3:编写函数代码,保存并运行

五、直接法的实验结果:
1、矩阵分解结果
ans =
1 0 0 0 0 0 0 0 0 0
1 1 0 0 0 0 0 0 0 0
0 1 1 0 0 0 0 0 0 0
0 0 1 1 0 0 0 0 0 0
0 0 0 2 5 2 0 0 0 0
0 0 0 0 2 5 2 0 0 0
0 0 0 0 0 2 5 2 0 0
0 0 0 0 0 0 2 5 2 0
0 0 0 0 0 0 0 2 5 2
0 0 0 0 0 0 0 0 2 5];
b=[5 18 27 36 45 54 63 72 81 68]';
相关文档
最新文档