数学建模报告——水的流出时间.doc

数学建模实验报告水的流出时间软件84周天牧08161102同组完成者 俞乐晨 参与讨论课讨论【问题描述】三个横截面积为常数A，高分别为H1，H2，H3的水池内都盛满了水，都由池底一横截面积为B的小孔放水。求在任意时刻的水面高度和将水放空所需的时间。【模型假设】由于我们建立的是理想化的模型，做如下假设：1. 题目中的水流从前一个水池流出后立即流入后一个水池化为水位的上升，在中间飞行的时间忽略不计。2. 题目中水池底部的出水孔很小，无论水位多么低都可以充满以保证流出速度，除非水位为0 有了这两个假设，我们可以开始建立模型。【模型建立】将此问题抽象成数学问题：求初值分别为H1，H2，H3的函数h1，h2，h3随时间变化的函数，以及他们变为0所需的时间（设水池1，2，3的流出速度为v1，v2，v3）。由机械能守恒定律可知，水从第i个水池小孔流出的速度为v（i）=sqrt（2*g*h（i）根据任何一个水池内水量在一个时间微元内的减少量等于流出量减去流入量可以得到如下关系：水池1：-Dh1*A= Ds1*B将D写成d，两边同除dt可得两边积分得把水放空所需时间为水池2：-Dh2*A= Ds2*B- Dh1*A注意到|Dh1*A|= |Ds1*B|,故两边写D为d并同除以dt可得化简并带入v1的函数表达式可得这是一个非线性常微分方程，难以得到解析解（并非不可求，但是在此处求出解析解并不是建模的重点，因为即使h2存在解析解，h3也不一定存在。后面给出了求此方程解析解的方法。），在这里我们采用数值方法进行分析。由于第二步的数值会对后续的水池产生影响，我们选用了精度较高的Runge-Kutta法（龙格-库塔法，以下简称RK法，具体方法参见数值计算方法第六章），matlab中给出了实现的函数：t,y=ODE45(odefun,tspan,y0)其中t，y用于输出广义时间与相空间向量，odefun为要求的微分方程，tspan为广义时间区间，y0为初值向量。定义函数：function y=funt(t,x)A=2;B=1;g=10;H=10;y=-1*B/A*(sqrt(2*g*x)-sqrt(2*g*H)+B/A*g*t)end再调用ode45函数，即可得到h2的数值解，即一组h2关于t的向量。需要注意的是，由于非线性微分方程解的稳定性问题，此迭代法是局部收敛而非全局收敛的，故当第二个水池流空了而第一个水池还未流空时会产生复数解，需舍去。水池3： -Dh3*A= Ds3*B- Dh2*A同理可得由于只有三步，此处我们采用简单一些的欧拉法(Euler法，具体方法见数值计算方法第六章)进行数值求解。采用此方法可以得到h3的数值解，即h3关于时间的离散函数。【代码实现】在这里给出欧拉法的代码实现j=1;h3(1)=H;arrayt3(1)=0;temp3=0while(ji-1)temp3=h3(j)+(h2(j)-h2(j+1)*(-1*B/A*(sqrt(2*g*h3(j)-sqrt(2*g*h2(j)h3=h3,temp3arrayt3=arrayt3,arrayt2(j+1)j=j+1end【模型分析】有了以上的方法，我们可以得到一张h1,h2,h3关于时间t的数值表格（前提是给定初值H1，H2，H3。其实最理想的情况其实是得到h1，h2，h3关于时间t的带参数H1，H2，H3的函数，但由于水平所限，无法用mathmatica软件实现RK方法和Euler方法，而matlab又不支持字母运算，故用数值代替之）。取H0=10，讨论以下初始情况下三个水池水位的关系：分别是三个水池初始量相等的情况，一多两少的情况和一少量多的情况H1H2H31 H0H0H02 H0H0/2H0/23 H0/2H0H0/24 H0/2H0/2H05 H0/2H0H06 H0H0/2H07H0H0H0/2情况一：H1=10,H2=10,H3=10 t1=2.83 t2=3.6492 t3=4.6333 情况二：H1=10,H2=5,H3=5 t1=2.83 t2=3.4458 t3=4.3738 情况三：H1=5,H2=10,H3=5 t1=2 t2=3.6 t3=4.6633 情况四：H1=5,H2=5,H3=10 t1=2 t2=3.4599 t3=4.4581 情况五：H1=5,H2=10,H3=10 t1=2 t2=3.6 t3=4.6751 情况六：H1=10,H2=5,H3=10 t1=2.83 t2=3.4458 t3=4.4707 情况七：H1=10,H2=10,H3=5 t1=2.83 t2=3.6492 t3=4.6233 【模型讨论】由以上的建模过程可以得到这样一些结论h1是单调减函数，并且其变化规律只与其初值有关系，给定初值的条件下流空的时间是一常数。由通式-dhi/dt=(v(i)-v(i-1)*B/A可知，任一时刻当第i个水池的水面高度大于第i-1个时，第i个水池的水面高度下降。反之当第i个水池的水面高度小于第i-1个时，第i个水池的水面高度上升。所以永远不会出现第i-1个水池流完之前，第i个水池就流完的了的情况。当外界有水以恒定速度流入的时候，经过足够长的时间，最终每个水池的高度会趋于一致，每个出水口的速度也会趋于一致当第2个和第3个水池的初始高度相同时，无论H1取何值，h3在相当长一段时间内都可以视为不变或者变化很小，可以理解为第2个水池为第3个水池做了缓冲。不难想象，当用n个相同初始高度的水池做缓冲系统时，第n个水池高度维持时间会更长，即效果会更好【模型改进】1.第二步到第三步的所采用的Euler法是最为简单的一种数值方法，但是精度较低，可以用其他方法改进。2.h3所有的分量均由h2参与迭代运算产生，故当h2为0之后，h3也停止了迭代。故在我们绘制的图中h3没有到0就停止了。可以改进为，当t时刻h2为0之后，h3采用h3(t)为初值，带入h1的表达式进行计算计算零点。3. 这里给出第二个方程的解析解法：由于是一个常数，我们将其简化为k。设h2=y2（sqrt(H2)=y）,则原式变为即两边积分计算得左边为h2-H2，右边为其中故解得h2关于t的表达式为：【问题拓展】当水池个数由3个变为N个时，跟以上第二个，第三个水池原理相同。在问题分析中所得到的结论也