第六讲:流体基本模型解析解——二维、稳定源1、非稳定源二维模型解析解(比照一维解析解的形式)2、基本模型稳定源解析解3、解析解的MATLAB实现方法第六讲/272一、非稳定源二维模型解析解根据物质平衡原理,流场中污染物浓度C(x,t)和瞬时投放污染物量M成正比。仿照一维模型的动态比例系数,可确定二维模型的动态比例系数:用tDDhtDtDhyxyx41441反映水深和弥散系数不同造成浓度随时间的变化。在一维纵向扩散中,考虑了横向(宽度B)和垂向(深度h)(其乘积为横截面积A)的影响,因为纵向的影响已经包括在后面的比例因子中。所以在二维扩散中,在该比例因子中考虑垂向(水深h)的影响,横向和纵向扩散的影响,包括在后面的比例因子中。用公式tDtuxxx4)(exp2×tDtuyyy4)(exp2反映由于介质在纵向和横向流动和弥散造成的浓度随时间的变化。同样用)exp(Kt反映由于污染物降解造成的浓度随时间变化。综合以上分析,可以得到瞬时投放的二维模型解析解:tDDhMtyxCyx4),,(tDtuxxx4)(exp2×tDtuyyy4)(exp2×)exp(Kt1、中心投放-边界无限对上述变换中的原函数乘以常数a/2,则根据拉普拉斯变换是线性变换性质可得:L[)4/exp(2)exp(2tattbta]=dpbpbpaap2)exp(=exp(-abp)所以可得到:L-1[J(x,p)=)4/exp(2)exp(20tattbtaCexp(ux/2Dx)=)exp()244exp(2220KtDutDxDtuttDxCxxxxxx=)exp(4)(exp220KttDtuxtDCuxxxx根据瞬时投放时初始浓度计算方法:C(x,t)=)exp(4)(exp42KttDtuxtDAMxxx一维瞬时源第六讲/273一、非稳定源二维模型解析解2、中心投放-边界有限如果污染物扩散在横向的扩散空间是有限的,则边界要反射。结果是造成环境介质中污染物浓度的增加。这种反射作用可以用虚拟的点源来考虑。把边界作为一个理想的反射镜面,在实源的对称位置设立一个源强与实源相等的虚源。虚源的设立可以把由于边界反射作用造成环境介质中污染物浓度的增加考虑进模型中。第六讲/274一、非稳定源二维模型解析解2、中心投放-边界有限河宽为B,中心瞬时投放:虚源1造成的浓度为:tDDhMtyxCyx4),,(1tDtuyBtDtuxyyxx4)(4)(exp22)exp(Kt虚源2造成的浓度为:tDDhMtyxCyx4),,(2tDtuyBtDtuxyyxx4)(4)(exp22)exp(KtBB/2B/2xy虚源2虚源1实源河宽为B中心排放情形·(x,y)第六讲/275一、非稳定源二维模型解析解2、中心投放-边界有限所以总浓度为:tDDhMtyxCyx4),,()exp(KttDtuyBtDtuxtDtuyBtDtuxtDtuytDtuxyyxxyyxxyyxx4)(4)(exp4)(4)(exp4)(4)(exp222222这就是中心投放边界有限时的解析解。第六讲/276一、非稳定源二维模型解析解3、边界投放-边界无限边界投放情况下,最多只有一个边界反射。在另一边界无限的情况下,相当于在上述问题中不考虑虚源2,并将虚源1和实源在边界处汇合。相应地对原来的模型,去掉虚源2造成的浓度影响项tDtuyBtDtuxyyxx4)(4)(exp22。同时,令B=0(虚源和实源在边界汇合),即得边界投放,另一边界无限时的解析解:tyDxDhMtyxC42),,(tDtuytDtuxyyxx4)(4)(exp22)exp(Ktxy边界源边界排放宽度无限根据物质平衡原理,流场中污染物浓度C(x,t)和瞬时投放污染物量M成正比。仿照一维模型的动态比例系数,可确定二维模型的动态比例系数:用tDDhtDtDhyxyx41441反映水深和弥散系数不同造成浓度随时间的变化。在一维纵向扩散中,考虑了横向(宽度B)和垂向(深度h)(其乘积为横截面积A)的影响,因为纵向的影响已经包括在后面的比例因子中。所以在二维扩散中,在该比例因子中考虑垂向(水深h)的影响,横向和纵向扩散的影响,包括在后面的比例因子中。用公式tDtuxxx4)(exp2×tDtuyyy4)(exp2反映由于介质在纵向和横向流动和弥散造成的浓度随时间的变化。同样用)exp(Kt反映由于污染物降解造成的浓度随时间变化。综合以上分析,可以得到瞬时投放的二维模型解析解:tDDhMtyxCyx4),,(tDtuxxx4)(exp2×tDtuyyy4)(exp2×)exp(Kt中心投放边界无限第六讲/277一、非稳定源二维模型解析解4、边界投放-边界有限对比上述模型:边界无限时边界排放的污染物浓度相当于中心排放时的2倍。边界有限情况下,只有一个边界反射,叠加实源和虚源的共同作用可得:tDDhMCyx4)exp(Kt×tDtuyBtDtuxtDtuytDtuxyyxxyyxx4)2(4)(4)(4)(exp2222河宽为B中心排放情形Bxy虚源边界实源Bexp第六讲/278二、非稳定源三维模型解析解同样类似于二维模型解析解的推求方法,可得均匀流场中,典型三维模型解析解为:zyxEEEtMtzyxC3)(8),,,(tEtuztEtuytEtuxzzyyxx4)(4)(4)(exp222×)exp(Kt式中Ex、Ey、Ez分别为x、y、z方向的湍流扩散系数。由于三维模型中考虑了流速的三维分量,没有流速的断面平均,因此不出现弥散项,仅出现湍流扩散项。三维模型在一般的水环境预测计算中应用不多,但是计算污染物在大气湍流扩散中的浓度分布应用较多。第六讲/279三、稳定源的解析解稳定源是指污染物的排放源是连续的,其源强大小不随时间变化。在流体环境介质中,污染物处于不断的流动状态,由于排放源的稳定性,环境介质中污染物的分布状况也是稳定的,这时污染物在某一空间位置浓度不会随时间变化,即dC/dt=0,所以模型的求解过程相对简单。上述稳定源排放下模型的解析解一般称为稳态解。由于稳态情况下模型的参数估计和求解都比较简单,因此有时为了研究问题的方便,如果污染物浓度在一个较长时间内平均值保持一定的稳定性,也可以把此类情况作为稳态来处理,一般这样就可以满足环境规划和管理的要求。第六讲/2710三、稳定源的解析解1.零维模型的稳态解稳态条件下的零维模型(dtdC=0):kCVQCQC0=0所以得:TKCKQVCKVQQCC11000式中T=Q/V称为理论停留时间。如果一个河流系列中考虑多个完全混合单元,并假设河水流过各个混合单元的时间相同,则上述解可变形得:TKCC101,TKCC112,。。。,TKCCnn11所以nnTKCC110。第六讲/2711三、稳定源的解析解2.一维模型的稳态解第六讲/2712三、稳定源的解析解3.二维模型:中心投放,无边界二维模型:22xCDx+22yCDy-xCux-yCuy-KC=0一般情况下纵向(x方向)推流迁移远远大于弥散作用,而在横向(y方向)弥散作用又远远大于推流迁移,所以可以忽略Dx和uy:22yCDy-xCux-KC=0该方程可以用傅立叶变换(或Laplace变换)求解,也可用一种简便方法求解:第六讲/2713三、稳定源的解析解3.二维模型:中心投放,无边界设污染源为点源,单位时间内释放的污染物质量为Q(即为源强,g/s)。将介质的流动看作是一系列平行条块组成,这些条块以速度u通过点源时,接纳扩散质团质量为Q△t,而△t=△x/ux,为该条块通过点源的时间。按瞬时源一维扩散(在y方向)的浓度分布模型可写出该条块运行t时刻后的浓度分布为:C(y,t)=)exp(4exp42KttDytDxhtQyy将t=x/ux,ux=△x/△t代入上式得到:C(x,y)=)exp(4exp/42xyxxyxuxKxDyuuxDhuQ这就是二维模型的稳态解。第六讲/2714三、稳定源的解析解4.二维模型:中心投放,有边界如果在污染源中心排放的河宽有限(B)条件下,污染物的扩散受到两岸边界的反射,如果两岸的反射是连续的,则反射就会成为连锁式:)exp()4)(exp()4)(exp()4exp(/412122xniyxniyxyxxyxuKxxDyiBuxDyiBuxDyuuxDhuQC一般在计算中取n=3或4,就能达到计算精度。第六讲/2715三、稳定源的解析解5.二维模型:边界投放,无反射在边界排放条件下,不考虑边界反射,就是假设边界无限宽:)exp(4exp/422xyxxyxuKxxDyuuxDhuQC可见:边界排放条件下,污染物浓度较大(2倍)第六讲/2716三、稳定源的解析解6.二维模型:边界投放,有反射如果设宽度为B,则考虑边界反射条件下的解为:)exp()4)2(exp()4)2(exp()4exp(/4212122xniyxniyxyxxyxuKxxDyiBuxDyiBuxDyuuxDhuQC第六讲/2717三、稳定源的解析解7.三维模型解析解类同于二维模型的推导过程,只考虑纵向流速、横向和垂向湍流扩散,忽略纵向湍流扩散、横向和垂向流速的情况下,得到三维模型的稳态解为:)exp(44exp422xzxyxzyuKxxEzuxEyuEExQC式中Ey、Ez分别为横向和垂向的湍流扩散系数。第六讲/2718四、解析解的MATLAB实现1、回顾MATLAB中,常微分方程(组)的解析解是由函数dsolve()实现的,调用格式为:y=dsolve(‘eq1’,’eq2’,…,’eqm’)或y=dsolve(‘eq1’,’eq2’,…,’eqm’,‘x’)其中eqi既可以表示微分方程,也可以表示初始边界条件;并用D3y这样的记号表示d3y/dx3,用D2y(3)=5这样的记号表示5322xdxyd;其默认自变量为t,如果不是t而是x,则可以第二种格式中指明。第六讲/2719四、解析解的MATLAB实现2、简单举例(1)dy/dx=1+y2的通解以及y(0)=1的特解(2)tteyxdtdyeyxdtdxdtxd434222y11=dsolve('D1y=1+y^2','x')y12=dsolve('D1y=1+y^2','y(0)=1','x')%求满足条件y(0)=1的特解[x3,y3]=dsolve('D2x+2*Dx=x+2*y-exp(-t)','Dy=4*x+3*y+4*exp(-t)','t')第六讲/2720四、解析解的MATLAB实现3、dsolve对环境流体模型求解举例(1)1)0(,0CtCCKVQCVQdtdCC=dsolve('Dc=Q/V*c0-(Q/V+K)*c','c(0)=c1','t')pretty(C)(Q+KV)texp(----------------)(-Qc0+c1Q+c1KV)Qc0VC=-------+------------