水塔流量问题.docx

1、水塔流量问题本科生课程设计报告实习课程数值分析学院名称管理科学学院专业名称学生姓名学生学号指导教师实验地点实验成绩 二 一 六 年 六 月 二 一 六 年 六 月估计水塔的水流量摘要水塔流量的估计是一个较为经典的数学建模问题，本问题最大的困难在于不知泵启动时水位的变化和向外水流的速度.解决该问题，先确定近似流速，利用中点数值求导公式计算出每个时间点出的流速，再利用插值与拟合计算出流速与时间的函数，对0到24小时积分可得总用水量，这是第一种方法.第二种方法，水泵没有开动时利用高度差计算用水量，水泵开动时利用积分，这样计算出的结果较为准确，2种方法比较，可得出误差.关键词：中点数值求导；插值与拟合

2、；积分第1章 前言1.1 内容及要求某地的用水管理机构要求各社区提供各个时刻的用水率以及每天所用的总用水量。但许多社区并没有测量流入或流出当地水塔的水量的设备，他们只能代之以每小时测量水塔中的水位，其误差不超过0.5%。更为重要的是，无论什么时候，只要水塔中的水位下降到最低水位L时，水泵就自动启动向水塔重新充水直到最高水位H时水泵自动停止，但也无法得到水泵的供水量的测量数据。因此，在水泵正在工作时，人们不容易建立水塔中水位与水泵工作时的用水量之间的关系。水泵每天向水塔充水一次或两次，每次约二小时。下表为某地一天中的真实的数据。表1某天水塔水位测量记录时刻t（秒）0331666351061913

3、9371792121240水位（0.01英尺）3175311030542994294728922850时刻t（秒）25223285433228435932393323943543318水位（0.01英尺）279527522697水泵启动水泵启动35503445时刻t（秒）46636499535393657254605746455468535水位（0.01英尺）3350326031673087301229272842时刻t（秒）71854750217925482649859688995393270水位（0.01英尺）27672697水泵启动水泵启动347533973340水塔是一个高40英尺、直

4、径57英尺的圆柱。按照设计，水塔水位降至约L=27英尺时，水泵自动启动加水；当水位升高到约H=35.5英尺米时，水泵自动停止工作。试估计在任何时刻（包括水泵正在供水时）水从水塔流出的流量f(t)，并估计一天的总用水量。1.2 研究思路及结构安排本问题最大的困难在于不知道水泵启动时水位的变化和向外水流的速度。用水量等于向外水流速度乘以时间。因此，如何确定流速是解决该问题的关键。第2章 模型建立与求解2.1模型假设（1）除了问题中特别说明的数据以外，其他给定的数据其测量误差不超过0.5%；（2）一天之中，任意从一个开始时刻，如从t=0起到t=24小时结束，一天开始时刻的不同不影响一天总水量；（3）

5、管理部门不考虑水流速度的瞬间值，感兴趣的是整个一天中的用水总趋势；（4）水泵抽水的速度是均匀的；（5）假设水流的速度是连续变化的，流速可以用一条光滑的曲线近似表示，即在给定的数据点上有连续的二阶导数。2.2确定近似流速水泵不供水时,水流速度不难用数值方法描述，但如何描述清楚水泵向水塔供水期间的流速?我们首先估计水泵不启动时,每个已知数据点上的流速。流速等于水量的变化除以时间的变化。假设为时刻水塔中水的容量，表示时刻水塔中水的容量为。使用中点数值求导公式，求处的导数的近似值。水泵不启动时段共有3段，对于每段的前两个点对于每段的后两个点对于每段的其余点如此，可得各时刻点上的流速。 (负号保证流速为

6、正数)2.3 确定水泵启动时的流量及总流量曲线在，根据假设（5），水流的速度可以表示为一条光滑的曲线。因此考虑使用已经求出的水泵未启动时的流速来近似模拟水泵启动时的流速，这种模拟有两种方式:插值和拟合。由于假设流速具有连续的二阶导数，因此采用三次样条插值多项式来模拟任意时刻的流速。用表示时刻 的流速。2.4确定总用水量前面已经求出水流的流速曲线,如何求总用水量呢? 如果画出流速曲线，就可以从流速曲线图看出,流速在0时刻和24时的流速近似相等,因此流速是近似地以一天为周期的。可以考虑以下两种方式：1）一种方式是使用流速在时刻0至24小时上积分 2)直接根据高度的变化求出水泵不启动时的用水量用积分

7、求出水泵启动时的用水量近似值,然后相加。选择024小时内的范围，08.9678小时内,高度从31.75下降至26.97英尺，用水量为在,高度从35.50降至2697英尺，用水量为然后用数值积分求出两次水泵启动期间的用水量总用水量为第3章 算法步骤3.1 中点数值求导函数步骤及流程图图 1 中点数值求导函数流程图3.2 三次样条插值函数步骤及流程图图2 三次样条插值函数流程图第4章 算法实现4.1 程序总体结构图 3 程序函数调用关系4.2 源程序清单#include#include#define N 24double Q1=0,Q2=0,q1=0,q2=0,nN=0,mN=0;/依次为法一积分

8、、法二积分、两次启动水泵期间用水量、每次积分void zhuigan2(double a,double b,double c,double d,double m) /第二类边界追赶法求解方程组的解 double eN,fN; e0=c1/b1; f0=d1/b1; for(int i=1;iN-3;i+) ei=ci+1/(bi+1-ai+1*ei-1); for(i=1;i0;i-) mi=fi-1-ei-1*mi+1;void chazhi2(double x,double y,double s) /第二类边界条件及自然插值 double hN,aN,bN,cN,dN,mN; for(in

9、t i=0;iN;i+) bi=2; for(i=0;iN-1;i+) /计算hN hi=xi+1-xi; for(i=1;iN-1;i+) /计算aN,cN ai=hi-1/(hi+hi-1); ci=1-ai; for(i=1;iN-1;i+) /计算dN di=(6/(hi+hi-1)*(yi+1-yi)/hi)-(yi-yi-1)/hi-1); d1=d1-a1*s0; dN-2=dN-2-cN-2*s1; m0=s0; mN-1=s1; zhuigan2(a,b,c,d,m); /追赶法求mN printf(n追赶法求出m:n); for(i=0;iN;i+) /输出mN print

10、f(m%d=%.2fn,i,mi); /printfcz(h,x,y,m);void median(double *t,double *f,double *v) v0=1/(t2-t0)*(-3*f0+4*f1-f2);/第一段数据中值求导 v1=1/(t2-t0)*(f2-f0); v8=1/(t9-t7)*(f9-f7); v9=1/(t9-t7)*(3*f9-4*f8+f7); for(int i=2;i8;i+) vi=1/(3*(ti+2-ti-2)*(fi-2-8*fi-1+8*fi+1-fi+2); v10=1/(t12-t10)*(-3*f10+4*f11-f12);/第二段数

11、据中值求导 v11=1/(t12-t10)*(f12-f10); v19=1/(t20-t18)*(f20-f18); v20=1/(t20-t18)*(3*f20-4*f19+f18); for(i=12;i19;i+) vi=1/(3*(ti+2-ti-2)*(fi-2-8*fi-1+8*fi+1-fi+2); v21=1/(t23-t21)*(-3*f21+4*f22-f23);/第三段数据中值求导 v22=1/(t23-t21)*(f23-f21); v23=1/(t23-t21)*(3*f23-4*f22+f21);void Integral(double h,double x,do

12、uble y) /积分法 int i; double k3,c,d; /k数组为三次样条插值函数各项系数 printf(*n); printf(计算一天总水量:n); printf(tt方法1：024小时函数直接积分计算总用水量n); printf(积分结果为:n); for(i=0;iN-1;i+) /计算每段函数的积分 k0=(mi+1-mi)/(6*hi); k1=mi/2; k2=(yi+1-yi)/hi-hi*(2*mi+mi+1)/6; if(i21) c=pow(xi+1-xi,4)*k0/4+pow(xi+1-xi,3)*k1/3+pow(xi+1-xi,2)*k2/2+yi*

13、xi+1; /代上限 d=yi*xi; /代下限 ni=c-d; else c=pow(24-xi,4)*k0/4+pow(24-xi,3)*k1/3+pow(24-xi,2)*k2/2+yi*24; /代上限 d=yi*xi; /代下限 ni=c-d; for(i=0;iN-1;i+) Q1+=ni; /每段积分求和,即0到24小时用水量 printf(%fn,Q1); printf(*n);void way2(double t,double f) /高度变化结合积分求水量 double PI=3.,DIA=57; q1=fabs(PI*(DIA/2.0)*(DIA/2.0)*(f9-f0)

14、*0.01; q2=fabs(PI*(DIA/2.0)*(DIA/2.0)*(f19-f10)*0.01; Q2+=q1+q2+n9+n20+n21; printf(tt方法2：积分与高度变化结合求用水量); printf(n第1次水泵启动时用水量为：%fn,n9); printf(第2次水泵启动时用水量为：%fn,n20); printf(n%f%f用水量为：%fn,t0,t9,q1); printf(%f%f用水量为：%fn,t10,t20,q2); printf(总用水量为：%fn,Q2); printf(*n);void CountError() /计算误差 double w; w=(

15、Q1-n9-n20-n21)-(Q2-n9-n20-n21)/(Q2-n9-n20-n21); printf(用水量误差为：%f%n,w*100);void main() double tN=0,0.9211,1.8431,2.9497,3.8714,4.9781,5.9000,7.0064,7.9286,8.9678,10.9256,12.0328,12.9544,13.8758,14.9822,15.9039,16.8261,17.9317,19.0375,19.9594,22.9581,23.8800,24.9869,25.9083; double fN=3175,3110,3054,2

16、994,2947,2892,2850,2795,2752,2697,3550,3445,3350,3260,3167,3087,3012,2927,2842,2767,2697,3475,3397,3340; double vN=0; double S=3.1415*28.5*28.5; for(int i=0;iN;i+) fi=fi*S*0.01; printf(数据已录入，现在计算n); median(t,f,v);/计算导数 for(i=0;iN;i+) vi=-vi; printf(中值求导计算完毕，各时间点流量导数为:n); for(i=0;iN;i+) printf(%10.2f

17、,vi); printf(现在开始三次样条插值近似出导数函数n); double s2=0; chazhi2(t,v,s);/计算拟合函数 for(i=0;iN;i+) fi=fi/S/0.01; Integral(f,t,v); way2(t,f); /方法2计算总水量 CountError(); /计算误差4.3 程序运行图4 水塔流量程序运行结果图5 水塔流量程序运行结果运行结果说明：通过积分的方法算出的用水量为37879，通过综合计算的方法计算出的用水量为36363，两次结果都与实验说明相差8000，原因应该是计算过程中，曲线的拟合处出现了较大误差；而两次计算的相对误差为4.7%，结果

18、相对比较接近实验说明的要求。所以，此次编程虽然有些许误差，但是实验结果仍然表明整个方法没有根本性错误。第5章 误差分析考虑前面两种求总用水量的结果，由数值积分求出，结合高度和数值积分求出。根据高度求出的水泵未启动时的用水量可以认为是精确的。由数值积分求出的水泵未启动时的用水量。根据高度求出的水泵未启动时的用水量，可以认为是精确值。我们使用水泵未启动时的用水量的误差近似替代总的误差。第6章 模型的评价和改进6.1 优点1）模型的灵活性较好，可以用到那些拥有垂直圆柱型水塔的社区或单位，模型中的数据只要是离散的时间和水位的高度；2）)模型所用的数学概念简单，并且也容易理解，一阶导数、三次样条插值和积

19、分都是被广泛使用的概念；3）只要提供了数据和计算机，模型很容易实现；4）模型为一天里的流速和总用水量提供了较为精确的估计。6.2 缺点1)模型受到水塔形状限制；2)光滑曲线不能模拟真实流速观测中的微小变化。6.3 模型的改进方向1) 可以增加一些参数,以使模型适用于不同形状的水塔；2) 也可以加入随机因素，并进行方差分析，使模型更适应于真实的情况。 参考文献/archive/2013/01/26/2878092.html学生学习心得通过对三次样条插值法的应用，更加深入地了解了其在实际生产生活中的应用；通过对中值求导法的应用，能够更加熟练地对实际数据进行处理，以得到我们想要的数据。在对水塔问题的处理中，充分体会到了一些简单的数值处理方法的有机组合能够达到的强大功能，能够处理的问题的复杂性。学生（签名）： 2016 年 6月 16日诚信承诺本人郑重声明所呈交的课程报告是本人在指导教师指导下进行的研究工作及取得的研究成果。据我所知，除了文中特别加以标注的地方外，论文中不包含其他人已经发表或撰写过的研究成果。与我一同工作的同学对本文研究所做的贡献均已在报告中作了明确的说明并表示谢意。学生（签名）：任课教师评语成绩评定：任课教师（签名）： 年 月 日