《数值分析》课程实验报告（5篇模版）

第一篇：《数值分析》课程实验报告

《数值分析》课程实验报告 姓 名：

学 号：

学 院：

机 电 学 院 日 期：

2015 年 X 月X 日 目 录 实验一 函数插值方法 1 实验二 函数逼近与曲线拟合 5 实验三 数值积分与数值微分 7 实验四 线方程组的直接解法 9 实验五 解线性方程组的迭代法 15 实验六 非线性方程求根 19 实验七 矩阵特征值问题计算 21 实验八 常微分方程初值问题数值解法 24 实验一 函数插值方法 一、问题提出 对于给定的一元函数的n+1个节点值。试用Lagrange公式求其插值多项式或分段二次Lagrange插值多项式。

数据如下：

（1）0.4 0.55 0.65 0.80 0.95 1.05 0.41075 0.57815 0.69675 0.90 1.00 1.25382 求五次Lagrange多项式，和分段三次插值多项式，计算, 的值。（提示：结果为,）（2）1 2 3 4 5 6 7 0.368 0.135 0.050 0.018 0.007 0.002 0.001 试构造Lagrange多项式，计算的,值。（提示：结果为,）二、要求 1、利用Lagrange插值公式 编写出插值多项式程序；

2、给出插值多项式或分段三次插值多项式的表达式；

3、根据节点选取原则，对问题（2）用三点插值或二点插值，其结果如何；

4、对此插值问题用Newton插值多项式其结果如何。Newton插值多项式如下：

其中：

三、目的和意义 1、学会常用的插值方法，求函数的近似表达式，以解决其它实际问题；

2、明确插值多项式和分段插值多项式各自的优缺点；

3、熟悉插值方法的程序编制；

4、如果绘出插值函数的曲线，观察其光滑性。

四、实验步骤（1）0.4 0.55 0.65 0.80 0.95 1.05 0.41075 0.57815 0.69675 0.90 1.00 1.25382 求五次Lagrange多项式，和分段三次插值多项式，计算, 的值。（提示：结果为,）第一步：先在matlab中定义lagran的M文件为拉格朗日函数代码为：

function[c,l]=lagran(x,y)w=length(x);n=w-1;l=zeros(w,w);for k=1:n+1 v=1;for j=1:n+1 if(k~=j)v=conv(v,poly(x(j)))/(x(k)-x(j));end end l(k,:)=v;end c=y*l;end 第二步：然后在matlab命令窗口输入：

>>>> x=[0.4 0.55 0.65 0.80,0.95 1.05];y=[0.41075 0.57815 0.69675 0.90 1.00 1.25382];>> lagran(x,y)回车得到：

ans =121.6264-422.7503 572.5667-377.2549 121.9718-15.0845 由此得出所求拉格朗日多项式为 p（x）=121.6264x5-422.7503x4+572.5667x3-377.2549x2+121.9718x-15.0845 第三步：在编辑窗口输入如下命令：

>> x=[0.4 0.55 0.65 0.80,0.95 1.05];>> y=121.6264*x.^5-422.7503*x.^4+572.5667*x.^3-377.2549*x.^2+121.9718*x-15.0845;>> plot(x,y)命令执行后得到如下图所示图形，然后 >> x=0.596;>> y=121.6264*x.^5-422.7503*x.^4+572.5667*x.^3-377.2549*x.^2+121.9718*x-15.084 y =0.6262 得到f（0.596）=0.6262 同理得到f（0.99）=1.0547（2）1 2 3 4 5 6 7 0.368 0.135 0.050 0.018 0.007 0.002 0.001 试构造Lagrange多项式，和分段三次插值多项式，计算的,值。（提示：结果为,）实验步骤：

第一步定义 function[c,l]=lagran(x,y)w=length(x);n=w-1;l=zeros(w,w);for k=1:n+1 v=1;for j=1:n+1 if(k~=j)v=conv(v,poly(x(j)))/(x(k)-x(j));end end l(k,:)=v;end c=y*l;end 定义完拉格朗日M文件 第二步：然后在matlab命令窗口输入：

>>>> x=[1 2 3 4 5 6 7];y=[0.368 0.135 0.050 0.018 0.007 0.002 0.001];>> lagran(x,y)回车得到：

ans =0.0001-0.0016 0.0186-0.1175 0.4419-0.9683 0.9950 由此得出所求拉格朗日多项式为 p（x）=0.0001x6-0.0016x5+0.0186x4-0.1175x3+0.4419x2-0.9683x+0.9950 第三步：在编辑窗口输入如下命令：

>> x=[1 2 3 4 5 6 7];>> y=0.0001*x.^6-0.0016*x.^5+0.0186*x.^4-0.1175*x.^3+0.4419*x.^2-0.9683*x+0.9950;>> plot(x,y)命令执行后得到如下图所示图形，然后 >> x=1.8;>> y=121.6264*x.^5-422.7503*x.^4+572.5667*x.^3-377.2549*x.^2+121.9718*x-15.084 y =0.1650 得到f（0.596）=0.6262 同理得到f（6.15）=2.3644 五、实验结论 插值是在离散数据的基础上补插连续函数，使得这条连续曲线通过全部给定的离散数据点，它是离散函数逼近的重要方法，利用它可通过函数在有限个点处的取值状况，估算出函数在其他点处的近似值。

实验二 函数逼近与曲线拟合 一、问题提出 从随机的数据中找出其规律性，给出其近似表达式的问题，在生产实践和科学实验中大量存在，通常利用数据的最小二乘法求得拟合曲线。

在某冶炼过程中，根据统计数据的含碳量与时间关系，试求含碳量与时间t的拟合曲线。

t(分)0 5 10 15 20 25 30 35 40 45 50 55 0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.58 4.02 4.64 二、要求 1、用最小二乘法进行曲线拟合；

2、近似解析表达式为；

3、打印出拟合函数，并打印出与的误差，；

4、另外选取一个近似表达式，尝试拟合效果的比较；

5、* 绘制出曲线拟合图。

三、目的和意义 1、掌握曲线拟合的最小二乘法；

2、最小二乘法亦可用于解超定线代数方程组；

3、探索拟合函数的选择与拟合精度间的关系 四、实验步骤：

第一步先写出线性最小二乘法的M文件 function c=lspoly(x,y,m)n=length(x);b=zeros(1:m+1);f=zeros(n,m+1);for k=1:m+1 f(:,k)=x.^(k-1);end a=f＇*f;b=f＇*y＇;c=a＼b;c=flipud(c);第二步在命令窗口输入：

>>lspoly([0,5,10,15,20,25,30,35,40,45,50,55],[0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.02,4.64],2)回车得到：

ans =-0.0024 0.2037 0.2305 即所求的拟合曲线为y=-0.0024x2+0.2037x+0.2305 在编辑窗口输入如下命令：

>> x=[0,5,10,15,20,25,30,35,40,45,50,55];>> y=-0.0024*x.^2+0.2037*x+0.2305;>> plot(x,y)命令执行得到如下图 五、实验结论  分析复杂实验数据时,常采用分段曲线拟合方法。利用此方法在段内可以实现最佳逼近,但在段边界上却可能不满足连续性和可导性。分段函数的光滑算法,给出了相应的误差分析.给出了该方法在分段曲线拟合中的应用方法以及凸轮实验数据自动分段拟合。

实验三 数值积分与数值微分 一、问题提出 选用复合梯形公式，复合Simpson公式，Romberg算法，计算（1）（2）（3）（4）二、要求 1、编制数值积分算法的程序；

2、分别用两种算法计算同一个积分，并比较其结果；

3、分别取不同步长，试比较计算结果（如n = 10, 20等）；

4、给定精度要求ε，试用变步长算法，确定最佳步长。

三、目的和意义 1、深刻认识数值积分法的意义；

2、明确数值积分精度与步长的关系；

3、根据定积分的计算方法，可以考虑二重积分的计算问题。

四、实验步骤 第一步：编写各种积分的程序 复合梯形程序如下：

function I=TX(x,y)n=length(x);m=length(y);if n~=m error(＇The lengths of X and Y must be equal＇);return;end h=(x(n)-x(1))/(n-1);a=[1 2*ones(1,n-2)1];I=h/2*sum(a.*y);复合Simpson程序如下：

function s = simpr1(f,a,b,n)h=(b-a)/(2*n);s1=0;s2=0;for k=1:10 x=a+h*(2*k-1);s1=s1+feval(f,x);end for k=1:(10-1)x=a+h*2*k;s2=s2+feval(f,x);end s=h*(feval(f,a)+feval(f,b)+4*s1+2*s2)/3;end Romberg程序如下：

function I = Romber_yang(fun,a,b,ep)if nargin<4 ep=1e-5;end;m=1;h=b-a;I=h/2*(feval(fun,a)+feval(fun,b));T(1,1)=I;while 1 N=2^(m-1);h=h/2;I=I/2;for i=1:N I=I+h*feval(fun,a+(2*i-1)*h);end T(m+1,1)=I;M=2*N;k=1;while M>1;T(m+1,k+1)=(4^k*T(m+1,k)-T(m,k))/(4^k-1);M=M/2;k=k+1;end if abs(T(k,k)-T(k-1,k-1))

2、对于积分Ι=01sin⁡XXdx，f(0)=1,梯形积分T=0.94607307，辛普森积分S=0.94607308，Romberg积分R=0.94607307。

3、对于积分Ι=01eX4+X2dx，梯形积分T=0.39081248，辛普森积分S=0.39081185，Romberg积分R=0.39081885。

4、对于积分Ι=01ln1+X1+X2dx，梯形积分T=0.27218912，辛普森积分S=0.27219844，Romberg积分R=0.27219827。

五、实验结论，通过本实验学会复合梯形公式，复合Simpson公式，Romberg公式的编程与应用，掌握MATLAB提供的计算积分的各种函数的使用方法。

实验四 线方程组的直接解法 一、问题提出 给出下列几个不同类型的线性方程组，请用适当算法计算其解。

1、设线性方程组  2、设对称正定阵系数阵线方程组  3、三对角形线性方程组 二、要求 1、对上述三个方程组分别利用Gauss顺序消去法与Gauss列主元消去法；

平方根法与改进平方根法；

追赶法求解（选择其一）；

2、应用结构程序设计编出通用程序；

3、比较计算结果，分析数值解误差的原因；

4、尽可能利用相应模块输出系数矩阵的三角分解式。

三、目的和意义 1、通过该课题的实验，体会模块化结构程序设计方法的优点；

2、运用所学的计算方法，解决各类线性方程组的直接算法；

3、提高分析和解决问题的能力，做到学以致用；

4、通过三对角形线性方程组的解法，体会稀疏线性方程组解法的特点。

四、实验步骤：

列主元高斯消去法的matlab的M文件程序 function [x,det,index]=Gauss(A,b)% 求线形方程组的列主元Gauss消去法，其中，% A为方程组的系数矩阵；

% b为方程组的右端项；

% x为方程组的解；

% det为系数矩阵A的行列式的值；

% index为指标变量，index=0表示计算失败，index=1表示计算成功。

[n,m]=size(A);nb=length(b);% 当方程组行与列的维数不相等时，停止计算，并输出出错信息。

if n~=m error(＇The rows and columns of matrix A must be equal!＇);return;end % 当方程组与右端项的维数不匹配时，停止计算，并输出出错信息 if m~=nb error(＇The columns of A must be equal the length of b!＇);return;end % 开始计算，先赋初值 index=1;det=1;x=zeros(n,1);for k=1:n-1 % 选主元 a_max=0;for i=k:n if abs(A(i,k))>a_max a_max=abs(A(i,k));r=i;end end if a_max<1e-10 index=0;return;end % 交换两行 if r>k for j=k:n z=A(k,j);A(k,j)=A(r,j);A(r,j)=z;end z=b(k);b(k)=b(r);b(r)=z;det=-det;end % 消元过程 for i=k+1:n m=A(i,k)/A(k,k);for j=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);% 回代过程 if abs(A(n,n))<1e-10 index=0;return;end for k=n:-1:1 for j=k+1:n b(k)=b(k)-A(k,j)*x(j);end x(k)=b(k)/A(k,k);end 然后在命令窗口输入 >> A=[4 2-3-1 2 1 0 0 0 0;8 6-5-3 6 5 0 1 0 0;4 2-2-1 3 2-1 0 3 1;0-2 1 5-1 3-1 1 9 4;-4 2 6-1 6 7-3 3 2 3;8 6-8 5 7 17 2 6-3 5;0 2-1 3-4 2 5 3 0 1;16 10-11-9 17 34 2-1 2 2;4 6 2-7 13 9 2 0 12 4;0 0-1 8-3-24-8 6 3-1];>> b=[5 12 3 2 3 46 13 38 19-21];>> gauss(A,b)ans = 1.0000-1.0000 0.0000 1.0000 2.0000 0.0000 3.0000 1.0000-1.0000 2.0000 高斯-约当消去法maltab的M文件程序 function [x,flag]=Gau_Jor(A,b)% 求线形方程组的列主元Gauss-约当法消去法，其中，% A为方程组的系数矩阵；

% b为方程组的右端项；

% x为方程组的解；

[n,m]=size(A);nb=length(b);% 当方程组行与列的维数不相等时，停止计算，并输出出错信息。

if n~=m error(＇The rows and columns of matrix A must be equal!＇);return;end % 当方程组与右端项的维数不匹配时，停止计算，并输出出错信息 if m~=nb error(＇The columns of A must be equal the length of b!＇);return;end % 开始计算，先赋初值 flag=＇ok＇;x=zeros(n,1);for k=1:n % 选主元 max1=0;for i=k:n if abs(A(i,k))>max1 max1=abs(A(i,k));r=i;end end if max1<1e-10 falg=＇failure＇;return;end % 交换两行 if r>k for j=k:n z=A(k,j);A(k,j)=A(r,j);A(r,j)=z;end z=b(k);b(k)=b(r);b(r)=z;end % 消元过程 b(k)=b(k)/A(k,k);for j=k+1:n A(k,j)=A(k,j)/A(k,k);end for i=1:n if i~=k for j=k+1:n A(i,j)=A(i,j)-A(i,k)*A(k,j);end b(i)=b(i)-A(i,k)*b(k);end end end % 输出x for i=1:n x(i)=b(i);end 然后保存后在命令窗口输入：

>> A=[4 2-4 0 2 4 0 0;2 2-1-2 1 3 2 0;-4-1 14 1-8-3 5 6;0-2 1 6-1-4-3 3;2 1-8-1 22 4-10-3;4 3-3-4 4 11 1-4;0 2 5-3-10 1 14 2;0 0 6 3-3-4 2 19];>> b=[0-6 20 23 9-22-15 45];>> Gau_Jor(A,b)ans = 121.1481-140.1127 29.7515-60.1528 10.9120-26.7963 5.4259-2.0185 五、实验结论 用LU法，调用matlab中的函数lu中，L往往不是一个下三角，但可以直接计算不用它的结果来计算，不用进行行变换。如果进行行变b也要变，这样会很麻烦。

实验五 解线性方程组的迭代法 一、问题提出 对实验四所列目的和意义的线性方程组，试分别选用Jacobi 迭代法，Gauss-Seidel迭代法和SOR方法计算其解。

二、要求 1、体会迭代法求解线性方程组，并能与消去法做以比较；

2、分别对不同精度要求，如由迭代次数体会该迭代法的收敛快慢；

3、对方程组2，3使用SOR方法时，选取松弛因子ω=0.8，0.9，1，1.1，1.2等，试看对算法收敛性的影响，并能找出你所选用的松弛因子的最佳者；

4、给出各种算法的设计程序和计算结果。

三、目的和意义 1、通过上机计算体会迭代法求解线性方程组的特点，并能和消去法比较；

2、运用所学的迭代法算法，解决各类线性方程组，编出算法程序；

3、体会上机计算时，终止步骤或k >（给予的迭代次数），对迭代法敛散性的意义；

4、体会初始解，松弛因子的选取，对计算结果的影响。

四、实验步骤 第一步编写实验所需的Jacobi迭代法，Gauss-Seidel迭代法，SOR迭代法的程序。

Jacobi迭代法：

function [x,k,index]=J(A,b,ep,itmax)if nargin<4 itmax=100;end if nargin<3 ep=1e-5;end n=length(A);k=0;x=zeros(n,1);y=zeros(n,1);index=1;while 1 for i=1:n y(i)=b(i);for j=1:n if j~=i y(i)=y(i)-A(i,j)*x(j);end end if abs(A(i,i))<1e-10|k==itmax index=0;return;end y(i)=y(i)/A(i,i);end if norm(y-x,inf)

function [x,k,index]=G(A,b,ep,itmax)if nargin<4 itmax=100;end if nargin<3 ep=1e-5;end n=length(A);k=0;x=zeros(n,1);y=zeros(n,1);index=1;while 1 y=x;for i=1:n z=b(i);for j=1:n if j~=i z=z-A(i,j)*x(j);end end if abs(A(i,i))<1e-10|k==itmax index=0;return;end z=z/A(i,i);x(i)=z;end if norm(y-x,inf)

function [x,k,index]=SOR(A,b,ep,w,itmax)if nargin<5 itmax=100;end if nargin<4 w=1;end if nargin<3 ep=1e-5;end n=length(A);k=0;x=zeros(n,1);y=zeros(n,1);index=1;while 1 y=x;for i=1:n z=b(i);for j=1:n if j~=i z=z-A(i,j)*x(j);end end if abs(A(i,i))<1e-10|k==itmax index=0;return;end z=z/A(i,i);x(i)=(1-w)*x(i)+w*z;end if norm(y-x,inf)

1、设线性方程组  2、设对称正定阵系数阵线方程组 3、三对角形线性方程组 五、实验结论 迭代法是解线性方程组的一个重要的实用方法，特别适用于求解在实际中大量出现的，系数矩阵为稀疏阵的大型线性方程组。通过此次实验学会了Jacobi迭代法，Gauss-Seidel迭代法，SOR迭代法的程序编写，并掌握了它们各自的优缺点及其适用条件。

实验六 非线性方程求根 一、问题提出 设方程有三个实根 现采用下面六种不同计算格式，求 f(x)=0的根或 1、2、3、4、5、6、二、要求 1、编制一个程序进行运算，最后打印出每种迭代格式的敛散情况；

2、用事后误差估计来控制迭代次数，并且打印出迭代的次数；

3、初始值的选取对迭代收敛有何影响；

4、分析迭代收敛和发散的原因。

三、目的和意义 1、通过实验进一步了解方程求根的算法；

2、认识选择计算格式的重要性；

3、掌握迭代算法和精度控制；

4、明确迭代收敛性与初值选取的关系。

四、实验步骤 第一步：编写实验所需的程序。

function [x_star,index,it]=DD(fun,x,ep,itmax)if nargin<4 itmax=100;end if nargin<3 ep=1e-5;end index=0;k=1;while k

1、，x1=0，x2=0。

2、，x1=无穷大，x2=-0.3473。

3、，x1=1.8794，x2=1.8794。

4、，x1=-0.3473，x2=-0.3473.。

5、，x1=1.8794，x2=1.8794。

6、，x1=1.8794，x2=-0.3473。

五、实验结论 对于非线性方程，求它的解析解有时候是很困难的，但采用数值方法可以很容易地求它的近似解。此次实验就是采用迭代法求非线性方程的根。对于一个非线性方程，选用不同的迭代形式，因为其收敛程度不一样，造成其效率与精确度有很大的差别。

实验七 矩阵特征值问题计算 一、问题提出 利用冪法或反冪法，求方阵的按模最大或按模最小特征值及其对应的特征向量。

设矩阵A的特征分布为：

且 试求下列矩阵之一（1）求，及 取 结果（2）求及 取 结果：

（3）求及 取 结果(4)取 这是一个收敛很慢的例子,迭代次才达到 结果（5）有一个近似特征值,试用幂法求对应的特征向量，并改进特征值（原点平移法）。

取 结果 二、要求 1、掌握冪法或反冪法求矩阵部分特征值的算法与程序设计；

2、会用原点平移法改进算法，加速收敛；

对矩阵B=A-PI取不同的P值，试求其效果；

3、试取不同的初始向量，观察对结果的影响；

4、对矩阵特征值的其它分布，如且如何计算。

三、目的和意义 1、求矩阵的部分特征值问题具有重要实际意义，如求矩阵谱半径，稳定性问题往往归于求矩阵按模最小特征值；

2、进一步掌握冪法、反冪法及原点平移加速法的程序设计技巧；

3、问题中的题（5），反应了利用原点平移的反冪法可求矩阵的任何特征值及其特征向量。

四、实验步骤 第一步：写出实验所需的幂法求最大特征值及反幂法求最小特征值的程序。

幂法程序：

function [m,u,index]=TZ(A,ep,itmax)if nargin<3 itmax=100;end if nargin<2 ep=1e-5;end n=length(A);u=ones(n,1);index=0;k=0;m1=0;while k<=itmax v=A*u;[vmax,i]=max(abs(v));m=v(i);u=v/m;if abs(m-m1)

function [m,u,index]=FTZ(A,ep,itmax)if nargin<3 itmax=100;end if nargin<2 ep=1e-5;end n=length(A);u=ones(n,1);index=0;k=0;m1=0;invA=inv(A);while k<=itmax v=invA*u;[vmax,i]=max(abs(v));m=v(i);u=v/m;if abs(m-m1)

λ3=3.4723，x3=（1.0000 0.5229 0.2422）T。，λ1=21.3053，X1=（0.8724 0.5401 0.9973 0.5644 0.4972 1.0000）T；

λ6=1.6214。

五、实验结论 求n阶方阵A的特征值和特征向量，也是实际中常常碰到的问题。通过此次实验掌握了用幂法和反幂法求一个方阵的最大特征值和特征向量，绝对值最小的特征值和特征向量。

实验八 常微分方程初值问题数值解法 一、问题提出 科学计算中经常遇到微分方程（组）初值问题，需要利用Euler法，改进Euler法，Rung-Kutta方法求其数值解，诸如以下问题：

（1）分别取h=0.1,0.2,0.4时数值解。

初值问题的精确解。

（2）用r=3的Adams显式和预-校式求解 取步长h=0.1，用四阶标准R-K方法求值。

（3）用改进Euler法或四阶标准R-K方法求解 取步长0.01,计算数值解，参考结果。

（4）利用四阶标准R-K方法求二阶方程初值问题的数值解（I）(II)(III)(IV) 二、要求 1、根据初值问题数值算法，分别选择二个初值问题编程计算；

2、试分别取不同步长，考察某节点处数值解的误差变化情况；

3、试用不同算法求解某初值问题，结果有何异常；

4、分析各个算法的优缺点。

三、目的和意义 1、熟悉各种初值问题的算法，编出算法程序；

2、明确各种算法的精度与所选步长有密切关系；

3、通过计算更加了解各种算法的优越性。

四、实验步骤 function [x,y]=euler(fun,x0,xfinal,y0,n);if nargin<5,n=50;end h=(xfinal-x0)/n;x(1)=x0;y(1)=y0;for i=1:n;x(i+1)=x(i)+h;y(i+1)=y(i)+h*feval(fun,x(i),y(i));end 实验程序及分析（Ⅰ）（1）、算法程序 function E =Euler_1(fun,x0,y0,xN,N)% Euler向前公式，其中 % fun为一阶微分方程的函数 % x0,y0为初始条件 % xN为取值范围的一个端点 % h为区间步长 % N为区间个数 % x为Xn构成的向量 % y为yn构成的向量 x=zeros(1,N+1);y=zeros(1,N+1);x(1)=x0;y(1)=y0;h=(xN-x0)/N;for n=1:N x(n+1)=x(n)+h;y(n+1)=y(n)+h*feval(fun,x(n),y(n));end T=[x＇,y＇] function z=f(x,y)z=4*x/y-x*y;(2)、运行程序 >> Euler_1(＇f＇,0,3,2,20)结果 ：

>> Euler_1(＇f＇,0,3,2,20)T = 0 3.0000 0.1000 2.9836 0.2000 2.9517 0.3000 2.9058 0.4000 2.8481 0.5000 2.7810 0.6000 2.7073 0.7000 2.6297 0.8000 2.5511 0.9000 2.4739 1.0000 2.4004 1.1000 2.3325 1.2000 2.2714 1.3000 2.2177 1.4000 2.1717 1.5000 2.1332 1.6000 2.1017 1.7000 2.0765 1.8000 2.0567 1.9000 2.0414 2.0000 2.0299 五、实验结论 很多科学技术和工程问题常用微分方程的形式建立数学模型，因此微分方程的求解是很有意义的，但对于绝大多数的微分方程问题很难或者根本不可能得到它的解析解，因此，研究微分方程的数值方法是非常有意义的。

第二篇：清华大学数值分析实验报告

数值分析实验报告

一、实验3.1

题目：

考虑线性方程组，，编制一个能自动选取主元，又能手动选取主元的求解线性代数方程组的Gauss消去过程。

（1）取矩阵，则方程有解。取计算矩阵的条件数。分别用顺序Gauss消元、列主元Gauss消元和完全选主元Gauss消元方法求解，结果如何？

（2）现选择程序中手动选取主元的功能，每步消去过程都选取模最小或按模尽可能小的元素作为主元进行消元，观察并记录计算结果，若每步消去过程总选取按模最大的元素作为主元，结果又如何？分析实验的结果。

（3）取矩阵阶数n=20或者更大，重复上述实验过程，观察记录并分析不同的问题及消去过程中选择不同的主元时计算结果的差异，说明主元素的选取在消去过程中的作用。

（4）选取其他你感兴趣的问题或者随机生成的矩阵，计算其条件数，重复上述实验，观察记录并分析实验的结果。

1.算法介绍

首先，分析各种算法消去过程的计算公式，顺序高斯消去法：

第k步消去中，设增广矩阵中的元素（若等于零则可以判定系数矩阵为奇异矩阵，停止计算），则对k行以下各行计算，分别用乘以增广矩阵的第行并加到第行，则可将增广矩阵中第列中以下的元素消为零；重复此方法，从第1步进行到第n-1步，则可以得到最终的增广矩阵，即；

列主元高斯消去法：

第k步消去中，在增广矩阵中的子方阵中，选取使得，当时，对中第行与第行交换，然后按照和顺序消去法相同的步骤进行。重复此方法，从第1步进行第n-1步，就可以得到最终的增广矩阵，即；

完全主元高斯消去法：

第k步消去中，在增广矩阵中对应的子方阵中，选取使得，若或，则对中第行与第行、第列与第列交换，然后按照和顺序消去法相同的步骤进行即可。重复此方法，从第1步进行到第n-1步，就可以得到最终的增广矩阵，即；

接下来，分析回代过程求解的公式，容易看出，对上述任一种消元法，均有以下计算公式：

2.实验程序的设计

一、输入实验要求及初始条件；

二、计算系数矩阵A的条件数及方程组的理论解；

三、对各不同方法编程计算，并输出最终计算结果。

3.计算结果及分析

（1）

先计算系数矩阵的条件数，结果如下，可知系数矩阵的条件数较大，故此问题属于病态问题，b或A的扰动都可能引起解的较大误差；

采用顺序高斯消去法，计算结果为：

最终解为x=(1.***,1.***,1.***,1.***,0.***,1.***,0.***,1.***,0.***,1.***)T

使用无穷范数衡量误差，得到=2.842***1e-14，可以发现，采用顺序高斯消元法求得的解与精确解之间误差较小。通过进一步观察，可以发现，按照顺序高斯消去法计算时，其选取的主元值和矩阵中其他元素大小相近，因此顺序高斯消去法方式并没有对结果造成特别大的影响。

若采用列主元高斯消元法，则结果为：

最终解为x=(1.***,1.***,1.***,1.***,1.***,1.***,1.***,1.***,1.***,1.***)T

同样使用无穷范数衡量误差，有=0；

若使用完全主元高斯消元法，则结果为

最终解x=(1.***,1.***,1.***,1.***,1.***,1.***,1.***,1.***,1.***,1.***)T

同样使用无穷范数衡量误差，有=0；

（2）

若每步都选取模最小或尽可能小的元素为主元，则计算结果为

最终解x=(1.***

1.***

1.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***)T

使用无穷范数衡量误差，有为2.842***1e-14；而完全主元消去法的误差为=0。

从（1）和（2）的实验结果可以发现，列主元消去法和完全主元消去法都得到了精确解，而顺序高斯消去法和以模尽量小的元素为主元的消去法没有得到精确解。在后两种消去法中，由于程序计算时的舍入误差，对最终结果产生了一定的影响，但由于方程组的维度较低，并且元素之间相差不大，所以误差仍比较小。

为进一步分析，计算上述4种方法每步选取的主元数值，并列表进行比较，结果如下：

第n次消元

顺序

列主元

完全主元

模最小

6.***

6.***

4.***

4.***

4.***

4.***

4.***3333

4.***3333

4.***

4.***

4.***

4.***

4.0***063

4.0***063

4.***

4.***

4.0039***

4.0039***

4.***

0.0***469

0.0***469

4.***

从上表可以发现，对这个方程组而言，顺序高斯消去选取的主元恰好事模尽量小的元素，而由于列主元和完全主元选取的元素为8，与4在数量级上差别小，所以计算过程中的累积误差也较小，最终4种方法的输出结果均较为精确。

在这里，具体解释一下顺序法与模最小法的计算结果完全一致的原因。该矩阵在消元过程中，每次选取主元的一列只有两个非零元素，对角线上的元素为4左右，而其正下方的元素为8，该列其余位置的元素均为0。在这样的情况下，默认的主元也就是该列最小的主元，因此两种方法所得到的计算结果是一致的。

理论上说，完全高斯消去法的误差最小，其次是列主元高斯消去法，而选取模最小的元素作为主元时的误差最大，但是由于方程组的特殊性（元素相差不大并且维度不高），这个理论现象在这里并没有充分体现出来。

（3）

时，重复上述实验过程，各种方法的计算结果如下所示，在这里，仍采用无穷范数衡量绝对误差。

顺序高斯消去法

列主元高斯消去

完全主元高斯消去

选取模最小或尽可能小元素作为主元消去

X

1.***

1.***

1.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

1.***

1.***

1.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

2.***e-11

0

0

2.***e-11

可以看出，此时列主元和完全主元的计算结果仍为精确值，而顺序高斯消去和模尽可能小方法仍然产生了一定的误差，并且两者的误差一致。与n=10时候的误差比相比，n=20时的误差增长了大约1000倍，这是由于计算过程中舍入误差的不断累积所致。所以，如果进一步增加矩阵的维数，应该可以看出更明显的现象。

（4）

不同矩阵维度下的误差如下，在这里，为方便起见，选取2-条件数对不同维度的系数矩阵进行比较。

维度

条件数

顺序消去

列主元

完全主元

模尽量小

1.7e+3

2.84e-14

0

0

2.84e-14

1.8e+6

2.91e-11

0

0

2.91e-11

5.7e+7

9.31e-10

0

0

9.31e-10

1.8e+9

2.98e-08

0

0

2.98e-08

1.9e+12

3.05e-05

0

0

3.05e-05

3.8e+16

3.28e+04

3.88e-12

3.88e-12

3.28e+04

8.5e+16

3.52e+13

4.2e-3

4.2e-3

3.52e+13

从上表可以看出，随着维度的增加，不同方法对计算误差的影响逐渐体现，并且增长较快，这是由于舍入误差逐步累计而造成的。不过，方法二与方法三在维度小于40的情况下都得到了精确解，这两种方法的累计误差远比方法一和方法四慢；同样地，出于与前面相同的原因，方法一与方法四的计算结果保持一致，方法二与方法三的计算结果保持一致。

4.结论

本文矩阵中的元素差别不大，模最大和模最小的元素并没有数量级上的差异，因此，不同的主元选取方式对计算结果的影响在维度较低的情况下并不明显，四种方法都足够精确。

对比四种方法，可以发现采用列主元高斯消去或者完全主元高斯消去法，可以尽量抑制误差，算法最为精确。不过，对于低阶的矩阵来说，四种方法求解出来的结果误差均较小。

另外，由于完全选主元方法在选主元的过程中计算量较大，而且可以发现列主元法已经可以达到很高的精确程度，因而在实际计算中可以选用列主元法进行计算。

附录：程序代码

clear

clc;

format

long;

%方法选择

n=input(＇矩阵A阶数:n=＇);

disp(＇选取求解方式＇);

disp(＇1

顺序Gauss消元法，2

列主元Gauss消元法，3

完全选主元Gauss消元法，4

模最小或近可能小的元素作为主元＇);

a=input(＇求解方式序号：＇);

%赋值A和b

A=zeros(n,n);

b=zeros(n,1);

for

i=1:n

A(i,i)=6;

if

i>1

A(i,i-1)=8;

end

if

i

A(i,i+1)=1;

end

end

for

i=1:n

for

j=1:n

b(i)=b(i)+A(i,j);

end

end

disp(＇给定系数矩阵为：＇);

A

disp(＇右端向量为：＇);

b

%求条件数及理论解

disp(＇线性方程组的精确解：＇);

X=(A＼b)＇

fprintf(＇矩阵A的1-条件数:

%f

＼n＇,cond(A,1));

fprintf(＇矩阵A的2-条件数:

%f

＼n＇,cond(A));

fprintf(＇矩阵A的无穷-条件数:

%f

＼n＇,cond(A,inf));

%顺序Gauss消元法

if

a==1

A1=A;b1=b;

for

k=1:n

if

A1(k,k)==0

disp(＇主元为零，顺序Gauss消元法无法进行＇);

break

end

fprintf(＇第%d次消元所选取的主元：%g＼n＇,k,A1(k,k))

%disp(＇此次消元后系数矩阵为：＇);

%A1

for

p=k+1:n

l=A1(p,k)/A1(k,k);

A1(p,k:n)=A1(p,k:n)-l*A1(k,k:n);

b1(p)=b1(p)-l*b1(k);

end

end

x1(n)=b1(n)/A1(n,n);

for

k=n-1:-1:1

for

w=k+1:n

b1(k)=b1(k)-A1(k,w)*x1(w);

end

x1(k)=b1(k)/A1(k,k);

end

disp(＇顺序Gauss消元法解为:＇);

disp(x1);

disp(＇所求解与精确解之差的无穷-范数为＇);

norm(x1-X,inf)

end

%列主元Gauss消元法

if

a==2

A2=A;b2=b;

for

k=1:n

[max_i,max_j]=find(A2(:,k)==max(abs(A2(k:n,k))));

if

max_i~=k

A2_change=A2(k,:);

A2(k,:)=A2(max_i,:);

A2(max_i,:)=A2_change;

b2_change=b2(k);

b2(k)=b2(max_i);

b2(max_i)=b2_change;

end

if

A2(k,k)==0

disp(＇主元为零，列主元Gauss消元法无法进行＇);

break

end

fprintf(＇第%d次消元所选取的主元：%g＼n＇,k,A2(k,k))

%disp(＇此次消元后系数矩阵为：＇);

%A2

for

p=k+1:n

l=A2(p,k)/A2(k,k);

A2(p,k:n)=A2(p,k:n)-l*A2(k,k:n);

b2(p)=b2(p)-l*b2(k);

end

end

x2(n)=b2(n)/A2(n,n);

for

k=n-1:-1:1

for

w=k+1:n

b2(k)=b2(k)-A2(k,w)*x2(w);

end

x2(k)=b2(k)/A2(k,k);

end

disp(＇列主元Gauss消元法解为:＇);

disp(x2);

disp(＇所求解与精确解之差的无穷-范数为＇);

norm(x2-X,inf)

end

%完全选主元Gauss消元法

if

a==3

A3=A;b3=b;

for

k=1:n

VV=eye(n);

[max_i,max_j]=find(A3(k:n,k:n)==max(max(abs(A3(k:n,k:n)))));

if

numel(max_i)==0

[max_i,max_j]=find(A3(k:n,k:n)==-max(max(abs(A3(k:n,k:n)))));

end

W=eye(n);

W(max_i(1)+k-1,max_i(1)+k-1)=0;

W(k,k)=0;

W(max_i(1)+k-1,k)=1;

W(k,max_i(1)+k-1)=1;

V=eye(n);

V(k,k)=0;

V(max_j(1)+k-1,max_j(1)+k-1)=0;

V(k,max_j(1)+k-1)=1;

V(max_j(1)+k-1,k)=1;

A3=W*A3*V;

b3=W*b3;

VV=VV*V;

if

A3(k,k)==0

disp(＇主元为零，完全选主元Gauss消元法无法进行＇);

break

end

fprintf(＇第%d次消元所选取的主元：%g＼n＇,k,A3(k,k))

%disp(＇此次消元后系数矩阵为：＇);

%A3

for

p=k+1:n

l=A3(p,k)/A3(k,k);

A3(p,k:n)=A3(p,k:n)-l*A3(k,k:n);

b3(p)=b3(p)-l*b3(k);

end

end

x3(n)=b3(n)/A3(n,n);

for

k=n-1:-1:1

for

w=k+1:n

b3(k)=b3(k)-A3(k,w)*x3(w);

end

x3(k)=b3(k)/A3(k,k);

end

disp(＇完全选主元Gauss消元法解为:＇);

disp(x3);

disp(＇所求解与精确解之差的无穷-范数为＇);

norm(x3-X,inf)

end

%模最小或近可能小的元素作为主元

if

a==4

A4=A;b4=b;

for

k=1:n

AA=A4;

AA(AA==0)=NaN;

[min_i,j]=find(AA(k:n,k)==min(abs(AA(k:n,k))));

if

numel(min_i)==0

[min_i,j]=find(AA(k:n,k)==-min(abs(AA(k:n,k:n))));

end

W=eye(n);

W(min_i(1)+k-1,min_i(1)+k-1)=0;

W(k,k)=0;

W(min_i(1)+k-1,k)=1;

W(k,min_i(1)+k-1)=1;

A4=W*A4;

b4=W*b4;

if

A4(k,k)==0

disp(＇主元为零，模最小Gauss消元法无法进行＇);

break

end

fprintf(＇第%d次消元所选取的主元：%g＼n＇,k,A4(k,k))

%A4

for

p=k+1:n

l=A4(p,k)/A4(k,k);

A4(p,k:n)=A4(p,k:n)-l*A4(k,k:n);

b4(p)=b4(p)-l*b4(k);

end

end

x4(n)=b4(n)/A4(n,n);

for

k=n-1:-1:1

for

w=k+1:n

b4(k)=b4(k)-A4(k,w)*x4(w);

end

x4(k)=b4(k)/A4(k,k);

end

disp(＇模最小Gauss消元法解为:＇);

disp(x4);

disp(＇所求解与精确解之差的无穷-范数为＇);

norm(x4-X,inf)

end

二、实验3.3

题目：

考虑方程组的解，其中系数矩阵H为Hilbert矩阵：

这是一个著名的病态问题。通过首先给定解（例如取为各个分量均为1）再计算出右端的办法给出确定的问题。

（1）选择问题的维数为6，分别用Gauss消去法（即LU分解）、J迭代法、GS迭代法和SOR迭代法求解方程组，其各自的结果如何？将计算结果与问题的解比较，结论如何。

（2）逐步增大问题的维数，仍用上述的方法来解它们，计算的结果如何？计算的结果说明的什么？

（3）讨论病态问题求解的算法。

1.算法设计

对任意线性方程组，分析各种方法的计算公式如下，（1）Gauss消去法：

首先对系数矩阵进行LU分解，有，则原方程转化为，令，则原方程可以分为两步回代求解：

具体方法这里不再赘述。

（2）J迭代法：

首先分解，再构造迭代矩阵，其中，进行迭代计算，直到误差满足要求。

（3）GS迭代法：

首先分解，再构造迭代矩阵，其中，进行迭代计算，直到误差满足要求。

（4）SOR迭代法：

首先分解，再构造迭代矩阵，其中，进行迭代计算，直到误差满足要求。

2.实验过程

一、根据维度n确定矩阵H的各个元素和b的各个分量值；

二、选择计算方法（Gauss消去法，J迭代法，GS迭代法，SOR迭代法），对迭代法设定初值，此外SOR方法还需要设定松弛因子；

三、进行计算，直至满足误差要求（对迭代法，设定相邻两次迭代结果之差的无穷范数小于0.0001；

3.计算结果及分析

（1）时，问题可以具体定义为

计算结果如下，Gauss消去法

第1次消元所选取的主元是：1

第2次消元所选取的主元是：0.0833333

第3次消元所选取的主元是：0.00555556

第4次消元所选取的主元是：0.000357143

第5次消元所选取的主元是：2.26757e-05

第6次消元所选取的主元是：1.43155e-06

解得X=(0.***

1.***

0.***

1.***

0.***

1.***)T

使用无穷范数衡量误差，可得=4.***e-10；

J迭代法

设定迭代初值为零，计算得到

J法的迭代矩阵B的谱半径为4.30853＞1，所以J法不收敛；

GS迭代法

设定迭代初值为零，计算得到GS法的迭代矩阵G的谱半径为：0.999998＜1，故GS法收敛，经过541次迭代计算后，结果为X=(1.001***6

0.***

0.***

1.***

1.***

0.***)T

使用无穷范数衡量误差，有=0.***；

SOR迭代法

设定迭代初值为零向量，并设定，计算得到SOR法迭代矩阵谱半径为0.***，经过100次迭代后的计算结果为

X=(1.***

0.***

1.03***59

1.06***81

1.***

0.9***527)T；

使用无穷范数衡量误差，有=0.***；

对SOR方法，可变，改变值，计算结果可以列表如下

迭代次数

迭代矩阵的谱半径

0.***

0.***

0.***

0.***

X

1.***

0.***

1.01***40

1.***

1.0***681

0.***

1.***

0.***

1.***

1.***

1.***

0.***

1.***

0.***

1.***

1.***

0.***

0.***

1.05***66

0.***

1.***

0.***

1.***

0.***

0.***

0.***

0.***

0.***

可以发现，松弛因子的取值对迭代速度造成了不同的影响，上述四种方法中，松弛因子=0.5时，收敛相对较快。

综上，四种算法的结果列表如下：

算法

Gauss消去法

Jacobi法

GS法

SOR法（取）

迭代次数

--

不收敛

541

迭代矩阵的谱半径

--

4.30853

0.999998

0.***

X

0.***

1.***

0.***

1.***

0.***

1.***

--

1.001***6

0.***

0.***

1.***

1.***

0.***

1.***

0.***

1.03***59

1.06***81

1.***

0.9***527

4.***e-10

--

0.***

0.***

计算可得，矩阵H的条件数为>>1，所以这是一个病态问题。由上表可以看出，四种方法的求解都存在一定的误差。下面分析误差的来源：

LU分解方法的误差存在主要是由于Hilbert矩阵各元素由分数形式转换为小数形式时，不能除尽情况下会出现舍入误差，在进行LU分解时也存在这个问题，所以最后得到的结果不是方程的精确解，但结果显示该方法的误差非常小；

Jacobi迭代矩阵的谱半径为4.30853，故此迭代法不收敛；

GS迭代法在迭代次数为541次时得到了方程的近似解，其误差约为0.05，比较大。GS迭代矩阵的谱半径为0.999998，很接近1，所以GS迭代法收敛速度较慢；

SOR迭代法在迭代次数为100次时误差约为0.08，误差较大。SOR迭代矩阵的谱半径为0.999999，也很接近1，所以时SOR迭代法收敛速度不是很快，但是相比于GS法，在迭代速度方面已经有了明显的提高；另外，对不同的，SOR方法的迭代速度会相应有变化，如果选用最佳松弛因子，可以实现更快的收敛；

（2）

考虑不同维度的情况，时，算法

Gauss消去

J法

GS法

SOR法（w=0.5）

计算结果

0.***

1.***

0.***

1.***

0.***

1.***

0.***

1.***

--

0.***

1.***

0.***

1.***

1.***

1.***

0.9968***

0.***

1.***

0.9397***

0.***

1.***

1.***

1.***

0.***

0.***

迭代次数

--

--

356

谱半径

--

6.04213

0.***

--

时，算法

Gauss消去法

Jacobi法

GS法

SOR法（w=0.5）

计算结果

0.***

1.***

0.***

1.000***1

0.***

1.***

0.***

1.***

0.***

1.***

0.***

--

0.***

1.***

0.***

0.***

0.***

1.02***91

1.***

1.***

1.***

0.***

0.947***7

1.0***572

0.***

0.***

0.***

1.***

1.***

1.***

1.***

0.***

0.***

0.***

迭代次数

--

--

1019

谱半径

--

8.64964

0.***

--

时

算法

Gauss消去法

Jacobi法

GS法

SOR法（w=0.5）

计算结果

0.***

1.***

0.***

0.***

1.***

0.***

2.***

-2.***

7.***

-7.***

7.***

-1.***

0.***

1.***

0.***

--

不收敛

1.***

1.***

0.907***9

0.***

0.***

1.***

1.09***64

1.***

1.***

1.***

1.0385***

0.***

0.942***3

0.***

0.***

迭代次数

--

--

262

谱半径

--

6.04213

＞1

1.***

8.***

--

--

0.***

分析以上结果可以发现，随着n值的增加，Gauss消去法误差逐渐增大，而且误差增大的速度很快，在维数小于等于10情况下，Gauss消去法得到的结果误差较小；但当维数达到15时，计算结果误差已经达到精确解的很多倍；

J法迭代不收敛，无论n如何取值，其谱半径始终大于1，因而J法不收敛，所以J迭代法不能用于Hilbert矩阵的求解；

对于GS迭代法和SOR迭代法，两种方法均收敛，GS迭代法是SOR迭代法松弛因子取值为1的特例，SOR方法受到取值的影响，会有不同的收敛情况。可以得出GS迭代矩阵的谱半径小于1但是很接近1，收敛速度很慢。虽然随着维数的增大，所需迭代的次数逐渐减少，但是当维数达到15的时候，GS法已经不再收敛。因此可以得出结论，GS迭代方法在Hilbert矩阵维数较低时，能够在一定程度上满足迭代求解的需求，不过迭代的速度很慢。另外，随着矩阵维数的增加，SOR法的误差水平基本稳定，而且误差在可以接受的范围之内。

经过比较可以得出结论，如果求解较低维度的Hibert矩阵问题，Gauss消去法、GS迭代法和SOR迭代法均可使用，且Gauss消去法的结果精确度较高；如果需要求解较高维度的Hibert矩阵问题，只有采用SOR迭代法。

（3）

系数矩阵的条件数较大时，为病态方程。由实验可知，Gauss法在解上述方程时，结果存在很大的误差。而对于收敛的迭代法，可以通过选取最优松弛因子的方法来求解，虽然迭代次数相对较多，但是结果较为精确。

总体来看，对于一般病态方程组的求解，可以采用以下方式：

1.低维度下采用Gauss消去法直接求解是可行的；

Jacobi迭代方法不适宜于求解病态问题；

GS迭代方法可以解决维数较低的病态问题，但其谱半径非常趋近于1，导致迭代算法收敛速度很慢，维数较大的时候，GS法也不再收敛；

SOR方法较适合于求解病态问题，特别是矩阵维数较高的时候，其优势更为明显。

2.采用高精度的运算，如选用双倍或更多倍字长的运算，可以提高收敛速度；

3.可以对原方程组作某些预处理，从而有效降低系数矩阵的条件数。

4.实验结论

（1）对Hibert矩阵问题，其条件数会随着维度的增加迅速增加，病态性会越来越明显；在维度较低的时候，Gauss消去法、GS迭代法和SOR迭代法均可使用，且可以优先使用Gauss消去法；如果需要求解较高维度的Hibert矩阵问题，只有SOR迭代法能够求解。

（2）SOR方法比较适合于求解病态问题，特别是矩阵维数较高的时候，其优点更为明显。从本次实验可以看出，随着矩阵维数的增大，SOR方法所需的迭代次数减少，而且误差基本稳定，是解决病态问题的适宜方法。

附录：程序代码

clear

all

clc;

format

long;

%矩阵赋值

n=input(＇矩阵H的阶数:n=＇);

for

i=1:n

for

j=1:n

H(i,j)=1/(i+j-1);

end

end

b=H*ones(n,1);

disp(＇H矩阵为：＇);

H

disp(＇向量b：＇);

b

%方法选择

disp(＇选取求解方式＇);

disp(＇1

Gauss消去法，2

J迭代法，3

GS迭代法，4

SOR迭代法＇);

a=input(＇求解方式序号：＇);

%Gauss消去法

if

a==1;

H1=H;b1=b;

for

k=1:n

if

H1(k,k)==0

disp(＇主元为零，Gauss消去法无法进行＇);

break

end

fprintf(＇第%d次消元所选取的主元是：%g＼n＇,k,H1(k,k))

for

p=k+1:n

m5=-H1(p,k)/H1(k,k);

H1(p,k:n)=H1(p,k:n)+m5*H1(k,k:n);

b1(p)=b1(p)+m5*b1(k);

end

end

x1(n)=b1(n)/H1(n,n);

for

k=n-1:-1:1

for

v=k+1:n

b1(k)=b1(k)-H1(k,v)*x1(v);

end

x1(k)=b1(k)/H1(k,k);

end

disp(＇Gauss消去法解为:＇);

disp(x1);

disp(＇解与精确解之差的无穷范数＇);

norm((x1-a),inf)

end

D=diag(diag(H));

L=-tril(H,-1);

U=-triu(H,1);

%J迭代法

if

a==2;

%给定初始x0

ini=input(＇初始值设定:x0=＇);

x0(:,1)=ini*diag(ones(n));

disp(＇初始解向量为:＇);

x0

xj(:,1)=x0(:,1);

B=(D^(-1))*(L+U);

f=(D^(-1))*b;

fprintf(＇(J法B矩阵谱半径为：%g＼n＇,vrho(B));

if

vrho(B)<1;

for

m2=1:5000

xj(:,m2+1)=B*xj(:,m2)+fj;

if

norm((xj(:,m2+1)-xj(:,m2)),inf)<0.0001

break

end

end

disp(＇J法计算结果为:＇);

xj(:,m2+1)

disp(＇解与精确解之差的无穷范数＇);

norm((xj(:,m2+1)-diag(ones(n))),inf)

disp(＇J迭代法迭代次数:＇);

m2

else

disp(＇由于B矩阵谱半径大于1，因而J法不收敛＇);

end

end

%GS迭代法

if

a==3;

%给定初始x0

ini=input(＇初始值设定:x0=＇);

x0(:,1)=ini*diag(ones(n));

disp(＇初始解向量为:＇);

x0

xG(:,1)=x0(:,1);

G=inv(D-L)*U;

fG=inv(D-L)*b;

fprintf(＇GS法G矩阵谱半径为：%g＼n＇,vrho(G));

if

vrho(G)<1

for

m3=1:5000

xG(:,m3+1)=G*xG(:,m3)+fG;

if

norm((xG(:,m3+1)-xG(:,m3)),inf)<0.0001

break;

end

end

disp(＇GS迭代法计算结果:＇);

xG(:,m3+1)

disp(＇解与精确解之差的无穷范数＇);

norm((xG(:,m3+1)-diag(ones(n))),inf)

disp(＇GS迭代法迭代次数:＇);

m3

else

disp(＇由于G矩阵谱半径大于1，因而GS法不收敛＇);

end

end

%SOR迭代法

if

a==4;

%给定初始x0

ini=input(＇初始值设定:x0=＇);

x0(:,1)=ini*diag(ones(n));

disp(＇初始解向量为:＇);

x0

A=H;

for

i=1:n

b(i)=sum(A(i,:));

end

x_star=ones(n,1);

format

long

w=input(＇松弛因子:w=＇);

Lw=inv(D-w*L)*((1-w)*D+w*U);

f=w*inv(D-w*L)*b;

disp(＇迭代矩阵的谱半径：＇)

p=vrho(Lw)

time_max=100;%迭代次数

x=zeros(n,1);%迭代初值

for

i=1:time_max

x=Lw*x+f;

end

disp(＇SOR迭代法得到的解为＇);

x

disp(＇解与精确解之差的无穷范数＇);

norm((x_star-x),inf)

end

pause

三、实验4.1

题目：

对牛顿法和拟牛顿法。进行非线性方程组的数值求解

（1）用上述两种方法，分别计算下面的两个例子。在达到精度相同的前提下，比较其迭代次数、CPU时间等。

（2）取其他初值，结果又如何？反复选取不同的初值，比较其结果。

（3）总结归纳你的实验结果，试说明各种方法适用的问题。

1.算法设计

对需要求解的非线性方程组而言，牛顿法和拟牛顿法的迭代公式如下，（1）牛顿法：

牛顿法为单步迭代法，需要取一个初值。

（2）拟牛顿法：（Broyden秩1法）

其中，拟牛顿法不需要求解的导数，因此节省了大量的运算时间，但需要给定矩阵的初值，取为。

2.实验过程

一、输入初值；

二、根据误差要求，按公式进行迭代计算；

三、输出数据；

3.计算结果及分析

（1）首先求解方程组（1），在这里，设定精度要求为，方法

牛顿法

拟牛顿法

初始值

计算结果X

x1

0.***

0.***

x2

1.***

1.0852***

x3

0.***

0.***

迭代次数

CPU计算时间/s

3.777815

2.739349

可以看出，在初始值相同情况下，牛顿法和拟牛顿法在达到同样计算精度情况下得到的结果基本相同，但牛顿法的迭代次数明显要少一些，但是，由于每次迭代都需要求解矩阵的逆，所以牛顿法每次迭代的CPU计算时间更长。

之后求解方程组（2），同样设定精度要求为

方法

牛顿法

拟牛顿法

初始值

计算结果X

x1

0.***

0.***

x2

0.***

0.***

x3

-0.***

-0.***

迭代次数

CPU计算时间/s

2.722437

3.920195

同样地，可以看出，在初始值相同情况下，牛顿法和拟牛顿法在达到同样计算精度情况下得到的结果是基本相同的，但牛顿法的迭代次数明显要少，但同样的，由于每次迭代中有求解矩阵的逆的运算，牛顿法每次迭代的CPU计算时间较长。

（2）对方程组（1），取其他初值，计算结果列表如下，同样设定精度要求为

初始值

方法

牛顿法

拟牛顿法

计算结果

0.***

1.***

0.***

9.21***94

-5.***

18.1***205

迭代次数

CPU计算时间/s

3.907164

4.818019

计算结果

0.***

1.***

0.***

9.21***91

-5.***

18.1***807

迭代次数

2735

CPU计算时间/s

8.127286

5.626023

计算结果

0.***

1.***

0.***

0.***

1.0852***

0.***

迭代次数

CPU计算时间/s

3.777815

2.739349

计算结果

0.***

1.***

0.***

0.***

1.***

0.***

迭代次数

188

CPU计算时间/s

3.835697

2.879070

计算结果

9.21***22

-5.***

18.1***605

Matlab警告矩阵接近奇异值，程序进入长期循环计算中

迭代次数

--

CPU计算时间/s

4.033868

--

计算结果

0.***

1.***

0.***

Matlab警告矩阵接近奇异值，程序进入长期循环计算中

迭代次数

--

CPU计算时间/s

12.243263

--

从上表可以发现，方程组（1）存在另一个在(9.2,-5.6,18.1)T附近的不动点，初值的选取会直接影响到牛顿法和拟牛顿法最后的收敛点。

总的来说，设定的初值离不动点越远，需要的迭代次数越多，因而初始值的选取非常重要，合适的初值可以更快地收敛，如果初始值偏离精确解较远，会出现迭代次数增加直至无法收敛的情况；

由于拟牛顿法是一种近似方法，拟牛顿法需要的的迭代次数明显更多，而且收敛情况不如牛顿法好（初值不够接近时，甚至会出现奇异矩阵的情况），但由于牛顿法的求解比较复杂，计算时间较长；

同样的，对方程组（2），取其他初值，计算结果列表如下，同样设定精度要求为

初始值

方法

牛顿法

拟牛顿法

计算结果

0.***

0.***

-0.***

0.***

0.***

-0.***

迭代次数

CPU计算时间/s

2.722437

3.920195

计算结果

0.***

0.***

-0.***

0.***

-0.***

76.***

迭代次数

CPU计算时间/s

5.047111

5.619752

计算结果

0.***

0.***

-0.***

1.0e+02

*

-0.***

-0.000***6

1.754***3

迭代次数

CPU计算时间/s

3.540668

3.387829

计算结果

0.***

0.***

-0.***

1.0e+04

*

0.***

-0.***

1.***

迭代次数

CPU计算时间/s

2.200571

2.640901

计算结果

0.***

0.***

-0.***

矩阵为奇异值，无法输出准确结果

迭代次数

--

CPU计算时间/s

1.719072

--

计算结果

0.***

0.***

-0.***

矩阵为奇异值，无法输出准确结果

迭代次数

149

--

CPU计算时间/s

2.797116

--

计算结果

矩阵为奇异值，无法输出准确结果

矩阵为奇异值，无法输出准确结果

迭代次数

--

--

CPU计算时间/s

--

--

在这里，与前文类似的发现不再赘述。

从这里看出，牛顿法可以在更大的区间上实现压缩映射原理，可以在更大的范围上选取初值并最终收敛到精确解附近；

在初始值较接近于不动点时，牛顿法和拟牛顿法计算所得到的结果是基本相同的，虽然迭代次数有所差别，但计算总的所需时间相近。

（3）

牛顿法在迭代过程中用到了矩阵的求逆，其迭代收敛的充分条件是迭代满足区间上的映内性，对于矩阵的求逆过程比较简单，所以在较大区间内满足映内性的问题适合应用牛顿法进行计算。一般而言，对于函数单调或者具有单值特性的函数适合应用牛顿法，其对初始值敏感程度较低，算法具有很好的收敛性。

另外，需要说明的是，每次计算给出的CPU时间与计算机当时的运行状态有关，同时，不同代码的运行时间也不一定一致，所以这个数据并不具有很大的参考价值。

4.实验结论

对牛顿法和拟牛顿法，都存在初始值越接近精确解，所需的迭代次数越小的现象；

在应用上，牛顿法和拟牛顿法各有优势。就迭代次数来说，牛顿法由于更加精确，所需的迭代次数更少；但就单次迭代来说，牛顿法由于计算步骤更多，且计算更加复杂，因而每次迭代所需的时间更长，而拟牛顿法由于采用了简化的近似公式，其每次迭代更加迅速。当非线性方程组求逆过程比较简单时，如方程组1的情况时，拟牛顿法不具有明显的优势；而当非线性方程组求逆过程比较复杂时，如方程组2的情况，拟牛顿法就可以体现出优势，虽然循环次数有所增加，但是CPU耗时反而更少。

另外，就方程组压缩映射区间来说，一般而言，对于在区间内函数呈现单调或者具有单值特性的函数适合应用牛顿法，其对初始值敏感程度较低，使算法具有很好的收敛性；而拟牛顿法由于不需要在迭代过程中对矩阵求逆，而是利用差商替代了对矩阵的求导，所以即使初始误差较大时，其倒数矩阵与差商偏差也较小，所以对初始值的敏感程度较小。

附录：程序代码

%方程1，牛顿法

tic;

format

long;

%%初值

disp(＇请输入初值＇);

a=input(＇第1个分量为：＇);

b=input(＇第2个分量为：＇);

c=input(＇第3个分量为：＇);

disp(＇所选定初值为＇);

x=[a;b;c]

%%误差要求

E=0.0001;

%%迭代

i=0;

e=2*E;

while

e>E

F=[12*x(1)-x(2)^2-4*x(3)-7;x(1)^2+10*x(2)-x(3)-11;x(2)^3+10*x(3)-8];

f=[12,-2*x(2),-4;2*x(1),10,-1;0,3*x(2)^2,10];

det_x=((f)^(-1))*(-F);

x=x+det_x;

e=max(norm(det_x));

i=i+1;

end

disp(＇迭代次数＇);

i

disp(＇迭代次数＇);

x

toc;

%方程1，拟牛顿法

tic;

format

long;

%%初值

%%初值

disp(＇请输入初值＇);

a=input(＇第1个分量为：＇);

b=input(＇第2个分量为：＇);

c=input(＇第3个分量为：＇);

disp(＇所选定初值为＇);

x0=[a;b;c]

%%误差要求

E=0.0001;

%%迭代

i=0;

e=2*E;

A0=eye(3);

while

e>E

F0=[12*x0(1)-x0(2)^2-4*x0(3)-7;x0(1)^2+10*x0(2)-x0(3)-11;x0(2)^3+10*x0(3)-8];

x1=x0-A0^(-1)*F0;

s=x1-x0;

F1=[12*x1(1)-x1(2)^2-4*x1(3)-7;x1(1)^2+10*x1(2)-x1(3)-11;x1(2)^3+10*x1(3)-8];

y=F1-F0;

A1=A0+(y-A0*s)*s＇/(s＇*s);

x0=x1;

A0=A1;

e=max(norm(s));

i=i+1;

end

disp(＇迭代次数＇);

i

disp(＇迭代次数＇);

x0

toc;

%方程2，牛顿法

tic;

format

long;

%%初值

disp(＇请输入初值＇);

a=input(＇第1个分量为：＇);

b=input(＇第2个分量为：＇);

c=input(＇第3个分量为：＇);

disp(＇所选定初值为＇);

x=[a;b;c]

%%误差要求

E=0.0001;

%%迭代

i=0;

e=2*E;

while

e>E

F=[3*x(1)-cos(x(2)*x(3))-0.5;x(1)^2-81*(x(2)+0.1)^2+sin(x(3))+1.06;exp(1)^(-x(1)*x(2))+20*x(3)+(10*pi-3)/3];

f=[3,x(3)*sin(x(2)*x(3)),x(2)*sin(x(2)*x(3));2*x(1),-162*x(2)-81/5,cos(x(3));-x(2)*exp(1)^(-x(1)*x(2)),-x(1)*exp(1)^(-x(1)*x(2)),20];

det_x=((f)^(-1))*(-F);

x=x+det_x;

e=max(norm(det_x));

i=i+1;

end

disp(＇迭代次数＇);

i

disp(＇迭代次数＇);

x

toc;

%方程2，拟牛顿法

tic;

format

long;

%%初值

%%初值

disp(＇请输入初值＇);

a=input(＇第1个分量为：＇);

b=input(＇第2个分量为：＇);

c=input(＇第3个分量为：＇);

disp(＇所选定初值为＇);

x0=[a;b;c]

%%误差要求

E=0.0001;

%%迭代

i=0;

e=2*E;

A0=eye(3);

while

e>E

F0=[3*x0(1)-cos(x0(2)*x0(3))-0.5;x0(1)^2-81*(x0(2)+0.1)^2+sin(x0(3))+1.06;exp(1)^(-x0(1)*x0(2))+20*x0(3)+(10*pi-3)/3];

x1=x0-A0^(-1)*F0;

s=x1-x0;

F1=[3*x1(1)-cos(x1(2)*x1(3))-0.5;x1(1)^2-81*(x1(2)+0.1)^2+sin(x1(3))+1.06;exp(1)^(-x1(1)*x1(2))+20*x1(3)+(10*pi-3)/3];

y=F1-F0;

A1=A0+(y-A0*s)*s＇/(s＇*s);

x0=x1;

A0=A1;

e=max(norm(s));

i=i+1;

end

disp(＇迭代次数＇);

i

disp(＇迭代次数＇);

x0

toc;

第三篇：数值分析实验报告写作范本

学号：

姓名：

实验二插值法

实验2.1（多项式插值的振荡现象）

问题提出：考虑一个固定的区间上用插值逼近一个函数。显然拉格朗日插值中使用的节点越多，插值多项式的次数就越高。我们自然关心插值多项式的次数增加时，L(x)是否也更加靠近被逼近的函数。龙格给出了一个极著名例子。设区间[-1，1]上函数

f(x)=1／(1+25x^2)

实验内容：考虑区间[-1，1]的一个等距划分，分点为：

x(i)=-1+2i/n,i=0,1,2…,n

则拉格朗日插值多项式为：

L(x)=∑l(i)(x)/(1+25x(j)^2)i=0,1,…n

其中l(i)(x), i=0,1,…n,n是n次拉格朗日插值基函数。

实验要求：

⑴ 选择不断增大的分点数目n=2,3…，画出f(x)及插值多项式函数L(x)在[-1，1]上的图象，比较分析实验结果。

（2）选择其它的函数，例如定义在区间[-5，5]上的函数

h(x)=x/(1+x^4),g(x)=arctanx

重复上述的实验看其结果如何。

（3）区间[a,b]上切比雪夫点的定义为：

xk=（b+a）/2+((b-a)/2)cos((2k-1)π/(2(n+1))),k=1,2,^,n+1

以x1,x2^x(n+1)为插值节点构造上述各函数的拉格朗日插值多项式，比较其结果。

实验过程：

程序：。。。。。。。。。。。。。。。。

数值实验结果及分析：。。。。。。。。。。。。。。

讨论。。。。。。。。。。。。。。。。。。。

实验总结：。。。。。。。。。。。。。。。

第四篇：“数值分析”课程教学改革浅谈

摘要：“数值分析”是计算机科学比较重要的基础课之一，从多方面就“数值分析”课程教学中存在的问题以及提高教学质量、学生兴趣的教学方法进行了探讨。

关键词：数值分析；教学方法；实践

作者简介：黄文芝（1978-），女，湖北武汉人，武汉工程大学计算机科学与工程学院，讲师；张蕾（1982-），女，湖北武汉人，武汉工程大学计算机科学与工程学院，讲师。（湖北 武汉 430073）

基金项目：本文系武汉工程大学青年科学基金项目（项目编号：q201107）的研究成果。

中图分类号：g642.0     文献标识码：a     文章编号：1007-0079（2012）05-0039-02

“数值分析”也称计算方法，它与计算工具发展密切相关。计算方法是数学的一个组成部分，很多方法都与当时的数学家名字相联系，如牛顿插值公式，方程求根的牛顿法，解线性方程组的高斯消去法，多项式求值的秦九韶算法，计算积分的辛普森公式等，这表明计算方法就是数学的一部分，它没有形成单独的学科分支。而计算机出现以后，计算方法迅速发展并形成数学科学的一个独立分支――计算数学。这说明了计算方法与计算机的密切联系，以及在计算机研究领域的重要性。并且数值分析在计算机相关领域应用比较广泛，比如在数学建模中，在图像处理中，在信号处理中等都会用到数值分析中相关的一些知识。这些都说明“数值分析”是计算机专业学生的一门核心专业基础课程。

“数值分析”课程的教学内容主要包括三部分，一部分是插值拟合，一部分是方程和方程组求解，另外一部分是常微分方程初值问题数值解。而数值积分也是在插值的基础进行，故笔者把它归为插值拟合部分。这些内容看上去都是以前学过的知识，积分是在高等数学里学过的，而方程和方程组求解更是中学就重点讲解过的知识，学生刚开始接触这门课的时候会和以前所学的纯数学学习的思想结合起来。通过“数值分析”课程的教学，培养学生用计算机解决问题的能力，并且为后续阶段的专业课程打下基础。

笔者是计算机科学与技术专业的一名老师，使用的教材是清华大学出版的李庆扬等编的《数值分析》，本文就当前“数值分析”课程在计算机科学与技术专业教学中存在的一些问题和教学方法、教学模式等方面进行讨论，其目的在于改进教学方法和手段，提高学生兴趣和教学效果。

一、“数值分析”课程教学中存在的问题

1.数学理论强，公式繁多冗长，学生学习兴趣不高

“数值分析”是数学的一部分，具有与其他数学课程一样的理论性强的特点，但“数值分析”又还有一些和以往学生所学各类数学课程不同的特点。首先，“数值分析”研究的是计算算法，用计算机来解决问题，以前学生学习数学课程大都是从理论学习到作业联系，涉及的知识逻辑推理的特性比较强，并且以往研究的大多数都是连续的，这种研究对象的差异使得学生不能很快接受，思想不能很快转变过来。其次，“数值分析”比以往所学的数学课程的公式更加繁多，更加冗长，比如解线性方程组，如果用以前的知识，学生都会解，但现在解线性方程组不仅仅是要得出结果，更重要的是解线性方程组的算法以及它的实现，这就涉及到至少4个公式，而我们要弄清楚了这些公式的来历才能通过编程实现这个算法，这也是学生不感兴趣的主要原因。

另外，由于学生对数学课程以及对数学公式的害怕，对“数值分析”这门课程的重要性认识不足，当学生学习遇到困难时，容易失去学习兴趣，从而放弃学习。虽然“数值分析”是计算机科学与技术专业的基础课，是大多数课程的基础，但学生还不能理会到“数值分析”这门课程对以后课程的重要性，对于大三的学生来说他们现在所学的课程还没能很好地得到应用，而对他们比较实际的用处――找工作也没有显现出比较重要的作用，因而学生会在潜意识里无视这门课，在课程学习遇到困难的情况下，他们往往会选择放弃学习。

2.知识点多，信息量大，掌握困难

这门课的知识点比较多，信息量比较大，对于理学的学生来说该课程学时比较多，但笔者承担的“数值分析”课程的学时是48学时，并且完全是讲授部分，然而相对于课程所包含的大量内容，这些学时数远远不够，比如函数逼近与快速傅里叶变换，它涉及到范数，赋范线性空间，欧氏空间，三角插值等许多概念，想让学生在规定的学时数内真正掌握这些概念比较困难，尤其是对计算机科学专业的学生而言。因为理学院的学生学过实变函数、泛函分析，所以理解这些概念就略显容易些。

3.重理论，轻实践

当前“数值分析”课程教学过程中，仍然存在理论与实践脱离的现象，虽然这门课实践比较重要，但鉴于课时的安排，大多数教师只能按书本知识来讲，学生听，学生没理解理论的用处，没能立刻就在实践中体现出来，因此使得很多学生只是为了考试而学习，为了学习而学习，不知道它的作用，考完就还给老师。这样他们也只获得了知识的皮毛，而没有抓住知识的精髓和实质。

二、“数值分析”课程教学方法浅谈

1.强调课程的重要性，提高学生的学习兴趣

为了让学生正确对待这门课，应该让学生充分认识到“数值分析”课程在计算机科学与技术专业中的重要性。在组织教学的过程中，可以安排一些有实践经验的学生介绍经验（这样学生更好理解，更容易相信，更实际），联系具体的研究方向，给出简单的例子，论述“数值分析”在计算机科学与技术专业方向中的应用，让学生切实感受到“数值分析”课程是后续课程学习的基础，应用比较广泛。另外，在教学中教师还必须联系实际，在课程中穿插一些有实际应用意义的例子，比如现在很多数学建模就用到“数值分析”的内容，可以就里面简单的例子引用一个。这样让学生了解到“数值分析”不是空洞抽象的理论，而是能够解决实际问题的工具，通过这些方法，使学生逐步树立“数值分析”比较有用，应该学好“数值分析”课程的观念。

然而仅有应该学好该课程的观念还不够，还应该从各个方面提高学生学习的兴趣，兴趣是最好的老师，只有有了兴趣，学生才会真正自主去学习，而不是被动的，为了考试而学习。如何让枯燥的课程变得生动有趣是值得研究的问题。在实际教学过程中，可以采用学生自己讲解，学生之间互相提问等方法，另外也可以编一些小程序，演示计算机解题的过程，这样让学生体会到虽然计算机的功能比较强大，还是需要人脑来控制，灵魂还是人。这样能使学生在整个课题中能主动思考，而不是被动接收。

2.合理取舍教学内容，把握全局，突出重点

“数值分析”课程所涉及的内容非常丰富，但现在课时有限，因此合理取舍教学内容非常重要，应该在有限的学时内，让学生掌握比较重要的理论方法，比如根据学生专业的特点，可以将主要的教学时间安排在讲解误差分析，插值，数值积分，方程和方程组的解法上面。在矩阵特征值计算方面，有时间的条件下可以简单介绍思想方法，而对于常微分方程初值问题的数值解可以舍去，因为本专业的学生没有学常微分方程，所以对常微分方程初值问题的数值解会无法理解。

3.合理使用多种教学方法和手段

传统的“黑板+粉笔”的教学模式对数学课程的教学非常重要，通过板书学生可以了解教师处理问题的思维过程，然而鉴于“数值分析”的特点，又不能完全用传统的教学模式，因为“数值分析”课程中有大量的矩阵和公式，如果单纯使用“黑板+粉笔”，黑板无法板书完整，如果擦掉原先板书的内容又无法把前后联系起来讲解，而使用多媒体就可以解决这一问题。另外，有条件的学校可以把上课安排到有投影的机房，在讲解算法时教师可以演示一些程序，学生学起来就不会觉得完全是在听数学课了。因为是计算机专业的学生，这样和他们的联系更紧密些，他们也可以通过编程来实现算法。

4.强调理论联系实践，培养解决问题的能力

“数值分析”这门课重点讲授的是算法，而学生如果没有很好的实践，对这些算法的应用只能停留在死记硬背上，这不是学习的目的。本来计算机专业也应该突出学生的动手能力，所以对讲授的每个算法都应尽可能让学生编程来实现，这样一来可以巩固学生学到的知识，二来也可以让学生明白这门课不是单纯的数学课，而是和实际联系比较紧密的一门课。当然要实现每个算法都编程，在所授课的学时内是无法完成的，这样就要鼓励学生自己主动去编程，可以采取一些奖励的措施，比如对编程完成比较好的学生可以适当提高平时成绩等。学生自己主动的学习有利于提高其学习兴趣，开发学生智力，培养学生解决问题的能力，从而提高学生的综合素质。

三、总结

随着计算机的广泛应用，“数值分析”课程作为计算机科学与技术的一门专业基础课程，在学生学习和工作中越来越重要，因此“数值分析”课程教学也应该不断更新知识结构，丰富教学内容，改进教学手段，以提高学生学习兴趣，提高教学质量，培养学习的能力，从而为后续课程的学习和将来的工作打下坚实的基础。

第五篇：《数值计算方法》课程教学大纲.

《数值计算方法》课程教学大纲

课程名称：数值计算方法/Mathods of Numerical Calculation 课程代码：0806004066 开课学期：4 学时/学分：56学时/3.5学分（课内教学 40 学时，实验上机 16 学时，课外 0 学时）先修课程：《高等代数》、《数学分析》、《常微分方程》、《C语言程序设计》 适用专业：信息与计算科学

开课院（系）：数学与计算机科学学院

一、课程的性质与任务

数值计算方法是数学与应用数学专业的核心课程之一。它是对一个数学问题通过计算机实现数值运算得到数值解答的方法及其理论的一门学科。本课程的任务是架设数学理论与计算机程序设计之间的桥梁，建立解决数学问题的有效算法，讨论其收敛性和数值稳定性并寻找误差估计式，培养学生数值计算的能力。

二、课程的教学内容、基本要求及学时分配

（一）误差分析

2学时 了解数值计算方法的主要研究内容。2 理解误差的概念和误差的分析方法。熟悉在数值计算中应遵循的一些基本原则。重点：数值计算中应遵循的基本原则。难点：数值算法的稳定性。

（二）非线性方程组的求根

8学时 理解方程求根的逐步搜索法的含义和思路 掌握方程求根的二分法、迭代法、牛顿法及简化牛顿法、非线性方程组求根的牛顿法 3 熟悉各种求根方法的算法步骤，并能编程上机调试和运行或能利用数学软件求非线性方程的近似根。

重点：迭代方法的收敛性、牛顿迭代方法。难点：迭代方法收敛的阶。

（三）线性方程组的解法

10学时 熟练掌握高斯消去法 熟练地实现矩阵的三角分解：Doolittle法、Crout法、Cholesky法、LDR方法。3 掌握线性方程组的直接解法：Doolittle法、Crout法、Cholesky法（平方根法）、改进平方根法、追赶法。

4能熟练地求向量和矩阵的1-范数、2-范数、-范数和条件数。5 理解迭代法的基本思想，掌握迭代收敛的基本定理。掌握解线性方程组的雅可比（Jacobi）迭代法、高斯-赛德尔(Gauss-Seidel)迭代法、逐次超松驰（SOR）迭代法。7能写出线性方程组的各种直接解法和间接解法的算法，并能编程上机运行或能利用数学软件求解线性方程组。

重点：矩阵的三角分解。

难点：线性方程组迭代解法的收敛问题。

（四）插值法

6学时

1.了解插值的一般概念和多项式插值的存在唯一性。

2.熟练掌握Lagrange插值、Newton插值、Hermite插值、分段低次插值及三次样条插值的求解。

3.熟悉曲线拟合的最小二乘法，能熟练地求矛盾方程组的最小二乘解。

4.能对Lagrange插值、Newton插值、Neville插值、Hermite插值、三次样条插值、线拟合的最小二乘法等编程上机调试和运行或借助数学软件求插值函数和曲线拟合。

重点：Lagrange插值、Newton插值、Hermite插值。难点：三次样条插值的求解。

（五）最佳逼近多项式的一般理论

5学时 了解最佳逼近的基本问题。掌握C[a,b]空间中最佳逼近的唯一性问题。3 了解切贝绍夫定理与Vallee-Poussin定理。

（六）数值微分与数值积分

5学时 了解数值积分的基本思想，能够熟练地确定具体求积公式的代数精度及确定求积公式的节点和系数。熟练地用Newton-cotes公式，Romberg公式，两点、三点Gauss公式等进行数值积分 重点：确定具体求积公式的代数精度及确定求积公式的节点和系数。难点：用待定系数法确定Gauss型求积公式的节点和系数。

（七）常微分方程的数值解

4学时 理解常微分方程的数值解的含义 掌握常微分方程的欧拉解法、R—K方法、亚当姆斯方法，理解其算法思想。重点：基于数值积分的方法。难点：R—K方法。

三、推荐教材及参考书

推荐教材：

1、张韵华等编著，数值计算方法与算法，科学出版社，2001。

2、冯天祥编著，数值计算方法，四川科技出版社，2003。参考书：

1、冯天祥编著，数值计算方法理论与实践研究，西南交通大学出版社，２００５。

2、李庆扬等著，数值分析，华中理工大学出版社，2000。

3、林成森著，数值计算方法，科学出版社出版，1999。

4、李庆扬等著，现代数值分析，高等教育出版社，1998。

5封建湖等，计算方法典型题分析解集，西北工业大学出版社，1999。

四、结合近几年的教学改革与研究，对教学大纲进行的新调整 增加了最佳逼近多项式的一般理论。

大纲制订者：冯玉明

大纲审定者：陈小春

制订日期：2008-11-15