数值代数实验2012
数值代数实验报告

数值代数实验报告数值代数实验报告引言:数值代数是一门研究数值计算方法和算法的学科,它在科学计算和工程应用中起着重要的作用。
本实验报告旨在通过实际的数值计算问题,探讨数值代数的应用和效果。
实验一:线性方程组求解线性方程组求解是数值代数中的一个重要问题。
在实验中,我们使用了高斯消元法和LU分解法两种求解线性方程组的方法,并对比了它们的效果。
首先,我们考虑一个3×3的线性方程组:2x + 3y - z = 54x - 2y + 2z = 1x + y + z = 3通过高斯消元法,我们将该方程组转化为上三角形式,并得到解x=1, y=2, z=0。
而通过LU分解法,我们将该方程组分解为LU两个矩阵的乘积,并得到相同的解。
接下来,我们考虑一个更大的线性方程组,例如10×10的方程组。
通过比较高斯消元法和LU分解法的运行时间,我们可以发现LU分解法在处理大规模方程组时更加高效。
实验二:特征值与特征向量计算特征值与特征向量计算是数值代数中的另一个重要问题。
在实验中,我们使用了幂法和QR方法两种求解特征值与特征向量的方法,并对比了它们的效果。
首先,我们考虑一个3×3的矩阵:1 2 34 5 67 8 9通过幂法,我们可以得到该矩阵的最大特征值为15.372,对应的特征向量为[0.384, 0.707, 0.577]。
而通过QR方法,我们也可以得到相同的结果。
接下来,我们考虑一个更大的矩阵,例如10×10的矩阵。
通过比较幂法和QR 方法的运行时间,我们可以发现QR方法在处理大规模矩阵时更加高效。
实验三:奇异值分解奇异值分解是数值代数中的一种重要技术,它可以将一个矩阵分解为三个矩阵的乘积,从而实现数据降维和信息提取的目的。
在实验中,我们使用了奇异值分解方法,并通过实际的数据集进行了验证。
我们选取了一个包含1000个样本和20个特征的数据集,通过奇异值分解,我们将该数据集分解为三个矩阵U、S和V的乘积。
数值计算方法实验报告

数值计算方法实验报告一、实验介绍本次实验是关于数值计算方法的实验,旨在通过计算机模拟的方法,实现对于数值计算方法的掌握。
本次实验主要涉及到的内容包括数值微积分、线性方程组的求解、插值与拟合、常微分方程的数值解等。
二、实验内容1. 数值微积分数值微积分是通过计算机模拟的方法,实现对于微积分中的积分运算的近似求解。
本次实验中,我们将会使用梯形公式和辛普森公式对于一定区间上的函数进行积分求解,并比较不同公式的计算误差。
2. 线性方程组的求解线性方程组求解是数值计算领域中的重要内容。
本次实验中,我们将会使用高斯消元法、LU分解法等方法对于给定的线性方程组进行求解,并通过比较不同方法的计算效率和精度,进一步了解不同方法的优缺点。
3. 插值与拟合插值与拟合是数值计算中的另一个重要内容。
本次实验中,我们将会使用拉格朗日插值法和牛顿插值法对于给定的数据进行插值求解,并使用最小二乘法对于给定的函数进行拟合求解。
4. 常微分方程的数值解常微分方程的数值解是数值计算中的难点之一。
本次实验中,我们将会使用欧拉法和龙格-库塔法等方法对于给定的常微分方程进行数值解的求解,并比较不同方法的计算精度和效率。
三、实验结果通过本次实验,我们进一步加深了对于数值计算方法的理解和掌握。
在数值微积分方面,我们发现梯形公式和辛普森公式都能够有效地求解积分,但是辛普森公式的计算精度更高。
在线性方程组求解方面,我们发现LU分解法相对于高斯消元法具有更高的计算效率和更好的数值精度。
在插值与拟合方面,我们发现拉格朗日插值法和牛顿插值法都能够有效地进行插值求解,而最小二乘法则可以更好地进行函数拟合求解。
在常微分方程的数值解方面,我们发现欧拉法和龙格-库塔法都能够有效地进行数值解的求解,但是龙格-库塔法的数值精度更高。
四、实验总结本次实验通过对于数值计算方法的模拟实现,进一步加深了我们对于数值计算方法的理解和掌握。
在实验过程中,我们了解了数值微积分、线性方程组的求解、插值与拟合、常微分方程的数值解等多个方面的内容,在实践中进一步明确了不同方法的特点和优缺点,并可以通过比较不同方法的计算效率和数值精度来选择合适的数值计算方法。
数值计算基础实验报告(3篇)

第1篇一、实验目的1. 理解数值计算的基本概念和常用算法;2. 掌握Python编程语言进行数值计算的基本操作;3. 熟悉科学计算库NumPy和SciPy的使用;4. 分析算法的数值稳定性和误差分析。
二、实验内容1. 实验环境操作系统:Windows 10编程语言:Python 3.8科学计算库:NumPy 1.19.2,SciPy 1.5.02. 实验步骤(1)Python编程基础1)变量与数据类型2)运算符与表达式3)控制流4)函数与模块(2)NumPy库1)数组的创建与操作2)数组运算3)矩阵运算(3)SciPy库1)求解线性方程组2)插值与拟合3)数值积分(4)误差分析1)舍入误差2)截断误差3)数值稳定性三、实验结果与分析1. 实验一:Python编程基础(1)变量与数据类型通过实验,掌握了Python中变量与数据类型的定义方法,包括整数、浮点数、字符串、列表、元组、字典和集合等。
(2)运算符与表达式实验验证了Python中的算术运算、关系运算、逻辑运算等运算符,并学习了如何使用表达式进行计算。
(3)控制流实验学习了if-else、for、while等控制流语句,掌握了条件判断、循环控制等编程技巧。
(4)函数与模块实验介绍了Python中函数的定义、调用、参数传递和返回值,并学习了如何使用模块进行代码复用。
2. 实验二:NumPy库(1)数组的创建与操作通过实验,掌握了NumPy数组的基本操作,包括创建数组、索引、切片、排序等。
(2)数组运算实验验证了NumPy数组在数学运算方面的优势,包括加、减、乘、除、幂运算等。
(3)矩阵运算实验学习了NumPy中矩阵的创建、操作和运算,包括矩阵乘法、求逆、行列式等。
3. 实验三:SciPy库(1)求解线性方程组实验使用了SciPy库中的线性代数模块,通过高斯消元法、LU分解等方法求解线性方程组。
(2)插值与拟合实验使用了SciPy库中的插值和拟合模块,实现了对数据的插值和拟合,并分析了拟合效果。
数值代数上机实验报告

数值代数课程设计实验报告姓名: 班级: 学号: 实验日期:一、实验名称 代数的数值解法 二、实验环境 MATLAB7.0实验一、平方根法与改进平方根法一、实验要求:用熟悉的计算机语言将不选主元和列主元Gasuss 消元法编写成通用的子程序,然后用编写的程序求解下列方程组⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡--⨯14151515157681681681681681612321n n n n n x x x x x x 用所编的程序分别求解40、84、120阶方程组的解。
二、算法描述及实验步骤GAuss 程序如下:(1)求A 的三角分解:LU A =;(2)求解b y =L 得y ; (3)求解y x =U 得x ;列主元Gasuss 消元法程序如下: 1求A 的列主元分解:LU PA =;2求解b y P L =得y ; 3求解y x =U 得x ;三、调试过程及实验结果:%----------------方程系数---------------->> A1=Sanduijiaozhen(8,6,1,40); >> A2=Sanduijiaozhen(8,6,1,84); >> A3=Sanduijiaozhen(8,6,1,120); >> b1(1)=7;b2(1)=7;b3(1)=7;>> for i=2:39b1(i)=15;end>> b1(40)=14;>> for i=2:83b2(i)=15;end>> b2(40)=14;>> for i=2:119b1(i)=15;end>> b3(120)=14;%----------------方程解---------------->> x11=GAuss(A1,b1')>> x12=GAuss Zhu(A1,b1')>> x21=GAuss(A2,b2')>> x22=GAuss Zhu(A3,b3')>> x31=GAuss(A3,b3')>> x32=GAuss Zhu(A3,b3')运行结果:(n=40)GAuss消元法的解即为x11 =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000列主元GAuss消元法的解即为x12 =1 1 1 1 1 1 1 1 1 1 111111111111111111111111111111六、源程序:function A=Sanduijiaozhen(a,b,c,n)%生成n阶以a,b,c为元素的三对角阵A=diag(b*ones(1,n),0)+diag(c*ones(1,n-1),1)+diag(a*ones(1,n-1),-1);function x=GAuss(A,b)n=length(b);x=b;%-------分解---------------for i=1:n-1for j=i+1:nmi=A(j,i)/A(i,i);b(j)=b(j)-mi*b(i);for k=i:nA(j,k)=A(j,k)-mi*A(i,k);endAB=[A,b];endend%-----------回代------------------x(n)=b(n)/A(n,n);for i=n-1:-1:1s=0;for j=i+1:ns=s+A(i,j)*x(j);endx(i)=(b(i)-s)/A(i,i);endfunction x=GAussZhu(A,b)n=length(b);x=b;%----------------------选主元---------------------for k=1:n-1a_max=0;for i=k:nif abs(A(i,k))>a_maxa_max=abs(A(i,k));r=i;endendif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;end%--------------消元-----------------for i=k+1:nm=A(i,k)/A(k,k);for j=k:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);endendif abs(A(n,n))==0return;endAbZhu=[A,b];%----------------回代-----------------------x(n)=b(n)/A(n,n);for i=n-1:-1:1for j=i+1:nb(i)=b(i)-A(i,j)*x(j);endx(i)=b(i)/A(i,i);end实验二、平方根法与改进平方根法一、实验要求:用计算机语言将平方根法和改进的平方根法编成通用的子程序,然后用编写的程序求解对称正定方程组100阶方程组AX=b,二、算法描述及实验步骤:平方根法函数程序如下:1、求A 的Cholesky 分解:L L A T=;2、求解b y =L 得y ;3、求解y x =TL 得x ; 改进平方根法函数程序如下:1、求A 的Cholesky 分解:T=LDL A ; 2、求解b y =L 得y ; 3、求解y x =TDL 得x ;三、调试过程及实验结果:clear;clc;%----------------方程系数---------------->> A=Sanduijiaozhen(1,10,1,100); >> b(1)=11; >> for i=2:99 b(i)=12; end>> b(100)=11;>> x1=Cholesky(A,b') >> x2=GJCholesky(A,b')运行结果:平方根法的解即为 x1 =1.0000 1.00001.0000 1.00001.0000 1.00001.0000 1.00001.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000改进平方根法解得的解即为x2 =1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00000.99991.00090.99081.09080.1010四、源程序:function x=Cholesky(A,b)n=size(A);n=n(1);% x=A^-1*b;% disp('Matlab自带解即为x');%-----------------Cholesky分解-------------------for k=1:nA(k,k)=sqrt(A(k,k));A(k+1:n,k)=A(k+1:n,k)/A(k,k);for j=k+1:n;A(j:n,j)=A(j:n,j)-A(j:n,k)*A(j,k);endend%------------------前代法求解Ly=b----------------------------for j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n);%-----------------回代法求解L'x=y-----------------------------A=A';for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1);disp('平方根法的解即为');function b=GJCholesky(A,b)n=size(A);n=n(1);v=zeros(n,1);%----------------------LDL'分解-----------------------------for j=1:nfor i=1:j-1v(i)=A(j,i)*A(i,i);endA(j,j)=A(j,j)-A(j,1:j-1)*v(1:j-1);A(j+1:n,j)=(A(j+1:n,j)-A(j+1:n,1:j-1)*v(1:j-1))/A(j,j);endB=diag(A);D=zeros(n);for i=1:nD(i,i)=B(i);A(i,i)=1;end%-------------------前代法---------------------------A=tril(A); %得到L和Dfor j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n);%-----------------回代法-----------------------------A=D*(A');for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1);disp('改进平方根法解得的解即为');实验三、二次多项式拟合一、实验要求:用计算机语言编制利用QR分解求解线性最小二乘问题的通用子程序,用编写的程序求解一个二次多项式使在残向量的范数最小的意义下拟合下面的数据t-1 -0.75 -0.5 0 0.25 0.5 0.75iy 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125i二、算法描述及实验步骤:QR分解求解程序如下:1、求A 的QR 分解;2、计算b c 11T =Q ;3、求解上三角方程1c x =R 得x ;三、调试过程及实验结果:>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75];>> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> plot(t,y,'r*');>> legend('实验数据(ti,yi)'); >> xlabel('t'), ylabel('y');>> title('二次多项式拟合的数据点(ti,yi)的散点图');运行后屏幕显示数据的散点图(略).(3)编写下列MATLAB 程序计算)(x f 在),(i i y x 处的函数值,即输入程序 >> syms a b c>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75]; >> fi=a.*t.^2+ b.*t+c%运行后屏幕显示关于 ,,a b c 的线性方程组fi =[a-b+c,9/16*a-3/4*b+c,1/4*a-1/2*b+c,c,1/16*a+1/4*b+c,1/4*a+1/2*b+c,9/16*a+3/4*b +c]编写构造残向量2范数的MATLAB 程序>> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> fy=fi-y; fy2=fy.^2; J=sum(fy.^2); 运行后屏幕显示误差平方和如下 J=(a-b+c-1)^2+(9/16*a-3/4*b+c-13/16)^2+(1/4*a-1/2*b+c-3/4)^2+(c-1)^2+(1/16*a+1/4*b+c-21/16)^2+(1/4*a+1/2*b+c-7/4)^2+(9/16*a+3/4*b+c-37/16)^2为求,,a b c 使J 达到最小,只需利用极值的必要条件0J a ∂=∂,0J b ∂=∂,0J c∂=∂,得到关于,,a b c 的线性方程组,这可以由下面的MATLAB 程序完成,即输入程序 >> Ja1=diff(J,a); Ja2=diff(J,b); Ja3=diff(J,c);>> Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3) 运行后屏幕显示J 分别对,,a b c 的偏导数如下 Ja11 =451/128*a-63/32*b+43/8*c-887/128 Ja21 =-63/32*a+43/8*b-3/2*c-61/32Ja31 =43/8*a-3/2*b+14*c-143/8解线性方程组112131000Ja Ja Ja ===,,,输入下列程序 >> A=[451/128, -63/32, -3/2 ;-63/32,43/8,-3/2;43/8,-3/2,14]; >> B=[887/128,61/32,143/8]; >> C=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f 及其系数C 如下 C =0.3081 0.8587 1.4018 f =924/2999*x^2+10301/11996*x+4204/2999 故所求的拟合曲线为2()0.30810.8581 1.4018f x x x =++四、源程序:>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75];>> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> plot(t,y,'r*');>> legend('实验数据(ti,yi)'); >> xlabel('t'), ylabel('y');>> title('二次多项式拟合的数据点(ti,yi)的散点图'); >> syms a b c>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75]; >> fi=a.*t.^2+ b.*t+c fi =[ a-b+c, 9/16*a-3/4*b+c, 1/4*a-1/2*b+c, c, 1/16*a+1/4*b+c, 1/4*a+1/2*b+c, 9/16*a+3/4*b+c]>> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> fy=fi-y; fy2=fy.^2; J=sum(fy.^2) J =(a-b+c-1)^2+(9/16*a-3/4*b+c-13/16)^2+(1/4*a-1/2*b+c-3/4)^2+(c-1)^2+(1/16*a+1/4*b+c-21/16)^2+(1/4*a+1/2*b+c-7/4)^2+(9/16*a+3/4*b+c-37/16)^2>> Ja1=diff(J,a); Ja2=diff(J,b); Ja3=diff(J,c);>> Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3)Ja11 =451/128*a-63/32*b+43/8*c-887/128Ja21 =-63/32*a+43/8*b-3/2*c-61/32Ja31 =43/8*a-3/2*b+14*c-143/8>> A=[451/128, -63/32, -3/2 ;-63/32,43/8,-3/2;43/8,-3/2,14]; >> B=[887/128,61/32,143/8];>> C=B/A, f=poly2sym(C)C =0.3081 0.8587 1.4018f =924/2999*x^2+10301/11996*x+4204/2999>>。
数值计算方法I实验报告

实验报告实验课程名称数值计算方法I开课实验室数学实验室学院理学院年级2012 专业班信息与计算科学2班学生姓名学号开课时间2012 至2013 学年第 2 学期实验一 误差分析试验1.1(病态问题)问题提出:考虑一个高次的代数多项式)1.1()()20()2)(1()(201∏=-=---=k k x x x x x p显然该多项式的全部根为1,2,…,20共计20个,且每个根都是单重的。
现考虑该多项式的一个扰动)2.1(0)(19=+x x p ε其中ε是一个非常小的数。
这相当于是对(1.1)中19x 的系数作一个小的扰动。
我们希望比较(1.1)和(1.2)根的差别,从而分析方程(1.1)的解对扰动的敏感性。
实验内容:为了实现方便,我们先介绍两个MA TLAB 函数:“roots ”和“poly ”。
roots(a)u =其中若变量a 存储n+1维的向量,则该函数的输出u 为一个n 维的向量。
设a 的元素依次为121,,,+n a a a ,则输出u 的各分量是多项式方程01121=+++++-n n n n a x a x a x a的全部根;而函数 poly(v)b =的输出b 是一个n+1维向量,它是以n 维向量v 的各分量为根的多项式的系数。
可见“roots ”和“poly ”是两个互逆的运算函数。
))20:1((;)2();21,1(;000000001.0ve poly roots ess ve zeros ve ess +===上述简单的MA TLAB 程序便得到(1.2)的全部根,程序中的“ess ”即是(1.2)中的ε。
实验要求:(1)选择充分小的ess ,反复进行上述实验,记录结果的变化并分析它们。
如果扰动项的系数ε很小,我们自然感觉(1.1)和(1.2)的解应当相差很小。
计算中你有什么出乎意料的发现?表明有些解关于如此的扰动敏感性如何?(2)将方程(1.2)中的扰动项改成18x ε或其它形式,实验中又有怎样的现象? (3)(选作部分)请从理论上分析产生这一问题的根源。
数值代数(第08周)实验报告()

A=sparse(i,j,s,n,n,3*n-2);
spy(A);
x0=A\b;
[x]=Zuiganfa(a,c,d,b);
ifnorm(x'-x0) > 1.0e-6
error('误差太大');
end
2.2.3实验(或测试)结果
得分:
??? Error using ==> Untitled2 at 14
(2)计算矩阵 的克拉特(Crout)分解
其中
【注:在实际计算中,将数组alpha和beta分别保存在数组b和c中.】
(3)解
(4)解
2.2程序清单
2.2.1写出追赶法求解三对角方程组的函数代码
得分:
function[x]=Zuiganfa(a,c,d,b);
n=length(a);
n1=length(c);
数值代数(第08周)实验报告
学号
20150
姓名
实验目的
追赶法实验
1、实验要求
给定一个 级主对角占优的三对角矩阵 ,和一个 的列向量 ,用追赶法求解线性方程组 ,其中
注:尝试用规模尽量大(即n尽量大)的矩阵来完成实验.
2、追赶法
2.1算法描述
得分:
(写出追赶法求解三对角方程组的计算过程)
(1)将矩阵 的三条对角线和常数向量 分别记为数组
%U(i,i+1)=q(i);
%end
end
2.2.2实验(或测试)代码
得分:
n=100000;
c=rand(1,n)*100;
d=rand(1,n)*100;
a = rand(1,n)*1000+c+d;
《数值计算方法》教学大纲

河北联合大学第2012-2013-1学期《数值计算方法》教学大纲依据我校章程,特制定了适合我校理工科各专业本科生的《数值计算方法》教学大纲。
一、课程计划课程名称:数值计算方法Numerical Calculation Methods开课单位:理学院课程类型:专业必修课开设学期:第五学期讲授学时:共15周,每周4学时,共60学时学时安排:课堂教学44学时+实验教学16学时适用专业:信科、数学、统计理科专业本科生教学方式:讲授(多媒体为主)+上机考核方式:闭卷40% +上机实验20%+课程报告20% +平时成绩10%学分:4学分与其它课程的联系预修课程:数学分析、高等代数、常微分方程、计算机高级语言等。
后继课程:偏微分方程数值解及其它专业课程。
二、课程介绍数值计算方法也称为数值分析,是研究用计算机求解各种数学问题的数值方法及其理论的一门学科。
随着计算科学与技术的进步和发展,科学计算已经与理论研究、科学实验并列成为进行科学活动的三大基本手段,作为一门综合性的新科学,科学计算已经成为了人们进行科学活动必不可少的科学方法和工具。
数值计算方法是科学计算的核心内容,它既有纯数学高度抽象性与严密科学性的特点,又有应用的广泛性与实际实验的高度技术性的特点,是一门与计算机使用密切结合的实用性很强的数学课程。
主要介绍数值计算的误差、插值法、函数逼近与曲线拟合、线性方程组迭代解法、数值积分与数值微分、非线性方程组解法、矩阵特征值与特征向量数值计算以及常微分方程数值解,并特别加强实验环节的训练以提高学生动手能力。
通过本课程的学习,不仅能使学生初步掌握数值计算方法的基本理论知识,了解算法设计及数学建模思想,而且能使学生具备一定的科学计算能力和分析与解决问题的能力,不仅为学习后继课程打下良好的理论基础,也为将来从事科学计算、计算机应用和科学研究等工作奠定必要的数学基础。
教学与实验教学课堂教学实验教学论文报告机动课内学时课外学时学时数44 16 8 2 60 10三、重点难点课程重点:理解各种常用数值计算方法的数学原理和理论分析过程,掌握各种数值计算方法的示范性上机程序,学会设计数值算法的基本思路、一般原理和各种数值算法的程序实现。
共轭梯度法优缺点

共轭梯度法优缺点数值代数实验报告一、实验名称:用共轭梯度法解线性方程组。
二、实验目的:进一步熟悉理解掌握共轭梯度法解法思路,提高matlab 编程能力三、实验要求:已知线性方程矩阵,应用共轭梯度法在相关软件编程求解线性方程组的解。
四、实验原理:1.共轭梯度法:考虑线性方程组A某b的求解问题,其中A是给定的n阶对称正定矩阵,b是给定的n维向量,某是待求解的n维向量.为此,定义二次泛函TT(某)某A某2b某.定理1设A对称正定,求方程组A某b的解,等价于求二次泛函(某)的极小值点.定理1表明,求解线性方程组问题就转化为求二次泛函(某)的极小值点问题.求解二次函数极小值问题,通常好像盲人下山那样,先给定一个初始向量某0,确定一个下山方向p0,沿着经过点某0而方向为p0的直线某某0p0找一个点某某p,1000使得对所有实数有某p某p,00000即在这条直线上某1使(某)达到极小.然后从某1出发,再确定一个下山的方向p1,沿着直线某某p再跨出一步,即找到1使得某在某2某11p1达到极小:11某p某p.11111重复此步骤,得到一串0,1,2,L和p0,p1,p2,L,称pk为搜索方向,k为步长.一般情况下,先在某k点找下山方向pk,再在直线某某p上确定步长k使kk某p某pkkkkk,最后求出某某p.然而对不同的搜索方向和步长,得到各种不同的算法.k1kkk1由此,先考虑如何确定k.设从某k出发,已经选定下山方向pk.令f某pkkTT某pA某p2b某pkkkkkk2T2TpAprp某,kkkkk其中rkbApk.由一元函数极值存在的必要条件有TTf2pAp2rp0kkkk所确定的即为所求步长,即k步长确定后,即可算出TrpkkkTpApkk.某某p.k1kkkT此时,只要rp0,就有kk即某某.k1k某某某p某2k1kkkkkTrpkk2T2T0pAprpkkkkkkTpApkk再考虑如何确定下山方向pk.易知负梯度方向是(某)减小最快的方向,但简单分析就会发现负梯度方向只是局部最佳的下山方向,而从整体来看并非最佳.故采用新的方法寻求更好的下山方向——共轭梯度法.下面给出共轭梯度法的具体计算过程:给定初始向量某0,第一步仍选用负梯度方向为下山方向,即p0r0,于是有Trr00,某某p,rbA某.0T100010pAp00对以后各步,例如第k+1步(k1),下山方向不再取rk,而是在过点由向量rk和pk1所张成的二维平面2{某|某某krkpk1,,R}内找出使函数下降最快的方向作为新的下山方向pk.考虑在2上的限制:,(某krkpk)1T(某rp)A(某rp)kkk1kkk1T2b(某rp).kkk1计算关于,的偏导得:TTT2rArrAprr,kkkk1kkTT2rAppAp,T其中最后一式用到了rp10,这可由rk的定义直接验证.令kk1k1k1kk0,即知在2内有唯一的极小值点某%某rp,k0k0k12其中0和0满足TTTrArrAprr0kk0kk1kk,TTrAppAp0kk10k1k10.由于rk0必有00,所以可取1p某%某rpkkkk100作为新的下山方向.显然,这是在平面2内可得的最佳下山方向.令0 k,则可1得Tkk1.k1TpApk1k1T注:这样确定的pk满足pAp10,即pk与pk1是相互共轭的.kk总结上面的讨论,可得如下的计算公式:Trpkk,某k1某kkpk,kTpApkkrbA某,k1k1Tk1k,pk1rk1kpk.kTpApkk在实际计算中,常将上述公式进一步简化,从而得到一个形式上更为简单而且对称的计算公式.首先来简化rk1的计算公式:r1bA某1bA(某p)rAp.kkkkkkkk因为Apk在计算k是已经求出,所以计算rk1时可以不必将某k1代入方程计算,而是从递推关系rk1bkApk得到.再来简化k和k的计算公式.此处需要用到关系式TTTrr1rp1r1p0,k1,2,K.kkkkkk从而可导出1rrr,,k1k1k111kTTTpApprrprkkkkk1kk1T1Tkkrrprr.kkk1k1kkkk由此可得TTrrrrkk11.kk,.kTkTpAprrkkkk从而有求解对称正定方程组的共轭梯度法算法如下:某初值rbA某;k000whilerk0kk1ifk1pr003eleTTkrkrkrkrk21122k1k1k2k2endTTkrkrkpkApk11111某某pkk1k1k1rrApkk1k1k1end某某k注:该算法每迭代一次仅需要使用系数矩阵A做一次矩阵向量积运算.定理2由共轭梯度法得到的向量组r和pi具有如下基本性质:iT(1)pr0,0ijk;T(2)rr0,ij,0i,jk;ijT(3)pAp0,ij,0i,jk;ij(4)pan{r,K,rk}pan{p,K,pk}(A,r,k1),000其中k(A,r,k1)pan{r,Ar,K,Ar},0000通常称之为Krylov子空间.下面给出共轭梯度法全局最优性定理:定理3用共轭梯度法计算得到的近似解某k满足某min某:某某(A,r,k)k00或某某某某某某Ark,某min某:0(,0,)kAA其中T某某A某,某某是方程组A某b的解,(A,r0,k)是由所定义的Krylov子空间.A定理2表明,向量组r0,K,rk和p0,K,pk分别是Krylov子空间(A,r0,k1)的正交基和共轭正交基.由此可知,共轭梯度法最多n步便可得到方程组的解某某.因此,理论上来讲,共轭梯度法是直接法.然而实际使用时,由于误差的出现,使rk之间的正交性很快损失,以致于其有限步终止性已不再成立.此外,在实际应用共轭梯度法时,由于一般n很大,以至于迭代On次所耗费的计算时间就已经使用户无法接受了.因此,实际上将共轭梯度法作为一种迭代法使用,而且通常是r是否已经很小及迭代次数是否已经达到最大允许的迭代次数kma某来终止迭代.从而得到解对称正定线性方程组的实用共轭梯度法,其算法如下:某初值4k0;rbA某;Trrwhilebandkkma某2kk1ifk1prele%;prpendTApp某某pTrr;%;rrend算法中,系数矩阵A的作用仅仅是用来由已知向量p产生向量Ap,这不仅可以充分利用A的稀疏性,而且对某些提供矩阵A较为困难而由已知向量p产生向量Ap又十分方便的应用问题是十分有益的。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数值代数实验
数值线性代数实验一
一、实验名称:矩阵的LU分解.
二、实验目的:用不选主元的LU分解和列主元LU分解求解线性方程组Ax=b, 并比较这
两种方法.
三、实验内容与要求
(1)用所熟悉的计算机语言将不选主元和列主元LU分解编成通用的子程序,然后用编写的程序求解下面的84阶方程组
将计算结果与方程组的精确解进行比较,并就此谈谈你对Gauss消去法的看法。
(2)写出追赶法求解三对角方程组的过程,并编写程序求该实验中的方程组
数值线性代数实验二
一、实验名称:实对称正定矩阵的A的Cholesky分解.
二、实验目的:用平方根法和改进的平方根方法求解线性方程组Ax=b.
三、实验内容与要求
用所熟悉的计算机语言将Cholesky分解和改进的Cholesky分解编成通用的子程序,然后用编写的程序求解对称正定方程组Ax=b,其中
(1)b随机的选取,系数矩阵为100阶矩阵
(2)系数矩阵为40阶Hilbert矩阵,即系数矩阵A的第i行第j列元素为
,向量b的第i个分量为
(3)用实验一的程序求解这两个方程组,并比较所有的计算结果,然后评价各个方法的优劣。
数值线性代数实验三
一、实验名称:矩阵A 的QR 分解
二、实验目的:应用改进的Gram —Schmidt 方法和Householder 变换的方法计算矩阵A 的
QR 分解. 其中),()(n m R a A n m ij ≥∈=⨯ rank A =n
三、实验内容与要求
输入:A 的各列),,1,],,[(,,,121n j a a a a a a T m j j j n ==
输出:Q 的各列元素(存放在A 的相应位置上)以及R 的元素
),,,,,1(n i j n i r ij ==
数值线性代数实验四
一、实验名称:用迭代法求解方程组及超松弛迭代和最佳松弛因子的确定.
二、实验目的:应用Jacobi 迭代法、Gauss -Seidel 迭代法和超松弛迭代方法求解线性
方程组 Ax=b, 并选择不同的松弛因子ω,观察松弛因子对松弛迭代法计算效果的影响.
三、实验内容与要求
将常微分方程22(0)0,(1)1d y dy a dx dx y y ε⎧+=⎪⎨⎪==⎩
(01a <<)离散化得到差分方程Ax=b ,取
10.5,100,a n h n
===,应用Jacobi 迭代法、Gauss -Seidel 迭代法和超松弛迭代方法求解线性方程组,分别取1,0.1,0.01,0.0001ε=,用SOR 迭代法计算对应的数值解,并与精确解进行比较. 写出这三种迭代法求解线性方程组的步骤,并对计算结果进行分析.
四、实验原理
将[0,1]区间n 等分,方程离散化得差分方程
211()(2)i i i h y h y y ah εεε+-+-++=
差分方程对应的系数矩阵和右端项分别为
20222(1)(1)(1)1(2)(2)(2)(2)()n n n n h h ah y h h ah A b h h ah h ah h y εεεεεεεεεεεε-⨯--⨯⎛⎫-++-⎛⎫ ⎪ ⎪-++ ⎪ ⎪ ⎪ ⎪==-+ ⎪ ⎪+ ⎪ ⎪ ⎪ ⎪-+-+⎝⎭⎝
⎭ 其中0(0),(1).n y y y y ==
SOR 迭代法的迭代矩阵和常数项分别为
11()[(1)],()B D L U D f D L b ωωωωωωω--=-+-=-
SOR 迭代法的迭代公式: 1k k y B y f ωω+=+
五、实验过程
利用Matlab 语言实现SOR 迭代公式的计算,分别对于1,0.1,0.01,0.0001ε=,松弛因子在(1,2)范围内,先按步长0.1搜索,在最少迭代步数附近按步长0.01再搜索. 算法终止准则:1410k k y y +--≤ 或者410k y y --≤或者迭代次数1000≥.
记录四种情况下的计算结果,并找出最佳松弛因子ω的值,并观察松弛因子ω与ε间的关系.。