贵州大学数值分析上机实验

合集下载

数值分析上机实验报告3

数值分析上机实验报告3

实验报告三题目:函数逼近——曲线拟合目的:掌握曲线拟合基本使用方法数学原理:[P,S]=polyfit(x,y,3)其中x,y为取样值,3为得出的结果的最高次数。

P为对应次数的系数,S为误差值向量,其中x,y是等长的向量,P是一个长度为m+1的向量。

结果分析和讨论:23.观察物体的运动,得出时间t与距离s的关系如表,求运动方程。

t=[0,0.9,1.9,3.0,3.9,5.0];s=[0,10,30,50,80,110];[P,S]=polyfit(t,s,5)P =-0.5432 6.4647 -26.5609 46.1436 -13.2601 -0.0000S =R: [6x6 double]df: 0normr: 1.2579e-012所以得到方程为:-13.2601x46.1436x-26.5609x6.4647x-0.5432x2345++24.在某化学反应堆里,根据实验所得分解物的质量分数y与时间t的关系,用最小拟合求y=F(t);>> x=0:5:55;y=[0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.62,4.64];>> [P,S]=polyfit(x,y,5)P =0.0000 -0.0000 0.0002 -0.0084 0.2851 0.0082S =R: [6x6 double]df: 6normr: 0.0487所以得到方程为:0082.02851.00084.00002.023++-xxx结论:在23题中计算的结果误差为4.5769,而在24中计算的结果误差为0.0487,说明对于曲线拟合来说,总会有误差,因为取样点并不是都过拟合的曲线的。

数值分析上机报告(1)

数值分析上机报告(1)

一.上机目的1. 通过上机编程,复习巩固以前所学程序设计语言及上机操作指令;2. 通过上机计算,了解舍入误差所引起的数值不稳定性;3. 熟悉并掌握拉格朗日插值多项式、牛顿插值多项式和分段低次插值,注意其不同特点;4. 了解最小二乘法的基本原理,能通过计算机解决实际问题。

二.上机环境MATLAB 软件等。

三.上机内容1.数值算法稳定性实验;2.插值法实验:拉格朗日插值、牛顿插值以及分段低次插值;3.曲线拟合实验:最小二乘法。

四.实验内容1、数值稳定性实验对n=0,1,2,…20计算定积分dx x x y n n ⎰+=105算法1 利用递推公式151--=n n y ny n=1,2,…,20 取182322.05ln 6ln 51100≈-=+-⎰dx x y 代码:y(1)=log(6)-log(5);for i=1:20y(i+1)=1/i-5*y(i);endk=ones(7,3);for i=1:7for j=1:3k(i,j)=y(3*(i-1)+j);endenddigits(6)vpa(k)结果:[ 0.182322, 0.0883922, 0.0580389][ 0.0431387, 0.0343063, 0.0284684][ 0.0243249, 0.0212326, 0.0188369][ 0.0169265, 0.0153676, 0.0140713][ 0.0129767, 0.0120398, 0.0112295][ 0.0105192, 0.00990388, 0.00930414][ 0.00903483, 0.00745741, 0.012713]算法2 利用递推公式n n y n y 51511-=- n=20,19,.…,1 注意到105151561126110200201020=≤+≤=⎰⎰⎰dx x dx x x dx x 取008730.0)12611051(2120≈+≈y 代码: y(21)=0.008730;for i=2:21j=22-i;y(j)=1/(5*j)-1/5*y(j+1);endk=ones(7,3);for i=1:7for j=1:3k(i,j)=y(3*(i-1)+j);endenddigits(6) ;vpa(k)结果:[ 0.182322, 0.0883922, 0.0580389][ 0.0431387, 0.0343063, 0.0284684][ 0.0243249, 0.0212326, 0.0188369][ 0.0169265, 0.0153676, 0.0140713][ 0.0129766, 0.0120399, 0.0112292][ 0.0105205, 0.0098975, 0.00933601][ 0.00887552, 0.008254, 0.00873]说明:从计算结果可以看出,算法1是不稳定的,而算法2是稳定的。

数值分析第七章实验报告7

数值分析第七章实验报告7

贵州师范大学数学与计算机科学学院学生实验报告课程名称: 数值分析 班级:实验日期:学 号: 点名序号:26姓名: 指导教师:实验成绩:一、实验名称实验六: 常微分方程初值问题数值解法二、实验目的及要求1. 让学生掌握用Euler 法, Runge-Kutta 法求解常微分方程初值问题.2. 培养Matlab 编程与上机调试能力.三、实验环境每人一台计算机,要求安装Windows XP 操作系统,Microsoft office2003、MATLAB6.5(或7.0).四、实验内容1. 取步长h=0.1,0.05,0.01, ,用Euler 法及经典4阶Runge-Kutta 法求解初值问题⎩⎨⎧=≤≤++-=1)0()10(2222'y t t t y y 要求:1) 画出准确解(准确解22t e y t +=-)的曲线,近似解折线;2) 把节点0.1和0.5上的精确解与近似解比较,观察误差变化情况.2. 用 Euler 法,隐式Euler 法和经典4阶R-K 法取不同步长解初值问题⎪⎩⎪⎨⎧=∈-=21)0(],1,0[,50'y x y y 并画出曲线观察稳定性.注:题1必须写实验报告五、算法描述及实验步骤4阶R-K 算法:功能输入 f(x),a,b,x0(x0=a),y0.输出 4阶R-K 解y.步1 m<=(b-a)/h,xn=a+n*h(n=1.2…m)步2 对n=0.1…m-1执行K1<=f(xn,yn)K2<=f(xn+0.5,y n+0.5*h*K1),K3<=f(xn+0.5,y n+0.5*h*K2).K4<f(x n+1,y n+h*K3).yn+1<=yn+(h/6)*(K1+2*K2+2*K3+K4);步3 输出y=(y1,y2,…,ym)’;结束Euler法算法功能解初值问题y’=f(x,y),y(x0)=y0;输入 f(x,y),a,b,h,x0(x0=a),y0;输出 Euler解y;步1 m<=(b-a)/h ,xn=a+n*h(n=1.2…m);步2 对n=0.1…m-1执行Yn+1<=yn+h*f(xn,yn);步3 输出y=(y1,y2,…,ym)’结束实验步骤:1.按要求在安装Windows2000或Windows XP操作系统,Matlab软件环境下编写源程序。

数值分析2024上机实验报告

数值分析2024上机实验报告

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

数值分析上机实验报告

数值分析上机实验报告

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

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

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

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

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

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

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

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

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

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

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

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

数值分析上机实验报告

数值分析上机实验报告

数值分析上机实验理学院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的值。

数值分析上机实践报告

数值分析上机实践报告

数值分析上机实践报告一、实验目的本次实验主要目的是通过上机操作,加深对数值分析算法的理解,并熟悉使用Matlab进行数值计算的基本方法。

在具体实验中,我们将实现三种常见的数值分析算法:二分法、牛顿法和追赶法,分别应用于解决非线性方程、方程组和线性方程组的求解问题。

二、实验原理与方法1.二分法二分法是一种常见的求解非线性方程的数值方法。

根据函数在给定区间端点处的函数值的符号,不断缩小区间的长度,直到满足精度要求。

2.牛顿法牛顿法是求解方程的一种迭代方法,通过构造方程的泰勒展开式进行近似求解。

根据泰勒展式可以得到迭代公式,利用迭代公式不断逼近方程的解。

3.追赶法追赶法是用于求解三对角线性方程组的一种直接求解方法。

通过构造追赶矩阵,采用较为简便的向前追赶和向后追赶的方法进行计算。

本次实验中,我们选择了一组非线性方程、方程组和线性方程组进行求解。

具体的实验步骤如下:1.调用二分法函数,通过输入给定区间的上下界、截止误差和最大迭代次数,得到非线性方程的数值解。

2.调用牛顿法函数,通过输入初始迭代点、截止误差和最大迭代次数,得到方程组的数值解。

3.调用追赶法函数,通过输入追赶矩阵的三个向量与结果向量,得到线性方程组的数值解。

三、实验结果与分析在进行实验过程中,我们分别给定了不同的参数,通过调用相应的函数得到了实验结果。

下面是实验结果的汇总及分析。

1.非线性方程的数值解我们通过使用二分法对非线性方程进行求解,给定了区间的上下界、截止误差和最大迭代次数。

实验结果显示,根据给定的输入,我们得到了方程的数值解。

通过与解析解进行比较,可以发现二分法得到的数值解与解析解的误差在可接受范围内,说明二分法是有效的。

2.方程组的数值解我们通过使用牛顿法对方程组进行求解,给定了初始迭代点、截止误差和最大迭代次数。

实验结果显示,根据给定的输入,我们得到了方程组的数值解。

与解析解进行比较,同样可以发现牛顿法得到的数值解与解析解的误差在可接受范围内,说明牛顿法是有效的。

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

数值分析上机实验报告课程名称:数值分析上机实验学院:机械工程学院专业:机械制造姓名:张法光学号:2012021691 年级:12级任课教师:代新敏老师2012年12月30日一.已知A 与b12.38412 2.115237 -1.061074 1.112336 -0.1135840.718719 1.742382 3.067813 -2.031743 2.11523719.141823-3.125432 -1.012345 2.1897361.563849-0.784165 1.112348 3.123124 -1.061074 -3.125A =43215.567914 3.123848 2.031454 1.836742-1.056781 0.336993 -1.010103 1.112336 -1.012345 3.12384827.108437 4.101011-3.741856 2.101023 -0.71828 -0.037585 -0.1135842.189736 2.031454 4.10101119.8979180.431637-3.111223 2.121314 1.784317 0.718719 1.563849 1.836742 -3.741856 0.4316379.789365-0.103458 -1.103456 0.238417 1.742382 -0.784165 -1.056781 2.101023-3.111223-0.10345814.7138465 3.123789 -2.213474 3.067813 1.112348 0.336993-0.71828 2.121314-1.103456 3.12378930.719334 4.446782 -2.031743 3.123124 -1.010103-0.037585 1.7843170.238417-2.213474 4.44678240.00001[ 2.1874369 33.992318 -25.173417 0.84671695 1.784317 -86.612343 1.1101230 4.719345 -5.6784392]TB ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦=(2)用超松弛法求解Bx=b (取松弛因子ω=1.4,x (0)=0,迭代9次)。

(3)用列主元素消去法求解Bx=b 。

解:(3)、用列主元素消去法求解Bx=b(一)、理论依据:其基本思想是选取绝对值尽量大的元素作为主元素,进行行与列的交换,再进行回代,求出方程的解。

将方阵A 和向量b 写成C=(A b )。

将C 的第1列中第1行的元素与其下面的此列的元素逐一进行比较,找到最大的元素1j c ,将第j 行的元素与第1行的元素进行交换,然后通过行变换,将第1列中第2到第n 个元素都消成0。

将变换后的矩阵(1)C 的第二列中第二行的元素与其下面的此列的元素逐一进行比较,找到最大的元素(1)2k c ,将第k 行的元素与第2行的元素进行交换,然后通过行变换,将第2列中第3到第n 个元素都消成0。

以此方法将矩阵的左下部分全都消成0。

(二)、计算程序: #include "math.h" #include "stdio.h" void main() { double u[9],x1[9],y[9],q[9],b1[9][10],x[9],a[9][9]={{12.38412,2.115237,-1.061074,1.112336,-0.113584,0.718719,1.742382,3.067813,-2.031743},{2.115237,19.141823,-3.125432,-1.012345,2.189736,1.563849,-0.784165,1.112348,3.123124},{-1.061074,-3.125432,15.567914,3.123848,2.031454,1.836742,-1.056781,0.336993,-1.010103}, {1.112336,-1.012345,3.123848,27.108437,4.101011,-3.741856,2.101023,-0.71828,-0.037585}, {-0.113584,2.189736,2.031454,4.101011,19.897918,0.431637,-3.111223,2.121314,1.784317}, {0.718719,1.563849,1.836742,-3.741856,0.431637,9.789365,-0.103458,-1.103456,0.238417}, {1.742382,-0.784165,-1.056781,2.101023,-3.111223,-0.103458,14.713846,3.123789,-2.213474}, {3.067813,1.112348,0.336993,-0.71828,2.121314,-1.103456,3.123789,30.719334,4.446782},{-2.031743,3.123124,-1.010103,-0.037585,1.784317,0.238417,-2.213474,4.446782,40.00001}};int sign(double x);double k,t,s,w,e,c,z;int i,j,n,r;doubleb[9]={2.1874369,33.992318,-25.173417,0.84671695,1.784317,-86.612343,1.1101230,4.719345,-5.6784392};for(r=0;r<=6;r++) /*Household 变换*/{e=0.0;for(i=r+1;i<=8;i++)e=e+a[i][r]*a[i][r];s=sqrt(e);t=s*s+fabs(a[r+1][r])*s;for(i=0;i<=r;i++) u[i]=0;c=a[r+1][r]; /*求u[i]的值*/u[r+1]=a[r+1][r]+s*sign(c);for(i=r+2;i<=8;i++)u[i]=a[i][r];for(i=0;i<=8;i++){y[i]=0;for(j=0;j<=8;j++)y[i]+=a[i][j]*u[j]/t;} /*求出y向量*/ k=0;for(i=0;i<9;++i)k+=0.5*(u[i]*y[i])/t;for(i=0;i<=8;i++)q[i]=y[i]-k*u[i];for(i=0;i<=8;i++)for(j=0;j<=8;j++)a[i][j]=a[i][j]-(q[j]*u[i]+u[j]*q[i]);} /*求结果*/printf("Household变换:\n");for(i=0;i<9;++i)for(j=0;j<9;++j) /*打印转化后的矩阵*/ {if (j%9==0)printf("\n");printf("%-9.5f",a[i][j]);}printf("\n");printf("超松弛变量法得解:\n");w=1.4; /*超松弛法*/for(i=0;i<9;i++)x1[i]=0;for(i=0;i<9;i++)for(j=0;j<9;j++){if(i==j)b1[i][j]=0;else b1[i][j]=-a[i][j]/a[i][i];} /*求出矩阵b1[9][9]和b1[i][9]的值*/ for(i=0;i<9;i++)b1[i][9]=b[i]/a[i][i];for(n=0;n<9;n++)for(i=0;i<9;i++){z=0;for(j=0;j<9;j++)z=z+b1[i][j]*x1[j]; /*执行本算法*/z=z+b1[i][9];x1[i]=x1[i]*(1-w)+w*z;}for(i=0;i<9;i++){if (i==5)printf("\n");printf("x%d=%-10.6f",i,x1[i]);}printf("\n");printf("列主元消去法得解:\n");u[0]=a[0][0]; /*以下是消去法*/y[0]=b[0];for(i=1;i<9;i++){q[i]=a[i][i-1]/u[i-1];u[i]=a[i][i]-q[i]*a[i-1][i];y[i]=b[i]-q[i]*y[i-1];}x[8]=y[8]/u[8]; /*执行本算法*/for(i=7;i>=0;i--)x[i]=(y[i]-a[i][i+1]*x[i+1])/u[i]; /*求出x的值*/ for(i=0;i<9;i++){if (i==5)printf("\n");printf(" x%d=%-10.6f",i,x[i]);}printf("\n");}int sign(double x){int z;z=(x>=(1e-6)?1:-1);return(z);}(三)、计算结果打印:(四)、问题讨论:a .由于选主元,使|lij|最小,这样去乘方程的每一系数时,系数的舍入误差不至扩大,并防止溢出与停机。

b . 由于选主元,回代时作除数主元aii(i)的绝对值也是最大,这样扩大的误差也是最小.c . 若detA ≠0,则选主元后的主元aii(i) ≠0,这是因为若aii(i)=0则必有aji(i)=0 (j>i)这样按行列式的Laplace 展开式就有detA=0而矛盾,因此用主元消去法中进行下去,不止中断停机.同时也由于选主元, aii(i)接近零的概率减少.运行结果基本与准确值无异,因为这种算法的无误差的算法。

三.试用三次样条插值求()4.563f 及(4.563)f '的近似值。

(一)方法由s(x i )=y i ,I=1,2,…N s’(x 0)=y 0’,s ’(x N )=y N ’可得S(x i )= ∑+-=11N j c j Ω3(x i -x j /h)=y iS ’(x 0)=1/h ∑+-=11N j c j Ω3’(x 0-x j /h)=y ’0S ’(x N )=1/h ∑+-=11N j c j Ω3’(x N -x j /h)=y ’N其中j c 必须的方程 :⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛---2.081551.1318338.1347667.12675460.11750557.106566268.93177664.85916738.6158839.4021011411411411411411411411411411411011098765432101c c c c c c c c c c c cq[0]=0 , u[0]=0 ,1,,2,1]),[][][(][][-⋅⋅⋅=⋅+-=n i i q i a i b i c i qn i i q i a i b i u i a i d i u ,,2,1]),1[][][(])[][][(}[⋅⋅⋅=-⋅+⋅-= x[9]=u[9]1,,2,1],[]1[][][⋅⋅⋅--=++⋅=n n i i u i x i q i x []333333)2()1(46)1(4)2(61)(+++++-+--++-+=Ωx x x x x x s ' (x)= ∑∑-=-=--Ω--+Ω10121012)21()21(j j j j j x c j x c(二)原程序及运行结果:#include<stdio.h> #include<math.h> #define N 9 void main() {double x[N+1]={1,2,3,4,5,6,7,8,9,10},y[N+1]={0,0.69314718,1.0986123,1.3862944,1.6094378,1.7917595,1.9459101,2.079445,2.1 972246,2.3025851},h[N+1],d[N+1],a[N+1],c[N+1],b[N+1]={2,2,2,2,2,2,2,2,2,2},s[N+1],t[N+1],l[N+1],M[N+1], f,f1;int i;for(i=1;i<=N;i++) /*使用对任意分化的三弯矩插值法*/h[i-1]=x[i]-x[i-1];d[0]=6/h[0]*((y[1]-y[0])/h[0]-1);d[N]=6/h[N-1]*(0.1-(y[N]-y[N-1])/h[N-1]);for(i=1;i<=N-1;i++){d[i]=6/(h[i-1]+h[i])*((y[i+1]-y[i])/h[i]-(y[i]-y[i-1])/h[i-1]);a[i]=h[i-1]/(h[i-1]+h[i]);c[i]=1-a[i];}c[0]=1;a[N]=1;s[0]=b[0];t[0]=d[0]; /*用追赶法求解三对角方程组*/for(i=0;i<=N-1;i++){l[i+1]=a[i+1]/s[i];s[i+1]=b[i+1]-l[i+1]*c[i];t[i+1]=d[i+1]-l[i+1]*t[i];}M[N]=t[N]/s[N];for(i=N-1;i>=0;i--)M[i]=(t[i]-c[i]*M[i+1])/s[i];f=M[3]*pow((x[4]-4.563),3)/6/h[3] /*求算4.563这点的函数值*/+M[4]*pow((4.563-x[3]),3)/6/h[3]+(y[3]-M[3]*h[3]*h[3]/6)*(x[4]-4.563)/h[3]+(y[4]-M[4]*h[3]*h[3]/6)*(4.563-x[3])/h[3];f1=-3*M[3]*pow((x[4]-4.563),2)/6/h[3] /*求算4.563这点的一阶导数值*/+3*M[4]*pow((4.563-x[3]),2)/6/h[3] -(y[3]-M[3]*h[3]*h[3]/6)/h[3] +(y[4]-M[4]*h[3]*h[3]/6)/h[3];printf("f(4.563)=%lf f'(4.563)=%lf\n",f,f1);}(三) 输出结果:(四) 问题讨论:其基本思想是对均匀分划的插值函数的构造,三次样条函数空间中不取1,x,,x 2,x 3,(x-x j )+3为基函数,而取B 样条函数Ω3(x-x j /h)为基函数.由于三次样条函数空间是N+3维的,故我们把分点扩大到X -1,X N+1,则任意三次样条函数可用Ω3(x-x j /h)线性组合来表示 S(x)= ∑+-=11N j c j Ω3(x-x j /h) 这样对不同插值问题,若能确定c j 由解的唯一性就能求得解S(x).同时样条插值效果比Lagrange 插值好,近似误差较小.没有Runge 现象。

相关文档
最新文档