一种基于悬浮泥沙沉降概率的磷输移与分布计算方法与流程

1.本发明涉及水环境模拟与调控领域，具体是一种基于悬浮泥沙沉降概率的磷输移与分布计算方法。


背景技术：

2.泥沙运动状态、颗粒级配和污染物吸附特征的变化直接关系到污染物的含量以及赋存形态。水中的磷可以分为可溶态和颗粒态，磷伴随着泥沙在水体中的运动而迁移。一方面，颗粒态磷随泥沙沉降到底床中，减小了赋存于水体中的污染物量，在一定条件下沉积在底床的污染物随着冲刷泥沙重新进入水体；另一方面，当水化学、水动力等外部条件发生改变，泥沙上的吸附相磷会与水体溶解相磷发生动态交换，水体中悬浮泥沙的运动、输移直接影响着磷在固-液两相间的分配关系和时空分布。泥沙是一种宽级配、多粒径组成的混合体，其粒径组成是影响沉降速率和污染物吸附解吸速率的主要因素之一。
3.河流泥沙磷输移过程的数学模型是解决水库工程管理问题的有效工具。最早模拟泥沙水质影响的是thomas模型，随后有velez-gannon模型、dobbins-comp模型，但它们将泥沙运动的各种过程用一个常量参数表述，不能反映泥沙输送对污染物迁移的影响(何用,李义天,郜会彩,王家生.泥沙污染水质模型研究[j].四川大学学报(工程科学版),2004,6:12-17.)。许多研究者提出并建立了含沙水体的水质模型，先后引入不同的参数来考虑泥沙运动对磷输移的影响。褚君达等(褚君达,徐惠慈.河流底泥冲刷沉降对水质影响的研究.水利学报,1994,11,42-47,69.)从污染对流扩散方程出发，引入颗粒相污染物，考虑了底泥冲刷、沉降及释放等因素的影响，推导得到一维水质模型；禹雪中建立了溶解相、吸附相污染物一维输移转化分相模型(禹雪中,杨志峰,钟德钰,彭期冬.河流泥沙与污染物相互作用数学模型[j].水利学报,2006,1:10-15.)。黄磊等引入河流动力学中床沙悬沙交换的理论，建立了河流水动力-泥沙-磷输移过程的分相模型(黄磊,方红卫,王靖宇,等.河流泥沙磷输移过程的数学模型研究[j].水利学报,2014,4:394-402.)，该模型从机理上考虑磷在底泥中的不同分配比例，重点关注冲淤变化下底泥层中有氧层与无氧层之间的磷交换，但是该模型磷输移转化方程中的沉降系数仅为一个固定经验常数(即群体稳恒沉速)，这种概化的方式忽略了水体层颗粒磷随不同粒径泥沙沉降受阻和多相相互作用，缺少沉降过程的基础理论，因此，这类模型更适用于外部扰动较小的大型湖泊型水库沉积物磷环境行为研究，但对于天然河流或者河道型水库而言，由于库区水动力条件年变幅较大、沿程空间变化明显，在人为水沙调控作用下不同粒径泥沙及其所携颗粒磷的输送沉降特性存在显著时空差异性，对沉降过程的简化处理可能给模拟结果带来重要影响。
[0004]
本发明人在实现本发明的过程中经过研究发现：现有水沙污染物耦合模型或技术存在较大不足：
[0005]
(1)现有泥沙-水质数学模型过度概化泥沙颗粒粒径这一关键影响因子，在模拟中均采用中值粒径作为泥沙颗粒粒径，这种简化处理，对于输送距离近、颗粒大小较为均匀的泥沙较为适用，但泥沙根据粒径大小通常可分为黏土、粉沙、极细沙、细沙和粗砂等，中值粒
径往往不能反映泥沙颗粒的真实颗粒组成，且对于长距离泥沙输送来说，随着水动力条件改变，不同粒径的泥沙其运动范围、沉降距离及运动状态均表现出巨大的差异性，比如粗砂通常沉降与库尾，细沙则淤积在坝前且细沙能够吸附富集更多污染物。
[0006]
(2)其次，泥沙颗粒运动过程中，磷等污染物在颗粒态与溶解相间进行再分配，分配系数的确定通常选取代表性粒径颗粒，通过恒温震荡吸附解吸实验，借助langmuir等吸附模型计算获得；泥沙粒径呈正态分布，震荡后的粒径与实验前已发生巨大变化，因此该分配系数也不能科学反映泥沙粒径或组成对颗粒态磷与溶解相磷相互转化的影响。


技术实现要素：

[0007]
针对现有技术存在的上述缺陷，本发明提供一种基于悬浮泥沙沉降概率的磷输移与分布计算方法，其考虑泥沙的全粒径组成，分析河流泥沙的粒径级配，基于泥沙沉降概率计算全粒径泥沙沉降速率及输送距离，基于室内实验获取不同粒径泥沙颗粒的磷分配系数，构建一种基于泥沙沉降概率的磷输移模型，可全面解析全粒径泥沙输送的环境效应。
[0008]
一种基于悬浮泥沙沉降概率的磷输移与分布计算方法，包括如下步骤：
[0009]
步骤1：选取典型断面，现场测量地形数据，同步测定悬浮泥沙浓度s、水温t、流速u，所述地形数据包括河宽b、水深h，采用原位激光粒度仪实时监测泥沙颗粒级配等变化；
[0010]
步骤2：以分组泥沙颗粒粒径尺寸大小初步判断流态，根据流态计算泥沙颗粒沉速，后通过计算雷诺数反证所判断流态是否正确，若不属于判别流层，则假定另一流态计算沉速，直到算出正确结果，根据悬浮泥沙浓度对泥沙群体沉速进行修正，并合成得到不同粒径分组泥沙沉速；
[0011]
步骤3：根据全粒径泥沙沉速，由泥沙沉降概率密度函数计算颗粒磷沉降概率；
[0012]
步骤4：根据颗粒磷沉降概率及合成后的全粒径颗粒沉降速率，计算颗粒磷沉降前的输移距离；
[0013]
步骤5：通过室内实验分别得到不同粒级分组下泥沙吸附速率k
p
及不同悬浮泥沙浓度下泥沙吸附速率ks，通过多项式拟合得到综合吸附速率k
α
；
[0014]
步骤6：采用磷分相输移模型分析溶解相磷及吸附相磷沿程变化，提取特定断面数据，得到磷分布特征，根据两种形态磷的比例得到溶解相磷及吸附相磷分配特征。
[0015]
进一步的，步骤2中采用公式(7)对泥沙群体沉速进行修正，采用公式(8)合成得到不同粒径分组泥沙沉速：
[0016][0017][0018]
式中，为泥沙群体沉速；pi为第i组粒径所占百分比；ωi为第i组粒径泥沙沉速；ω
s0
为泥沙颗粒在实验条件下的沉速，sf为实验条件下的悬沙浓度，α是一个与雷诺数有关的系数。
[0019]
进一步的，步骤3具体包括：
[0020]
模型假定泥沙颗粒沉降到底部时会以一定的概率沉降下来，其沉降概率在0～1之
间变化，泥沙颗粒运动速度服从gauss正态分布，分布密度函数为：
[0021][0022]
式中x为颗粒速度；u为均值；σ为标准差，颗粒磷与泥沙相互依存，将gauss分布密度函数积分，得颗粒态磷沉降概率函数：
[0023][0024]
式中，ps为颗粒态磷沉降概率；x为泥沙颗粒速度；ω为颗粒磷沉速，即全粒径泥沙沉速；σ
ν
为垂向紊动强度标准差。
[0025]
进一步的，步骤(4)中通过式(11)计算颗粒磷沉降前的输移距离：
[0026]
泥沙颗粒在水流垂线平均流速及沉速ω之下，由距河床为h处沉降至河底时纵向输移距离为：
[0027][0028]
式中lc为泥沙落淤纵向输移距离；ω为沉速；为水流垂线平均流速；β为泥沙颗粒与水流紊动跟随性的参数。
[0029]
进一步的，步骤5中吸附速率k
α
采用室内吸附动力学实验得到，动力学方程采用颗粒内扩散模型、一级动力学、二级动力学等分别拟合，由可决系数选取最佳拟合参数，
[0030]
其中：
[0031]kp
＝f(d,pi)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)
[0032]ks
＝f(s)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(13)
[0033]kα
＝f(k
p
,ks)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(14)。
[0034]
进一步的，步骤6中吸附相磷的表达为：
[0035][0036][0037]
式中cs为泥沙吸附相磷含量，n为泥沙粒径分组，q
si
为第i组泥沙吸附相磷含量，pi为第i组泥沙质量百分数；
[0038]
磷输移的控制方程采用分相模型，包括水体中溶解相磷cw、悬沙上吸附相磷cs，具体表达式如下：
[0039][0040][0041]
式中，a为控制断面面积；cb、cs、cw分别为床沙、悬浮泥沙及水体中污染物含量；s为
断面平均含沙量；s
*
为断面平均挟沙力；d
x
为纵向离散系数；α为恢复饱和系数；ω为颗粒沉速，即包含了泥沙粒径分组信息的颗粒沉速；b为河宽；k
α
为吸附速率，即单位时间单位质量泥沙吸附量的变化率，k
α
》0表示泥沙吸附，反之表示解吸；u为断面平均流速；dcw/dt为生化反应项；σ为单位河长污染物的释放速率；r为排放源强；
[0042]
挟沙力采用张瑞瑾公式，引入粒径级配校正系数：
[0043][0044]s*
为挟沙力，单位为kg/m3；u为断面平均流速，单位为m/s；g为重力加速度，单位为m/s2；h为平均水深，单位为m；ω为平均沉速，单位为m/s；m为挟沙力指数，一般取0.92；k为挟沙力系数；pi为第i组泥沙所占质量百分数；n为泥沙粒径分组数；
[0045]
不同泥沙组分中吸附相磷的含量关系如下：
[0046][0047]
式中cs为泥沙吸附相磷含量，n为泥沙粒径分组，q
si
为第i组泥沙吸附相磷含量，pi为第i组泥沙质量百分数。
[0048]
由于采用了上述方案，本发明可模拟磷随泥沙输移的过程，具有如下有益效果：
[0049]
(1)结合悬沙粒径水下原位监测技术和颗粒沉速公式，计算出随一个包含粒径级配的颗粒沉速，用以获取真实的全组分泥沙沉降速率变化过程；采用实时观测的粒径数据，解决了传统水沙营养盐模型中泥沙沉降速率给定固定常数的局限，实现了全粒径泥沙和生源要素在干流水体中的沉降输移过程的精准模拟，更符合实际情况；
[0050]
(2)模拟了不同粒径颗粒磷的沉降轨迹，计算全组分颗粒磷的输移距离，分析颗粒磷的沉降概率及拦截效率。
附图说明
[0051]
图1是本发明构建的基于悬浮泥沙沉降概率的磷输移模型的框架示意图；
[0052]
图2是不同紊动速度下颗粒磷沉降概率分布图；
[0053]
图3是颗粒磷运动轨迹图；
[0054]
图4是汛期各断面表层颗粒磷输送轨迹；
[0055]
图5是本发明实施例颗粒态磷模拟值与实测值对比图；
[0056]
图6是本发明基于悬浮泥沙沉降概率的磷输移与分布计算方法其中一个实施例的流程图。
具体实施方式
[0057]
下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整的描述。
[0058]
请参阅图1，本发明建的基于悬浮泥沙沉降概率的磷输移模型包括水动力模块、悬浮泥沙泥沙输移模块和磷输送分配模块3个部分，其中水动力模块和磷输送分配模块均考
虑了泥沙级配及含沙量的影响，颗粒磷沉降速率通过gauss方程表示。
[0059]
其中水动力模块包括水流连续方程和水流运动方程，具体如下：
[0060]
水流连续方程：
[0061][0062]
水流运动方程：
[0063][0064]
悬浮泥沙输移模块中的泥沙输移采用非饱和输沙方法计算，考虑悬移质泥沙与底泥之间的交换及泥沙粒径组成。泥沙输移的控制方程为连续方程，具体如下：
[0065][0066][0067]
式中，ω为泥沙沉速，角标i为不同粒径组泥沙分组号；q为流量；a为过水面积；t为时间；x为沿流程坐标；z为水位；k为断面流量模数；s为含沙量；s
*
为水流挟沙力；ρ为淤积物干容重；b为断面宽度；g为重力加速度；u为平均流速；α为恢复饱和系数；ad为悬移质河床冲淤面积；d为粒径；u为流速；ab为推移质河床冲淤面积。
[0068]
磷输送分配模块采用分相模型描述磷的输移过程，包括水体中溶解相磷(cw)、悬沙上吸附相磷(cs)，具体表达式如下：
[0069][0070][0071]
式中，a为控制断面面积；cb、cs、cw分别为床沙、悬浮泥沙及水体中污染物含量；s为断面平均含沙量；s
*
为断面平均挟沙力；d
x
为纵向离散系数；α为恢复饱和系数；ω为颗粒沉速；b为河宽；k
α
为吸附速率，即单位时间单位质量泥沙吸附量的变化率，k
α
》0表示泥沙吸附，反之表示解吸；u为断面平均流速；dcw/dt为生化反应项；σ为单位河长污染物的释放速率；r为排放源强。
[0072]
请继续参阅图6，本发明实施例提供一种基于悬浮泥沙沉降概率的磷输移与分布计算方法，包括如下步骤：
[0073]
步骤1：选取典型断面，现场测量河宽b、水深h等地形数据，同步测定悬浮泥沙浓度s、水温t、流速u等，采用原位激光粒度仪实时监测泥沙颗粒级配等变化；
[0074]
步骤2：以分组泥沙颗粒粒径尺寸大小初步判断流态，根据流态计算泥沙颗粒沉速，后通过计算雷诺数反证所判断流态是否正确，若不属于判别流层，则假定另一流态计算沉速，直到算出正确结果；根据悬浮泥沙浓度对泥沙群体沉速进行修正公式(式7)，并合成得到不同粒径分组泥沙沉速(式8)：
[0075][0076][0077]
式中，为泥沙群体沉速；pi为第i组粒径所占百分比；ωi为第i组粒径泥沙沉速；ω
s0
为泥沙颗粒在实验条件下的沉速，sf为实验条件下的悬沙浓度，α是一个与雷诺数有关的系数；
[0078]
步骤3：根据全粒径泥沙沉速(式8)，由泥沙沉降概率密度函数(式10)计算颗粒磷沉降概率。
[0079]
模型假定泥沙颗粒沉降到底部时会以一定的概率沉降下来，其沉降概率在0～1之间变化。大量实验观测证实泥沙颗粒运动速度服从gauss正态分布，分布密度函数为：
[0080][0081]
式中x为颗粒速度；u为均值；σ为标准差。颗粒磷与泥沙相互依存，将gauss分布密度函数积分，可得颗粒态磷沉降概率函数：
[0082][0083]
式中，ps为颗粒态磷沉降概率；x为泥沙颗粒速度；ω为颗粒磷沉速(即全粒径泥沙沉速)；σ
ν
为垂向紊动强度标准差。
[0084]
当vs取不同值的时，相应的泥沙(或颗粒磷)沉降概率如图2所示。
[0085]
步骤4：根据颗粒磷沉降概率(式10)及合成后的全粒径颗粒沉降速率(式8)，通过式(11)计算颗粒磷沉降前的输移距离。
[0086]
泥沙颗粒在水流垂线平均流速及沉速ω之下，由距河床为h处沉降至河底时纵向输移距离为：
[0087][0088]
式中lc为泥沙落淤纵向输移距离；ω为沉速；为水流垂线平均流速；β为泥沙颗粒与水流紊动跟随性的参数。
[0089]
步骤5：通过室内实验分别得到不同粒级分组下泥沙吸附速率k
p
(式12)及不同悬浮泥沙浓度下泥沙吸附速率ks(式13)，通过多项式拟合得到综合吸附速率k
α
(式14)。吸附速率k
α
可采用室内吸附动力学实验得到。动力学方程采用颗粒内扩散模型、一级动力学、二级动力学等分别拟合，由可决系数选取最佳拟合参数。
[0090]
其中：
[0091]kp
＝f(d,pi)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)
[0092]ks
＝f(s)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(13)
[0093]kα
＝f(k
p
,ks)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(14)
[0094]
步骤6：采用磷分相输移模型分析溶解相磷及吸附相磷沿程变化，提取特定断面数
据，得到磷分布特征。根据两种形态磷的比例得到溶解相磷及吸附相磷分配特征。
[0095]
吸附相磷的表达为：
[0096][0097][0098]
式中cs为泥沙吸附相磷含量，n为泥沙粒径分组，q
si
为第i组泥沙吸附相磷含量，pi为第i组泥沙质量百分数。
[0099]
磷输移的控制方程采用分相模型，包括水体中溶解相磷(cw)、悬沙上吸附相磷(cs)、具体表达式如下：
[0100][0101][0102]
式中，a为控制断面面积；cb、cs、cw分别为床沙、悬浮泥沙及水体中污染物含量；s为断面平均含沙量；s
*
为断面平均挟沙力；d
x
为纵向离散系数；α为恢复饱和系数；ω为颗粒沉速；b为河宽；k
α
为吸附速率，即单位时间单位质量泥沙吸附量的变化率，k
α
》0表示泥沙吸附，反之表示解吸；u为断面平均流速；dcw/dt为生化反应项；σ为单位河长污染物的释放速率；r为排放源强。
[0103]
挟沙力采用张瑞瑾公式，引入粒径级配校正系数：
[0104][0105]s*
为挟沙力，kg/m3；u为断面平均流速，m/s；g为重力加速度，m/s2；h为平均水深，m；ω为平均沉速，m/s；m为挟沙力指数，一般取0.92；k为挟沙力系数；pi为第i组泥沙所占质量百分数；n为泥沙粒径分组数。
[0106]
不同泥沙组分中吸附相磷的含量关系如下：
[0107][0108]
式中cs为泥沙吸附相磷含量，n为泥沙粒径分组，q
si
为第i组泥沙吸附相磷含量，pi为第i组泥沙质量百分数。
[0109]
引入相应的边界条件和初始条件，导入步骤1测定的地形及水动力参数，采用有限差分法对上述方程进行数值求解，由于悬浮泥沙浓度与紊动参数之间、溶解相磷含量与吸附相磷含量之间均相互影响和制约，因此，在每一个时间步均需进行迭代计算，当两次迭代结果最大绝对误差《10-5
时方可进入下一个时间步的计算。计算过程中，先通过水流运动方
程和悬浮泥沙运动方程求得悬浮泥沙浓度，代入磷输移方程中即可求得溶解相磷和吸附相磷的时空分布。
[0110]
下面以一个具体实例进行说明：
[0111]
步骤1：选取典型河道为研究范围，根据需要划分典型监测断面。采用水下激光原位粒子尺寸分析仪lisst-200x开展悬移质泥沙级配现场在线监测。现场测量地形数据，同步测定悬浮泥沙浓度s、水温t、流速u，所述地形数据包括河宽b、水深h。
[0112]
步骤2：根据泥沙级配分组，考虑现场温度、分组泥沙浓度、泥沙容重等参数，选用适当的计算公式，计算得到各组平均泥沙沉降速率，通过每组泥沙所占比例，合成得到全粒径分组颗粒磷沉降速度。
[0113]
步骤3：根据颗粒磷运动轨迹计算原理，考虑沿水流和水深两个方向的流速，分层分段计算颗粒磷的运动轨迹(如图3所示)。从第1点运动到第2点，沿水流方向的运动距离为l
x
，沿水深方向的运动距离为lz；t为颗粒磷从第1点运动至第2点的时长。
[0114]
当ps＝50％时，泥沙输移距离可由下式计算：
[0115][0116]
式中lc为泥沙落淤纵向距离；ω为颗粒沉速；为水流垂线平均流速。
[0117]
当ps＝15.85％时，取则泥沙颗粒输移距离为：
[0118][0119]
j为水力坡降。
[0120]
步骤4：以现场监测数据为基础，根据上游水沙特征和水质监测结果设置水动力模型、泥沙输移模型及磷输移模型边界条件，确定模型中主要参数的取值。
[0121]
步骤5：引入相应的边界条件和初始条件(表1)，采用有限差分法对上述方程进行数值求解。图4为经过计算后4个断面表层水体颗粒磷输送轨迹，图中这个四个点泥沙颗粒输送距离可达200～350km。图5为由分相模型计算得到的水体颗粒相磷浓度与实测颗粒相磷浓度的对比情况，模型计算结果与真实至相差一般小于20％。
[0122]
表1模型初始条件
[0123][0124]
最后应当说明的是：以上实施例仅用以说明本发明的技术方案而非对其限制，尽管参照上述实施例对本发明进行了详细的说明，所属领域的普通技术人员应当理解：依然可以对本发明的具体实施方式进行修改或者等同替换，而未脱离本发明精神和范围的任何
修改或者等同替换，其均应涵盖在本发明的权利要求保护范围之内。