《教学项目4数值分析法模型》由会员分享,可在线阅读,更多相关《教学项目4数值分析法模型(63页珍藏版)》请在金锄头文库上搜索。
1、1 教学项目四教学项目四 数值分析法模型数值分析法模型24.1 4.1 插值法建模插值法建模 拉格朗日插值拉格朗日插值分段线性插值分段线性插值三次样条插值三次样条插值一、插值的定义一、插值的定义二、插值的方法二、插值的方法三、用三、用Matlab解插值问题解插值问题3已知已知 n+1个节点个节点其中其中互不相同,不妨设互不相同,不妨设求任一插值点求任一插值点处的插值处的插值节点可视为由节点可视为由产生,产生,表达式复杂表达式复杂,或无封闭形式或无封闭形式,或未知。或未知。4.1.1 4.1.1 插值法的定义插值法的定义 4 构造一个构造一个(相对简单的相对简单的)函数函数通过全部节点通过全部节
2、点, 即即再用再用计算插值,即计算插值,即54.1.2 4.1.2 插值的方法插值的方法 已知函数f(x)在n+1个点x0,x1,xn处的函数值为 y0,y1,yn 。求一n次多项式函数Pn(x),使其满足:拉格朗日拉格朗日(Lagrange)插值插值6 称为拉格朗日插值基函数拉格朗日插值基函数。 解决此问题的拉格朗日插值多项式公式如下其中Li(x) 为n次多项式:7特殊情况特殊情况两点一次两点一次( (线性线性) )插值多项式插值多项式: :三点二次三点二次(抛物抛物)插值多项式插值多项式:8 采用拉格朗日多项式插值:选取不同插值节点个数n+1,其中n为插值多项式的次数,当n取10时,绘出插
3、值结果图形.例例9分段线性插值分段线性插值计算量与n无关;n越大,误差越小.xjxj-1xj+1x0xnxoy10比分段线性插值更光滑。比分段线性插值更光滑。xyxi-1 xiab 在数学上,光滑程度的定量描述是:函数(曲线)的k阶导数存在且连续,则称该曲线具有k阶光滑性。 光滑性的阶次越高,则越光滑。是否存在较低次的分段多项式达到较高阶光滑性的方法?三次样条插值就是一个很好的例子。三次样条插值三次样条插值11 三次样条插值g g( (x x) )为为被插值函数被插值函数。12用用MATLABMATLAB作插值计算作插值计算一维插值函数:一维插值函数:yi=interp1(x,y,xi,met
4、hod)插值方法插值方法被插值点被插值点插值节点插值节点xixi处的插处的插值结果值结果nearest :最邻近插值:最邻近插值linear : 线性插值;线性插值;spline : 三次样条插三次样条插值;值;cubic : 立方插值。立方插值。缺省时:缺省时: 分段线性插值。分段线性插值。注意:所有的插值方法都要求注意:所有的插值方法都要求x x是单调的,并且是单调的,并且xi不能够不能够超过超过x的范围。的范围。13 例:在例:在1-121-12的的1111小时内,每隔小时内,每隔1 1小时测量一次温小时测量一次温度,测得的温度依次为:度,测得的温度依次为:5 5,8 8,9 9,151
5、5,2525,2929,3131,3030,2222,2525,2727,2424。试估计每隔。试估计每隔1/101/10小时的温小时的温度值。度值。hours=1:12;temps=5 8 9 15 25 29 31 30 22 25 27 24;h=1:0.1:12;t=interp1(hours,temps,h,spline); (直接输出数据将是很多的)plot(hours,temps,+,h,t,hours,temps,r:) %作图xlabel(Hour),ylabel(Degrees Celsius)14xy实验习题实验习题 已知飞机下轮廓线上数据如下,求已知飞机下轮廓线上数据如
6、下,求x每改变每改变0.1时的时的y值。值。机翼下轮廓线154.2 4.2 拟合法简介拟合法简介 2.2.拟合的基本原理拟合的基本原理1. 拟合问题引例拟合问题引例16拟拟 合合 问问 题题 引引 例例 1 1温度温度t(0C) 20.5 32.7 51.0 73.0 95.7电阻电阻R( ) 765 826 873 942 1032已知热敏电阻数据:已知热敏电阻数据:求求60600C时的电阻时的电阻R。 设设 R=at+ba,b为待定系数为待定系数17拟拟 合合 问问 题题 引引 例例 2 2 t (h) 0.25 0.5 1 1.5 2 3 4 6 8c ( g/ml) 19.21 18.
7、15 15.36 14.10 12.89 9.32 7.45 5.24 3.01已知一室模型快速静脉注射下的血药浓度数据已知一室模型快速静脉注射下的血药浓度数据(t=0注射注射300mg)求血药浓度随时间的变化规律求血药浓度随时间的变化规律c(t).作半对数坐标系作半对数坐标系(semilogy)下的图形下的图形18曲曲 线线 拟拟 合合 问问 题题 的的 提提 法法已知一组(二维)数据,即平面上已知一组(二维)数据,即平面上 n个点个点(xi,yi) i=1,n, 寻求一个函数(曲线)寻求一个函数(曲线)y=f(x), 使使 f(x) 在某种准则下与所有在某种准则下与所有数据点最为接近,即曲
8、线拟合得最好。数据点最为接近,即曲线拟合得最好。 +xyy=f(x)(xi,yi)i i 为点为点(xi,yi) 与与曲线曲线 y=f(x) 的距离的距离19拟合与插值的关系拟合与插值的关系 函数插值与曲线拟合都是要根据一组数据构造一个函数作函数插值与曲线拟合都是要根据一组数据构造一个函数作为近似,由于近似的要求不同,二者的数学方法上是完全不同为近似,由于近似的要求不同,二者的数学方法上是完全不同的。的。 实例:实例:下面数据是某次实验所得,希望得到X和 f之间的关系?问题:问题:给定一批数据点,需确定满足特定要求的曲线或曲面解决方案:解决方案:若不要求曲线(面)通过所有数据点,而是要求它反映
9、对象整体的变化趋势,这就是数据拟合数据拟合,又称曲线拟合或曲面拟合。若要求所求曲线(面)通过所给所有数据点,就是插值问题插值问题;20最临近插值、线性插值、样条插值与曲线拟合结果:最临近插值、线性插值、样条插值与曲线拟合结果:21曲线拟合问题最常用的解法曲线拟合问题最常用的解法线性最小二乘法的基本思路线性最小二乘法的基本思路第一步: :先选定一组函数先选定一组函数 r1(x), r2(x), rm(x), m0)k(0)模型假设模型假设1. 1. 机体看作一个房室,室内血药浓度均匀机体看作一个房室,室内血药浓度均匀一室模型一室模型模型建立模型建立 在此,在此,d=300mg,t及及c(t)在某
10、些点处的值见前表,)在某些点处的值见前表,需经拟合求出参数需经拟合求出参数k、v45用线性最小二乘拟合用线性最小二乘拟合c(t)计算结果:计算结果:d=300;t=0.25 0.5 1 1.5 2 3 4 6 8;c=19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01;y=log(c);a=polyfit(t,y,1)k=-a(1)v=d/exp(a(2)程序:程序:46给药方案给药方案 设计设计cc2c10t 设每次注射剂量D, 间隔时间 血药浓度c(t) 应c1 c(t) c2 初次剂量D0 应加大给药方案记为:给药方案记为:2、1、计算结果
11、:计算结果:给药方案:给药方案:c1=10,c2=25k=0.2347v=15.0247故可制定给药方案:故可制定给药方案:即即: 首次注射首次注射375mg, 其余每次注射其余每次注射225mg, 注射的间隔时间为注射的间隔时间为4小时。小时。48估计水塔的流量估计水塔的流量2、解题思路、解题思路3、算法设计与编程、算法设计与编程1、问题、问题49 某居民区有一供居民用水的园柱形水塔,一般可以通过测量其水位来估计水的流量,但面临的困难是,当水塔水位下降到设定的最低水位时,水泵自动启动向水塔供水,到设定的最高水位时停止供水,这段时间无法测量水塔的水位和水泵的供水量通常水泵每天供水一两次,每次约
12、两小时.水塔是一个高12.2米,直径17.4米的正园柱按照设计,水塔水位降至约8.2米时,水泵自动启动,水位升到约10.8米时水泵停止工作表1 是某一天的水位测量记录,试估计任何时刻(包括水泵正供水时)从水塔流出的水流量,及一天的总用水量50流量估计的解题思路流量估计的解题思路拟合水位拟合水位时间函数时间函数确定流量确定流量时间函数时间函数估计一天总用水量估计一天总用水量51 拟合水位拟合水位时间函数时间函数 测量记录看,一天有两个供水时段(以下称第1供水时段和第2供水时段),和3个水泵不工作时段(以下称第1时段t=0到t=8.97,第2次时段t=10.95到t=20.84和第3时段t=23以
13、后)对第1、2时段的测量数据直接分别作多项式拟合,得到水位函数为使拟合曲线比较光滑,多项式次数不要太高,一般在36由于第3时段只有3个测量记录,无法对这一时段的水位作出较好的拟合52 2、确定流量确定流量时间函数时间函数 对于第1、2时段只需将水位函数求导数即可,对于两个供水时段的流量,则用供水时段前后(水泵不工作时段)的流量拟合得到,并且将拟合得到的第2供水时段流量外推,将第3时段流量包含在第2供水时段内533、一天总用水量的估计一天总用水量的估计 总用水量等于两个水泵不工作时段和两个供水时段用水量之和,它们都可以由流量对时间的积分得到。54算法设计与编程算法设计与编程1、拟合第拟合第1、2
14、时段的水位,并导出流量时段的水位,并导出流量2、拟合供水时段的流量拟合供水时段的流量3、估计一天总用水量估计一天总用水量4、流量及总用水量的检验、流量及总用水量的检验55 1、拟合第拟合第1时段的水位,并导出流量时段的水位,并导出流量 设t,h为已输入的时刻和水位测量记录(水泵启动的4个时刻不输入),第第1时段时段各时刻的流量可如下得:1) c1=polyfit(t(1:10),),h(1:10),),3);); %用3次多项式拟合第1时段水位,c1输出3次多项式的系数2)a1=polyder(c1);); % a1输出多项式(系数为c1)导数的系数 3)tp1=0:0.1:9; x1=-po
15、lyval(a1,tp1); % x1输出多项式(系数为a1)在tp1点的函数值(取负后边为正值),即tp1时刻的流量 4)流量函数为:流量函数为:56 2、拟合第拟合第2时段的水位,并导出流量时段的水位,并导出流量 设t,h为已输入的时刻和水位测量记录(水泵启动的4个时刻不输入),第第2时段时段各时刻的流量可如下得:1) c2=polyfit(t(10.9:21),h(10.9:21),3); %用3次多项式拟合第2时段水位,c2输出3次多项式的系数2) a2=polyder(c2); % a2输出多项式(系数为c2)导数的系数 3)tp2=10.9:0.1:21; x2=-polyval(
16、a2,tp2); % x2输出多项式(系数为a2)在tp2点的函数值(取负后边为正值),即tp2时刻的流量4)流量函数为:流量函数为:57 3、拟合供水时段的流量拟合供水时段的流量 在第1供水时段(t=911)之前(即第1时段)和之后(即第2时段)各取几点,其流量已经得到,用它们拟合第1供水时段的流量为使流量函数在t=9和t=11连续,我们简单地只取4个点,拟合3次多项式(即曲线必过这4个点),实现如下: xx1=-polyval(a1,8 9); %取第1时段在t=8,9的流量 xx2=-polyval(a2,11 12); %取第2时段在t=11,12的流量 xx12=xx1 xx2; c
17、12=polyfit(8 9 11 12,xx12,3); %拟合3次多项式 tp12=9:0.1:11; x12=polyval(c12,tp12); % x12输出第1供水时段 各时刻的流量拟合的流量函数为:拟合的流量函数为:58 在第2供水时段之前取t=20,20.8两点的流水量,在该时刻之后(第3时段)仅有3个水位记录,我们用差分得到流量,然后用这4个数值拟合第2供水时段的流量如下: dt3=diff(t(22:24)); %最后3个时刻的两两之差 dh3=diff(h(22:24)); %最后3个水位的两两之差 dht3=-dh3./dt3; %t(22)和t(23)的流量 t3=2
18、0 20.8 t(22) t(23); xx3=-polyval(a2,t3(1:2),dht3); %取t3各时刻的流量 c3=polyfit(t3,xx3,3); %拟合3次多项式 t3=20.8:0.1:24; x3=polyval(c3,tp3);% x3输出第2供水时段 (外推至t=24)各时刻的流量拟合的流量函数为:拟合的流量函数为:59 3、一天总用水量的估计一天总用水量的估计 第1、2时段和第1、2供水时段流量的积分之和,就是一天总用水量虽然诸时段的流量已表为多项式函数,积分可以解析地算出,这里仍用数值积分计算如下: y1=0.1*trapz(x1); %第1时段用水量(仍按高
19、 度计),0.1为积分步长 y2=0.1*trapz(x2); %第2时段用水量 y12=0.1*trapz(x12); %第1供水时段用水量 y3=0.1*trapz(x3); %第2供水时段用水量 y=(y1+y2+y12+y3)*237.8*0.01; %一天总用水量( ) 计算结果:计算结果:y1=146.2, y2=266.8, y12=47.4, y3=77.3,y=1250.460 4、流量及总用水量的检验流量及总用水量的检验 计算出的各各时时刻刻的的流流量量可用水位记录的数值微分来检验用水量y1可用第1时段水位测量记录中下降高度968-822=146来检验,类似地,y2用108
20、2-822=260检验供供水水时时段段流流量量的一种检检验验方方法法如下:供水时段的用水量加上水位上升值260是该时段泵入的水量,除以时段长度得到水泵的功率(单位时间泵入的水量),而两个供水时段水泵的功率应大致相等第1、2时段水泵的功率可计算如下: p1=(y12+260)/2; %第1供水时段水泵的功率 (水量仍以高度计) tp4=20.8:0.1:23; xp2=polyval(c3,tp4); % xp2输出第2供水时段 各时刻的流量 p2=(0.1*trapz(xp2)+260)/2.2; %第2供水时段水泵的功率 (水量仍以高度计)计算结果计算结果:p1=154.5 ,p2=140.161计算结果计算结果流量函数为:流量函数为:62流量曲线见图流量曲线见图n=(3,4)n=(5,6)63