第五章水塔用水量的估计一插值 §5-1水塔用水量问题 本题中所使用的长度单位为E(=30.24cm):容积单 位为G(=3.785L(升) 某些州的用水管理机构需估计公众的用水速度 (单位是G/h)和每天总用水量的数据。许多地方没有测 量流入或流出水箱流量的设备,而只能测量水箱中的 水位(误差不超过5%)当水箱水位低于某最低水位L 时,水泵抽水,灌入水箱内直至水位达到某最高水位 H为止。但是也无法测量水泵的流量,因此在水泵启 动时无法立即将水箱中的水位和水量联系起来。水泵 天灌水1~2次,每次约2h.试估计在任意时刻(包括 水泵灌水期间)流出水箱的流量f(t),并估计一天的总 用水量。 下表给出了某镇中某一天的真实用水数据,表中 测量时间以秒为单位,水位以00E为单位。例如在 3316s时,水箱中的水深为31.10E 时间)水位00D时间()水位00E 3175 46636 3350 3316 3110 49953 3260 6635 3054 59363167 10619 3087 13937 2947 60574 3012 17921 2892 64554 2927 21240 2850 68535 2842 25223 2795 71854 2767 75021 2697 32284 2697 79254 泵水 35932 82649 泵水 39332 泵水 85968 3475 39435 89953 3397 43318 3445 93270 3340 已知水箱是直径为57E、高为40E的正圆柱体 当水位落到27E以下时,水泵自动启动把水泵入水箱; 当水位回升至355E时,水泵停止工作。 §5-3水塔用水量的计算 5.3.1问题分析 1.假设 28次采集数据,得到24对数据(时刻:秒。水位: 英尺),要估计水流出水箱的速度ft(即:任意时刻t 时的水流量)。先作几条假设。 (1)只考虑通常情况下的正常用水,不考虑突发 事变,如水管破裂、自然灾害等 (2)全天24小时能保证用户的正常用水,不会发 生某个用户想用水时打开水龙头而没水流出这样的情 况。即:水箱中水的流出量ft完全由用户的需求量决 定,与水箱自身的情况如水位高度、出口口径等因素 无关。 (3)题中所给“水位约2700E时,水泵启动:水
第五章 水塔用水量的估计----插值 §5—1 水塔用水量问题 本题中所使用的长度单位为 E(=30.24cm)；容积单 位为 G(=3.785L(升)). 某些州的用水管理机构需估计公众的用水速度 (单位是 G/h)和每天总用水量的数据。许多地方没有测 量流入或流出水箱流量的设备，而只能测量水箱中的 水位(误差不超过 5%). 当水箱水位低于某最低水位 L 时，水泵抽水，灌入水箱内直至水位达到某最高水位 H 为止。但是也无法测量水泵的流量，因此在水泵启 动时无法立即将水箱中的水位和水量联系起来。水泵 一天灌水 1~2 次，每次约 2h. 试估计在任意时刻(包括 水泵灌水期间)t 流出水箱的流量 f(t)，并估计一天的总 用水量。 下表给出了某镇中某一天的真实用水数据，表中 测量时间以秒为单位，水位以 0.01E 为单位。例如在 3316s 时，水箱中的水深为 31.10E . 时间(s) 水位(0.01E) 时间(s) 水位(0.01E) 0 3175 46636 3350 3316 3110 49953 3260 6635 3054 53936 3167 10619 2994 57254 3087 13937 2947 60574 3012 17921 2892 64554 2927 21240 2850 68535 2842 25223 2795 71854 2767 28543 2752 75021 2697 32284 2697 79254 泵水 35932 泵水 82649 泵水 39332 泵水 85968 3475 39435 3550 89953 3397 43318 3445 93270 3340 已知水箱是直径为 57E、高为 40E 的正圆柱体。 当水位落到 27E 以下时，水泵自动启动把水泵入水箱； 当水位回升至 35.5E 时，水泵停止工作。 §5—3 水塔用水量的计算 5.3.1 问题分析 1．假设 28 次采集数据，得到 24 对数据(时刻：秒。水位： 英尺)，要估计水流出水箱的速度 f(t)(即：任意时刻 t 时的水流量)。先作几条假设。 (1). 只考虑通常情况下的正常用水，不考虑突发 事变，如水管破裂、自然灾害等。 (2). 全天 24 小时能保证用户的正常用水，不会发 生某个用户想用水时打开水龙头而没水流出这样的情 况。即：水箱中水的流出量 f(t)完全由用户的需求量决 定，与水箱自身的情况如水位高度、出口口径等因素 无关。 (3). 题中所给“水位约 27.00E 时，水泵启动；水
位约3550时,水泵停止;每次工作约2小时”这样 的已知信息说的比较含糊,故这里将其作为一条假设 同时也假设水泵的灌水速度为常数、水泵不会损坏或 不需要维修。 4)水箱的出水流量f(t),前面已假设是由用户决 定的。单个用户的用水流量是时间t的突变函数(即: 突然很大,突然为0),当大量的用户随机用水时(本题 的实际背景正是这样),总的用水流量函数f(t)就可以 近似当作是关于时间t连续的、光滑的,故假设水流 量函数f(t)阶连续可导 本题目标是计算出水箱的出水流量函数f(t)、及 天24小时总的用水量。已知条件为24个时刻的水位, 我们的想法是:分析确定水泵的两次启动、停止工作 的时刻,用数值计算的方法确定测量点时刻的水流量, 以这些测量点作为节点做三次样条插值,得出水流量 函数f(t 2.水泵的工作情况 为了便于计算,将时间单位变成小时h,长度单位 变成英尺E,见下表 时间)水位(E)时间h)[水位(E) 0 31.75 12954 33.50 0.921 31.10 13.876 32.60 1.843 30.54 3167 950 29.94 15904 30.87 3.871 2947 16826 30.12 4.978 28.92 17932 29.27 5.900 28.50 19.038 28.42 7.006 279 19959 2767 97 8968 26.97 22015 泵水 9981 22.958 泵水 10926 泵水 23.880 35.50 24.987 33.97 12.033 3445 25.908 33.40 为了方便分析,从上述28对信息中提取出一小部 分来: 时间(h)水位E)时间(h)水位(E) 8.968 26.97 20.839 26.97 9981 22015 泵水 10926 泵水 22.958 泵水 10954 35.50 23.880 根据已知,水位降低到约27E时,水泵开始启动 在t=8968h时刻,水位是26.97E,低于27E,水泵应 该工作了,但实际情况是此时刻水泵没工作。显然, 过了t=8968h后很短的时间内,水泵将开始工作,故 近似确定,t=8.968h是水泵第一次工作的启动时刻。 同理,近似确定t20.839h是水泵第二次工作的启动 时刻。 根据已知,水位上升到约355E时,水泵停止工 作。在t=10.954h时刻,水位是3550E,且水泵已停
位约 35.50 时，水泵停止；每次工作约 2 小时”这样 的已知信息说的比较含糊，故这里将其作为一条假设， 同时也假设水泵的灌水速度为常数、水泵不会损坏或 不需要维修。 (4). 水箱的出水流量 f(t)，前面已假设是由用户决 定的。单个用户的用水流量是时间 t 的突变函数(即： 突然很大，突然为 0)，当大量的用户随机用水时(本题 的实际背景正是这样)，总的用水流量函数 f(t)就可以 近似当作是关于时间 t 连续的、光滑的，故假设水流 量函数 f(t)一阶连续可导。 本题目标是计算出水箱的出水流量函数 f(t)、及一 天 24 小时总的用水量。已知条件为 24 个时刻的水位， 我们的想法是：分析确定水泵的两次启动、停止工作 的时刻，用数值计算的方法确定测量点时刻的水流量， 以这些测量点作为节点做三次样条插值，得出水流量 函数 f(t) . 2. 水泵的工作情况 为了便于计算，将时间单位变成小时 h，长度单位 变成英尺 E，见下表： 时间(h) 水位(E) 时间(h) 水位(E) 0 31.75 12.954 33.50 0.921 31.10 13.876 32.60 1.843 30.54 14.982 31.67 2.950 29.94 15.904 30.87 3.871 29.47 16.826 30.12 4.978 28.92 17.932 29.27 5.900 28.50 19.038 28.42 7.006 27.95 19.959 27.67 7.929 27.52 20.839 26.97 8.968 26.97 22.015 泵水 9.981 泵水 22.958 泵水 10.926 泵水 23.880 34.75 10.954 35.50 24.987 33.97 12.033 34.45 25.908 33.40 为了方便分析，从上述 28 对信息中提取出一小部 分来： 时间(h) 水位(E) 时间(h) 水位(E) 8.968 26.97 20.839 26.97 9.981 泵水 22.015 泵水 10.926 泵水 22.958 泵水 10.954 35.50 23.880 34.75 根据已知，水位降低到约 27E 时，水泵开始启动。 在 t=8.968h 时刻，水位是 26.97E，低于 27E，水泵应 该工作了，但实际情况是此时刻水泵没工作。显然， 过了 t=8.968h 后很短的时间内，水泵将开始工作，故 近似确定，t=8.968h 是水泵第一次工作的启动时刻。 同理，近似确定 t=20.839h 是水泵第二次工作的启动 时刻。 根据已知，水位上升到约 35.5E 时，水泵停止工 作。在 t=10.954h 时刻，水位是 35.50E，且水泵已停
止工作。再考虑到已知水泵每次工作约2小时,而 10954-8968=19862,实际情况是, 在22958时刻水泵还再工作,本次工作时长已经超过 2小时了,故可推断在此刻之后的很短时间内水泵就 要停止工作了。 理由(2):(考察邻近时间段的用水情况)从19959 至20.839时段,对应水位从2767降到2697,平均 降速为0.796E/h;从2388至24987时段,对应水位 从34.75降到33.97,平均降速为0.705E/h;实际中 T>22.958,从T至2388时段,对应水位从355降到 3475,平均降速为 35.5-34.7535.5-34.75 23.88-T23.88 229580.813Eh 注意到T离22958越远,水位降速会越比0.813大, 而0.813已明显比0.796、0.705大,该时段内的水位 降速不应该再大了,故推断T离22958很近。 基于上述两个理由,我们近似确定t=22958h是水 泵第二次工作的停止时刻,且此刻水位是35.50E (注:原来已知24对数据,现在又多了一对已知数据) 测量点时刻的水流量 t时刻的水流量为f(t),设水箱中水的体积为vt) 则,当水泵不工作的时候必有 d()(t) 把25个测量点时刻从小到大排列,记为 h1,2 t时刻的水位记为h,i=1,2,…,25 L时刻水箱中水的体积V(t1)记为V,则有 H=x{(2)角,1=12-235 1时刻的水流量f(1)记为f 根据已知数据(1,V)i=12,…,25,我们用“数值微分 的方法来计算∫: 数学知识:数值微分 问题1:y=f(x)连续可导,不知到表达式
止工作。再考虑到已知水泵每次工作约 2 小时，而 10.954-8.968=1.9862，实际情况是， 在 22.958 时刻水泵还再工作，本次工作时长已经超过 2 小时了，故可推断在此刻之后的很短时间内水泵就 要停止工作了。 理由(2)：（考察邻近时间段的用水情况）从 19.959 至 20.839 时段，对应水位从 27.67 降到 26.97，平均 降速为 0.796E/h；从 23.88 至 24.987 时段，对应水位 从 34.75 降到 33.97，平均降速为 0.705E/h；实际中 T>22.958，从 T 至 23.88 时段，对应水位从 35.5 降到 34.75，平均降速为 0.813 23.88 22.958 35.5 34.75 23.88 35.5 34.75 − − − − T E/h 注意到 T 离 22.958 越远，水位降速会越比 0.813 大， 而 0.813 已明显比 0.796、0.705 大，该时段内的水位 降速不应该再大了，故推断 T 离 22.958 很近。 基于上述两个理由，我们近似确定 t=22.958h 是水 泵第二次工作的停止时刻，且此刻水位是 35.50E . (注：原来已知 24 对数据，现在又多了一对已知数据) 3．测量点时刻的水流量 t 时刻的水流量为 f(t)，设水箱中水的体积为 V(t)， 则，当水泵不工作的时候必有 dt V t dt dV t f t ( ) ( ) ( ) = = − . 把 25 个 测 量 点 时 刻 从 小 到 大 排 列 ， 记 为 1 2 25 t ,t , ,t . i t 时刻的水位记为 h , i =1,2, ,25. i i t 时刻水箱中水的体积 ( ) i V t 记为 Vi ，则有 , 1,2, ,25 2 57 2 = V = h i i i i t 时刻的水流量 ( ) i f t 记为 i f . 根据已知数据 (t i ,Vi ), i =1,2, ,25 ，我们用“数值微分” 的方法来计算 i f ： 数学知识：数值微分 问题 1：y=f(x)连续可导，不知到表达式
已知 x1<x2 xn 求y=f(x)在每个节点x处的一阶导数 分析:学过插值法(如: Lagrange插值、三次样条插 值)和最小二乘法的人,很容易想到用这样的 方法:先拟合出y=f(x)的近似表达式,再求导 函数,再计算导函数在每个节点处的值。此法 不好,应该直接用已知数据(x,y)计算节点x 处的一阶导数 数值微分法:节点x处的一阶导数记为∫ i=1,2时 f -3+4y-y12 3≤≤n-2时:f=-M+2+8-8+y-2 2(x+1-x) 1,n时:f y 2( 问题2:y=f(x)表达式已知,且可导,但求不出导函数 f(x),要求任意一点x处的一阶导数值 数值微分法:取小步长h f(x0)= f(x0+2h)+8f(x+h)-8f(x0-h)+f(x0-2h 12h 补充知识完毕 将25个节点分成3个组 第一组 第二组 第三组 1,12,……,l t1,412,…,l2 l2,l23,l24,l25 每组的最前两点,用公式 ,= 3:+4V.,-V 2 每组的最末两点,用公式: f,= 2(1-l1-1) 每组的其余点,用公式: ,= ,+8 12(t1 5.3.2模型建立与求解 前面已经得到测量点时刻的水流量
已知 n n y y y y x x x x 1 2 1 2 : : ， 求 y=f(x) 在每个节点 i x 处的一阶导数。 分析：学过插值法（如：Lagrange 插值、三次样条插 值）和最小二乘法的人，很容易想到用这样的 方法：先拟合出 y=f(x)的近似表达式，再求导 函数，再计算导函数在每个节点处的值。此法 不好，应该直接用已知数据 ( , ) i i x y 计算节点 i x 处的一阶导数。 数值微分法：节点 i x 处的一阶导数记为 i f i=1,2 时： 2( ) 3 4 1 1 2 i i i i i i x x y y y f − − + − = + + + 3 i n−2 时： 12( ) 8 8 1 2 1 1 2 i i i i i i i x x y y y y f − − + − + = + + + − − i = n −1, n 时： 2( ) 3 4 1 1 2 − − − − − + = i i i i i i x x y y y f . 问题 2：y=f(x)表达式已知，且可导，但求不出导函数 f (x) ，要求任意一点 0 x 处的一阶导数值。 数值微分法：取小步长 h h f x h f x h f x h f x h f x 12 ( 2 ) 8 ( ) 8 ( ) ( 2 ) ( ) 0 0 0 0 0 − + + + − − + − = 补充知识完毕. 将 25 个节点分成 3 个组 第一组 1 2 10 t ,t , ,t 第二组 11 12 21 t ,t , ,t 第三组 22 23 24 25 t ,t ,t ,t 每组的最前两点，用公式： 2( ) 3 4 1 1 2 i i i i i i t t V V V f − − + − = + + + 每组的最末两点，用公式： 2( ) 3 4 1 1 2 − − − − − + = i i i i i i t t V V V f 每组的其余点，用公式： 12( ) 8 8 1 2 1 1 2 i i i i i i i t t V V V V f − − + − + = + + + − − 5.3.2 模型建立与求解 前面已经得到测量点时刻的水流量
用三次样条插值可得水流量函数f=f() 对f(t)积分,可得任何一个时间段的用水总量 下面程序计算并画出f=f(1)图象(程序中的数据t0,v0 取自书P67之表54) clear tO=0,0.921,1,8432.949,3.871,4.978,5.9,7.006,7928,89 67,10.954,12.032,12954,13875,14.982,15.903,16826,1 7931,19.037,19959,20.839,22958,23.88,24986,25908] v0=[2294,2247,2206,2163,2129,2089,2059,2020,1988,19 48,2564,2489,2420,2355,2288,2230,2176,2114,2053,19 9,1948,2564,2510,2454,2413 for F1: 25 vitch case{1,2,112,22,23} fo(i)=abs(-3*vO()+4*vo(计+1)v0(+2)(2*tO(i+1)tO(i) ) case{9,10,20,21,24.25} f0()=abs(3*vO()-4*v0(-1)+v0(-2)/(2*(tO()-t0(-)) otherwise f0()=abs(V0(计+2)+8*v0(计+1)-8*v0(-1)+vO(-2)/(12*(t 0(+1)to() end t=0: 0.05: 26, f-spline(to, f o, t): plot(t, f), grid 上面程序中,时间t以小步长005遍历计算得到 了每个点处的f值,用“数值积分中的复化梯形法” 计算区间[0,24]上定积分就得一天总用水量。程序如 下 (续上面程序内容) zyl=0; e2ys=ysl+005*f0+0+1)2 sI 计算结果为:zs=1.2582e+003 5.3.3模型的检验 从几个侧面,分别对模型进行检验 1.以不同时刻作为起点,计算总用水量 下面程序实现书P70表56 (续上面程序内容) for qd=0.250.25:1.7 j0=qd/0.0
(t , f ), i =1,2,...,25. i i 用三次样条插值可得水流量函数 f = f (t). 对 f (t) 积分，可得任何一个时间段的用水总量。 下面程序计算并画出 f = f (t). 图象（程序中的数据 t0, v0 取自书 P67 之表 5.4）： clear t0=[0,0.921,1.843,2.949,3.871,4.978,5.9,7.006,7.928,8.9 67,10.954,12.032,12.954,13.875,14.982,15.903,16.826,1 7.931,19.037,19.959,20.839,22.958,23.88,24.986,25.908] ; v0=[2294,2247,2206,2163,2129,2089,2059,2020,1988,19 48,2564,2489,2420,2355,2288,2230,2176,2114,2053,199 9,1948,2564,2510,2454,2413]; for i=1:25 switch i case {1,2,11,12,22,23} f0(i)=abs((-3*v0(i)+4*v0(i+1)-v0(i+2))/(2*(t0(i+1)-t0(i)) )); case {9,10,20,21,24,25} f0(i)=abs((3*v0(i)-4*v0(i-1)+v0(i-2))/(2*(t0(i)-t0(i-1)))); otherwise f0(i)=abs((-v0(i+2)+8*v0(i+1)-8*v0(i-1)+v0(i-2))/(12*(t 0(i+1)-t0(i)))); end end t=0:0.05:26;f=spline(t0,f0,t);plot(t,f),grid 上面程序中，时间 t 以小步长 0.05 遍历计算得到 了每个点处的 f 值，用“数值积分中的复化梯形法” 计算区间[0,24]上定积分就得一天总用水量。程序如 下： （续上面程序内容） zysl=0; for j=1:480 zysl=zysl+0.05*(f(j)+f(j+1))/2; end zysl 计算结果为：zysl= 1.2582e+003 5.3.3 模型的检验 从几个侧面，分别对模型进行检验： 1．以不同时刻作为起点，计算总用水量 下面程序实现书 P70 表 5.6 （续上面程序内容） for qd=0.25:0.25:1.75 j0=qd/0.05;
zsl=0 forj=j0+1j0+480 qd, ysl end 计算结果为 qd=0.2500zysl=1.2593e+003 qd=0.5000zysl=1.2608e+003 qd=0.7500zys=1.2627e+003 qd=1.2500zysl=1.267le+003 zysl=1.2688e+003 qd=1.7500zysl=1.2698c003 2.分三段,实际用水量与模型计算量比较。 略算,看书P70表5.7 3.估计水泵的功率 略算,书中得 水泵第一次工作:充水量=733446,工作时长=1987 水泵第一次工作:充水量=736.705,工作时长=2.119 733446 第一次功率= 369 1.987 第一次功_736705 2.119 两次功率的相对误差为 369-348 348 §5-3二维插值的应用 (黑板:简略介绍“二维插值问题”) 已知节点(x,y)及节点处的函数值==f(x,y),求 二元函数=f(x,y),并且画曲面图形。 关键: interp2( 下面程序实现书P72之图: clear yf=1:12;wd=5:10:85 qX=[24,18.7,20.8,22.1,373,48.2,256,53,0.3;1.6,214,18 5,20.1,288,366,24.2,5.3,0; 2.4,16.2,18.2,20.5,27.8,35.5,25.5,54,03.2,9.2,16.6,25.1 37.2,40,24.6,4.90.3; 1.0,2.8,12929.2,403,37.6,21.1,4.9,0:0.5,1.7,10.1,32.64 7,354,22.2,7.1,0 04,14,83,33.0,46.2,35,20.2,5.3,0.102,24,11.2,310,39 9,347,21.2,7.3,0.2;
zysl=0; for j=j0+1:j0+480 zysl=zysl+0.05*(f(j)+f(j+1))/2; end qd,zysl end 计算结果为： qd =0.2500 zysl =1.2593e+003 qd =0.5000 zysl =1.2608e+003 qd =0.7500 zysl =1.2627e+003 qd =1 zysl =1.2649e+003 qd =1.2500 zysl =1.2671e+003 qd =1.5000 zysl =1.2688e+003 qd =1.7500 zysl =1.2698e+003 2．分三段，实际用水量与模型计算量比较。 略算，看书 P70 表 5.7 3．估计水泵的功率 略算，书中得： 水泵第一次工作：充水量=733.446，工作时长=1.987 水泵第一次工作：充水量=736.705，工作时长=2.119 第一次功率 369 1.987 733.446 = = 第一次功率 348 2.119 736.705 = = 两次功率的相对误差为 6% 348 369 348 = − §5—3 二维插值的应用 （黑板：简略介绍“二维插值问题”） 已知节点 ( , ) i j x y 及节点处的函数值 ( , ) ij i j z = f x y ，求 二元函数 z = f (x, y) ，并且画曲面图形。 关键：interp2(…………) 下面程序实现书 P72 之图： clear yf=1:12;wd=5:10:85; qx=[2.4,18.7,20.8,22.1,37.3,48.2,25.6,5.3,0.3;1.6,21.4,18 .5,20.1,28.8,36.6,24.2,5.3,0; 2.4,16.2,18.2,20.5,27.8,35.5,25.5,5.4,0;3.2,9.2,16.6,25.1, 37.2,40,24.6,4.9,0.3; 1.0,2.8,12.9,29.2,40.3,37.6,21.1,4.9,0;0.5,1.7,10.1,32.6,4 1.7,35.4,22.2,7.1,0; 0.4,1.4,8.3,33.0,46.2,35,20.2,5.3,0.1;0.2,2.4,11.2,31.0,39. 9,34.7,21.2,7.3,0.2;
0.5,58,12.5,28.6,35.9,357,22.6,7,0.3:0.8,9.2,21.1,32.0,4 03,39.5,28.5,86,0; 24,103,23.9,28.1,38.2,40,253,63,0.1;3.6,16,25.5,256,4 34,41.9,24.366,0.3] [xy]= meshgrid(1:1:12,5:585) ZF=interp(yf, wd, x, x,y, cubic) mesh(x,y, z)
0.5,5.8,12.5,28.6,35.9,35.7,22.6,7,0.3;0.8,9.2,21.1,32.0,4 0.3,39.5,28.5,8.6,0; 2.4,10.3,23.9,28.1,38.2,40,25.3,6.3,0.1;3.6,16,25.5,25.6,4 3.4,41.9,24.3,6.6,0.3]; [x,y]=meshgrid(1:1:12,5:5:85); z=interp2(yf,wd,qx',x,y,'cubic'); mesh(x,y,z)