屋檐水槽构建模型555

1、 组别：第六小组 姓名：梁伟仙、叶佳盛、冼深发屋檐水槽模型摘要 对于倾斜屋顶采用何种结构进行有组织排水，是民用建筑一个普通而有实际的问题。在民用建筑中一般是在房顶的边缘安装一个檐槽和一个竖立的排水管来排水，这种排水方式的可行性探讨实际上是解决水槽的容量在单位时间内能否足以排除雨水的问题。构建这种排水方式的微分方程模型，运用Maple对设计方案进行了可行性论证，提出了屋檐水槽模型的优化方案。关键词: 速度平衡原理；微分方程；数值解；maple一、问题的重述 为了雨天出入方便，房屋管理部门想在房顶边缘安装一个檐槽。现在有一个公司想承担这项业务，他们承诺：提供一种新型的可持久的檐槽，这种檐槽不管天气

2、如何都能排掉房顶的雨水。房屋部门希望检验公司的承诺能否实现。简单来说, 从屋脊到屋檐的房顶可以看成是一个12米长,6米宽的矩形平面.房顶与水平方向的倾斜角度要视具体的房屋而定,一般来说,这个角度通常在20°50°之间。现在有一个公司想承接这项业务,他们允诺:提供一种新型的可持久的檐槽,它包括一个横截面为半圆形(半径为7.5厘米)的水槽和一个竖直的排水管(直径为10厘米),并且不管天气情况如何,这种檐槽都能排掉房顶的雨水。如图1所示。 二、模型的假设2.1、合理假设 (1)、降雨分布均匀并以垂直下落，并且直接落在房顶上；(2)、所有落在房顶上的雨水迅速流入水槽中；(3)、雨水

3、不断从水槽中溅出；(4)、排水管道顺畅，没有任何障碍或阻塞；(5)、假设雨开始下时槽内有雨水，且水深高度为0.01m。2.2、符号说明 ： 如下图表所示 各因素与符号说明有关的因素因素类型 符号 单位降水强调输入变量rm/s时间变量ts房顶的倾斜度输入参数 弧度房顶长度输入参数dm房顶的宽度输入参数bm水槽的半径输入参数am水槽中水的深度输入参数hm水槽中水的容量变量Vm流入水槽的流速变量Qm/s流出水槽的流速变量Qm/s排水管横截面积参数Am重力加速度常数gm/s 表1三、模型的建立参照图1, 根据速度平衡的原理, 对于房顶排水系统有:水槽中水的流量的变化率= 雨水的流入流量- 排水流出的流

4、量。即Vc(t)= Q- Q,这里Q、Q 分别是单位时间流入水槽和从水槽流出的雨水流量。房顶雨水的流动情况, 如图2 所示。房顶的面积是bd, 由于房顶是倾斜的, 根据合理假设(1), 实际受雨的水平面积应为bdcos, 房顶上雨水的流量就是r(t)bdcos,雨水的流动是沿倾斜的房顶向下的, 从而流入水槽的流量应该是它在铅垂方向的分量 , 即Q= r( t)cossin (1) 而根据能量守恒定律 得出， （2） 所以，Vc(t)= Q- Q水槽中水的深度h<a时的状况如图3所示，槽中水的体积为： v(t)= （3）由图3有： cos （4） 由（4）式得： （5） sin2=2(a-

5、h) (6)把(5)、（6）式代入(1)式得： V(t)=adarcos()- (7)下面运用maple对（7）式中t求导【用v(t)表示d(v)/d(t)】：>V(t):=a2*d*(arcos(a-h(t)/a)-(a-h(t)*sqrt(2*a*h(t)-h(t)2)/a2)定义v(t);V(t):=adarccos->V1：=diff(v(t),t);#求v（t） >v2:=a2*d(diff(h(t),t)/(2*a*h(t)-h(t)2)(1/2)+diff(h(t),t)*(2*a*h(t)-h(t)2)(1/2)/a2-1/2*(a-h(t)/(2*a*h(t

6、)-h(t)2)(1/2)/a2*(2*a*diff(h(t),t)-2*h(t)*diff(h(t),t);>v(t):=simplify(%);#对上式化简 v(t):=- （8）所以根据（1）、（2）、（8）这样就得到模型： （9）即 （10）： 四、 模型的求解与分析设定一组数据：a=0.075m,b=6m,d=12m,g=9.8m/s,A=0.0025m,h(0)=0.01m（假设槽内有些积水）。下面，用Maple 求（10）数值解。 >restart;eq:=diff(h(t),t)=(r*b*d*sin(alpha)*cos(alpha)-A*sqrt(2*g*h(t

7、)/(2*d*sqrt(2*a*h(t)-(h(t)2);#定义（10）式；>eq:=subs(a=0.075,b=6,d=12,g=9.8,A=0.0025*Pi,alpha=Pi/6,%);把这组数据代到上式eq := > eq:=evalf(eq);化简针对房屋管理部门的要求，所以我们可能遇到两种情况：1、r(t)为常数；2、r(t)为周期函数：（1） 当r(t)为常数时。考虑水槽的深度趋于一个低于0.075m的稳定值，即h(t)=0> eq1:=0 = .4166666667e-1*(31.17691454*r-.3477105892e-1*h(t)(1/2)/(.1

8、50*h(t)-1.*h(t)2)(1/2);eq1:=0=> solve(eq1,h(t);8.03952576210r> eq2:=803952.5762*r2>=0.075; eq2:=0.075<=8.03952576210r> solve(eq2,r);r<=-0.0003054326190,0.0003054326190<=r 根据以上数据，我们可以得知，当r=0.00305时，水不溢出；当r>0.000305时，水溢出。 所以我们可以令r=0.00021; r=0.000305; r=0.00031 分别代入eq中求出其数值解，并作

9、出h(t)的图如图（4）所示；>eq3:=diff(h(t),t)=0.4166666667e-1*(31.19458600*r-0.3477105892e-1*h(t)(1/2)/(0.150*h(t)-1*h(t)2)(1/2);#定义方程h(t)= >eq3:=subs(r=0.00021,eq)=0； #将r=0.00021代入eq3并令其为0 =）>initvals:=h(0)=0.01>s1:=dsolve(eq4,initvals,h(t),type=numeric)；#求方程eq3的数值解 h(0)=0.01 proc(x_rkf45)end proc&

10、gt;sl(5)；#可以看出s1的第2项为好h(t)>plot(rhs(s1(t)2,t=0.110):所以根据图（4）知道，当r=0.00021m/s时，h为t的增函数，h的最大值在0.035m附近；而根据图（5）知道当r=0.000305时，h为t的增函数，h的最大值为0.074m左右；而根据图（6）知道当r=0.00031m/s时，h为t的增函数，h的最大值超过0.075m。 由于各地具体情况不同,各地气象预报部门对于当地各类降水的标准也有些自己的规定。一般而言, 当地气象部门规定24 h降水量在60 mm 以上的雨为特大暴雨。对于r(t)= 常数这种情形, r> 0.000

11、 305 m/s 的强降雨机率几乎为0,因此,这个公司的承诺是能兑现的。（2）当r(t)为周期函数时，并且设为正弦函数： r(t)= ，h（t）=0.01这表明下雨过程是在60s内发生的一个短促的强震雨行为，最大的降雨强度是0.000001m/s，由方程eq3得到如下的微分方程：而且h(0)=0.01这样我们可以运用上面的数值解法，得到h(t)的图（7）。 而且我们可以从图（7）中可以看出，h(t)的最大值才只有0.01m，也就是水槽的高度，而上图的最大值为0.075m，所以不会造成溢出现象，对于第2种情况，水槽的水也不会出现溢出的情况，所以这个公司的承诺可以兑现。五、 对模型的优化并且改进：

12、如果是长时间下暴雨，可能会造成溢出现象，因此有两种优化方案；a):增大排水管的横截面积A,即增大排水管的半径。当排水管半径 增大为0. 056 m,对于模型的第1 种情形,r 0. 000 36 m/ s, 水槽的水不会出现溢出的情况。图8 是h( t)随时间t 变化的形。h( t)的最大值不会超过0.07 m。对于模型的第2 种情形, 水槽的水不会溢出。图9 是h( t)随时间t 变化的图形。h( t)的最大值不会超过0.01 m。b):改变水槽的连接方式, 如图10, 让水槽往屋檐倾斜一定角度, 这相当于增加水槽的容水高度。图9 改进后的h随t 变化示意参考文献1刘斌中国三农问题报告M北京

13、：中国发展出版社，2O04 2黎诣远微 经济分析M北京：清华大学出版社，20(13 3葛志华wro与当代中国农民M南京：江苏人民出版社，2001 4鲜祖德，农民收入增长问题研究M北京：中国经济出版社，2002附件1restart;eq:=diff(h(t),t)=(r*b*d*sin(alpha)*cos(alpha)-A*sqrt(2*g*h(t)/(2*d*sqrt(2*a*h(t)-(h(t)2);> eq:=subs(a=0.075,b=6,d=12,g=9.8,A=0.0025*Pi,alpha=Pi/6,%);> eq:=evalf(eq);> eq1:=subs

14、(diff(h(t),t)=0,eq);> h(t):=solve(eq1,h(t);> eq2:=803952.5762*r2>=0.075; solve(eq2,r);> restart;eq := diff(h(t),t) = .4166666667e-1*(31.17691454*r-.3477105892e-1*h(t)(1/2)/(.150*h(t)-1.*h(t)2)(1/2);> eq3:=subs(r=1/2000*sin(Pi/40*t),eq);eq4:=evalf(eq3);> s1:=dsolve(eq4,h(0)=0.01,h(t),numeric);> s1(5);#可以看出是（1）第二项为h(t)> plot('rhs(s1(t)2)',t=0.50);> eq4:=subs(r=0.000305,eq);>