哈工大数值分析报告上机实验报告

合集下载

数值分析2024上机实验报告

数值分析2024上机实验报告

数值分析2024上机实验报告数值分析是计算数学的一个重要分支,它研究如何用数值方法来解决数学问题。

在数值分析的学习过程中,学生需要通过上机实验来巩固理论知识,并学会使用相应的数值方法来解决实际问题。

本篇报告将详细介绍2024年度数值分析上机实验的内容和结果。

一、实验内容2024年度数值分析上机实验分为四个部分,分别是:方程求根、插值与拟合、数值积分和常微分方程的数值解。

1.方程求根这部分实验要求使用数值方法求解给定的非线性方程的根。

常见的数值方法有二分法、牛顿法、割线法等。

在实验过程中,我们需要熟悉这些数值方法的原理和实现步骤,并对不同方法的收敛性进行分析和比较。

2.插值与拟合这部分实验要求使用插值和拟合方法对给定的一组数据进行拟合。

插值方法包括拉格朗日插值、牛顿插值等;拟合方法包括最小二乘拟合、多项式拟合等。

在实验中,我们需要熟悉插值和拟合方法的原理和实现步骤,并对不同方法的精度和稳定性进行比较。

3.数值积分这部分实验要求使用数值方法计算给定函数的积分。

常见的数值积分方法有梯形法则、辛普森法则、龙贝格积分等。

在实验过程中,我们需要熟悉这些数值积分方法的原理和实现步骤,并对不同方法的精度和效率进行比较。

4.常微分方程的数值解这部分实验要求使用数值方法求解给定的常微分方程初值问题。

常见的数值方法有欧拉法、改进的欧拉法、四阶龙格-库塔法等。

在实验中,我们需要熟悉这些数值解方法的原理和实现步骤,并对不同方法的精度和稳定性进行比较。

二、实验结果在完成2024年度数值分析上机实验后,我们得到了以下实验结果:1.方程求根我们实现了二分法、牛顿法和割线法,并对比了它们的收敛速度和稳定性。

结果表明,割线法的收敛速度最快,但在一些情况下可能会出现振荡;二分法和牛顿法的收敛速度相对较慢,但稳定性较好。

2.插值与拟合我们实现了拉格朗日插值和最小二乘拟合,并对比了它们的拟合效果和精度。

结果表明,拉格朗日插值在小区间上拟合效果较好,但在大区间上可能出现振荡;最小二乘拟合在整体上拟合效果较好,但可能出现过拟合。

数值分析上机实习报告

数值分析上机实习报告

指导教师:姓名:学号:专业:联系电话:上海交通大学目录序言 (3)实验课题(一) 雅可比迭代法和高斯-塞得尔迭代法的收敛性和收敛速度 (4)数值分析 (6)实验课题(二) 松弛因子对SOR法收敛速度的影响 (6)数值分析 (12)总结 (13)附录(程序清单) (14)1.雅可比迭代法和高斯-塞得尔迭代法的收敛性和收敛速度 (14)雅可比迭代法: (14)高斯-塞得尔迭代法: (16)2.松弛因子对SOR法收敛速度的影响 (18)松弛法(SOR) (18)序言随着科学技术的发展,提出了大量复杂的数值计算问题,在实际解决这些计算问题的长期过程中,形成了计算方法这门学科,专门研究各种数学问题的数值解法(近似解法),包括方法的构造和求解过程的误差分析,是一门内容丰富,有自身理论体系的实用性很强的学科。

解决工程问题,往往需要处理很多数学模型,这就要花费大量的人力和时间,但是还有不少数学模型无法用解析法得到解。

使用数值方法并利用计算机,就可以克服这些困难。

事实上,科学计算已经与理论分析、科学实验成为平行的研究和解决科技问题的科学手段,经常被科技工作者所采用。

作为科学计算的核心内容——数值分析(数值计算方法),已逐渐成为广大科技工作者必备的基本知识并越来越被人重视。

由于数值方法是解数值问题的系列计算公式,所以数值方法是否有效,不但与方法本身的好坏有关,而且与数值问题本身的好坏也有关,因此,研究数值方法时,不但需要研究数值方法的好坏,即数值稳定性问题,而且还需要研究数值问题本身的好坏,即数值问题的性态,以及它们的判别问题。

数值计算的绝大部分方法都具有近似性,而其理论又具有严密的科学性,方法的近似值正是建立在理论的严密性基础上,根据计算方法的这一特点。

因此不仅要求掌握和使用算法,还要重视必要的误差分析,以保证计算结果的可靠性。

数值计算还具有应用性强的特点,计算方法的绝大部分方法如求微分方程近似解,求积分近似值,求解超越方程,解线性方程组等都具有较强的实用性,而插值法,最小二乘法,样条函数等也都是工程技术领域中常用的,有实际应用价值的方法。

数值分析上机实验报告(插值)

数值分析上机实验报告(插值)

数值分析第一次上机练习实验报告——Lagrange 插值与三次样条插值一、 问题的描述设()2119f x x =+, []1,1x ∈-,取15iix =-+,0,1,2,...,10i =.试求出10次Lagrange 插值多项式()10L x 和三次样条插值函数()S x (采用自然边界条件),并用图画出()f x ,()10L x ,()S x .二、 方法描述——Lagrange 插值与三次样条插值我们取15i ix =-+,0,1,2,...,10i =,通过在i x 点的函数值()2119i i f x x =+来对原函数进行插值,我们记插值函数为()g x ,要求它满足如下条件:()()21,0,1,2,...,1019i i ig x f x i x ===+ (1) 我们在此处要分别通过Lagrange 插值(即多项式插值)与三次样条插值的方法对原函数()2119f x x=+进行插值,看两种方法的插值结果,并进行结果的比较。

10次的Lagrange 插值多项式为:()()10100i i i L x y l x ==∑ (2)其中:()21,0,1,2,...,1019i i i y f x i x ===+以及()()()()()()()()()011011......,0,1,2,...,10......i i n i i i i i i i n x x x x x x x x l x i x x x x x x x x -+-+----==----我们根据(2)进行程序的编写,我们可以通过几个循环很容易实现函数的Lagrange 插值。

理论上我们根据区间[]1,1-上给出的节点做出的插值多项式()n L x 近似于()f x ,而多项式()n L x 的次数n 越高逼近()f x 的精度就越好。

但实际上并非如此,而是对任意的插值节点,当n →+∞的时候()n L x 不一定收敛到()f x ;而是有时会在插值区间的两端点附近会出现严重的()n L x 偏离()f x 的现象,即所谓的Runge 现象。

数值分析第一次上机练习实验报告

数值分析第一次上机练习实验报告

数值分析第一次上机练习实验报告一、实验目的本次实验旨在通过上机练习,加深对数值分析方法的理解,并掌握实际应用中的数值计算方法。

二、实验内容1. 数值计算的基本概念和方法在本次实验中,我们首先回顾了数值计算的基本概念和方法。

数值计算是一种通过计算机进行数值近似的方法,其包括近似解的计算、误差分析和稳定性分析等内容。

2. 方程求解的数值方法接下来,我们学习了方程求解的数值方法。

方程求解是数值分析中非常重要的一部分,其目的是找到方程的实数或复数解。

我们学习了二分法、牛顿法和割线法等常用的数值求解方法,并对它们的原理和步骤进行了理论学习。

3. 插值和拟合插值和拟合是数值分析中常用的数值逼近方法。

在本次实验中,我们学习了插值和拟合的基本原理,并介绍了常见的插值方法,如拉格朗日插值和牛顿插值。

我们还学习了最小二乘拟合方法,如线性拟合和多项式拟合方法。

4. 数值积分和数值微分数值积分和数值微分是数值分析中的两个重要内容。

在本次实验中,我们学习了数值积分和数值微分的基本原理,并介绍了常用的数值积分方法,如梯形法和辛卜生公式。

我们还学习了数值微分的数值方法,如差商法和牛顿插值法。

5. 常微分方程的数值解法常微分方程是物理和工程问题中常见的数学模型,在本次实验中,我们学习了常微分方程的数值解法,包括欧拉法和四阶龙格-库塔法。

我们学习了这些方法的步骤和原理,并通过具体的实例进行了演示。

三、实验结果及分析通过本次实验,我们深入理解了数值分析的基本原理和方法。

我们通过实际操作,掌握了方程求解、插值和拟合、数值积分和数值微分以及常微分方程的数值解法等数值计算方法。

实验结果表明,在使用数值计算方法时,我们要注意误差的控制和结果的稳定性。

根据实验结果,我们可以对计算结果进行误差分析,并选择适当的数值方法和参数来提高计算的精度和稳定性。

此外,在实际应用中,我们还需要根据具体问题的特点和条件选择合适的数值方法和算法。

四、实验总结通过本次实验,我们对数值分析的基本原理和方法有了更加深入的了解。

(完整版)哈工大-数值分析上机实验报告

(完整版)哈工大-数值分析上机实验报告

实验报告一题目:非线性方程求解摘要:非线性方程的解析解通常很难给出,因此线性方程的数值解法就尤为重要。

本实验采用两种常见的求解方法二分法和Newton法及改进的Newton法。

前言:(目的和意义)掌握二分法与Newton法的基本原理和应用。

数学原理:对于一个非线性方程的数值解法很多。

在此介绍两种最常见的方法:二分法和Newton法。

对于二分法,其数学实质就是说对于给定的待求解的方程f(x),其在[a,b]上连续,f(a)f(b)<0,且f(x)在[a,b]内仅有一个实根x*,取区间中点c,若,则c恰为其根,否则根据f(a)f(c)<0是否成立判断根在区间[a,c]和[c,b]中的哪一个,从而得出新区间,仍称为[a,b]。

重复运行计算,直至满足精度为止。

这就是二分法的计算思想。

Newton法通常预先要给出一个猜测初值x0,然后根据其迭代公式产生逼近解x*的迭代数列{x k},这就是Newton法的思想。

当x0接近x*时收敛很快,但是当x0选择不好时,可能会发散,因此初值的选取很重要。

另外,若将该迭代公式改进为其中r为要求的方程的根的重数,这就是改进的Newton法,当求解已知重数的方程的根时,在同种条件下其收敛速度要比Newton法快的多。

程序设计:本实验采用Matlab的M文件编写。

其中待求解的方程写成function的方式,如下function y=f(x);y=-x*x-sin(x);写成如上形式即可,下面给出主程序。

二分法源程序:clear%%%给定求解区间b=1.5;a=0;%%%误差R=1;k=0;%迭代次数初值while (R>5e-6) ;c=(a+b)/2;if f12(a)*f12(c)>0;a=c;elseb=c;endR=b-a;%求出误差k=k+1;endx=c%给出解Newton法及改进的Newton法源程序:clear%%%% 输入函数f=input('请输入需要求解函数>>','s')%%%求解f(x)的导数df=diff(f);%%%改进常数或重根数miu=2;%%%初始值x0x0=input('input initial value x0>>');k=0;%迭代次数max=100;%最大迭代次数R=eval(subs(f,'x0','x'));%求解f(x0),以确定初值x0时否就是解while (abs(R)>1e-8)x1=x0-miu*eval(subs(f,'x0','x'))/eval(subs(df,'x0','x'));R=x1-x0;x0=x1;k=k+1;if (eval(subs(f,'x0','x'))<1e-10);breakendif k>max;%如果迭代次数大于给定值,认为迭代不收敛,重新输入初值ss=input('maybe result is error,choose a new x0,y/n?>>','s');if strcmp(ss,'y')x0=input('input initial value x0>>');k=0;elsebreakendendendk;%给出迭代次数x=x0;%给出解结果分析和讨论:1.用二分法计算方程在[1,2]内的根。

数值分析上机实验报告

数值分析上机实验报告

数值分析上机实验报告摘要:本报告是对数值分析课程上机实验的总结和分析,涵盖了多种算法和数据处理方法,通过对实验结果的分析,探究了数值计算的一般过程和计算的稳定性。

1. 引言数值计算是数学的一个重要分支,广泛应用于物理、金融、工程等领域。

本次实验是对数值分析课程知识的实际应用,通过上机实现算法,探究数值计算的可靠性和误差分析。

2. 实验方法本次实验中,我们实现了多种算法,包括:(1)牛顿迭代法求方程的根;(2)高斯消元法求线性方程组的解;(3)最小二乘法拟合数据点;(4)拉格朗日插值法估计函数值;(5)梯形公式和辛普森公式求积分近似值。

对于每个算法,我们都进行了多组数值和不同参数的实验,并记录了相关数据和误差。

在实验过程中,我们着重考虑了算法的可靠性和计算的稳定性。

3. 实验结果与分析在实验中,我们得到了大量的实验数据和误差分析,通过对数据的展示和分析,我们得到了以下结论:(1)牛顿迭代法求解非线性方程的根能够对算法的初始值和迭代次数进行适当的调整,从而达到更高的稳定性和可靠性。

(2)高斯消元法求解线性方程组的解需要注意到矩阵的奇异性和精度的影响,从而保证计算的准确性。

(3)最小二乘法拟合数据点需要考虑到拟合的函数形式和数据的误差范围,采取适当的数据预处理和拟合函数的选择能够提高计算的准确性。

(4)拉格朗日插值法估计函数值需要考虑到插值点的选择和插值函数的阶数,防止出现龙格现象和插值误差过大的情况。

(5)梯形公式和辛普森公式求积分近似值需要考虑到采样密度和拟合函数的选择,从而保证计算的稳定性和收敛速度。

4. 结论通过本次实验的分析和总结,我们得到了深入的认识和理解数值计算的一般过程和算法的稳定性和可靠性,对于以后的数值计算应用也提供了一定的指导和参考。

优质文档精选——数值分析上机实验报告

优质文档精选——数值分析上机实验报告

数值分析上机实验报告《数值分析》上机实验报告1.用Newton 法求方程 X 7-X 4+14=0在(0.1,1.9)中的近似根(初始近似值取为区间端点,迭代6次或误差小于0.00001)。

1.1 理论依据:设函数在有限区间[a ,b]上二阶导数存在,且满足条件{}αϕ上的惟一解在区间平方收敛于方程所生的迭代序列迭代过程由则对任意初始近似值达到的一个中使是其中上不变号在区间],[0)(3,2,1,0,)(')()(],,[x |))(),((|,|,)(||)(|.4;0)(.3],[)(.20)()(.110......b a x f x k x f x f x x x Newton b a b f a f mir b a c x f ab c f x f b a x f b f x f k k k k k k ==-==∈≤-≠>+令)9.1()9.1(0)8(4233642)(0)16(71127)(0)9.1(,0)1.0(,1428)(3225333647>⋅''<-=-=''<-=-='<>+-=f f x x x x x f x x x x x f f f x x x f故以1.9为起点⎪⎩⎪⎨⎧='-=+9.1)()(01x x f x f x x k k k k 如此一次一次的迭代,逼近x 的真实根。

当前后两个的差<=ε时,就认为求出了近似的根。

本程序用Newton 法求代数方程(最高次数不大于10)在(a,b )区间的根。

1.2 C 语言程序原代码:#include<stdio.h>#include<math.h> main(){double x2,f,f1;double x1=1.9; //取初值为 1.9 do{x2=x1;f=pow(x2,7)-28*pow(x2,4)+14; f1=7*pow(x2,6)-4*28*pow(x2,3); x1=x2-f/f1;}while(fabs(x1-x2)>=0.00001||x1<0.1); //限制循环次数 printf("计算结果:x=%f\n",x1);}1.3 运行结果:1.4 MATLAB 上机程序function y=Newton(f,df,x0,eps,M) d=0;for k=1:Mif feval(df,x0)==0d=2;breakelsex1=x0-feval(f,x0)/feval(df,x0);ende=abs(x1-x0);x0=x1;if e<=eps&&abs(feval(f,x1))<=epsd=1;breakendendif d==1y=x1;elseif d==0y='迭代M次失败';elsey= '奇异'endfunction y=df(x)y=7*x^6-28*4*x^3;Endfunction y=f(x)y=x^7-28*x^4+14;End>> x0=1.9;>> eps=0.00001;>> M=100;>> x=Newton('f','df',x0,eps,M);>> vpa(x,7)1.5 问题讨论:1.使用此方法求方解,用误差来控制循环迭代次数,可以在误差允许的范围内得到比较理想的计算结果。

数学实验报告-哈工大

数学实验报告-哈工大

实验成绩实验报告
课程名称:数学实验
实验项目:常微分方程数值解
所在院系:
学生姓名:
学生学号:
授课学期: 2020秋
完成时间:
当1)0(=y ,取区间为[0,10],用改进欧拉法得到的结果与ode23得到的结果如下图所示:
(2)0cos =+''x y y 1)0(=y ,0)0(='y 取区间[0,10],用改进欧拉法得到的结果与ode23得到的结果如下图所示:
四、实验结果或结论分析
由实验结果可知,改进欧拉法与用ode23输出的结果十分相近,欧拉法在求解微分方程中具有极其重要的作用。

若改进式欧拉法取点越多,迭代次数越多,所得结果越接近真实值。

五、Matlab源程序
问题1的程序:
n=101;
y=zeros(n,1);
x=zeros(n,1);
y(1)=1;
x(1)=0;
long=10;
b=0;
m=0;
for i=1:100
m=y(i+1);
x(i+1)=(long/n)*i;
K1=8;
while k1>0.1
m=b;
b=y(i)+ long/n/2*(fun(x(i+1),m)+fun(x(i),y(i)));
%b=2*y(i+1)-y(i)+(long/n)^2/2*(fun(x(i+1),m)+fun(x(i),y(i))); K1=abs(m-b);
end
y(i+1)=b;
end
[m,n]=ode23(@(x,y)(x^2-y^2),[0,10],0.10);
plot(x,y);。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

实验报告一题目: Gauss 列主元消去法摘要:求解线性方程组地方法很多,主要分为直接法和间接法.本实验运用直接法地Guass 消去法,并采用选主元地方法对方程组进行求解.前言:(目地和意义)1. 学习Gauss 消去法地原理.2. 了解列主元地意义.3. 确定什么时候系数阵要选主元数学原理:由于一般线性方程在使用Gauss 消去法求解时,从求解地过程中可以看到,若)1(-k kk a =0,则必须进行行交换,才能使消去过程进行下去.有地时候即使≠-)1(k kk a 0,但是其绝对值非常小,由于机器舍入误差地影响,消去过程也会出现不稳定得现象,导致结果不正确.因此有必要进行列主元技术,以最大可能地消除这种现象.这一技术要寻找行r ,使得)1()1(max ||->-=k ik ki k rk a a 并将第r 行和第k 行地元素进行交换,以使得当前地)1(-k kk a 地数值比0要大地多.这种列主元地消去法地主要步骤如下:1. 消元过程对k =1,2,…,n -1,进行如下步骤.1) 选主元,记ik ki rk a a >=max || 若||rk a 很小,这说明方程地系数矩阵严重病态,给出警告,提示结果可能不对.2) 交换增广阵A 地r ,k 两行地元素.kj rj a a ↔ (j=k,…,n +1)3) 计算消元kk kj ik ij ij a a a a a /-= (i=k+1,…,n ; j =k +1,……,n +1)2. 回代过程对k = n , n -1,…,1,进行如下计算)/(11,∑-=+-=nk j kk j kj n k k a x a a x至此,完成了整个方程组地求解.程序设计:本实验采用Matlab地M文件编写.Gauss消去法源程序:cleara=input('输入系数阵:>>\n')b=input('输入列阵b:>>\n')n=length(b);A=[a b]x=zeros(n,1);%%%函数主体for k=1:n-1;%%%是否进行主元选取if abs(A(k,k))<yipusilong;%事先给定地认为有必要选主元地小数yzhuyuan=1;else yzhuyuan=0;endif yzhuyuan;%%%%选主元t=A(k,k);for r=k+1:n;if abs(A(r,k))>abs(t)p=r;else p=k;endend%%%交换元素if p~=k;for q=k:n+1;s=A(k,q);A(k,q)=A(p,q);A(p,q)=s;endendend%%%判断系数矩阵是否奇异或病态非常严重if abs(A(k,k))< yipusilongdisp(‘矩阵奇异,解可能不正确’)end%%%%计算消元,得三角阵for r=k+1:n;m=A(r,k)/A(k,k);for q=k:n+1;A(r,q)=A(r,q)-A(k,q)*m;endendend%%%%求解xx(n)=A(n,n+1)/A(n,n);for k=n-1:-1:1;s=0;for r=k+1:n;s=s+A(k,r)*x(r);endt=(A(k,n+1)-s)x(k)=(A(k,n+1)-s)/A(k,k)end结果分析和讨论:例:求解方程⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡321643.5072.12-623.4712.31-32108-z y x .求解地结果为:x =[]367257386.0,05088607.0-49105822.0-, 例:求解方程⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡73109104-10172-42-4z y x 求得地结果为:x =[]857142857.1,89285714.0-196428571.0, 结论:采用Gauss 消去法时,如果在消元时对角线上地元素始终较大(假如大于10-5),那么本方法不需要进行列主元计算,计算结果一般就可以达到要求,否则必须进行列主元这一步,以减少机器误差带来地影响,使方法得出地结果正确.实验报告二题目: Rung 现象产生和克服摘要:由于高次多项式插值不收敛,会产生Runge 现象,本实验在给出具体地实例后,采用分段线性插值和三次样条插值地方法有效地克服了这一现象,而且还取地很好地插值效果.前言:(目地和意义)1. 深刻认识多项式插值地缺点.2. 明确插值地不收敛性怎样克服.3. 明确精度与节点和插值方法地关系.数学原理:在给定n+1个节点和相应地函数值以后构造n 次地Lagrange 插值多项式,实验结果表明(见后面地图)这种多项式并不是随着次数地升高对函数地逼近越来越好,这种现象就是Rung 现象.解决Rung 现象地方法通常有分段线性插值、三次样条插值等方法.分段线性插值:设在区间[a, b ]上,给定n+1个插值节点a=x 0<x 1<…<x n =b和相应地函数值y 0,y 1,…,y n ,,求作一个插值函数)(x φ,具有如下性质:1) j j y x =)(φ,j=0,1,…,n .2) )(x φ在每个区间[x i , x j ]上是线性连续函数.则插值函数)(x φ称为区间[a, b ]上对应n 个数据点地分段线性插值函数.三次样条插值:给定区间[a, b ]一个分划⊿:a=x 0<x 1<…<x N =b若函数S(x)满足下列条件:1) S(x)在每个区间[x i , x j ]上是不高于3次地多项式.2) S(x)及其2阶导数在[a, b ]上连续.则称S(x)使关于分划⊿地三次样条函数. 程序设计流程:本实验采用Matlab 地M 文件编写.其中待插值地方程写成function 地方式,如下function y=f(x);y=1/(1+25*x*x );写成如上形式即可,下面给出主程序Lagrange 插值源程序:n=input('将区间分为地等份数输入:\n');s=[-1+2/n*[0:n]];%%%给定地定点,Rf为给定地函数x=-1:0.01:1;f=0;for q=1:n+1;l=1;%求插值基函数for k=1:n+1;if k~=q;l=l.*(x-s(k))./(s(q)-s(k));elsel=l;endendf=f+Rf(s(q))*l;%求插值函数endplot(x,f,'r')%作出插值函数曲线grid onhold on分段线性插值源程序clearn=input('将区间分为地等份数输入:\n');s=[-1+2/n*[0:n]];%%%给定地定点,Rf为给定地函数m=0;hh=0.001;for x=-1:hh:1;ff=0;for k=1:n+1;%%%求插值基函数switch kcase 1if x<=s(2);l=(x-s(2))./(s(1)-s(2));elsel=0;endcase n+1if x>s(n);l=(x-s(n))./(s(n+1)-s(n));elsel=0;endotherwiseif x>=s(k-1)&x<=s(k);l=(x-s(k-1))./(s(k)-s(k-1));else if x>=s(k)&x<=s(k+1);l=(x-s(k+1))./(s(k)-s(k+1));elsel=0;endendendff=ff+Rf(s(k))*l;%%求插值函数值endm=m+1;f(m)=ff;end%%%作出曲线x=-1:hh:1;plot(x,f,'r');grid onhold on三次样条插值源程序:(采用第一边界条件)clearn=input('将区间分为地等份数输入:\n');%%%插值区间a=-1;b=1;hh=0.001;%画图地步长s=[a+(b-a)/n*[0:n]];%%%给定地定点,Rf为给定地函数%%%%第一边界条件Rf"(-1),Rf"(1)v=5000*1/(1+25*a*a)^3-50/(1+25*a*a)^4;for k=1:n;%取出节点间距h(k)=s(k+1)-s(k);endfor k=1:n-1;%求出系数向量lamuda,miula(k)=h(k+1)/(h(k+1)+h(k));miu(k)=1-la(k);end%%%%赋值系数矩阵Afor k=1:n-1;for p=1:n-1;switch pcase kA(k,p)=2;case k-1A(k,p)=miu(p+1);case k+1A(k,p)=la(p-1);otherwiseA(k,p)=0;endendend%%%%求出d阵for k=1:n-1;switch kcase 1d(k)=6*f2c([s(k) s(k+1) s(k+2)])-miu(k)*v;case n-1d(k)=6*f2c([s(k) s(k+1) s(k+2)])-la(k)*v;otherwised(k)=6*f2c([s(k) s(k+1) s(k+2)]);endend%%%%求解M阵M=A\d';M=[v;M;v];%%%%m=0;f=0;for x=a:hh:b;if x==a;p=1;elsep=ceil((x-s(1))/((b-a)/n));endff1=0;ff2=0;ff3=0;ff4=0;m=m+1;ff1=1/h(p)*(s(p+1)-x)^3*M(p)/6;ff2=1/h(p)*(x-s(p))^3*M(p+1)/6;ff3=((Rf(s(p+1))-Rf(s(p)))/h(p)-h(p)*(M(p+1)-M(p))/6)*(x-s(p));ff4=Rf(s(p))-M(p)*h(p)*h(p)/6;f(m)=ff1+ff2+ff3+ff4 ;end%%%作出插值图形x=a:hh:b;plot(x,f,'k')hold ongrid on结果分析和讨论: 本实验采用函数22511)(xx f +=进行数值插值,插值区间为[-1,1],给定节点为 x j =-1+jh ,h=0.1,j =0,…,n .下面分别给出Lagrang e 插值,三次样条插值,线性插值地函数曲线和数据表.图中只标出Lagrang e 插值地十次多项式地曲线,其它曲线没有标出,从数据表中可以看出具体地误差.表中,L10(x)为Lagrang e插值地10次多项式,S10(x),S40(x)分别代表n=10,40地三次样条插值函数,X10(x),X40(x)分别代表n=10,40地线性分段插值函数.x f(x)L10(x)S10(x) S40(x) X10(x) X40(x) -1.00000000000000 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 -0.95000000000000 0.04244031830239 1.92363114971920 0.04240833151040 0.04244031830239 0.04355203619910 0.04244031830239 -0.90000000000000 0.04705882352941 1.57872099034926 0.04709697585458 0.04705882352941 0.04864253393665 0.04705882352941 -0.85000000000000 0.05245901639344 0.71945912837982 0.05255839923979 0.05245901639344 0.05373303167421 0.05245901639344 -0.80000000000000 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 -0.75000000000000 0.06639004149378 -0.23146174989674 0.06603986172744 0.06639004149378 0.06911764705882 0.06639004149378 -0.70000000000000 0.07547169811321 -0.22619628906250 0.07482116198866 0.07547169811321 0.07941176470588 0.07547169811321 -0.65000000000000 0.08648648648649 -0.07260420322418 0.08589776360849 0.08648648648649 0.08970588235294 0.08648648648649 -0.60000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 -0.55000000000000 0.11678832116788 0.21559187891257 0.11783833017713 0.11678832116788 0.12500000000000 0.11678832116788 -0.50000000000000 0.13793103448276 0.25375545726103 0.14004371555730 0.13793103448276 0.15000000000000 0.13793103448276 -0.45000000000000 0.16494845360825 0.23496854305267 0.16722724315883 0.16494845360825 0.17500000000000 0.16494845360825 -0.40000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 -0.35000000000000 0.24615384615385 0.19058046675376 0.24054799403464 0.24615384615385 0.27500000000000 0.24615384615385 -0.30000000000000 0.30769230769231 0.23534659131080 0.29735691695860 0.30769230769231 0.35000000000000 0.30769230769231 -0.25000000000000 0.39024390243902 0.34264123439789 0.38048738140327 0.39024390243902 0.42500000000000 0.39024390243902 -0.20000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 -0.15000000000000 0.64000000000000 0.67898957729340 0.65746969368431 0.64000000000000 0.62500000000000 0.64000000000000 -0.10000000000000 0.80000000000000 0.84340742982890 0.82052861660828 0.80000000000000 0.75000000000000 0.80000000000000 -0.05000000000000 0.94117647058824 0.95862704866073 0.94832323122810 0.94117647058824 0.87500000000000 0.941176470588240 1.00000000000000 1.00000000000000 1.00000000000000 1.00000000000000 1.000000000000001.000000000000000.05000000000000 0.94117647058824 0.95862704866073 0.94832323122810 0.94117647058824 0.87500000000000 0.941176470588240.10000000000000 0.80000000000000 0.84340742982890 0.82052861660828 0.80000000000000 0.75000000000000 0.800000000000000.15000000000000 0.64000000000000 0.67898957729340 0.65746969368431 0.64000000000000 0.62500000000000 0.640000000000000.20000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.500000000000000.25000000000000 0.39024390243902 0.34264123439789 0.38048738140327 0.39024390243902 0.42500000000000 0.390243902439020.30000000000000 0.30769230769231 0.23534659131080 0.29735691695860 0.30769230769231 0.350000000000000.307692307692310.35000000000000 0.24615384615385 0.19058046675376 0.24054799403464 0.24615384615385 0.27500000000000 0.246153846153850.40000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.200000000000000.45000000000000 0.16494845360825 0.23496854305267 0.16722724315883 0.16494845360825 0.17500000000000 0.164948453608250.50000000000000 0.13793103448276 0.25375545726103 0.14004371555730 0.13793103448276 0.15000000000000 0.137931034482760.55000000000000 0.11678832116788 0.21559187891257 0.11783833017713 0.11678832116788 0.12500000000000 0.116788321167880.60000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.100000000000000.65000000000000 0.08648648648649 -0.07260420322418 0.08589776360849 0.08648648648649 0.08970588235294 0.086486486486490.70000000000000 0.07547169811321 -0.22619628906250 0.07482116198866 0.07547169811321 0.07941176470588 0.075471698113210.75000000000000 0.06639004149378 -0.23146174989674 0.06603986172744 0.06639004149378 0.06911764705882 0.066390041493780.80000000000000 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.058823529411760.85000000000000 0.05245901639344 0.71945912837982 0.05255839923979 0.05245901639344 0.05373303167421 0.052459016393440.90000000000000 0.04705882352941 1.57872099034926 0.04709697585458 0.04705882352941 0.04864253393665 0.047058823529410.95000000000000 0.04244031830239 1.92363114971920 0.04240833151040 0.04244031830239 0.04355203619910 0.042440318302391.00000000000000 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154从以上结果可以看到,用三次样条插值和线性分段插值,不会出现多项式插值是出现地Runge现象,插值效果明显提高.进一步说,为了提高插值精度,用三次样条插值和线性分段插值是可以增加插值节点地办法来满足要求,而用多项式插值函数时,节点数地增加必然会使多项式地次数增加,这样会引起数值不稳定,所以说这两种插值要比多项式插值好地多.而且在给定节点数地条件下,三次样条插值地精度要优于线性分段插值,曲线地光滑性也要好一些.实验报告三题目: 多项式最小二乘法摘要:对于具体实验时,通常不是先给出函数地解析式,再进行实验,而是通过实验地观察和测量给出离散地一些点,再来求出具体地函数解析式.又因为测量误差地存在,实际真实地解析式曲线并不一定通过测量给出地所有点.最小二乘法是求解这一问题地很好地方法,本实验运用这一方法实现对给定数据地拟合. 前言:(目地和意义)1. 学习使用最小二成法地原理2. 了解法方程地特性 数学原理:对于给定地测量数据(x i ,f i )(i=1,2,…,n ),设函数分布为∑==mj j j x a x y 0)()(ϕ特别地,取)(x j ϕ为多项式j j x x =)(ϕ (j=0, 1,…,m )则根据最小二乘法原理,可以构造泛函∑∑==-=n i mj i j j i m x a f a a a H 110))((),,,(ϕ令0=∂∂ka H(k=0, 1,…,m ) 则可以得到法方程⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡),(),(),(),(),(),(),(),(),(),(),(),(1010101111000100m m m m m m m m f f f a a a ϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕ求该解方程组,则可以得到解m a a a ,,,10 ,因此可得到数据地最小二乘解∑=≈mj j j x a x f 0)()(ϕ程序设计:本实验采用Matlab 地M 文件编写.其中多项式函数j j x =ϕ写成function 地方式,如下function y=fai(x,j)y=1;for i=1:jy=x.*y;end写成如上形式即可,下面给出主程序.多项式最小二乘法源程序clear%%%给定测量数据点(s,f)s=[3 4 5 6 7 8 9];f=[2.01 2.98 3.50 5.02 5.47 6.02 7.05];%%%计算给定地数据点地数目n=length(f);%%%给定需要拟合地数据地最高次多项式地次数m=10;%%%程序主体for k=0:m;g=zeros(1,m+1);for j=0:m;t=0;for i=1:n;%计算内积(fai(si),fai(si))t=t+fai(s(i),j)*fai(s(i),k);endg(j+1)=t;endA(k+1,:)=g;%法方程地系数矩阵t=0;for i=1:n;%计算内积(f(si),fai(si))t=t+f(i)*fai(s(i),k);endb(k+1,1)=t;enda=A\b%求出多项式系数x=[s(1):0.01:s(n)]';y=0;for i=0:m;y=y+a(i+1)*fai(x,i);endplot(x,y)%作出拟合成地多项式地曲线grid onhold onplot(s,f,'rx') %在上图中标记给定地点结果分析和讨论:例用最小二乘法处理下面地实验数据.并作出)f地近似分布图.(x分别采用一次,二次和五次多项式来拟合数据得到相应地拟合多项式为:y1=-0.38643+0.82750x;y2=-1.03024+1.06893x-0.02012x2;y5=-50.75309+51.53527x-19.65947x2+3.66585x3-0.32886x4+0.01137x5;分别作出它们地曲线图,图中点划线为y1曲线,实线为y2曲线,虚线为y5曲线.’x’为给定地数据点.从图中可以看出并不是多项式次数越高越好,次数高了,曲线越能给定点处和实际吻合,但别地地方就很差了.因此,本例选用一次和两次地多项式拟合应该就可以了.实验报告四题目: Romberg 积分法摘要:对于实际地工程积分问题,很难应用Newton-Leibnitz 公式去求解.因此应用数值方法进行求解积分问题已经有着很广泛地应用,本文基于Romberg 积分法来解决一类积分问题.前言:(目地和意义)1. 理解和掌握Romberg 积分法地原理;2. 学会使用Romberg 积分法;3. 明确Romberg 积分法地收敛速度及应用时容易出现地问题. 数学原理:考虑积分⎰=ba dx x f f I )()(,欲求其近似值,通常有复化地梯形公式、Simpsion公式和Cotes 公式.但是给定一个精度,这些公式达到要求地速度很缓慢.如何提高收敛速度,自然是人们极为关心地课题.为此,记T 1,k 为将区间[a,b ]进行2k 等分地复化地梯形公式计算结果,记T 2,k 为将区间[a,b ]进行2k 等分地复化地Simpsion 公式计算结果,记T 3,k 为将区间[a,b ]进行2k 等分地复化地Cotes 公式计算结果.根据Richardson 外推加速方法,可以得到收敛速度较快地Romberg 积分法.其具体地计算公式为: 1. 准备初值,计算)]()([21,1b f a f ba T +-=2. 按梯形公式地递推关系,计算∑-=-+-+-+-+=1201,11,11))5.0(2(221k i k k k k i ab a f a b T T 3. 按Romberg 积分公式计算加速值1441,11,11,--=----+---m mk m m k m m m k m T T T m=2,…,k4. 精度控制.对给定地精度R ,若R T T m m <--1,11,则终止计算,并取1,m T 为所求结果;否则返回2重复计算,直至满足要求地精度为止. 程序设计:本实验采用Matlab 地M 文件编写.其中待积分地函数写成function 地方式,例如如下function yy=f(x,y); yy=x.^3;写成如上形式即可,下面给出主程序Romberg 积分法源程序%%% Romberg 积分法 clear%%%积分区间 b=3; a=1;%%%精度要求 R=1e-5;%%%应用梯形公式准备初值 T(1,1)=(b-a)*(f(b)+f(a))/2; T(1,2)=T(1,1)/2+(b-a)/2*f((b+a)/2); T(2,1)=(4*T(1,2)-T(1,1))/(4-1); j=2; m=2;%%%主程序体%%%while(abs(T(m,1)-T(m-1,1))>R);%%%精度控制 j=j+1; s=0;for p=1:2^(j-2);s=s+f(a+(2*p-1)*h/(2^(j-1))); endT(1,j)=T(1,j-1)/2+h*s/(2^(j-1)); %%%梯形公式应用 for m=2:j; k=(j-m+1);T(m,k)=((4^(m-1))*T(m-1,k+1)-T(m-1,k))/(4^(m-1)-1); end end%%%给出 Romberg 积分法地函数表 I=T(m,1)结果分析和讨论: 1. 求积分dx x10063.精确解I= 24999676.运行程序得Romberg 积分法地函数表为1.0e+007 *4.70101520000000 3.05022950000000 2.63753307500000 2.49996760000000 2.49996760000000 0 2.49996760000000 0 0由函数表知Romberg 积分给出地结果为2.4999676*10^7,与精确没有误差,精度很高.2. 求积分dx xx⎰10sin . 直接按前面方法进行积分,会发现系统报错,出现了0为除数地现象.出现这种情况地原因就是当x=0时,被积函数分母出现0,如果用一个适当地小数ε(最好不要小于程序给定地最小误差值,但不能小于机器地最大精度)来代替可以避免这个问题.本实验取R =ε,可得函数表为:0.92073548319659 0.93979327500190 0.94451351171417 0.94569085359489 0.94598501993743 0.94614587227034 0.94608692395160 0.94608330088846 0.94608307538495 0 0.94608299406368 0.94608305935092 0.94608306035138 0 0 0.94608306038722 0.94608306036726 0 0 0 0.94608306036718 0 0 0 0故该函数地积分为0.94608306036718,取8位有效数字.3. 求积分dx x ⎰12sin本题地解析解很难给出,但运用Romberg 积分可以很容易给出近似解,函数表为:0.42073549240395 0.33406972582924 0.31597536075922 0.31168023948094 0.31062036680949 0.31035626065456 0.30518113697100 0.30994390573588 0.31024853238818 0.31026707591900 0.31026822526959 0 0.31026142365354 0.31026884083167 0.31026831215439 0.31026830189296 0 0 0.31026895856465 0.31026830376269 0.31026830173008 0 0 0 0.31026830119484 0.31026830172211 0 0 0 0 0.31026830172262 0 0 0 0 0故该函数地积分为0.31026830172262,取8位有效数字.结论:Romberg 积分通常要求被积函数在积分区间上没有奇点.如有奇点,且奇点为第一间断点,那么采用例3地方法,还是能够求出来地,否则,必须采用其它地积分方法.当然,Romberg 积分地收敛速度还是比较快地.。

相关文档
最新文档