实验报告实验项目名称估计水塔的水流量实验项目类型验证 演示 综合 设计 其他指导教师 成绩一、实验目的（1）学会对实际问题的分析方法（2）学会利用所学的知识解决实际问题（3）设计出相应的算法，编制相应的应用程序二、实验内容某居民区，其自来水是有一个圆柱形水塔提供，水塔高12.2m ，塔的直径为17.4m ，水塔是由水泵根据水塔中的水位自动加水，一般水泵每天工作两次。
按照设计，当水塔中的水位降低至最低水位，约8.2m 时，水泵自动启动加水。
当水位升至最高水位，约10.8m 时，水泵停止工作。
下表给出了某一天的测量记录，测量了28个时刻的数据，但由于水泵正向水塔供水，由3个时刻无法测量到水位（表中为—）。
试建立数学模型，计算居民的用水速度和日总用水量。
三、实验原理、方法（算法）、步骤时刻 0 0.921 1.8432.9493.8714.9785.900 水位 9.677 9.479 9.308 9.125 8.982 8.814. 8.686 时刻 7.006 7.928 8.967 9.981 10.925 10.945 12.032 水位 8.525 8.388 8.220 — — 10.820 10.500 时刻 12.954 13.875 14.982 15.903 16.826 17.931 19.037 水位 10.210 9.936 9.653 9.409 9.1808.921 8.662 时刻 19.959 20.839 22.015 22.958 23.880 24.986 25.908 水位 8.433 8.220 —10.82010.59110.35410.180（1）水泵不供水时水流速度的计算居民的用水速度即单位时间流出的水的体积，称之为流量。
由于水塔是一个圆柱体，体积24V D hπ=可以通过水位高度h 计算出来，这样在水泵不工作的时间段，水流速度就可以从体积对时间的导数计算出来，由于没有水的体积关于时间的函数表达式，而只能利用问题中给定的原始数据表公式hD V 24π=，计算出离散的在测量时刻的体积V ，并用差商代替微商，对各离散点求导。
设()t y v t =为t 时刻水塔中水的容量，（,i i t v ）表示i t 时刻水塔中水的容量为i f对于每段的前两个点利用公式：'1123311()(34)v t v v v t t =-+-- 231311'()()v t f f t t =--对于每段的后两个点利用公式：'1221()()n n n n n v t v v t t ---=--'1221()(34)n n n n n n v t v v v t t ---=-+-对于每段的其余点：'2112121()(88)3()i i i i i i i v t v v v v t t --+++-=-+--（2）水泵供水时的水流速度的计算水泵供水时段的流量不能直接由数据计算得到，因此考虑使用已经求出的水泵未启动时的流速来近似模拟水泵启动时的流速。
对第一、二段供水时段分别选取供水前、后的两个时间进行拟合。
选取多个多项式函数进行求解对比检验，得到最佳的多项式函数。
（3）日总的用水量的计算为了求得总的用水量，可以对已求得的流速在时刻0至24小时上积分24()Q V t dx=⎰四、实验环境（所用软件、硬件等）及实验数据文件实验环境：M atlab数据文件：ch1.m,nihe.m,niusunihe.m,jifen.m五、实验结果及实例分析(1)计算未供水时的流速t1=[0,0.921,1.843,2.949,3.871,4.978,5.900,7.006,7.928,8.967,10.945,12. 032,12.954,13.875,14.982,15.903,16.826,17.931,19.037,19.959,20.839,22. 958,23.880,24.986,25.908];t=0.1*t1;h=[9.677,9.479,9.308,9.125,8.982,8.814,8.686,8.525,8.388,8.220,10.820, 10.500,10.210,9.936,9.653,9.409,9.180,8.921,8.662,8.433,8.220,10.820,1 0.591,10.354,10.180];v=(pi*17.4^2/4)*hi=1;f(i)=-(-v(i+2)+4*v(i+1)-3*v(i))/(t(i+2)-t(i));i=2;f(i)=-(v(i+1)-v(i-1))/(t(i+1)-t(i-1));for i=3:8f(i)=-(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/ (3*(t(i+2)-t(i-2)));endi=9;f(i)=-(v(i+1)-v(i-1))/(t(i+1)-t(i-1));i=10;f(i)=-(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));i=11;f(i)=-(-v(i+2)+4*v(i+1)-3*v(i))/(t(i+2)-t(i));i=12;f(i)=-(v(i+1)-v(i-1))/(t(i+1)-t(i-1));for i=13:19f(i)=-(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(3*(t(i+2)-t(i-2)));endi=20;f(i)=-(v(i+1)-v(i-1))/(t(i+1)-t(i-1));i=21;f(i)=-(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));i=22;f(i)=-(-v(i+2)+4*v(i+1)-3*v(i))/(t(i+2)-t(i));i=23;f(i)=-(v(i+1)-v(i-1))/(t(i+1)-t(i-1));i=24;f(i)=-(v(i+1)-v(i-1))/(t(i+1)-t(i-1));i=25;f(i)=-(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));disp('水塔流速')fplot(t,f,'b*')title('流速散点图');xlabel('时间（小时）');ylabel('流速（立方米/小时）')结果：水塔流速f =Columns 1 through 18545.7621 476.0904 437.5709 379.6076 364.5640 345.4175 335.6489 355.6773 369.8372 419.9609 793.0184 722.0018 656.7567 689.0155 652.1951 541.1182 569.3991 619.0908Columns 19 through 25581.4813 583.2515 553.9360 527.6342 546.3945 481.9059 367.512700.51 1.52 2.53300400500600700800流速散点图时间（小时）流速（立方米/小时）（2）拟合水泵供水时的水流速度首先对第一、二段未供水时段进行拟合：t=[0 0.921 1.843 2.949 3.871 4.978 5.900 7.006 7.928 8.967 10.954 12.032 12.954 13.875 14.982 15.903 16.826 17.931 19.037 19.959 20.839 22.958 23.880 24.986 25.908];h=[9.677 9.479 9.308 9.125 8.982 8.814 8.686 8.525 8.388 8.220 10.820 10.500 10.210 9.936 9.653 9.409 9.180 8.921 8.662 8.433 8.220 10.820 10.591 10.354 10.180]; c1=polyfit(t(1:10),h(1:10),3); tp1=0:0.1:8.9; x1=polyval(c1,tp1); plot(tp1,x1)c1=polyfit(t(1:10),h(1:10),3); c2=polyfit(t(11:21),h(11:21),3); a1=polyder(c1); a2=polyder(c2);tp1=0:0.01:8.967;tp2=10.954:0.01:20.839; x13=-polyval(a1,tp1);x113=-polyval(a1,[0:0.01: 8.967]); wgsysl1=100*trapz(tp1,x113); x14=-polyval(a1,[ 7.928, 8.967]);x23=-polyval(a2,tp2);x114=-polyval(a2,[ 10.954:0.01: 20.839]); wgsysl2=100*trapz(tp2,x114); x24=-polyval(a2,[ 10.954, 12.032]);x25=-polyval(a2,[ 19.959,20.839]); subplot(1,2,1)plot(tp1,x13*100) subplot(1,2,2)plot(tp2,x23*100)结果：51014161820222410152025232425262728293031进行五次多项式拟合，得到的各系数为c1 =-0.0000 0.0004 -0.0041 0.0250 -0.2346 9.6769 c2 =0.0000 -0.0019 0.0626 -0.9892 7.3067 -8.7188对第一供水及混合时段进行拟合：对第一段供水时段进行拟合dygsdsj=[ 7.928,8.967, 10.954,12.032]; dygsdls=[x14, x24];nhjg=polyfit(dygsdsj, dygsdls,3); nhsj=7.928:0.1:12.032;nhlsjg=polyval(nhjg ,nhsj); gssj1=8.967:0.01:10.954;gs1=polyval(nhjg,[8.967:0.01:10.954]);gsysl1=100*trapz(gssj1,gs1) */该语句计算第1供水时段的总用水量*/ plot(nhsj, 100*nhlsjg)对混合时段进行拟合t3=[19.037,19.959,20.839,t(22),t(23)]; ls3=[x25*100,24.8373,21.4286,18.8720];nhhddxsxs=polyfit(t3,ls3,3); tp3=19.96:0.01:25.91;xx3=polyval(nhhddxsxs, tp3);gssj2=20.839:0.01:23;gs2=polyval(nhhddxsxs,[ 20.839:0.01:23]); gsysl2=trapz(gssj2,gs2); plot(tp3,xx3)结果：1820222426152025303540nhjg =-0.0001 0.0019 -0.0190 0.0640 0 0 gsysl1 = 45.7842 nhhddxsxs =1.0e+003 *0.0000 -0.0004 0.0143 -0.2048 1.0970 0（3）用水总量s1=100*trapz(tp1,x13) s2=100*trapz(tp2,x23) s3=100*trapz(gssj1,gs1) gsysl2=100*trapz(gssj2,gs2)结果：s1 =145.5843 S2 =789101112131520253035260.0198s3 =45.7842gsysl2 =51.5332故总的用水量为上面四项之和为502.9215结果分析：从拟合得到的各时间段的流速可以看出，其与根据实际数据进行数值微分得到的图形的趋势大致相同，认为结果较为合理。