MATLAB精品教程课件第13章 数值计算方法实际应用案例.ppt

13.1 水塔水流量的估计 根据题目的要求，关键在于确定水流量函数。流量是指单位时间内流出水的体积。由于水塔是正圆柱形，横截面积是常数，所以在水泵不工作时段，流量很容易根据水位相对时间的变化算出。问题的难点在于如何估计水泵供水时段的流量。 水泵供水时段的流量只能靠供水时段前后的流量经插值或拟合得到。作为用于插值或拟合的原始数据，我们希望水泵不工作时段的流量越准确越好。这此流量大体上可由两种方法计算，一是直接对书本表13-1中的水量用数值微分算出各时段的流量，用它们拟合其它时刻或连续时间的流量；二是先用表中数据拟合水位一时间函数，求导数即可得到连续时间的流量。 有了任何时刻的流量，就不难计算一天的总用水量。其实，水泵不工作时段的用水量可以由测量记录直接得到，由表13-1中下降水位乘以水搭的截面积就是这一时段的用水量。这个数值可以用来检验数据插值或拟合的结果。 为了表示方便，我们将表13-1中的数据全部化为国际标准单位，时间用小时(h)，高度用米(m)。为了书写的方便，首先将表13-1中的数据保存为data.mat文件，具体程序如下： t=[0,3316,6635,10619,13937,17921,21240,25223,28543,32284,... 35932,39332,39435,43318,46636,49953,53936,57254,60574,64554,... 68535,71854,75021,79254,82649,85968,89953,93270]; water=[3175,3110,3054,2994,2947,2892,2850,2795,2752,2697,... NaN,NaN,3550,3445,3350,3260,3167,3087,3012,2927,... 2842,2767,2697,NaN,NaN,3475,3397,3340]; save data t water % 保存为data.mat文件 单位标准化： load data t water % 读取data.mat文件 T=t/3600 h=0.01*water*0.3048 流量估计方法 首先我们需要确定水塔充水时间：由表13-1可知，对于第一次充水时间，当t=32284s时，水位为26.97ft，略低于最低水位27.0ft，因此可作为第一次开始充水时刻，当t=39435s时，水位为35.50ft，恰好为最高水位，因此可作为第一次充水结束时刻，故第一次充水时间为dt=(39435-32284)s=1.9864h；对于第二次充水时间，当t=82649s时，水泵在工作，但充水时间dt=(82649-75021)s=2.1189h，且下一时刻t=85968s时水塔水位为34.75ft，低于最高充水水位，因此可将t=82649s作为第二次充水的结束时刻，且该时刻的水位为最大充水高度35.50ft。 根据上面的介绍以及表13-2所给数据，利用水塔的截面积是常数，作出时间-水塔中水的体积图如图所示。 k=find(t==82649); % 寻找t=82649s的位置 h(k)=35.50*0.3048; % 将t=82649s时刻的值设为最高水位 T(isnan(h))=[]; % 将水泵启动的时刻置空 h(isnan(h))=[]; % 将水泵启动的时刻的水位置空 V=pi*(57*0.3048/2)^2*h; % 水塔中水的体积 plot(T,V,'K*') % 绘制图形 xlabel('时间/h') % 标注x轴 ylabel('水的体积/m^3') % 标注y轴 下面来计算水箱流量与时间的关系。水流量应该是水塔中水的体积对时间的导数，即 但此处由于没有水的体积关于时间的函数表达式，只有一个离散的函数值表，因此考虑数值微分。这里首先选用均差近似导数， % 利用均值近似表示导数 t1=T(1:10); % 第一段时间向量 t2=T(11:21); % 第二段时间向量 t3=T(22:end); % 第三段时间向量 V1=V(1:10); % 第一段水的体积 V2=V(11:21); % 第二段水的体积 V3=V(22:end); % 第三段水的体积 dV=abs([gradient(V1,t1),gradient(V2,t2),gradient(V3,t3)]); % 用均差近似导数 ti=linspace(min(T),max(T)); fi=spline(T,dV,ti); % fi=interp1(t,dV,ti,'spline') % 三次样条插值 plot(T,dV,'*',ti,fi,'k') xlabel('时间/h') ylabel('流速/m