数值分析上机实验报告
昆明理工大学数值分析上机报告3

数值分析实验报告(数值积分)姓名:学号:2006231011专业:材料学学院:云南省新材料制备与加工重点实验室授课教师:昆明理工大学06工科硕士 《数值分析》上机实验报告专业: 材料物化 姓名: 学号: 2006231011 任课教师: 作业完成实验室:实验内容:1.题目/要求:1、 利用Lagrange 插值公式()k n ki i ik in k n y x x x x x L ⎪⎪⎪⎭⎫ ⎝⎛--∑=∏≠==00 编写出插值多项式程序; 2、 给出插值多项式或分段三次插值多项式的表达式;3、 根据节点选取原则,对问题(2)用三点插值或二点插值,其结果如何; 对此插值问题用Newton 插值多项式其结果如何2.作业环境(包括选用的程序语言、运行环境) Visual C++ 6.03.数学(理论背景)描述1. Lagrange 插值多项式定义:若n 次多项式上满足条件个节点在,1),1,0)((10n j x x x n n j x l ⋅⋅⋅<+⋅⋅⋅=: 当0)(1)(=≠==k j k j x l j k x l j k 时;当时, 其中n k j ⋅⋅⋅=1,0, ,就称这1+n个次多项式)(),(),(10x l x l x l n ⋅⋅⋅为节点n x x x ⋅⋅⋅10,上的n 次插值基函数。
插值多项式可 表示为:)()(0x ly x L knk k n ∑==,称为Lagrange 插值多项式。
2.分段线形插值就是通过插值点用折线段连接起来逼近).(x f 设已知节点b x x x a n =<⋅⋅⋅<=10上的函数值n f f f ⋅⋅⋅,,21,记k kk k k h h x x h ma x ,1=-=+求一折线函数)(x I h 满足:[][]为分段线形插值函数上是线形函数,则称在每个区间记)(,)(3),,2,1,0()(2,,)(11000x I x x x I k f x I b a C x I h k k h k h h +==∈4.数值计算公式Lagrange 插值多项式()k nki i ik in k n y x x x x x L ⎪⎪⎪⎭⎫ ⎝⎛--∑=∏≠==00; 分段线形插函数,在每个小区间[]1,+k k x x 上可表示为:1111)(++++--+--=k kk kk k k k h f x x x x f x x x x h I (1+≤≤k k x x x )在整个区间[]b a ,上为:)()(0x lf x I jnj j h ∑==其中基函数),1,0,()(n k j x l jk k j ⋅⋅⋅==δ,其形式是:[][],.,,,,0);(,);0(,)(11111111⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧∉∈=≤≤--=≤≤--=+-+++---j j j j j j j j j j j j j x x x b a x n j x x x x x x x j x x x x x x x x l5.算法程序流程图Lagrange插值算法程序流程图分段低次插值算法程序流程图6.程序结构(程序中的函数调用关系图6.实验数据和实验结果(打印或用屏幕图形拷屏表示,可加为附页)(1)Lagrange 插值算法,62573238.0)596.0(=f 05422977.0)99.0(=f 分段低次插值算法,63270600.0)596.0(=f 10152800.1)99.0(=f (2Lagrange 插值算法16476189.0)8.1(=f 00126583.0)15.6(=f分段低次插值算法18160000.0)8.1(=f 00185000.0)15.6(=f7.讨论(包括题目要求的讨论和方法的适用性讨论)对于插值多项式L n (x),当∞→n 时,L n (x)不一定收敛到f(x),此时需用分段线性插值比L n (x)逼近f(x)好得多,题(2)就是一例。
数值分析2024上机实验报告

数值分析2024上机实验报告数值分析是计算数学的一个重要分支,它研究如何用数值方法来解决数学问题。
在数值分析的学习过程中,学生需要通过上机实验来巩固理论知识,并学会使用相应的数值方法来解决实际问题。
本篇报告将详细介绍2024年度数值分析上机实验的内容和结果。
一、实验内容2024年度数值分析上机实验分为四个部分,分别是:方程求根、插值与拟合、数值积分和常微分方程的数值解。
1.方程求根这部分实验要求使用数值方法求解给定的非线性方程的根。
常见的数值方法有二分法、牛顿法、割线法等。
在实验过程中,我们需要熟悉这些数值方法的原理和实现步骤,并对不同方法的收敛性进行分析和比较。
2.插值与拟合这部分实验要求使用插值和拟合方法对给定的一组数据进行拟合。
插值方法包括拉格朗日插值、牛顿插值等;拟合方法包括最小二乘拟合、多项式拟合等。
在实验中,我们需要熟悉插值和拟合方法的原理和实现步骤,并对不同方法的精度和稳定性进行比较。
3.数值积分这部分实验要求使用数值方法计算给定函数的积分。
常见的数值积分方法有梯形法则、辛普森法则、龙贝格积分等。
在实验过程中,我们需要熟悉这些数值积分方法的原理和实现步骤,并对不同方法的精度和效率进行比较。
4.常微分方程的数值解这部分实验要求使用数值方法求解给定的常微分方程初值问题。
常见的数值方法有欧拉法、改进的欧拉法、四阶龙格-库塔法等。
在实验中,我们需要熟悉这些数值解方法的原理和实现步骤,并对不同方法的精度和稳定性进行比较。
二、实验结果在完成2024年度数值分析上机实验后,我们得到了以下实验结果:1.方程求根我们实现了二分法、牛顿法和割线法,并对比了它们的收敛速度和稳定性。
结果表明,割线法的收敛速度最快,但在一些情况下可能会出现振荡;二分法和牛顿法的收敛速度相对较慢,但稳定性较好。
2.插值与拟合我们实现了拉格朗日插值和最小二乘拟合,并对比了它们的拟合效果和精度。
结果表明,拉格朗日插值在小区间上拟合效果较好,但在大区间上可能出现振荡;最小二乘拟合在整体上拟合效果较好,但可能出现过拟合。
数值分析上机实验报告

数值分析上机实验理学院11级统计01班41108030125鲁庆实验报告一一.实验名称误差与误差估计二.实验目的掌握数值运算的误差估计方法三.数学原理 1.绝对误差(*)e x设某一量的准确值为x ,近似值为x*,则x*与x 之差叫做近似值x*的绝对误差(简称误差),记为*(*)*e e x x x ==- 2.绝对误差限适当小的正数,使|(*)||*|*e x x x ε=-≤则称*ε为近似值 x * 的绝对误差限。
(有时用*x x ε*=±表示近似值x *的精度或准确值的所在范围。
3.相对误差(*)r e x绝对误差与准确值之比*(*)*(*),0r r e x x xe e x x x x-===≠称为x *的相对 误差。
4.相对误差限(*)r x ε若指定一个适当小的正数 (*)r x ε,使|(*)||(*)|(*)||r r e x e x x x ε=≤则称(*)r x ε为近似值 x *的相对误差限。
5.有效数字若近似值x*的绝对误差限是某一位的半个单位,该位到x*的第一位非零数字一共有n 位,则称近似值x*有n 位有效数字,或说x*精确到该位。
6.绝对误差的运算:)()()(2121x x x x εεε+=± )()()(122121x x x x x x εεε+≈22122121+=x x x x x x x )()()(εεε (f(x))()(x)f x εε'≈四.实验内容1. 计算I n=e 1-⎰10nxe x 2dx (n=0,1,...)并估计误差。
解: >> I0 = exp(-1)*quad('(x.^0).*exp(x.^2)',0,1,10^(-10));>> vpa(I0,10) ans =.5380795069>> I1= exp(-1)*quad('(x.^1).*exp(x.^2)',0,1,10^(-10)); >> vpa(I1,10) ans =.3160602794>> I2 = exp(-1)*quad('(x.^2).*exp(x.^2)',0,1,10^(-10)); >> vpa(I2,10) ans =.2309602465>> I3 = exp(-1)*quad('(x.^3).*exp(x.^2)',0,1,10^(-10)); >> vpa(I3,10) ans =.1839397206>> I4 = exp(-1)*quad('(x.^4).*exp(x.^2)',0,1,10^(-10)); >> vpa(I4,10) ans =.1535596302>> I5 = exp(-1)*quad('(x.^5).*exp(x.^2)',0,1,10^(-10)); >> vpa(I5,10) ans =.1321205588>> I6 = exp(-1)*quad('(x.^6).*exp(x.^2)',0,1,10^(-10)); >> vpa(I6,10) ans =.1161009245>> I7 = exp(-1)*quad('(x.^7).*exp(x.^2)',0,1,10^(-10)); >> vpa(I7,10) ans =.1036383235>> I8 = exp(-1)*quad('(x.^8).*exp(x.^2)',0,1,10^(-10)); >> vpa(I8,10) ans =.9364676413e-1>> I9 = exp(-1)*quad('(x.^9).*exp(x.^2)',0,1,10^(-10)); >> vpa(I9,10) ans =.8544670595e-1 2.计算x255的值。
数值分析上机实验报告(插值)

数值分析第一次上机练习实验报告——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 i g 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 iy 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 现象。
数值分析实验报告5篇

1.69376699767424 0.92310666706964 0.08471614569741 0.40804026409411
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
讨论:
利用这种方法进行这类实验,可以很精确的扰动敏感性的一般规律。即 当对扰动项的系数越来越小时,对其多项式扰动的结果也就越来越小, 即扰动敏感性与扰动项的系数成正比,扰动项的系数越大,对其根的扰 动敏感性就越明显,当扰动的系数一定时,扰动敏感性与扰动的项的幂 数成正比,扰动的项的幂数越高,对其根的扰动敏感性就越明显。
解线性方程组的直接方法
实验 (主元的选取与算法的稳定性) 问题提出:Gauss消去法是我们在线性代数中已经熟悉的。但由于计算 机的数值运算是在一个有限的浮点数集合上进行的,如何才能确保 Gauss消去法作为数值算法的稳定性呢?Gauss消去法从理论算法到数值 算法,其关键是主元的选择。主元的选择从数学理论上看起来平凡,它 却是数值分析中十分典型的问题。 实验内容:考虑线性方程组 编制一个能自动选取主元,又能手动选取主元的求解线性方程组的 Gauss消去过程。 实验要求: (1)取矩阵,则方程有解。取n=10计算矩阵的条件数。让程序自动选 取主元,结果如何? (2)现选择程序中手动选取主元的功能。每步消去过程总选取按模最 小或按模尽可能小的元素作为主元,观察并记录计算结果。若每步消去 过程总选取按模最大的元素作为主元,结果又如何?分析实验的结果。 (3)取矩阵阶数n=20或者更大,重复上述实验过程,观察记录并分析 不同的问题及消去过程中选择不同的主元时计算结果的差异,说明主元
东南大学数值分析上机报告完整版

数值分析上机实验报告目录1.chapter1舍入误差及有效数 (1)2.chapter2Newton迭代法 (3)3.chapter3线性代数方程组数值解法-列主元Gauss消去法 (7)4.chapter3线性代数方程组数值解法-逐次超松弛迭代法 (8)5.chapter4多项式插值与函数最佳逼近 (10)1.chapter1舍入误差及有效数1.1题目设S N =∑1j 2−1N j=2,其精确值为)11123(21+--N N 。
(1)编制按从大到小的顺序11131121222-+⋯⋯+-+-=N S N ,计算S N 的通用程序。
(2)编制按从小到大的顺序1211)1(111222-+⋯⋯+--+-=N N S N ,计算S N 的通用程序。
(3)按两种顺序分别计算64210,10,10S S S ,并指出有效位数。
(编制程序时用单精度)(4)通过本次上机题,你明白了什么? 1.2编写相应的matlab 程序 clear;N=input('please input N:'); AValue=((3/2-1/N-1/(N+1))/2); sn1=single(0); sn2=single(0); for i=2:Nsn1=sn1+1/(i*i-1); %从大到小相加的通用程序% endep1=abs(sn1-AValue); for j=N:-1:2sn2=sn2+1/(j*j-1); %从小到大相加的通用程序% endep2=abs(sn2-AValue);fprintf('精确值为:%f\n',AValue);fprintf('从大到小的顺序累加得sn=%f\n',sn1); fprintf('从大到小相加的误差ep1=%f\n',ep1); fprintf('从小到大的顺序累加得sn=%f\n',sn2); fprintf('从小到大相加的误差ep2=%f\n',ep2); disp('================================='); 1.3matlab 运行程序结果 >> chaper1please input N:100 精确值为:0.740050从大到小的顺序累加得sn=0.740049 从大到小相加的误差ep1=0.000001 从小到大的顺序累加得sn=0.740050 从小到大相加的误差ep2=0.000000 >> chaper1please input N:10000 精确值为:0.749900从大到小的顺序累加得sn=0.749852 从大到小相加的误差ep1=0.000048 从小到大的顺序累加得sn=0.749900 从小到大相加的误差ep2=0.000000please input N:1000000精确值为:0.749999从大到小的顺序累加得sn=0.749852 从大到小相加的误差ep1=0.000147 从小到大的顺序累加得sn=0.749999 从小到大相加的误差ep2=0.0000001.4结果分析以及感悟按照从大到小顺序相加的有效位数为:5,4,3。
数值分析上机实验
目录1 绪论 (1)2 实验题目(一) (2)2.1 题目要求 (2)2.2 NEWTON插值多项式 (3)2.3 数据分析 (4)2.3.1 NEWTON插值多项式数据分析 (4)2.3.2 NEWTON插值多项式数据分析 (6)2.4 问答题 (6)2.5 总结 (7)3 实验题目(二) (8)3.1 题目要求 (8)3.2 高斯-塞德尔迭代法 (8)3.3 高斯-塞德尔改进法—松弛法 (9)3.4 松弛法的程序设计与分析 (9)3.4.1 算法实现 (9)3.4.2 运算结果 (9)3.4.3 数据分析 (11)4 实验题目(三) (13)4.1 题目要求 (13)4.2 RUNGE-KUTTA 4阶算法 (13)4.3 RUNGE-KUTTA 4阶算法运算结果及数值分析 (14)总结 (16)附录A (17)1绪论数值分析是计算数学的一个主要部分,它主要研究各类数学问题的数值解法,以及分析所用数值解法在理论上的合理性。
实际工程中的数学问题非常复杂,所以往往需要借助计算机进行计算。
运用数值分析解决问题的过程:分析实际问题,构建数学模型,运用数值计算方法,进行程序设计,最后上机计算求出结果。
数值分析这门学科具有面向计算机、可靠的理论分析、好的计算复杂性、数值实验、对算法进行误差分析等特点。
本学期开设了数值分析课程,该课程讲授了数值分析绪论、非线性方程的求解、线性方程组的直接接法、线性方程组的迭代法、插值法、函数逼近与曲线拟合、数值积分和数值微分、常微分方程初值问题的数值解法等内容。
其为我们解决实际数学问题提供了理论基础,同时我们也发现课程中很多问题的求解必须借助计算机运算,人工计算量太大甚至无法操作。
所以学好数值分析的关键是要加强上机操作,即利用计算机程序语言实现数值分析的算法。
本报告就是基于此目的完成的。
本上机实验是通过用计算机来解答数值分析问题的过程,所用的计算工具是比较成熟的数学软件MATLAB。
数值分析上机实践报告
数值分析上机实践报告一、实验目的本次实验主要目的是通过上机操作,加深对数值分析算法的理解,并熟悉使用Matlab进行数值计算的基本方法。
在具体实验中,我们将实现三种常见的数值分析算法:二分法、牛顿法和追赶法,分别应用于解决非线性方程、方程组和线性方程组的求解问题。
二、实验原理与方法1.二分法二分法是一种常见的求解非线性方程的数值方法。
根据函数在给定区间端点处的函数值的符号,不断缩小区间的长度,直到满足精度要求。
2.牛顿法牛顿法是求解方程的一种迭代方法,通过构造方程的泰勒展开式进行近似求解。
根据泰勒展式可以得到迭代公式,利用迭代公式不断逼近方程的解。
3.追赶法追赶法是用于求解三对角线性方程组的一种直接求解方法。
通过构造追赶矩阵,采用较为简便的向前追赶和向后追赶的方法进行计算。
本次实验中,我们选择了一组非线性方程、方程组和线性方程组进行求解。
具体的实验步骤如下:1.调用二分法函数,通过输入给定区间的上下界、截止误差和最大迭代次数,得到非线性方程的数值解。
2.调用牛顿法函数,通过输入初始迭代点、截止误差和最大迭代次数,得到方程组的数值解。
3.调用追赶法函数,通过输入追赶矩阵的三个向量与结果向量,得到线性方程组的数值解。
三、实验结果与分析在进行实验过程中,我们分别给定了不同的参数,通过调用相应的函数得到了实验结果。
下面是实验结果的汇总及分析。
1.非线性方程的数值解我们通过使用二分法对非线性方程进行求解,给定了区间的上下界、截止误差和最大迭代次数。
实验结果显示,根据给定的输入,我们得到了方程的数值解。
通过与解析解进行比较,可以发现二分法得到的数值解与解析解的误差在可接受范围内,说明二分法是有效的。
2.方程组的数值解我们通过使用牛顿法对方程组进行求解,给定了初始迭代点、截止误差和最大迭代次数。
实验结果显示,根据给定的输入,我们得到了方程组的数值解。
与解析解进行比较,同样可以发现牛顿法得到的数值解与解析解的误差在可接受范围内,说明牛顿法是有效的。
哈工大数值分析报告上机实验报告
实验报告一题目: 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 积分地收敛速度还是比较快地.。
数值分析上机实验
数值分析上机实验⽬录1 绪论 (1)2 实验题⽬(⼀) (2)2.1 题⽬要求 (2)2.2 NEWTON插值多项式 (3)2.3 数据分析 (4)2.3.1 NEWTON插值多项式数据分析 (4)2.3.2 NEWTON插值多项式数据分析 (6)2.4 问答题 (6)2.5 总结 (7)3 实验题⽬(⼆) (8)3.1 题⽬要求 (8)3.2 ⾼斯-塞德尔迭代法 (8)3.3 ⾼斯-塞德尔改进法—松弛法 (9)3.4 松弛法的程序设计与分析 (9)3.4.1 算法实现 (9)3.4.2 运算结果 (9)3.4.3 数据分析 (11)4 实验题⽬(三) (13)4.1 题⽬要求 (13)4.2 RUNGE-KUTTA 4阶算法 (13)4.3 RUNGE-KUTTA 4阶算法运算结果及数值分析 (14)总结 (16)附录A (17)1绪论数值分析是计算数学的⼀个主要部分,它主要研究各类数学问题的数值解法,以及分析所⽤数值解法在理论上的合理性。
实际⼯程中的数学问题⾮常复杂,所以往往需要借助计算机进⾏计算。
运⽤数值分析解决问题的过程:分析实际问题,构建数学模型,运⽤数值计算⽅法,进⾏程序设计,最后上机计算求出结果。
数值分析这门学科具有⾯向计算机、可靠的理论分析、好的计算复杂性、数值实验、对算法进⾏误差分析等特点。
本学期开设了数值分析课程,该课程讲授了数值分析绪论、⾮线性⽅程的求解、线性⽅程组的直接接法、线性⽅程组的迭代法、插值法、函数逼近与曲线拟合、数值积分和数值微分、常微分⽅程初值问题的数值解法等内容。
其为我们解决实际数学问题提供了理论基础,同时我们也发现课程中很多问题的求解必须借助计算机运算,⼈⼯计算量太⼤甚⾄⽆法操作。
所以学好数值分析的关键是要加强上机操作,即利⽤计算机程序语⾔实现数值分析的算法。
本报告就是基于此⽬的完成的。
本上机实验是通过⽤计算机来解答数值分析问题的过程,所⽤的计算⼯具是⽐较成熟的数学软件MATLAB。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数值分析实验报告学院:理学院班级:信计1001班 姓名:兰文利 学号:100350131实验一、用Jacobi 迭代法计算迭代次数一、实验目的:1、学习使用matlab 编写数值计算程序。
2、了解Jacobi 迭代法的迭代原理和编程方法。
3、根据Jacobi 迭代法的迭代原理编写matlab 程序,并运行的出结果。
二、实验原理:设n 阶非奇异矩阵A的主对角元nn a a a ,,,⋯⋯2211全不为0,记D=diag ( nn a a a ,,,⋯⋯2211)是非奇异对角阵,做A的一个分裂:)A D (D A --=。
记b D g A,D I )A D (D B -111=-=-=--。
迭代过程式为g Bx x k 1k +=+。
这种迭代方法称为Jacobi 迭代法。
其迭代矩阵记为:J B 右端向量记为:J g⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋯⋯⋯⋯⋯⋯⋯⋯⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋯-=----0a a a a 0a a a a 0a a a a 0a aa a B n3n2n13n 32312n 23211n 13121nn 133122111J = -⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⋯⋯⋯⋯⋯⋯⋯⋯0a a a a a a a a 0a a a a a a a a 0a a a a a a a a 0nnn3nnn2nnn1333n 33323331222n 22232221111n 11131112记Tk n k k k k x x x x x )()()(3)(2)(1,,,,⋯=,Jacobi 迭代过程为:⎪⎪⎪⎭⎪⎪⎪⎬⎫⎪⎪⎪⎩⎪⎪⎪⎨⎧⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋯+⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋯⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋯⋯⋯⋯⋯⋯⋯⋯=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋯-++++n k n k k k n n n n n n k n k k k b b b b x x x x a a a a a a a a a a a a D x x x x 321)()(3)(2)(13213323122321113121)()()()(10000113121 对照等号两边,得到计算)1(+k i x (i=1,2,…,n)的公式为⎥⎥⎦⎤⎢⎢⎣⎡+⎪⎪⎭⎫ ⎝⎛+-=∑∑-=+=+i i j ni j k jijk j ij iik ib xax a a x111)()()1(1,或⎪⎪⎭⎫ ⎝⎛+-=∑=+)(1)()1(1k i n j k j ij i ii k i x x a b a x ,i=1,2,…,n 三、实验内容与步骤:1、实验内容:依照实验原理编写Jacobi 迭代法的程序。
2、实验步骤:首先,在电脑上安装matlab ,然后,启动matlab ,新建一个M 文件。
1、实验送的初值为:A =[-2 1 0; 1 -2 1 ; 0 1 -2]; b=[-2 0 -3]';2、实验结果为:五、实验分析:1、在本次实验中所输的初值所对应的矩阵为收敛的,所以会产生结果,如果当输入一个Jacobi 矩阵为发散的初值时,则matlab 会提示出错。
同时也可能会导致死机。
2、在本次实验中,很容易把迭代矩阵中各分量的代码的写错,因为编程技术还不够熟练,所以在使用数学指令上还不够全面,因此常翻阅书籍找指令代码,这在以后的实验中是 应该加强训练的。
3、通过本次实验也增强了我对matlab 这个软件的了解。
实验二、用最小二乘法计算一、实验目的:1、学习使用matlab 编写数值计算程序。
2、了解最小二乘法的拟合原理和编程方法。
3、根据最小二乘法的拟合原理编写matlab 程序,并运行的出结果。
二、实验原理:设A 是N M ⨯阶矩阵)(N M >,称线性方程组:b Ax = (1)为超定方程组,这里mm R b R x ∈∈,。
如果A的秩r(A)=n ,称A为列满秩矩阵。
记残向量r=b-Ax ,考虑确定一个向量x,使2222rAx b -=达到最小的问题称为线性最小二乘问题,这样的x称为方程组(1)的最小二乘解。
线性最小二乘问题与欧式空间nR 的正交性等有密切联系。
设超定方程组(1)的系数矩阵A 是列满秩矩阵,由定理(当m>n 时,超定方程组(1)的最小二乘解总是存在。
最小二乘解惟一的充分必要条件是0)(=A null )知,方程组(1)的最小二乘解总是存在而且惟一的。
设x 是最小二乘解,21b b b +=。
由定理(当m>n 时,超定方程组(1)的最小二乘解总是存在。
最小二乘解惟一的充分必要条件是0)(=A null ),残向量)()(221TA N b b Ax b Ax b r ∈=+-=-=,因此,0)(=-=Ax b A r A TT(2)也就是说。
X 是方程组(1)的最小二乘解等价于x 是方程组(2)的解。
方程组(2)称为最小二乘问题的正则方程组或法方程组,由于A 是列满秩的,因此,方程组(2)是正定方程组。
用正则方程组求最小二乘解的计算步骤为如下:①形成b A g A A B T T ==,;②计算B 的Cholesky 分解,记T LL B =;③求解下三角方程组g Ly =; ④求解上三角方程组y x L T=。
三、实验内容与步骤:1、实验内容:依照实验原理编写最小二乘法的程序。
2、实验步骤:首先,在电脑上安装matlab ,然后,启动matlab ,新建一个M 文件。
1、实验送的初值为: x=[-100 -75 -50 -25 0 25 50 75 100];y=[4.06 6.78 9.49 16.27 40.67 97.62 146.43 151.85 162.7]; m=4; 拟合波形为:2、实验送的初值为:x=[-100 -75 -50 -25 0 25 50 75 100];y=[ 3.855 7.555 25.554 36.245 45.547 65.558 98.587 129.548 256.3525]; m=4; 实验结果为:拟合波形为:3、实验送的初值为:x=[-100 -75 -50 -25 0 25 50 75 100];y=[25.6 65.4 75.58 50.63 17.58 47.57 12.55 35.25 12.35]; m=4; 实验结果为:五、实验分析:1、在本次实验中所输的初值x与y对应的不太好时,就会出现像初值3那样的,拟合点和 拟合曲线之间有很大的差距,但基本上拟合点还是均匀分布在曲线周围的,初值1和初 值2的拟合曲线和拟合点均对应的很好。
基本在曲线上,或是曲线附近很小的范围内。
2、这种算法只适用于拟合函数是多项式的情形。
用正则方程组求最小二乘解得方法,有两个注意的问题:①形成正则方程组时,计算出来的A A T 可能不是正定的,也可能是奇异的,以致 后面不能用Cholesky 分解来求解正则方程组。
比如说:⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=--3310001011A 采用六位有限数位计算,得⎪⎪⎭⎫ ⎝⎛=1111A A T ,这时A A T是奇异的。
②矩阵A的状态数与A A T 的状态数密切相关,对列满秩的)(n m n m >⨯阶矩阵A, 有[])()(222A A cond A cond T =,因此,用正则方程组求最小二乘解时,会增加误差的敏感性实验三、迭代法求线性方程组一、 实验概述:【实验目的】1、了解线性方程组Ax=b 迭代法的一些基本理论思想;2、通过三对角形线性方程组的解法,体会稀疏线性方程组解法的特点;3、学习Gauss-Seidel 迭代法和松弛迭代法的使用方法;4、掌握迭代法的格式及收敛条件,判别迭代法收敛的基本方法及运用。
【实验原理】1、Gauss –Seidel 迭代法 迭代格式:b x U L D b Ax =--⇔=)((若L D -可逆)b L D Ux L D x 11)()(---+-=⇔ C Gx x k k +=+)()1(, ,2,1,0=kU L D G 1)(--=,b L D C 1)(--=,)0(x 为初始向量(矩阵形式)⎪⎭⎫ ⎝⎛+--=∑∑-=+=++111)()1()1(1i j ni j i k j ij k j ij ii k ib x a x a a x,),,2,1(n i =(分量形式)2、弛迭代法(简称SOR)迭代格式:在Seidel 迭代公式中,⎪⎭⎫ ⎝⎛+--=∑∑-=+=++111)()1()1(1i j ni j i k j ij k j ij ii k ib x a x a a x=)()(k i k i x x ∆+,)(k i x ∆⎪⎭⎫ ⎝⎛+--=∑∑-==+11)()1(1i j ni j i k j ij k j ij ii b x a x a a 。
若在修正项)(k i x ∆前加上一个参数w ,便得到松弛法的迭代公式:)()()1(k i k i k ix x x ∆+=+ω,n i ,,2,1 =,k =0, 1, 2, … 将上式写成矩阵形式得b L D x x k k 1)()1()(-+-+Φ=ωωω, k=0, 1, 2, … ])1[()(1U D L D ωωωω+--=Φ-其中w 叫做松弛因子, 当w >1时叫做超松弛, w < 1时叫做低松弛。
注:(1)w =1时就是Seidel 迭代法;(2)适当选取松弛因子w , 可加快收敛速度;(3)关于迭代格式(3),引进收敛速度概念:)(ln G R ρ-= 【实验环境】Microsoft visual C++6.0 二、 实验内容: 【实验方案】1、明确实验要求任务:(1) 用Gauss-Seidel 迭代法求解下列方程组,当满足4)()1(10-∞+<-k k x x 时结束迭代。
⎪⎩⎪⎨⎧=+-=++--=++1105210241225321321321x x x x x x x x x(2)用松弛迭代法(取25.1=ω)求解下列方程组,当满足5)()1(10-∞+<-k k x x 时结束迭代,说明迭代的收敛性。
⎪⎩⎪⎨⎧=+--=-+=-+054224521222321321321x x x x x x x x x并讨论初值的变化对收敛的影响。
2、确定具体解决算法;3、选择一门可行计算机语言进行编程;4、对程序进行调试、改进和验证。
【实验过程】(实验步骤、记录、数据、分析) 1、 打开Microsoft visual C++6.0编程软件; 2、 编写Gauss-Seidel 迭代法的程序对第一个方程组进行运算,并显示出迭代结果; 3、 编写松弛因子ω为1.26的SOR 方法对第二个方程组进行运算,并显示出迭代结果; 4、 分别对所编写的程序进行保存、调试、运行和改进,并保存结果图;5、 分析结果数据; 【实验结论】(结果)1、由Gauss-Seidel 迭代法求解的程序运行结果可知初值的选取影响迭代格式的收敛快慢。