数值分析(计算方法)实验一
数值分析计算方法实验报告
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
数值分析实验报告一
实验报告课程名称:数值分析课题名称:比较算法专业:勘查技术与工程姓名:韩鹏洋班级:061132班完成日期:2015 年10 月11 日实验报告一、实验名称比较两种算法收敛性及复杂度二、实验目的(1)培养编程与上机调试能力(2)观察不同算法的差异(3)评估各算法稳定性三、实验要求利用matlab计算算法,并绘图观察收敛性。
四、实验原理利用泰勒展开式逼近函数值五、实验题目求ln 2的近似值六、实验步骤(1)写出ln(1+x)展开式(2)利用Matlab编程计算(3)最后结果分析七、实验整体流程图或算法八、程序及其运行结果clear all;ticn=1:100;s=0;for i=1:100s1=(-1).^(i-1)/i;s=s+s1;y(i)=s;endplot(n,y,'ro');tocclear all;ticn=1:50;s=0;for i=1:50s1=2*(1/3).^(2*i-1)/(2*i-1);s=s+s1;y(i)=s;endhold on;plot(n,y,'b-');toc运行结果:方法1时间已过0.369496 秒。
方法2时间已过0.025252 秒。
九、实验结果分析方法一趋近速度慢,复杂度100+100+(1+2+…+99)=5150 方法二趋近快,复杂度150+3+5+7+…+99=2499选用第二种方法更好十、实验体会。
数值分析实验报告——方程求根
《数值分析》实验报告实验一方程求根一、实验目的:掌握二分法、Newton法、不动点迭代法、弦截法求方程的根的各种计算方法、并实施程序调试和运行,学习应用这些算法于实际问题。
二、实验内容:二分法、Newton法、不动点迭代法、弦截法求方程的根、程序的调试和运行,给出实例的计算结果。
观察初值对收敛性的影响。
三、实验步骤:①、二分法:定义:对于区间[a,b]上连续不断且f(a)·f(b)<0的函数y=f(x),通过不断地把函数f(x)的零点所在的区间一分为二,使区间的两个端点逐步逼近零点,进而得到零点近似值的方法叫二分法。
实现方法:首先我们设一方程400*(x^4)-300*(x^3)+200*(x^2)-10*x-1=0,并求其在区间[0.1,1]上的根,误差限为e=10^-4。
PS:本方法应用的软件为matlab。
disp('二分法')a=0.1;b=1;tol=0.0001;n0=100;fa=400*(a.^4)-300*(a.^3)+200*(a.^2)-10*a-1;for i=1:n0 p=(a+b)/2;fp=400*(p.^4)-300*(p.^3)+200*(p.^2)-10*p-1;if fp==0||(abs((b-a)/2)<tol)disp('用二分法求得方程的根p=')disp(p)disp('二分迭代次数为:')disp(i)break;end;if fa*fp>0 a=p;else b=p;end;end;if i==n0&&~(fp==0||(abs((b-a)/2)<tol)) disp(n0) disp('次二分迭代后没有求出方程的根')end;程序调试:运行结果:用二分法求得方程的根p=0.1108二分迭代次数为:14②Newton法定义:取定初值x0,找到函数对应的点,然后通过该点作函数切线,交x轴,得到新的横坐标值,然后找函数对应的点,做切线,得到新的横坐标值,重复上述步骤,多次迭代,直到收敛到需要的精度。
计算方法与计算 实验一误差分析
% 输出的量--每次迭代次数k和迭代值xk,
%
--每次迭代的绝对误差juecha和相对误差xiangcha,
误差分析
误差问题是数值分析的基础,又是数值分析中一个困难的课题。在实际计算 中,如果选用了不同的算法,由于舍入误差的影响,将会得到截然不同的结果。 因此,选取算法时注重分析舍入误差的影响,在实际计算中是十分重要的。同时, 由于在数值求解过程中用有限的过程代替无限的过程会产生截断误差,因此算法 的好坏会影响到数值结果的精度。 一、实验目的
因为运行后输出结果为: y 1.370 762 168 154 49, yˆ =1.370 744 664 189
38, R 1.750 396 510 491 47e-005, WU= 1.782 679 830 970 664e-005 104 . 所
以, yˆ 的绝对误差为 10 4 ,故 y
③ 运行后输出计算结果列入表 1–1 和表 1-2 中。
④ 将算法 2 的 MATLAB 调用函数程序的函数分别用 y1=15-2*x^2 和
y1=x-(2*x^2+x-15)/(4*x+1)代替,得到算法 1 和算法 3 的调用函数程序,将其保
存,运行后将三种算法的前 8 个迭代值 x1, x2 ,, x8 列在一起(见表 1-1),进行
的精确解 x* 2.5 比较,观察误差的传播.
算法 1 将已知方程化为同解方程 x 15 2x2 .取初值 x0 2 ,按迭代公式
xk1 15 2xk2
数值分析实验报告模板
数值分析实验报告模板篇一:数值分析实验报告(一)(完整)数值分析实验报告12345篇二:数值分析实验报告实验报告一题目:非线性方程求解摘要:非线性方程的解析解通常很难给出,因此线性方程的数值解法就尤为重要。
本实验采用两种常见的求解方法二分法和Newton法及改进的Newton法。
利用二分法求解给定非线性方程的根,在给定的范围内,假设f(x,y)在[a,b]上连续,f(a)xf(b) 直接影响迭代的次数甚至迭代的收敛与发散。
即若x0 偏离所求根较远,Newton法可能发散的结论。
并且本实验中还利用利用改进的Newton法求解同样的方程,且将结果与Newton法的结果比较分析。
前言:(目的和意义)掌握二分法与Newton法的基本原理和应用。
掌握二分法的原理,验证二分法,在选对有根区间的前提下,必是收敛,但精度不够。
熟悉Matlab语言编程,学习编程要点。
体会Newton使用时的优点,和局部收敛性,而在初值选取不当时,会发散。
数学原理:对于一个非线性方程的数值解法很多。
在此介绍两种最常见的方法:二分法和Newton法。
对于二分法,其数学实质就是说对于给定的待求解的方程f(x),其在[a,b]上连续,f(a)f(b) Newton法通常预先要给出一个猜测初值x0,然后根据其迭代公式xk?1?xk?f(xk) f'(xk)产生逼近解x*的迭代数列{xk},这就是Newton法的思想。
当x0接近x*时收敛很快,但是当x0选择不好时,可能会发散,因此初值的选取很重要。
另外,若将该迭代公式改进为xk?1?xk?rf(xk) 'f(xk)其中r为要求的方程的根的重数,这就是改进的Newton 法,当求解已知重数的方程的根时,在同种条件下其收敛速度要比Newton法快的多。
程序设计:本实验采用Matlab的M文件编写。
其中待求解的方程写成function的方式,如下function y=f(x);y=-x*x-sin(x);写成如上形式即可,下面给出主程序。
数值分析实验报告1
显然该多项式的全部根为1,2,…,20共计20个,且每个根都是单重的。现考虑该多项式的一个扰动
其中是一个非常小的数。这相当于是对(1。1)中的系数作一个小的扰动.我们希望比较(1。1)和(1。2)根的差别,从而分析方程(1。1)的解对扰动的敏感性.
实验内容:为了实现方便,我们先介绍两个Matlab函数:“roots”和“poly”。
(3)取矩阵阶数n=20或者更大,重复上述实验过程,观察记录并分析不同的问题及消去过程中选择不同的主元时计算结果的差异,说明主元素的选取在消去过程中的作用。
(4)选取其他你感兴趣的问题或者随机生成矩阵,计算其条件数。重复上述实验,观察记录并分析实验结果。
实验过程:
程序:
建立M文件:
function x=gauss(n,r)
实验总结:
利用MATLAB来进行病态问题的实验,虽然其得出的结果是有误差的,但是可以很容易的得出对一个多次的代数多项式的其中某一项进行很小的扰动,对其多项式的根会有一定的扰动的,所以对于这类病态问题可以借助于MATLAB来进行问题的分析。
学号:06450210
姓名:万轩
实验二插值法
实验2.1(多项式插值的振荡现象)
学号:06450210
姓名:万轩
实验五解线性方程组的直接方法
实验5。1(主元的选取与算法的稳定性)
问题提出:Gauss消去法是我们在线性代数中已经熟悉的.但由于计算机的数值运算是在一个有限的浮点数集合上进行的,如何才能确保Gauss消去法作为数值算法的稳定性呢?Gauss消去法从理论算法到数值算法,其关键是主元的选择。主元的选择从数学理论上看起来平凡,它却是数值分析中十分典型的问题.
其中若变量a存储n+1维的向量,则该函数的输出u为一个n维的向量。设a的元素依次为,则输出u的各分量是多项式方程
数值分析实验报告1
p
得到m=(00)T
即M0=0 ;M1=;M2=;M3=;M4=0
则根据三次样条函数定义,可得:
S(x)=
接着,在Command Window里输入画图的程序代码,
下面是画牛顿插值以及三次样条插值图形的程序:
x=[ ];
y=[ ];
plot(x,y)
hold on
for i=1:1:5
y(i)= *(x(i)*(x(i)*(x(i)*(x(i)*(x(i)*(x(i)*(x(i)
Pn=f(x0)+f[x0,x1](x-x0)+ f[x0,x1,x2](x-x0) (x-x1)+···+ f[x0,x1,···xn](x-x0) ···(x-xn-1)
我们要知道牛顿插值多项式的系数,即均差表中得部分均差。
在MATLAB的Editor中输入程序代码,计算牛顿插值中多项式系数的程序如下:
【实验原理】
《数值分析》第二章“插值法”的相关内容,包括:牛顿多项式插值,三次样条插值,拉格朗日插值的相应算法和相关性质。
【实验环境】(使用的软硬件)
软件:
MATLAB 2012a
硬件:
电脑型号:联想 Lenovo 昭阳E46A笔记本电脑
操作系统:Windows 8 专业版
处理器:Intel(R)Core(TM)i3 CPU M 350 @
实验内容:
【实验方案设计】
第一步,将书上关于三种插值方法的内容转化成程序语言,用MATLAB实现;第二步,分别用牛顿多项式插值,三次样条插值,拉格朗日插值求解不同的问题。
【实验过程】(实验步骤、记录、数据、分析)
实验的主要步骤是:首先分析问题,根据分析设计MATLAB程序,利用程序算出问题答案,分析所得答案结果,再得出最后结论。
数值分析计算方法
void main(void)
{
double fa=fc(1),fb=fc(3),a=1,b=3,f,x0;
int k=0;
for(;b>a&&b-a>=pow(10,-4)*0.5;)
{
仁fc((a+b)/2);
if(仁=0)
{ x0=(a+b)/2; break;
}Байду номын сангаас
else if(fa*f<0)
{
double y;
y=pow((3*x+1),1.0/3);
return y;
}
double Derivative1(double x)
{
double y;
y=pow((3*x+1),-2.0/3);
return y;
} double Iterate2(double x) {
double y;
y=(1-x*x*x)/3.0;
1.19 817
1.23223
作五次插值,并求x1=0.46,x2=0.55,x3=0.60时的函数近似值.
1■实验目的:通过拉格朗日插值和牛顿插值的实例,了解两种求解方法,并分析 各自的优缺点。
2.算法描述:
3.源程序:
拉格朗日插值:
#include<stdio.h>
#define k2
void main()
三•实验三:分别用复化梯形公式和复化辛卜生公式计算f(x)=s in( x)/x的积分, 并与准确值比较判断精度。
1■实验目的:通过实例体会各种算法的精度。熟练掌握复化梯形,复化辛普森, 复化柯特斯求积方法的程序。
实验报告数值分析
《数值分析》实验报告姓 名: 学 号: 专 业:指导教师: 刘 建 生 教 授 日 期: 2015年12月25日实验一 Lagrange/newton 插值一:对于给定的一元函数)(x f y =的n+1个节点值(),0,1,,j j y f x j n ==L 。
试用Lagrange 公式求其插值多项式或分段二次Lagrange 插值多项式。
数据如下: 求计算(0.596)f ,(0.99)f 的值(提示:结果为(0.596)0.625732f ≈,(0.99) 1.05423f ≈ )试构造Lagrange 多项式6L ()x ,计算的(1.8)f ,(6.15)f 值。
(提示:结果为(1.8)0.164762f ≈, (6.15)0.001266f ≈ )二:实验程序及注释MATLAB 程序:function f=lagrange(x0,y0,x )n=length(x0); m=length(y0); format long s=0.0; for k=1:n p=1.0; for j=1:n if j~=kp=p*(x-x0(j))/(x0(k)-x0(j));endends=s+y0(k)*p;Endf=s;end结果运行:结果与提示值完全吻合,说明Lagrange插值多项式的精度是很高的;同时,若采用三点插值和两点插值的方法,用三点插值的精度更高。
若同时采用两点插值,选取的节点距离x越近,精度越高。
三:采用newton插值进行计算算法程序如下:format long;x0=[0.4 0.55 0.65 0.80 0.95 1.05 ];y0=[0.41075 0.57815 0.69675 0.90 1.00 1.25382 ];x=0.596;n=max(size(x0));y=y0(1);%disp(y);s=1;dx=y0;for i=1:n-1dx0=dx; for j=1:n-idx(j)=(dx0(j+1)-dx0(j))/(x0(i+j)-x0(j)); end df=dx(1); s=s*(x-x0(i));y=y+s*df; %计算 %%disp(y); end disp(y)运行结果:绘制出曲线图:与结果相吻合。
数值分析实验 实验报告
数值分析实验实验报告数值分析实验实验报告引言在现代科学与工程领域,数值分析是一项重要的技术手段。
通过数值方法,我们可以利用计算机模拟和解决各种实际问题,如物理、化学、生物、经济等领域中的方程求解、优化问题、数据拟合等。
本实验旨在通过实际案例,探讨数值分析的应用和效果。
实验一:方程求解首先,我们考虑一个简单的方程求解问题。
假设我们需要求解方程f(x) = 0的根,其中f(x)是一个在给定区间[a, b]上连续且单调的函数。
为了实现这个目标,我们可以采用二分法、牛顿法、弦截法等数值方法。
在本实验中,我们选择使用二分法来求解方程f(x) = 0。
这种方法的基本思想是通过不断缩小区间[a, b]的范围,直到找到一个近似的根。
我们首先选取一个中间点c,计算f(c)的值,然后根据f(c)与0的关系,将区间[a, b]分成两部分。
重复这个过程,直到找到满足精度要求的根。
实验二:数据拟合接下来,我们考虑一个数据拟合的问题。
假设我们有一组离散的数据点,我们希望找到一个函数,使得该函数与这些数据点的拟合误差最小。
为了实现这个目标,我们可以采用最小二乘法等数值方法。
在本实验中,我们选择使用最小二乘法来进行数据拟合。
这种方法的基本思想是通过最小化数据点与拟合函数之间的误差平方和,来确定拟合函数的参数。
我们首先选择一个拟合函数的形式,如线性函数、多项式函数等。
然后,通过最小化误差平方和的方法,计算出拟合函数的参数。
实验三:优化问题最后,我们考虑一个优化问题。
假设我们需要在给定的约束条件下,找到一个使得目标函数取得最大或最小值的变量。
为了实现这个目标,我们可以采用梯度下降法、遗传算法等数值方法。
在本实验中,我们选择使用梯度下降法来解决优化问题。
这种方法的基本思想是通过迭代的方式,不断调整变量的取值,直到找到一个满足约束条件的最优解。
我们首先计算目标函数关于变量的梯度,然后根据梯度的方向和大小,更新变量的取值。
通过不断迭代,我们可以逐步接近最优解。
数值分析—实验报告1
wc=yi-jqz
wc = 0.00487093023460 数据分析 :
表:ln(0.6)
内插
外插
x
0.4
0.5
0.7
0.7
0.8
0.9
f (x) -0.916291 -0.693147 -0.356675 -0.356675 -0.223144 -0.105361
yi
-0.50660833333333
31
左上方框里填写学号后两位,学习委员按此顺号(报告展开排序)交给老师
数值分析实验报告
专业 信息与计算科学 班级 15 级 1 班 组别 指导教师 汪玉霞
姓名
史博强
同组人
实验时间 2017
年 10 月 20 日
实验地点 k7 —403
实验名称
插值函数与数据拟合
实验目的:
(1)由函数 f ( x) 的 n 1个节点处函数值得出 n 次 Lagrange 插值函数; (2)由函数 f ( x) 的 n 1个节点处函数值得出 n 次 Newton插值函数; (3)由函数 f ( x) 的个 n 1节点处函数值得出 Hermite 插值函数或分段三次 Hermite 函数; (4)由未知函数的离散数据 f ( xi ), i 1,2, ,n 得出最小二乘拟合函数。
for i=1:n-k Y(i,k+1)=(Y(i+1,k)-Y(i,k))/(x(i+k)-x(i));
end end ndcz=0; for i=1:n
z=1; for k=1:i-1
z=z*(xin-x(k)); end ndcz=ndcz+Y(1,i)*z; End
题目 :设 f (x)=lnx 并已知 f(x) 的数据如下表:
数值计算方法实验报告
数值分析实验报告实验一、解线性方程组的直接方法——梯形电阻电路问题利用追赶法求解三对角方程组的方法,解决梯形电阻电路问题:电路中的各个电流{1i ,2i ,…,8i }须满足下列线性方程组:R V i i =- 22 210 252321=-+-i i i 0 252 432=-+-i i i 0 252 543=-+-i i i 0 252 654=-+-i i i 0 252 765=-+-i i i 0 252 876=-+-i i i 052 87=+-i i设V 220=V ,Ω=27R ,运用追赶法,求各段电路的电流量。
问题分析:上述方程组可用矩阵表示为:⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡--------------00000001481.8522520000002520000002520000002520000002520000002520000002287654321i i i i i i i i问题转化为求解A x b =,8阶方阵A 满足顺序主子式(1,2...7)0i A i =≠,因此矩阵A存在唯一的Doolittle 分解,可以采用解三对角矩阵的追赶法!追赶法a=[0 -2 -2 -2 -2 -2 -2 -2]; b=[2 5 5 5 5 5 5 5];c=[-2 -2 -2 -2 -2 -2 -2 0]; d=[220/27 0 0 0 0 0 0 0];Matlab 程序function x= zhuiganfa( a,b,c,d )%追赶法实现要求:|b1|>|C1|>0,|bi|>=|ai|+|ci| n=length(b); u=ones(1,n); L=ones(1,n); y=ones(1,n); u(1)=b(1); y(1)=d(1); for i=2:nL(i)=a(i)/u(i-1);u(i)=b(i)-c(i-1)*L(i); y(i)=d(i)-y(i-1)*L(i); endx(n)=y(n)/u(n); for k=n-1:-1:1x(k)=(y(k)-c(k)*x(k+1))/u(k); end endMATLAB 命令窗口输入:a=[0 -2 -2 -2 -2 -2 -2 -2]; b=[2 5 5 5 5 5 5 5];c=[-2 -2 -2 -2 -2 -2 -2 0] d=[220/27 0 0 0 0 0 0 0];x= zhuiganfa(a,b,c,d )运行结果为:x =8.1478 4.0737 2.0365 1.0175 0.5073 0.2506 0.1194 0.0477存在问题根据电路分析中的所讲到的回路电流法,可以列出8个以回路电流为独立变量的方程,课本上给出的第八个回路电流方程存在问题,正确的应该是78240i i -+=;或者可以根据电路并联分流的知识,同样可以确定78240i i -+=。
实验一实验报告1
贵州师范大学数学与计算机科学学院学生实验报告课程名称:数值分析班级数学与应用数学2班实验日期:年月日学号: 110701020060 姓名周明祥指导教师:杨一都实验成绩:一、实验名称实验一:递推法的稳定性,秦九韶算法二、实验目的及要求1. 熟悉数值稳定的概念, 通过上机计算,了解舍入误差所引起的数值不稳定性.2. 培养Matlab编程与上机调试能力.三、实验环境每人一台计算机,要求安装Windows XP操作系统,Microsoft office2003、MATLAB6.5(或7.0).四、实验内容1.教材例1.13中,取15位数字计算,并分析、比较计算结果.2.设100999832,用秦九韶算法编程计=+++++++f x x x x x x x()101100994321算()f x在1,2,3,4x=上的值.五、算法描述及实验步骤1、输入函数f(x)=1/(x+5);输出 I(14);步1 I(1)=int(f,0,1);步2 用for循环语句执行递推关系I(n+1)=1/n-5*I(n);步3 输出I(14);结束。
2、输入多项式系数a(0),a(2),...,a(n),给定点x。
(x)在x处的值y。
输出多项式Pn步1 v<==a(n)。
步2 对k=1,2,...,n执行v<==v*x+a(n-k)。
步3 y<==v。
步4 输出多项式的值y。
结束。
六、调试过程及实验结果1、 y=ditui(14)ans =[ .18232155679396, .8839221603023e-1, .580389198489e-1, .43138734088e-1, .3430632956e-1, .2846835222e-1, .243249055e-1, .21232615e-1, .18836925e-1, .1692649e-1, .153675e-1, .140713e-1,.12978e-1, .1204e-1, .1123e-1]2、 a=[101:-1:1];x=[1,2,3,4];Qinjiu(a,100,x)ans =1.0e+062 *0.0000 0.0000 0.0000 2.1569vpa(ans,8)ans =[ 5151., .25353012e33, .77693161e50, .21568680e63]七、总结八、附录(源程序清单)1、 function y=ditui(n)syms xf=1./(x+5);I(1)=int(f,0,1);for i=1:nI(i+1)=1/i-5*I(i);endvpa(I,15)2、 function y=Qinjiu(a,n,x) v=a(1);for k=1:nv=v.*x+a(k+1);endy=v;。
数值分析实验报告
《数值分析》实验报告学院:计算机科学与软件学院姓名:XXX班级:计算机XX班学号:XXXXXX实验一:舍入误差与数值稳定性实验目的:1、 通过上机编程,复习巩固以前所学程序设计语言;2、 通过上机计算,了解舍入误差所引起的数值不稳定性。
3、 通过上机计算,了解运算次序对计算结果的影响,从而尽量避免大数吃小数的现象。
实验内容:用两种不同的顺序计算644834.11000012≈∑=-n n ,分析其误差的变化。
实验流程图:实验源程序:#include <stdio.h>#include <math.h>void main(){ int i;float s1=0,s2=0,d1,d2;for (i=1;i<=10000;i++)s1=s1+1.0f/(i*i);for (i=10000;i>=1;i--)s2=s2+1.0f/(i*i);d1=(float)(fabs(1.644834-s1));d2=(float)(fabs(1.644834-s2));printf("正向求和结果为%f\n 误差为%f\n\n",s1,d1);printf("反向求和结果为%f\n 误差为%f\n\n",s2,d2);if(d1<d2)printf("正向求和误差小于负向求和误差\n");else if(d1==d2)printf("正向求和误差等于负向求和误差\n"); elseprintf("正向求和误差大于负向求和误差\n");}实验结果:实验分析:第一次做数值实验,又一次使用C语言编程,没有了刚学习C语言的艰难,能够将实验步骤转换成流程图并编写出完整的实验代码,在经过多次调试、改正后得到正确的程序和结果。
这个实验较简单,计算误差时如果输入数据有误差,而在计算过程中舍入误差不增长,则称此算法是稳定的,否则称此算法是数值不稳定的,减少运算次数可以减小舍入误差。
数值分析实验指导
实验一 算法设计及误差分析1.1 实验目的○1根据给定的问题,给出稳定的算法,并分析误差的控制。
○2了解舍入误差所引起的数值不稳定性 1.2 实验内容 1、递推计算的稳定性 实验题目,1dx xa xI nn ⎰+=n=0,1, (10)其中a 为参数,分别对a=0.05及a=15按下列两种方法计算,列出结果,并对其可靠性进行分析比较,说明原因。
方案1:用递推公式,11n aII n n +-=- n=1,2,…,10递推初值可由积分直接得.1ln0aa I +=方案2:用递推公式),1(11nI aI n n +-=- n=1,...1,-N N根据估计式,)1(1)1)(1(1+<<++n a I n a n 当1+≥n n a或,1)1)(1(1n I n a n ≤<++ 当10+<≤n n a取递推初值])1(1)1)(1(1[21++++≈N a N a I N=,)1)(1(212_NIN a a a ∆+++ 当1+≥N N a或,]1)1)(1(1[21_NN INN a I ∆+++=当10+≤≤N N a2、不同方案收敛速度的比较实验题目三种求ln2的算法比较。
构造逼近ln2的数列,求出ln2的近似值,要求精度51021-⨯=ε,观察比较3种计算方案的收敛速度。
方案1:利用级数∑∞=--=+-+-=11)1( (4)1312112ln k k k设∑=--=nk k n kS 11)1(,则n S ≈2ln .方案2:对上述∑=--=nk k n kS 11)1(,按21212)(ˆ---+---=n n n n n nn S S S S S S S ,...,4,3=n 生成新数列n S ˆ,则n S ˆ2ln ≈方案3 利用级数∑∞==+⨯+⨯+⨯+⨯=143221 (2)412312212112ln k kk设∑==nk kn k S 121,则nS ≈2ln .实验二 插值方法2.1 实验目的① 掌握牛顿插值法的基本思路和步骤; ② 培养编程与上机调试能力。
数值分析实验报告1
数值分析实验报告1数值分析上机实验报告(注:本实验报告中所有程序均为MATLAB语⾔程序)班级:姓名:学号:⼀章1、利⽤数值积分计算n I =21n x ex e -?dx (n=0,1,2,……). ⽬的:定积分数值求解原理:梯形公式法程序:clearformat long ;k=input('k=');m=input('m=');for n=1:kh=1/m;x=0:h:1;f=x.^n.*exp(x.^2);for i=1:ms(i)=(f(i)+f(i+1))*h/2;ends=sum(s);I(n)=exp(-1)*s;endI 运⾏结果:k=9m=1000I =Columns 1 through 60.3160604988 0.2309605799 0.1839401373 0.1535601302 0.1321211422 0.1161015912Columns 7 through 90.1036390735 0.0936475974 0.08544762262、利⽤秦九韶算法计算当0a =5,n a =21n a -+3;n=100,x=0.5;n=150,x=13多项式n p (x )=n a n x +…11n n a x --…1a x +0a 的值。
⽬的:通过调整程序以简化计算步骤,减少运算次数原理:秦久韶算法程序:n=input('n=');x=input('x=');a(1)=5;for k=1:n;a(k+1)=2.*a(k)+3;ends(n+1)=a(n+1);for i=n:-1:1s(i)=x.*s(i+1)+a(i);endPnx=s(1)运⾏结果:n=100x=0.5Pnx =802.0000000n=150x=13Pnx =1.4659714820e+2133、设0Y =28,按递推公式n Y = 1n Y -100Y ,500Y ≈27.982(五位有效数字),试问计算100Y 、500Y 将有多⼤的误差。
《数值分析》_实验
( )
( )
4、另外选取一个近似表达式,尝试拟合效果的比较; 5、* 绘制出曲线拟合图。 三、目的和意义 1、掌握曲线拟合的最小二乘法; 2、最小二乘法亦可用于解超定线代数方程组; 3、探索拟合函数的选择与拟合精度间的关系
第5页
数值分析实验,宋伦继
6
实验四 数值积分与数值微分
一、问题提出 选用复合梯形公式,复合 Simpson 公式,Romberg 算法,计算
x
*
⎡ 0 ⎤ ⎢ −6 ⎥ ⎥ ⎢ ⎢ 20 ⎥ ⎥ ⎢ ⎢ 23 ⎥ ⎢ 9 ⎥ ⎥ ⎢ ⎢− 22⎥ ⎢ − 15 ⎥ ⎥ ⎢ ⎥ ⎢ 45 ⎦ ⎣
=
( 1,
-1,
0,
2,
1,
-1,
0,
2)
T
第7页
数值分析实验,宋伦继
8
3、三对角形线性方程组
⎡ 4 − 1 0 0 0 0 0 0 0 0 ⎤ ⎡ x1 ⎤ ⎡ 7 ⎤ ⎢− 1 4 − 1 0 0 0 0 0 0 0 ⎥ ⎢ x ⎥ ⎢ 5 ⎥ ⎥ ⎢ 2⎥ ⎥ ⎢ ⎢ ⎢ 0 − 1 4 − 1 0 0 0 0 0 0 ⎥ ⎢ x3 ⎥ ⎢ − 13⎥ ⎥⎢ ⎥ ⎥ ⎢ ⎢ ⎢ 0 0 − 1 4 − 1 0 0 0 0 0 ⎥ ⎢ x4 ⎥ ⎢ 2 ⎥ ⎢ 0 0 0 − 1 4 − 1 0 0 0 0 ⎥ ⎢ x5 ⎥ ⎢ 6 ⎥ ⎥⎢ ⎥ = ⎢ ⎥ ⎢ ⎢ 0 0 0 0 − 1 4 − 1 0 0 0 ⎥ ⎢ x6 ⎥ ⎢− 12⎥ ⎢ 0 0 0 0 0 −1 4 −1 0 0 ⎥ ⎢ x ⎥ ⎢ 14 ⎥ ⎥ ⎢ 7⎥ ⎥ ⎢ ⎢ ⎢ 0 0 0 0 0 0 − 1 4 − 1 0 ⎥ ⎢ x8 ⎥ ⎢−4⎥ ⎢ 0 0 0 0 0 0 0 − 1 4 − 1⎥ ⎢ x ⎥ ⎢ 5 ⎥ ⎥ ⎢ 9⎥ ⎥ ⎢ ⎢ ⎥ ⎥⎢ ⎥ ⎢ 0 0 0 0 0 0 0 0 −1 4 ⎦ ⎢ −5⎦ ⎣ ⎣ ⎣ x10 ⎦
数值计算方法实验指导(Matlab版)
《数值计算方法》实验指导(Matlab版)学院数学与统计学学院计算方法课程组《数值计算方法》实验1报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验1 算法设计原则验证(之相近数相减、大数吃小数和简化计算步骤) 2. 实验题目(1) 取1610=z ,计算z z -+1和)1/(1z z ++,验证两个相近的数相减会造成有效数字的损失.(2) 按不同顺序求一个较大的数(123)与1000个较小的数(15310-⨯)的和,验证大数吃小数的现象.(3) 分别用直接法和九韶算法计算多项式n n n n a x a x a x a x P ++++=--1110)(在x =1.00037处的值.验证简化计算步骤能减少运算时间.对于第(3)题中的多项式P (x ),直接逐项计算需要2112)1(+=+++-+n n n 次乘法和n 次加法,使用九韶算法n n a x a x a x a x a x P ++++=-)))((()(1210则只需要n 次乘法和n 次加法. 3. 实验目的验证数值算法需遵循的若干规则. 4. 基础理论设计数值算法时,应避免两个相近的数相减、防止大数吃小数、简化计算步骤减少运算次数以减少运算时间并降低舍入误差的积累.两相近的数相减会损失有效数字的个数,用一个大数依次加小数,小数会被大数吃掉,乘法运算次数太多会增加运算时间. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab6. 实验过程(1) 直接计算并比较;(2) 法1:大数逐个加1000个小数,法2:先把1000个小数相加再与大数加; (3) 将由高次项到低次项的系数保存到数组A[n]中,其中n 为多项式次数.7. 结果与分析 (1) 计算的z z -+1= ,)1/(1z z ++.分析:(2) 123逐次加1000个6310-⨯的和是 ,先将1000个6310-⨯相加,再用这个和与123相加得.分析:(3) 计算次的多项式:直接计算的结果是,用时;用九韶算法计算的结果是,用时.分析:8. 附录:程序清单(1) 两个相近的数相减.%*************************************************************%* 程序名:ex1_1.m *%* 程序功能:验证两个相近的数相减会损失有效数字个数 *%*************************************************************z=1e16;x,y======================================================================(2) 大数吃小数%*************************************************************%* 程序名:ex1_2.m *%* 程序功能:验证大数吃小数的现象. *%*************************************************************clc; % 清屏clear all; % 释放所有存变量format long; % 按双精度显示浮点数z=123; % 大数t=3e-15; % 小数x=z; % 大数依次加小数% 重复1000次给x中加上ty=0; % 先累加小数% 重复1000次给y中加上ty=z + y; % 再加到大数x,y======================================================================(3) 九韶算法%*************************************************************%* 程序名:ex1_3.m *%* 程序功能:验证九韶算法可节省运行时间. *%*************************************************************clc; % 清屏clear all; % 释放所有存变量format long; % 按双精度显示浮点数A=[8,4,-1,-3,6,5,3,2,1,3,2,-1,4,3,1,-2,4,6,8,9,50,-80,12,35,7,-6,42,5,6,23,74,6 5,55,80,78,77,98,56];A(10001)=0; % 扩展到10001项,后面的都是分量0% A为多项式系数,从高次项到低次项x=1.00037;n=9000; % n为多项式次数% 直接计算begintime=clock; % 开始执行的时间 % 求x的i次幂% 累加多项式的i次项endtime=clock; % 完毕执行的时间time1=etime(endtime,begintime); % 运行时间disp('直接计算');disp(['p(',num2str(x),')=',num2str(p)]);disp([' 运行时间: ',num2str(time1),'秒']);% 九韶算法计算begintime=clock; % 开始执行的时间% 累加九韶算法中的一项endtime=clock; % 完毕执行的时间time2=etime(endtime,begintime); % 运行时间disp(' ');disp('九韶算法计算');disp(['p(',num2str(x),')=',num2str(p)]);disp([' 运行时间: ',num2str(time2),'秒']);《数值计算方法》实验1报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验1 算法设计原则验证(之数值稳定性) 2. 实验题目 计算定积分⎰==-1110,1,0,d n x e xI x nn ,分别用教材例1-7推导出的算法A 和B ,其中:算法A :⎩⎨⎧≈-=-6321.0101I nI I n n 算法B :⎪⎩⎪⎨⎧≈-=-0)1(1101I I nI n n 验证算法不稳定时误差会扩大.3. 实验目的验证数值算法需遵循的若干规则. 4. 基础理论设计数值算法时,应采用数值稳定性好的算法.数值稳定的算法,误差不会放大,甚至会缩小;而数值不稳定的算法会放大误差. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab6. 实验过程分别用数组IA[ ]和IB[ ]保存两种算法计算的结果. 7. 结果与分析 运行结果:(或拷屏)8. 附录:程序清单%*************************************************************%* 程序名:ex1_4.m *%* 程序功能:验证数值稳定性算法可控制误差. *%*************************************************************clc; % 清屏clear all; % 释放所有存变量format long; % 按双精度显示浮点数I=[0.856, 0.144, 0.712, 0.865, ...0.538, 0.308, 0.154, 0.938, ...0.492, 0.662, 0.843];% 保留14位小数的精确值, …是Matlab中的续行符% 算法AIA(1) = 0.6321; % Matlab下标从1开始,所以要用IA(n+1)表示原问题中的I(n)% 算法Bdisp('n 算法A 算法B 精确值');for n=1:11fprintf('%2d %14.6f %14.6f %14.6f\n',n-1,IA(n),IB(n),I(n));end% n显示为2位整数, 其它显示为14位其中小数点后显示6位的小数《数值计算方法》实验1报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验1 算法设计原则(除数绝对值不能太小) 2. 实验题目将线性方程组增广矩阵利用初等行变换可化为⎪⎪⎭⎫⎝⎛→-⎪⎪⎭⎫ ⎝⎛→-⎪⎪⎭⎫ ⎝⎛''0'0''02221112'12221121112222211121122121121b a b a r r b a b a a r r b a a b a a a a a a由此可解得'/',/'22221111a b x a b x ==.分别解增广矩阵为161011212-⎛⎫ ⎪⎝⎭和162121011-⎛⎫⎪⎝⎭的方程组,验证除数绝对值远小于被除数绝对值的除法会导致结果失真. 3. 实验目的验证数值算法需遵循的若干规则. 4. 基础理论设计数值算法时,应避免除数绝对值远小于被除数绝对值的除法,否则绝对误差会被放大,使结果失真. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab6. 实验过程用二维数组A 和B 存放方程组的增广矩阵,利用题目所给初等行变换求解方程组. 7. 结果与分析第1种顺序的方程组的解为x =,y =;第2种顺序的方程组的解为x =,y =. 分析:8. 附录:程序清单%************************************************************* %* 程 序 名:ex1_5.m * %* 程序功能:验证除数的绝对值太小可能会放大误差. * %*************************************************************clc;A=[1e-16, 1, 1; 2, 1, 2];B=[2, 1, 2; 1e-16, 1, 1]; % 增广矩阵% 方程组A% m = - a_{21}/a_{11} 是第2行加第1行的倍数% 消去a_{21}% m = - a_{12}/a_{22} 是第1行加第2行的倍数% 消去a_{12}, 系数矩阵成对角线% 未知数x1的值% 未知数x2的值disp(['方程组A的解: x1=',num2str(A(1,3)),', x2=',num2str(A(2,3))]); disp(' ');% 方程组B% m = - b_{21}/b_{11} 是第2行加第1行的倍数% 消去b_{21}% m = - b_{12}/b_{22} 是第1行加第2行的倍数% 消去b_{12}, 系数矩阵成对角线% 未知数x1的值% 未知数x2的值disp(['方程组B的解: x1=',num2str(B(1,3)),', x2=',num2str(B(2,3))]);《数值计算方法》实验2报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验2 非线性方程的迭代解法(之简单迭代法) 2. 实验题目用简单迭代法求方程010423=-+x x 在区间[1,2]的一个实根,取绝对误差限为410-.3. 实验目的掌握非线性方程的简单迭代法. 4. 基础理论简单迭代法:将方程0)(=x f 改写成等价形式)(x x ϕ=,从初值0x 开始,使用迭代公式)(1k k x x ϕ=+可以得到一个数列,若该数列收敛,则其极限即为原方程的解.取数列中适当的项可作为近似解. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程7. 结果与分析8. 附录:程序清单《数值计算方法》实验2报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验2 非线性方程的迭代解法(之Newton 迭代法) 2. 实验题目用Newton 迭代法求方程010423=-+x x 在区间[1,2]的一个实根,取绝对误差限为410-.3. 实验目的掌握求解非线性方程的Newton 迭代法. 4. 基础理论Newton 迭代法:解方程0)(=x f 的Newton 迭代公式为)(')(1k k k k x f x f x x -=+.5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程7. 结果与分析8. 附录:程序清单《数值计算方法》实验2报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验2 非线性方程的迭代解法(之对分区间法) 2. 实验题目用对分区间法求方程310x x --=在区间[1, 1.5]的一个实根,取绝对误差限为410-. 3. 实验目的掌握求解非线性方程的对分区间法. 4. 基础理论对分区间法:取[a ,b ]的中点p ,若f (p ) ≈ 0或b – a < ε,则p 为方程0)(=x f 的近似解;若f (a ) f (p ) < 0,则说明根在区间取[a ,p ]中;否则,根在区间取[p ,b ]中.将新的有根区间记为 [a 1,b 1],对该区间不断重复上述步骤,即可得到方程的近似根. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程用宏定义函数f (x );为了循环方便,得到的新的有根区间始终用[a ,b ]表示;由于新的有根区间可能仍以a 为左端点,这样会反复使用函数值f (a ),为减少运算次数,将这个函数值保存在一个变量fa 中;同样在判断新的有根区间时用到函数值f (p ),若新的有根区间以p 为左端点,则下一次用到的f (a )实际上就是现在的f (p ),为减少运算次数,将这个函数值保存在一个变量fp 中.算法的伪代码描述:Input :区间端点a ,b ;精度要求(即误差限)ε;函数f (x );最大对分次数N Output :近似解或失败信息7. 结果与分析8. 附录:程序清单说明: 源程序中带有数字的空行,对应着算法描述中的行号%**********************************************************%* 程序名:Bisection.m *%* 程序功能:使用二分法求解非线性方程. *%**********************************************************f=inline('x^3-x-1'); % 定义函数f(x)a=input('有根区间左端点: a=');b=input('右端点:b=');epsilon=input('误差限:epsilona=');N=input('最大对分次数: N=');1 % 对分次数计数器n置12 % 左端点的函数值给变量fafprintf('\n k p f(p) a(k) f(a(k))'); fprintf(' b(k) b-a\n');% 显示表头fprintf('%2d%36.6f%12.6f%12.6f%12.6f\n',0,a,fa,b,b-a);% 占2位其中0位小数显示步数0, 共12位其中小数6位显示各值3% while n≤ N 4 % 取区间中点p5% 求p 点函数值给变量fpfprintf('%2d%12.6f%12.6f',n,p,fp); % 输出迭代过程中的中点信息p 和f(p)6 % 如果f(p)=0或b-a 的一半小于误差限εfprintf('\n\n 近似解为:%f\n',p);% 则输出近似根p (7)return;% 并完毕程序 (7)89 % 计数器加110% 若f(a)与f(p)同号11% 则取右半区间为新的求根区间, 即a 取作p 12 % 保存新区间左端点的函数值 13% 否则14 % 左半区间为新的求根区间, 即b 取作p 15fprintf('%12.6f%12.6f%12.6f%12.6f\n',a,fa,b,b-a); %显示新区间端点与左端函数值、区间长度 16fprintf('\n\n 经过%d 次迭代后未达到精度要求.\n',N); % 输出错误信息(行17)《数值计算方法》实验2报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验2 非线性方程的迭代解法(之Aitken-Steffensen 加速法) 2. 实验题目用Aitken-Steffensen 加速法求方程010423=-+x x 在区间[1,2]的一个实根,取绝对误差限为410-.3. 实验目的熟悉求解非线性方程的Aitken-Steffensen 加速法. 4. 基础理论将方程0)(=x f 改写成等价形式)(x x ϕ=,得到从初值0x 开始的迭代公式)(1k k x x ϕ=+后,基于迭代公式)(1k k x x ϕ=+的Aitken-Steffensen 加速法是通过“迭代-再迭代-加速”完成迭代的,具体过程为kk k k k k k k k k k x y z z y x x y z x y +---===+2)(),(),(21ϕϕ. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程为了验证Aitken-Steffensen 加速法可以把一些不收敛的迭代加速成迭代收敛,我们使用将方程组变形为31021x x -=,取迭代函数31021)(x x -=ϕ,并利用宏定义出迭代函数.由于不用保存迭代过程,所以用x0表示初值同时也存放前一步迭代的值,y 和z 是迭代过程中产生的y k 和z k ,x 存放新迭代的结果.算法的伪代码描述:Input :初值x 0;精度要求(即误差限)ε;迭代函数φ(x );最大迭代次数N7. 结果与分析8. 附录:程序清单%************************************************************* %* 程 序 名:Aitken_Steffensen.m * %* 程序功能:用Aitken-Steffensen 加速法求方程. * %************************************************************* clc;clear all;phi=inline('0.5 * sqrt( 10 - x^3)'); % 迭代函数x0=input('初值: x0 = ');epsilon=input('误差限: epsilon='); N=input('最大迭代次数: N=');disp(' n 迭代中间值y(n-1) 再迭代结构z(n-1) 加速后的近似值x(n)'); fprintf('%2d%54.6f\n',0,x0);% 占2位整数显示步数0, 为了对齐, 占54位小数6位显示x01 % n 是计数器2 % while n<=Ny= 3 ; % 迭代 z= 3 ; % 再迭代 x= 3 ; % 加速% x0初值与前一步的近似值, y 和z 是中间变量, x 是下一步的近似值fprintf('%2d%18.6f%18.6f%18.6f\n',n,y,z,x);%显示中间值和迭代近似值6 % 如果与上一步近似解差的绝对值不超过误差限 fprintf('\n\n 近似解 x≈x(%d)≈%f \n',n,x);% 则输出近似根 (7), 可简略为: fprintf('\n\n 近似解 x=%f',x); return; % 并完毕程序(7) 8 % 相当于endif9 % 计数器加110 % 新近似值x 作为下一次迭代的初值 11fprintf('\n 迭代%d 次还不满足误差要求.\n\n',N); %输出错误信息(12)《数值计算方法》实验2报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验2 非线性方程的迭代解法(之Newton 下山法) 2. 实验题目用Newton 下山法求方程010423=-+x x 在区间[1,2]的一个实根,取绝对误差限为410-.3. 实验目的熟悉非线性方程的Newton 下山法. 4. 基础理论Newton 下山法:Newton 下山法公式为)(')(1k k kk k x f x f x x λ-=+,使|)(||)(|1k k x f x f <+,其中10≤<k λ.5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程定义函数f(x)和df(x),其中df(x)是f(x)的导函数.每步迭代时先取下山因子为1,尝试迭代,判断尝试结果是否满足下山因子,若满足则作为这步的迭代结果;否则将下山因子减半,然后再尝试.为防止当前的x k 是极小值点,附近不会有满足下述条件的其它点,使尝试陷入死循环,同时计算机中能表示出的浮点数也有下界,因此我们设置了最大尝试次数.当超过最大尝试次数时,不再进行下山尝试.由于反复尝试迭代且要判断下山条件,所以f (x 0)和f ‘(x 0)会反复使用,为避免重复计算浪费运行时间,将这两个值分别保存在变量fx0和dfx0.而尝试产生的节点,判断下山条件时要用到它的函数值,若尝试成功,这个点会作为下一步的初值再使用,所以把该点的函数值也保存在变量fx 中.算法的伪代码描述:Input :初值x 0;精度要求(即误差限)ε;函数与其导函数f (x )和f’(x);最大迭代次数N ;K 下山尝试最大次数Output :近似解或失败信息7. 结果与分析8. 附录:程序清单%*************************************************************%* 程序名:NewtonDownhill.m *%* 程序功能:用Newton下山法求解非线性方程. *%*************************************************************clc;clear all;f=inline('x^3-x-1'); % 函数f(x)df=inline('3*x^2-1'); % 函数f(x)的导函数x0=input('初值: x0 = ');epsilon=input('误差限: epsilon=');N=input('最大迭代次数: N=');K=input('最大下山尝试次数: K=');1 % 迭代次数计数器2 % 存x0点函数值fprintf('\n\n n x(n) f(x(n))\n'); % 显示表头fprintf('%2d%14.6f%14.6f\n',0,x0,fx0); % 2位整数显示0, 共14位小数6位显示x0和fx03 % while n≤ Ndisp(''); % 换行显示下山尝试过程的表头disp(' 下山因子尝试x(n) 对应f(x(n)) 满足下山条件');disp('');4 % 存x0点导数值, 每次下山尝试不用重新计算ifdfx0==0 % 导数为0不能迭代disp(‘无法进行Newton迭代’);return;endlambda=1.0; % 下山因子从1开始尝试k=1; % k下山尝试次数计数器while k<=K % 下山最多尝试K次% 下山公式fx=f(x); % 函数值fprintf('%22.6f%14.6f%14.6f',lambda,x,fx); % 显示尝试结果if (abs(fx)<abs(fx0)) % 判断是否满足下山条件fprintf(' 满足\n');break; % 是, 则退出下山尝试的循环elsefprintf(' 不满足\n');endlambda=lambda/2; % 不是, 则下山因子减半k=k+1; % 计数器加1endif k>Kfprintf('\n 下山条件无法满足, 迭代失败.\n\n');return;endfprintf('%2d%14.6f%14.6f\n',n,x,fx);% 2位整数显示步数n, 共14位小数6位显示下步迭代结果22 % 达到精度要求否fprintf('\n\n 方程的近似解为: x≈%f\n\n',x); % (23)return; % 达到, 则显示结果并完毕程序(23) end % (24)% 用x0,fx0存放前一步的近似值和它的函数值, 进行循环迭代25262728fprintf('\n 迭代%d次还不满足误差要求.\n\n',N);《数值计算方法》实验2报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验2 非线性方程的迭代解法(之弦截法) 2. 实验题目用弦截法求方程010423=-+x x 在区间[1,2]的一个实根,取绝对误差限为410-. 3. 实验目的熟悉非线性方程的弦截法. 4. 基础理论将Newton 迭代法中的导数用差商代替,得到弦截法(或叫正割法)公式)()()(111k k k k k k k x f x f x f x x x x --+---=.5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程不保存迭代过程,所以始终以x 0和x 1分别存放x k -1和x k ,而x 存放新产生的迭代值x k +1,这样,下一次迭代时需要把上一步的x 1(即x k )赋值于x 0(做新的x k -1).这些点的函数值会重复用到,在迭代公式中也要用到,上一步的x 1作为下一步的x 0也会再一次用它的函数值,为减少重新计算该点函数值的运行时间,将x 1点的函数值保存在变量fx1中.算法的伪代码描述:Input :初值x 0,x 1;精度要求(即误差限)ε;函数f (x );最大迭代次数N7. 结果与分析8. 附录:程序清单%*************************************************************%* 程序名:SecantMethod.m *%* 程序功能:用弦截法求解非线性方程. *%*************************************************************clc;clear all;f=inline('2*x^3-5*x-1'); % 函数f(x)x0=input('第一初值: x0 = ');x1=input('第二初值: x1 = ');epsilon=input('误差限: epsilon=');N=input('最大迭代次数: N=');fprintf('\n n x(n)\n'); % 显示表头fprintf('%2d%14.6f\n', 0, x0); % 占2位显示步数0, 共14位其中小数6位显示x0fprintf('%2d%14.6f\n', 1, x1); % 占2位显示步数1, 共14位其中小数6位显示x11 % 存x0点函数值2 % 存x1点函数值3 % 迭代计数器4 % while n≤ N% 弦截法公式fprintf('%2d%14.6f\n', n, x); %显示迭代过程6 % 达到精度要求否fprintf('\n\n 方程的近似解为: x≈%f\n\n', x);return; % 达到, 则显示结果并完毕程序89 % 原x1做x0为前两步的近似值10 % 现x做x1为一两步的近似值11 % x0点函数值12 % 计算x1点函数值, 为下一次循环13 % 计数器加1 14fprintf('\n 迭代%d 次还不满足误差要求.\n\n',N);《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验3 解线性方程组的直接法(之Gauss 消去法) 2. 实验题目用Gauss 消去法求解线性方程组⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--000.3000.2000.1643.5072.1000.2623.4712.3000.1000.3000.2001.0321x x x . 3. 实验目的掌握解线性方程组的Gauss 消去法. 4. 基础理论Gauss 消去法是通过对增广矩阵的初等行变换,将方程组变成上三角方程组,然后通过回代,从后到前依次求出各未知数.Gauss 消去法的第k 步(1≤k≤n -1)消元:若0≠kk a ,则依次将增广矩阵第k 行的kk ik a a /-倍加到第i 行(k+1≤i≤n),将第k 列对角线下的元素都化成0.5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程7. 结果与分析8. 附录:程序清单《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验3 解线性方程组的直接法(之Gauss 列主元消去法) 2. 实验题目用Gauss 列主元消去法求解线性方程组⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--000.3000.2000.1643.5072.1000.2623.4712.3000.1000.3000.2001.0321x x x . 3. 实验目的掌握解线性方程组的Gauss 列主元消去法. 4. 基础理论Gauss 列主元消去法也是通过对增广矩阵的初等行变换,将方程组变成上三角方程组,然后通过回代,从后到前依次求出各未知数.Gauss 列主元消去法的第k 步(1≤k≤n -1)消元:先在nk k k kk a a a ,,,,1 +中找绝对值最大的,将它所在的行与第k 行交换,然后将第k 行的kk ik a a /-倍加到第i 行(k+1≤i≤n),将第k 列对角线下的元素都化成0. 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程7. 结果与分析8. 附录:程序清单《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验3 解线性方程组的直接法(之Doolittle 分解) 2. 实验题目对矩阵A 进行Doolittle 分解,其中⎪⎪⎪⎪⎪⎭⎫⎝⎛----=3101141101421126A .3. 实验目的掌握矩阵的Doolittle 分解. 4. 基础理论矩阵的Doolittle 分解是指将矩阵n n ij a A ⨯=)(可以分解为一个单位下三角矩阵和一个上三角矩阵的乘积.若设⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=nn n n n n n n u u u u u u u u u u U l l ll l l L000000,1010010001333223221131211321323121则可依如下顺序公式计算⎪⎪⎩⎪⎪⎨⎧++=-=+=-=∑∑-=-=1111,,2,1,/)(,,1,,k t kk tk it ik ik k r rj kr kj kj nk k i u u l a l nk k j u l a u其中k = 1,2,…,n .5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程(1)按计算公式依次计算一行u 同时计算一列l ;(2)因为计算完u ij (或l ij )后,a ij 就不再使用,为节省存储空间,将计算的u ij (和l ij )仍存放在矩阵A 中的相应位置;(3)使用L 矩阵和U 矩阵时需要根据元素所在位置取固定值或A 中相应位置的值.L 对角线上的元素为1,上三角部分为0,下三角部分为A 中对应的元素;U 的下三角部分为0,上三角部分为A 中对应的元素.算法的伪代码描述: Input :阶数n ;矩阵A7. 结果与分析8. 附录:程序清单%****************************************************% 程序名: Doolittle.m *% 程序功能: 矩阵LU分解中的Doolittle分解. *%****************************************************clc;clear all;n=4; % 矩阵阶数A=[6 2 1 -1;2 4 1 0; 1 1 4 -1; -1 0 -1 3]disp('A=');disp(A);% LU分解(Doolittle分解)for k=1:n% 计算矩阵U的元素u_{kj}% (可参照下面l_{ik}的公式填写)% 计算矩阵L的元素l_{ik}% L 在A 下三角, U 在上三角(对角线为1) enddisp('分解结果:'); disp('L='); for i=1:n for j=1:nif i>j % 在下三角部分, 则取A 对于的元素显示 fprintf(' %8.4f',A(i,j));elseif i==j % 在对角线上, 则显示1 fprintf(' %8d',1);else % 在上三角部分, 则显示0 fprintf(' %8d',0); end endfprintf('\n'); % 换行 enddisp('U='); for i=1:n for j=1:nif i<=j % 在上三角部分或对角线上, 则取A 对于的元素显示 fprintf(' %8.4f',A(i,j));else % 在下三角部分, 则显示0 fprintf(' %8d',0); end endfprintf('\n'); % 换行 end《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验3 解线性方程组的直接法(之LU 分解法) 2. 实验题目用LU 分解(Doolittle 分解)法求解线性方程组⎪⎩⎪⎨⎧=++=++=++104615631552162321321321x x x x x x x x x 3. 实验目的熟悉解线性方程组LU 分解法.4. 基础理论若将矩阵A 进行了Doolittle 分解,A = LU ,则解方程组b x A=可以分解求解两个三角方程组b y L=和y x U =.它们都可直接代入求解,其中b y L=的代入公式为∑-==-=11,,2,1,k j j kj k k n k y l b y而y x U=的代入公式为∑+=-=-=nk j kk j kjk k n n k u x uy x 11,,1,,/)( .5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程(1)Doolittle 分解过程依次计算一行u 同时计算一列l 完成,并将计算的u ij (和l ij )仍存放在矩阵A 中的相应位置;(2)求解方程组的代入公式中用到的u ij 和l ij 都直接在A 的相应位置取值即可. 算法的伪代码描述:Input :阶数n ;矩阵A ;常数项向量b7. 结果与分析8. 附录:程序清单%**************************************************** % 程序名: LinearSystemByLU.m *% 程序功能: 利用LU分解(Doolittle分解)解方程组. *%****************************************************clc;clear all;n=3; % 矩阵阶数A=[1 2 6; 2 5 15; 6 15 46];b=[1;3;10];% LU分解(Doolittle分解)for k=1:n% 计算矩阵U的元素u_{kj}% (可参照下面l_{ik}的公式填写)% 计算矩阵L的元素l_{ik}% L在A下三角, U在上三角(对角线为1) endfor k=1:n % 用代入法求解下三角方程组Ly=by(k)=b(k);3 %∑-==-=11,,2,1,kjj kjk knkylby33enddisp('方程组Ly=b的解:y=');disp(y');for k=n:-1:1 % 回代求解上三角方程组Ux=y x(k)=y(k);6 %∑+=-=-=nkjj kjk knnkxuyx11,,1,,666 enddisp('原方程组的解:x='); disp(x');《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X成绩:1. 实验名称实验3 解线性方程组的直接法(之Cholesky 分解) 2. 实验题目对矩阵A 进行Cholesky 分解,其中⎪⎪⎪⎪⎪⎭⎫⎝⎛----=3101141101421126A . 3. 实验目的理解矩阵的Cholesky 分解. 4. 基础理论矩阵的Cholesky 分解是指将矩阵n n ij a A ⨯=)(可以分解为一个下三角矩阵L 和L 转置的乘积,即A =LL T,其中L 各元素可依如下顺序公式计算⎪⎪⎩⎪⎪⎨⎧++=-=-=∑∑-=-=11112,,2,1,/)(k t kktk it ik ik k r kr kk kk nk k i l l l a l l a l其中k = 1,2,…,n .5. 实验环境操作系统:Windows xp ; 程序设计语言:VC++ 6. 实验过程(1)按计算公式依次先计算一列对角线上的元素l kk ,再计算这列其他元素l ik ,且对称位置的元素也取同一个值;(2)因为计算完l ij 后,a ij 就不再使用,为节省存储空间,将计算的l ij 仍存放在矩阵A 中的相应位置;(3)使用L 矩阵时需要根据元素所在位置取固定值或A 中相应位置的值.L 上三角部分为0,对角线和下三角部分为A 中对应的元素.算法的伪代码描述:Input :阶数n ;矩阵AOutput :矩阵L (合并存储在数组A 中)行号 伪代码注释1 for k ← 1 to n2∑-=-=112k r krkk kk l a l3 for i ← k to n4 ∑-=-=11/)(k t kk tk it ik ik l l l a l计算结果存放在a ij5 endfor6 endfor7return L输出L7. 结果与分析8. 附录:程序清单%************************************************************* %* 程 序 名:Cholesky.m * %* 程序功能:对称正定矩阵的Cholesky 分解. * %*************************************************************n=4; % 矩阵阶数 A=[6,2,1,-1; 2,4,1,0; 1,1,4,-1; -1,0,-1,3];disp('A ='); for i=1:n for j=1:nfprintf('%10.4f',A(i,j)); % 共占14位endfprintf('\n');% 一行完毕换行end% Cholesky 分解 for k=1:n % 计算对角线上的l _{kk}% 计算其他的l _{ik} % 和l _{ki}end % L 在A 下三角, L^T 在上三角disp('分解结果:'); disp('L='); for i=1:n for j=1:n if i>=j % 在下三角部分或对角线上, 则取A 对于的元素显示fprintf('%10.4f',A(i,j));else % 在上三角部分, 则显示0 fprintf('%10d',0); end endfprintf('\n'); % 换行 end《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X成绩:1. 实验名称实验3 解线性方程组的直接法(之改进的Cholesky 分解) 2. 实验题目对矩阵A 进行改进的Cholesky 分解,其中⎪⎪⎪⎪⎪⎭⎫⎝⎛----=3101141101421126A .3. 实验目的理解矩阵改进的Cholesky 分解. 4. 基础理论矩阵的改进的Cholesky 分解是指将矩阵n n ij a A ⨯=)(可以分解为一个单位下三角矩阵L 和对角矩阵D 与L 转置的乘积,即A =LDL T,其中L 和D 各元素可依如下顺序公式计算⎪⎪⎩⎪⎪⎨⎧++=-=-=∑∑-=-=11112,,2,1,/)(k t k kt it t ik ik k r kr r kk k nk k i d l l d a l l d a d其中k = 1,2,…,n .5. 实验环境操作系统:Windows xp ; 程序设计语言:VC++ 6. 实验过程(1)按计算公式依次先计算D 的一个元素d k ,再计算L 中这列的元素l ik ,且对称位置的元素也取同一个值;(2)因为计算完d k 和l ij 后,a kk 或a ij 就不再使用,为节省存储空间,将计算的a kk 或l ij 仍存放在矩阵A 中的相应位置;(3)使用L 矩阵时需要根据元素所在位置取固定值或A 中相应位置的值.L 对角线和上三角部分为0,下三角部分为A 中对应的元素;D 对角线为A 中对应的元素,其余都是0.算法的伪代码描述: Input :阶数n ;矩阵AOutput :矩阵L (合并存储在数组A 中)7. 结果与分析8. 附录:程序清单%************************************************************* %* 程 序 名:ImprovedCholesky.m * %* 程序功能:对称正定矩阵的改进的Cholesky 分解. * %*************************************************************n=4; % 矩阵阶数A=[6,2,1,-1; 2,4,1,0; 1,1,4,-1; -1,0,-1,3];disp('A =');for i=1:nfor j=1:nfprintf('%10.4f',A(i,j)); % 共占14位endfprintf('\n'); % 一行完毕换行end% Cholesky分解for k=1:n% 计算D对角线上的u_{kk}% 计算L的元素l_{ik}% 和L转置的元素l_{ki} end % L在A下三角, D在对角线disp('分解结果:');disp('L=');for i=1:nfor j=1:nif i>j % 在下三角部分, 则取A对于的元素显示fprintf('%10.4f',A(i,j));elseif i==j % 在对角线上, 则显示1fprintf('%10d',1);else % 在上三角部分, 则显示0fprintf('%10d',0);endendfprintf('\n'); % 换行enddisp('D='); for i=1:n for j=1:n if i==j % 在对角线上, 则取A 对于的元素显示fprintf('%10.4f',A(i,j));else % 其余显示0fprintf('%10d',0); end endfprintf('\n'); % 换行 end《数值计算方法》实验3报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验3 解线性方程组的直接法(之追赶法) 2. 实验题目用追赶法求解线性方程组⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-----101053001210023100124321x x x x 3. 实验目的熟悉解线性方程组的追赶法. 4. 基础理论对于系数矩阵为三对角矩阵的方程组,其Crout 分解可分解为⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=------11111211122111122211n n nn n n nn n n t t t s a s a s a s b a c b a c b a c b A这样,解方程组可以由如下2步完成:“追”:,,,3,2,/)(,,/,/,1111111111n i s y a f y t a b s s c t s f y b s i i i i i i i i i i i i =-=-====-----其中:Tn f f ),,(1 为方程组的常数项,n t 没用;“赶”:.1,,2,1,,1 --=-==+n n i x t y x y x i i i i n n5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程在“追”的过程中,向量s 和y 都有n 个元素,t 只有n -1个元素,又1s 和1y 的计算公式与其它i s 和i y 不同,所以先单独计算1s 和1y ,然后在一个n -1次循环中,求其它i s 和i y 以与i t .由于在“追”的过程中,i b ,i c 和i f 在分别计算完对应的i s ,i t 和i y 后就不再使用,所以借用数组b ,c 和f 存储向量s ,t 和y ;同样在“赶”的过程中,i y 在计算完对应的i x 后就不再使用,所以再一次借用数组f 存储向量x .追赶法算法的伪代码描述:Input :阶数n ;三对角矩阵的三条对角线向量a ,b ,c ,常数项向量f Output :方程组的解x改进的追赶法算法的伪代码描述:Input :阶数n ;三对角矩阵的三条对角线向量a ,b ,c ,常数项向量f Output :方程组的解x7. 结果与分析8. 附录:程序清单%*************************************************************%* 程序名:ChaseAfter.m *%* 程序功能:用追赶法求解三对角线性方程组. *%*************************************************************clc;clear all;n=4;a=[0,-1,-1,-3];b=[2, 3, 2, 5];c=[-1, -2, -1, 0];f=[0, 1, 0, 1];% "追"s(1) = b(1);y(1) = f(1); % 先单独求s_1和y_1 for k = 1 : n-1% 再求t_i(i=1,2,…,n-1)% s_i(i=2,3,…,n)% y_i(i=2,3,…,n)end% "赶"x(n) = y(n); % 先单独求x_nfor k = n-1 : -1 : 1% 再求x_i(i=n-1,n-2, (1)endx=x' % 输出解向量-------------------------------------------------------------------------------------------------------------------改进的程序:%*************************************************************%* 程序名:ChaseAfter.m *%* 程序功能:用追赶法求解三对角线性方程组. *%*************************************************************clc;clear all;n=4;a=[0,-1,-1,-3];b=[2, 3, 2, 5];c=[-1, -2, -1, 0];f=[0, 1, 0, 1];% "追"% b(1)=b(1); % s_1仍在b_1中,不用重新计算y(1)=f(1)/b(1); % 先单独y_1for k=1:n-1% 再求t_i(i=1,2,…,n-1)% s_i(i=2,3,…,n)% y_i(i=2,3,…,n)end% "赶"% f(n)=f(n); % x_n等于y_n仍在f_n中for k=n-1:-1:1% 再求x_i(i=n-1,n-2, (1)endx=f' % 输出解向量《数值计算方法》实验4报告班级:20##级####x班学号:20##2409####:##X 成绩:1. 实验名称实验4 解线性方程组的迭代法(之Jacobi迭代)2. 实验题目用Jacobi迭代法求解线性方程组1231231232251223x x x x x x x x x +-=⎧⎪++=⎪⎨++=⎪⎪⎩任取3. 实验目的掌握解线性方程组的Jacobi 迭代法. 4. 基础理论将第i (n i ≤≤1)个方程i n in i i b x a x a x a =+++ 2211移项后得到等价方程ii n in i i i i i i i i i a x a x a x a x a b x /)(11,11,11------=++--便可构造出Jacobi 迭代公式,1,0,/)()()(11,)(11,)(11)1(=------=++--+k a x a x a x a x a b x ii k n in k i i i k i i i k i i k i . 5. 实验环境操作系统:Windows xp ; 程序设计语言:Matlab 6. 实验过程7. 结果与分析8. 附录:程序清单《数值计算方法》实验4报告班级: 20##级####x 班 学号: 20##2409#### : ##X 成绩:1. 实验名称实验4 解线性方程组的迭代法(之Gauss-Seidel 迭代) 2. 实验题目用Gauss-Seidel 迭代法求解线性方程组。
数值分析的实验报告
数值分析的实验报告数值分析的实验报告导言数值分析是一门研究数值计算方法和数值计算误差的学科,它在科学计算、工程技术和社会经济等领域具有广泛的应用。
本实验旨在通过对数值分析方法的实际应用,验证其有效性和可靠性。
实验一:方程求根方程求根是数值分析中的基础问题之一。
我们选取了一个非线性方程进行求解。
首先,我们使用二分法进行求解。
通过多次迭代,我们得到了方程的一个近似解。
然后,我们使用牛顿法进行求解。
与二分法相比,牛顿法的收敛速度更快,但需要选择一个初始点。
通过比较两种方法的结果,我们验证了牛顿法的高效性。
实验二:插值与拟合插值与拟合是数值分析中常用的数据处理方法。
我们选取了一组实验数据,通过拉格朗日插值法和最小二乘法进行插值和拟合。
通过对比两种方法的拟合效果,我们验证了最小二乘法在处理含有噪声数据时的优势。
同时,我们还讨论了插值和拟合的精度与样本点数量之间的关系。
实验三:数值积分数值积分是数值分析中的重要内容之一。
我们选取了一个定积分进行计算。
首先,我们使用复化梯形公式进行积分计算。
通过增加分割区间的数量,我们得到了更精确的结果。
然后,我们使用复化辛普森公式进行积分计算。
与复化梯形公式相比,复化辛普森公式具有更高的精度。
通过比较两种方法的结果,我们验证了复化辛普森公式的优越性。
实验四:常微分方程数值解常微分方程数值解是数值分析中的重要应用之一。
我们选取了一个常微分方程进行数值解的计算。
首先,我们使用欧拉方法进行数值解的计算。
然后,我们使用改进的欧拉方法进行数值解的计算。
通过比较两种方法的结果,我们验证了改进的欧拉方法的更高精度和更好的稳定性。
实验五:线性方程组的数值解法线性方程组的数值解法是数值分析中的重要内容之一。
我们选取了一个线性方程组进行数值解的计算。
首先,我们使用高斯消元法进行数值解的计算。
然后,我们使用追赶法进行数值解的计算。
通过比较两种方法的结果,我们验证了追赶法在求解三对角线性方程组时的高效性。
《数值分析》课程实验报告范文
《数值分析》课程实验报告范文《数值分析》课程实验报告姓名:学号:学院:机电学院日期:2022年某月某日目录实验一函数插值方法1实验二函数逼近与曲线拟合5实验三数值积分与数值微分7实验四线方程组的直接解法9实验五解线性方程组的迭代法15实验六非线性方程求根19实验七矩阵特征值问题计算21实验八常微分方程初值问题数值解法24实验一函数插值方法一、问题提出对于给定的一元函数的n+1个节点值。
试用Lagrange公式求其插值多项式或分段二次Lagrange插值多项式。
实验二函数逼近与曲线拟合一、问题提出从随机的数据中找出其规律性,给出其近似表达式的问题,在生产实践和科学实验中大量存在,通常利用数据的最小二乘法求得拟合曲线。
在某冶炼过程中,根据统计数据的含碳量与时间关系,试求含碳量与时间t的拟合曲线。
t(分)051015202530354045505501.272.162.863.443.874.154.374.51 4.584.024.64二、要求1、用最小二乘法进行曲线拟合;2、近似解析表达式为;3、打印出拟合函数,并打印出与的误差,;4、另外选取一个近似表达式,尝试拟合效果的比较;5、某绘制出曲线拟合图。
三、目的和意义1、掌握曲线拟合的最小二乘法;2、最小二乘法亦可用于解超定线代数方程组;3、探索拟合函数的选择与拟合精度间的关系四、实验步骤:第一步先写出线性最小二乘法的M文件functionc=lpoly(某,y,m)n=length(某);b=zero(1:m+1);f=zero(n,m+1); fork=1:m+1f(:,k)=某.^(k-1);enda=f'某f;b=f'某y';c=a\b;c=flipud(c);第二步在命令窗口输入:>>lpoly([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)回车得到:an=-0.00240.20370.2305即所求的拟合曲线为y=-0.0024某2+0.2037某+0.2305在编辑窗口输入如下命令:>>某=[0,5,10,15,20,25,30,35,40,45,50,55];>>y=-0.0024某某.^2+0.2037某某+0.2305;>>plot(某,y)命令执行得到如下图五、实验结论分析复杂实验数据时,常采用分段曲线拟合方法。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
《数值分析》 课程实验指导书
实验一 函数插值方法
一、问题提出
对于给定的一元函数)(x f y =的n+1个节点值(),0,1,
,j j y f x j n ==。
试用
Lagrange 公式求其插值多项式或分段二次Lagrange 插值多项式。
数据如下:
(1)
求五次Lagrange 多项式5L ()x ,和分段三次插值多项式,计算(0.596)f ,(0.99)f 的值。
(提示:结果为(0.596)0.625732f ≈, (0.99) 1.05423f ≈ )
(2)
试构造Lagrange 多项式6L ()x ,计算的(1.8)f ,(6.15)f 值。
(提示:结果为(1.8)0.164762f ≈, (6.15)0.001266f ≈ )
二、要求
1、 利用Lagrange 插值公式
00,()n n
i n k k i i k k i x x L x y x x ==≠⎛⎫-= ⎪-⎝⎭∑∏编写出插值多项式程序; 2、 给出插值多项式或分段三次插值多项式的表达式;
3、 根据节点选取原则,对问题(2)用三点插值或二点插值,其结果如何;
4、 对此插值问题用Newton 插值多项式其结果如何。
四、实验分析:
Lagrange 插值多项式的表达式:
1,,2,1,)()()(,
)()(1111+=--==∏∑+≠=+=n i x x x x x l x l y x L n i j j j i j i n i i i 。
其中)(x l i 被称为插值基函数,实际上是一个n 次多项式。
)(x l i 的这种表示具有较好的对称性。
公式具有两大优点:(1)求插值多项式,不需要求解线性方程组,当已知数据点较多时,此公式更能显示出优越性。
(2)函数值可以用符号形式表示,数据点未确定的纵坐标可用多项式表示。
Newton 插值多项式如下:
1
0010,()()[,,]()k n n j k k j j k N x f x f x x x x -==≠=+•-∑∏ 其中: 00,0()()[,,]k i k i i j j j i k f x x x f x x ==≠-=∑∏
Newton 插值多项式的优点是:当每增加一个节点时,只增加一项多项式。
三、实验程序及注释
1、m 程序:
function [c,l]=lagran(x,y)
% x 为n 个节点的横坐标组成的向量,y 为纵坐标所组成的向量
% c 为所得插值函数的系数所组成的向量
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;
function fi=Lagran_(x,f,xi)
fi=zeros(size(xi));
n=length(f);
for i=1:n
z=ones(size(xi));
for j=1:n
if i~=j,z=z.*(xi-x(j))/(x(i)-x(j));
end
end
fi=fi+z*f(i);
end
Return
2、主程序:
对于数据(1),在主窗口输入:
x=[0.4 0.55 0.65 0.80 0.95 1.05]; y=[0.41075 0.57815 0.68675 0.90 1.00 1.25382];
lagran (x,y)
ans =
143.8486 -506.0836 694.2890 -463.5132 151.5174 -18.9859
format long
xi=[0.596,0.990]; yi=Lagran_(x,y,xi)
yi =
0.61916350788547 1.05159915747938
xx=0.3:0.01:1.1
[c,l]=lagran(x,y)
yy= polyval(c,xx)
plot(xx,yy, x,y,'o');
得出在拉格朗日的插值多项式曲线图,如图1-1所示:
图1-1
由此得到所求的lagrange 多项式为:
9859.185174.1515132.4632890.6940836.5068486.143L 23455
-+-+-=x x x x x x )( 对于数据(2),在主窗口输入:
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
xi=[1.8,6.15]; yi=Lagran_(x,y,xi)
yi =
0.16476 0.0012658
xx=0:7
[c,l]=lagran(x,y)
yy= polyval(c,xx)
plot(xx,yy, x,y,'o');
得出在拉格朗日的插值多项式曲线图,如图1-2所示:
图1-2
所以得到所求的lagrange 多项式为:
9950.09683.04419.01175.00186.00016.0001.0L 234566--+-+-=x x x x x x x )
(对于要求(3),在主窗口输入:
x=[0.4,0.65,0.95];
y=[0.41075,0.69675,1.00];
[c,l]=lagran(x,y)
c =
-0.24212121212121 1.39822727272727 -0.10980151515152
xx=04:0.01:1.0
yy= polyval(c,xx)
plot(xx,yy, x,y,'o');
得出在拉格朗日的插值多项式曲线图,如图1-3所示:
图1-3
xi=[0.596,0.990]; yi=Lagran_(x,y,xi)
yi =
0.63753661090909 1.03714048484848
4.1基于Newton插值多项式对插值问题的分析。
function[c,d]=newploy(x,y)
% 这里x为n个节点的横坐标所组成的向量,y为纵坐标所组成的向量% c为所求的Newton插值多项式的系数构成的向量
n=length(x);
d=zeros(n,n);
d(:,1)=y';
for j=2:n
for k=j:n
d(k,j)=(d(k,j-1)-d(k-1,j-1))/(x(k)-x(k-j+1));
end
end
c=d(n,n);
for k=(n-1):-1:1
c=conv(c,poly(x(k)));
m=length(c);
c(m)=c(m)+d(k,k);
end
x=[0.4 0.55 0.65 0.80 0.95 1.05];
y=[0.41075 0.57815 0.68675 0.90 1.00 1.25382];
>> newploy(x,y)
ans =
143.8486 -506.0836 694.2890 -463.5132 151.5174 -18.9859
三、目的和意义
1、学会常用的插值方法,求函数的近似表达式,以解决其它实际问题;
2、明确插值多项式和分段插值多项式各自的优缺点;
3、熟悉插值方法的程序编制;
4、如果绘出插值函数的曲线,观察其光滑性。
五、实验步骤:
1.进入C或matlab开发环境;
2.根据实验内容和要求编写程序;
3.调试程序;
4.运行程序;
5.撰写报告,讨论分析实验结果.。