首页 > 分享 > 一种自适应时间步长的场地地表污染扩散方法

一种自适应时间步长的场地地表污染扩散方法

花匠小妙招
2024-12-29 09:44

1.本发明属于污染扩散模拟领域，具体为一种自适应时间步长的场地地表污染扩散方法。

背景技术：

2.污染事故发生时，污染物首先会到达地表，从地势高处向地势地处进行流动，同时根据土壤性质，向下渗入土壤层，污染土壤和地下水，给生态造成持久性危害。一般来说，污染物在地表的扩散范围越大，其在土壤层和地下水中的蔓延范围也越广。因此，准确地模拟污染物在地表的扩散过程，并确定污染物在污染事件发生后的扩散范围是管理决策者首先要解决的问题。
3.现有对地表污染扩散的研究主要集中在地表流域和地表径流方面。对于污染物在流域中的扩散模拟，多采用swat模型、hspf模型、sparrow模型，但是这些模型结构复杂，主要面向较大尺度范围的模拟，且受空间不均匀性影响，带有一定的不确定性，并不适用于小范围详细场地模拟。而对于污染物在地表径流中的扩散模拟，多是基于降雨条件，采用城市暴雨模型和地表水污染模型来估算水污染负荷，然而这些模型多用于水污染评价，无法很好地模拟污染物在降雨下的时空动态扩散过程；也有采用流体动力学模型(sph)模拟地表污染扩散，这种方法能较好地模拟污染物在复杂地表的扩散过程，但是，描述流体扩散的sph 方程的相关参数难以得到，空间边界难以确定，因此这种方法比较适合理想实验环境，普适性较弱。
4.元胞自动机模型是一种时间空间都离散的动力学模型，它能很好地简化原有模型，能通过局部推演整体来描述复杂的地理情况，已经有大量的研究证实了元胞自动机在计算机模拟复杂世界方面的优势，如在城市扩展模拟、土地利用变迁、雨洪模拟、森林火灾模拟等领域，元胞自动机模型都显示了强大的功能。
5.元胞自动机模型的时间步长(即元胞状态根据转化规则更新一次代表的时间)对准确模拟污染扩散过程至关重要。当时间步长设置的过小时，需要多次更新元胞状态才能模拟污染扩散较为缓慢的过程，造成计算资源的浪费；当时间步长设置的过大时，对扩散过程较快的情况，模拟过程会漏掉许多扩散细节，影响管理者进行决策。由于液态污染物的扩散过程受地表随重力、地表粗糙度以及外界降雨和水流条件等多重因素的影响，其扩散过程是一个非均匀的过程，而现有的元胞自动机大都以固定的时间步长来模拟，无法准确地模拟出污染物的扩散过程。因此，此方法对传统的元胞自动机模型进行改进，根据相邻元胞之间扩散的时间差异、污染物扩散速度随质量衰减的特征，设计一种自适应调整时间步长的元胞动态时间演变算法，考虑地表污染物受不同坡度和粗糙度影响下扩散速度不均匀的情况，分别设计无降雨和降雨两种模拟情形，对污染物的时空扩散过程进行仿真。

技术实现要素：

6.针对现有污染扩散模拟领域存在的问题，本发明以bing地图影像和高程数据、苯污染数据和残差平方和小于0.0005的一组v-g模型参数为基础，基于元胞自动机模型，引入地表高差和粗糙度影响因子，建立地表污染物扩散模型，设计降雨和非降雨条件，并计算两种条件下的污染物流速，设计引入自适应时间步长演变算法，探索更贴近客观世界的污染扩散模拟能力，寻求更准确地模拟污染物在地表的扩散过程和结果分布。
7.具体的技术方案为：
8.一种自适应时间步长的场地地表污染扩散方法，包括以下步骤：
9.(1)数据预处理：对影像和高程数据进行预处理，包括高精度影像的切片和发布；测定污染物和土壤属性数据，根据实验确定van genuchten(v-g)模型、horton方程的参数；测定场地边界特征点坐标，基于最小包围盒思想，确定模型边界条件和模型空间，然后确定元胞最小单元大小，最后划分元胞空间；
10.(2)构建一种基于元胞自动机的场地地表污染扩散模型：基于元胞自动机模型，考虑地表高差和粗糙度的影响，构建一种地表污染扩散模型，对地表污染物的迁移过程进行模拟；
11.所述的基于元胞自动机的场地地表污染扩散模型为：
[0012][0013]
式中：m为污染物质量；i、j分别为元胞行列号；分别为垂向和斜向扩散系数，影响污染物在垂向和斜向的扩散速率；t为当前时刻；δt为演变时间步长；k为影响系数；usδt为污染物的损失项，k通过高差影响系数kd和粗糙度影响系数kr计算得到，计算公式为： k＝kd×
kr+1；usδt为当前时间步长内污染物损失的质量，其中u为污染物损失系数，s为污染扩散方向上通过的横截面积。
[0014]
(3)降雨非降雨下地表污染物流速计算：考虑降雨和非降雨两种外界条件，分别计算降雨条件下污染物随地表扩散时的流速和非降雨条件下受自身重力扩散的流速，增强模型的适用性；
[0015]
(4)将自适应时间步长算法应用到扩散模型中：根据污染物在不同地表位置的扩散速度设计一种随污染物质量和坡度衰减的元胞自适应时间步长演变算法，将其应用到污染扩散模型之中，以动态确定一次演变时间的间隔，同时动态表示一次演变中的污染先后过程，从而更真实地模拟污染物在地表随时间的扩散过程；
[0016]
(5)模拟降雨非降雨地表污染物的扩散过程：模拟非降雨条件下两种下垫面时的污染扩散结果；分别模拟降雨条件下不同降雨时长、降雨强度下的扩散结果；对比降雨条件下污染同一位置的浓度变化过程；模拟降雨条件下运用自适应时间步长算法的污染扩散过程；
[0017]
(6)结果分析与验证：分析非降雨条件下不同下垫面的污染扩散结果差异、降雨对污染物扩散的影响、自适应时间步长算法模拟结果；以降雨条件为例，根据坡度图，验证污染扩散结果。
[0018]
实验结果表明，不同下垫面对污染物扩散速度有很大的影响，污染物在水泥地表扩散速度远大于一般性土壤地表；降雨强度和时长能够加快污染物的扩散，且扩散速度随着降雨曲线变化而改变，此外，污染物会向地势更低的地方进行扩散，服从坡度分布特征，且高污染区域范围慢慢变小，一定时间之后，浓度变化渐渐趋于平稳；自适应步长算法能够较好地体现一次演变过程中污染物扩散在不同邻域元胞方向上的差异，能更好地模拟地表污染物在降雨条件下的时空不均匀的扩散过程。本发明旨在探索一种针对传统扩散模型难以动态模拟地表污染物时空不均匀扩散过程问题的解决方式，以期更好更准确地模拟地表污染物在降雨和非降雨条件下的扩散过程。
附图说明
[0019]
图1a为实施例实验区地理位置及边界条件；
[0020]
图1b为实施例实验区的元胞空间；
[0021]
图1c为实施例实验区的空间划分网格；
[0022]
图2为本发明元胞自适应时间演变算法流程图；
[0023]
图3为本发明的方法流程图；
[0024]
图4a-1为非降雨条件下下垫面为一般性土壤时污染物的扩散结果图；
[0025]
图4a-2为非降雨条件下下垫面为一般性土壤时污染物的扩散结果坐标图；
[0026]
图4b-1为非降雨条件下下垫面为水泥路面时污染物的扩散结果图；
[0027]
图4b-2为非降雨条件下下垫面为水泥路面时污染物的扩散结果坐标图；
[0028]
图5a-1为实施例降雨条件下受降雨时长影响时污染物的扩散结果图之一；
[0029]
图5a-2为实施例降雨条件下受降雨时长影响时污染物的扩散结果图之二；
[0030]
图5a-3为实施例降雨条件下受降雨时长影响时污染物的扩散结果图之三；
[0031]
图5a-4为实施例降雨条件下受降雨时长影响时污染物的扩散结果图之四；
[0032]
图5b-1为实施例降雨条件下受降雨强度影响时污染物的扩散结果图之一；
[0033]
图5b-2为实施例降雨条件下受降雨强度影响时污染物的扩散结果图之二；
[0034]
图5b-3为实施例降雨条件下受降雨强度影响时污染物的扩散结果图之三；
[0035]
图5b-4为实施例降雨条件下受降雨强度影响时污染物的扩散结果图之四；
[0036]
图5b-5为实施例降雨条件下受降雨强度影响时污染物的扩散结果图之五；
[0037]
图6为实施例不同地点的污染物浓度随时间变化的示意图；
[0038]
图7为实施例一次演变中污染物不同方向上的扩散差异示意图；
[0039]
图8-1为实施例污染扩散结果与实验区坡度分布特征的对比图之一；
[0040]
图8-2为实施例污染扩散结果与实验区坡度分布特征的对比图之二；
25～35mm/h。根据不同降雨强度，确定表层土壤渗透系数，计算对应的土壤含水率，再计算当前雨强的稳定水流入渗量。
[0057]
降雨参数选取：：
[0058]
恒定35mm/h雨强降雨1小时；最大降雨量为60mm/h，最低降雨量为0，降雨时长为30min，抛物线方程为r＝t(30-t)。
[0059]
实验区坡度分布图：
[0060]
分析实验区的坡度坡向，如图8-1到图8-3所示。
[0061]
3基于元胞自动机的场地地表污染扩散模型：
[0062]
污染物在地表的流动和在静水中的流动过程不同，还会受到其他因素的影响。其中最主要的影响因素为重力和坡度，污染物会因重力作用从地势高处向地势低处进行扩散，此外，不同的地表粗糙度也会影响污染物的扩散，地面越光滑，污染物的扩散速度就越快，反之则越慢。本发明考虑相邻元胞间的高差与地表粗糙度对污染扩散的影响以及污染物扩散的主要损耗，建立了场地地表污染扩散元胞自动机模型：
[0063][0064]
式中：m为污染物质量；i、j分别为元胞行列号；分别为垂向和斜向扩散系数，影响污染物在垂向和斜向的扩散速率；t为当前时刻；δt为演变时间步长；k为影响系数；usδt为污染物的损失项，k通过高差影响系数kd和粗糙度影响系数kr计算得到，计算公式为： k＝kd×
kr+1；usδt为当前时间步长内污染物损失的质量，其中u为污染物损失系数，s为污染扩散方向上通过的横截面积。
[0065]
第i个元胞的高差影响系数的计算公式为：
[0066][0067]
式中：i、j都分别表示高度大于中央元胞高度的邻域元胞索引；hi、hj分别为第i、j个元胞的高程；h0为中央元胞的高程；∑(h
j-h0)是中央元胞与比其高程小的邻域元胞之间的高差之和。
[0068]
第i个元胞的粗糙度影响系数的计算公式为：
[0069]
[0070]
式中：i、j分别表示邻域元胞的索引；ri、rj表示粗糙度；∑(1/rj)表示与中央元胞存在质量输出的元胞的粗糙度倒数之和，其中粗糙度的计算公式为：
[0071][0072]
式中：slope为当前元胞与中央元胞之间的坡度角。
[0073]
4降雨非降雨下地表污染物流速计算：
[0074]
采用曼宁公式计算地表污染物的平均流速：
[0075][0076]
式中：k是转换常数，值为1；n是糙率，又称粗糙系数，反映影响流体速度的综合阻力因素，无量纲，糙率越大，流体速度越小；rh是水力半径，是流体截面积与湿周长的比值，根据经验可取流体平均深度来代替；s是坡度，为地表两点之间的高差与距离的比值。
[0077]
引入雷诺数的概念，常采用雷诺数来区分水流流态，雷诺数越小，意味着流体的粘性影响显著，越大则自身的惯性成为主导：
[0078]
re＝ρvd/μ
[0079]
式中：re表示雷诺数；ρ为液体的密度；v为液体流速；d为一固定长度，和液体流动环境有关；μ为液体粘性系数。
[0080]
液态污染物与水流在流动过程中有一定的相似，可以通过雷诺数的大小对两者进行区分。基于此，本发明对经典曼宁公式进行进一步处理，在原方程的基础上添加雷诺数转换因子，得到一种更普适的液体均匀流流速计算公式。
[0081][0082]
式中：f
re
为转换因子，作用是将一般水流转换为特定的液态污染物，其计算公式为：
[0083][0084]
式中：是转换目标流体的雷诺数；是水的雷诺数。转换因子的计算规则是假定外界条件如温度、地表粗糙度等因素都相同的情况下，只受自身密度和粘性系数影响下的比值。
[0085]
(2)计算非降雨条件下污染物的平均流速：
[0086]
根据下一时刻污染物的扩散质量以及网格的大小来计算污染物在每个网格内的平均高度，即水力半径：
[0087][0088]
式中：m是当前网格内的污染物质量；ρ是污染物的密度；s是网格的底面面积。
[0089]
在非饱和土壤中，流体在土壤的渗透系数可以表示其向下渗流的强弱，表达式为：
[0090]
k＝pρg/η
[0091]
式中：k为渗透系数；p为土壤渗透率，和孔隙度有关；ρ为流体密度；η为粘性系数。水的渗透系数一般比较好获取，实验测定了场地非饱和土壤水的渗透系数，通过公式将其
转换为污染物在土壤中的渗透系数。
[0092]
由于影响渗透系数的因素主要为流体密度和粘性系数，因此可以采用雷诺数比值f
re
进行转换，转换公式为：
[0093]kt
＝f
rek[0094]
式中：k
t
为目标流体的渗透系数；k为水的渗透系数。
[0095]
考虑污染物向下渗透的质量流失，根据污染物的渗透系数及一次演变的向下渗入时间 δt，计算下渗的污染物质量：
[0096][0097]
得到了水力半径之后，根据曼宁公式计算得到流速。
[0098]
(2)计算降雨条件下污染物平均流速：
[0099]
本发明考虑降雨作用和下渗作用，计算了不同时刻的净雨量，采用在地表形成的径流深度，来代替平均水力半径，从而计算污染物在降雨过程中的流速，下式是产流的计算公式。
[0100]
r＝p-i
[0101]
式中：r是产流量；p是降雨量；i是下渗量；单位都为mm。
[0102]
采用改进horton方程计算某个时刻的下渗量：
[0103][0104]
式中：f(t)为t时刻的下渗率；fc为稳定下渗率；f0为初始下渗率；u为衰减系数。θ0和 θc分别是土壤的初始含水率和饱和含水率；θ为t时刻的土壤水分；k为改进系数。
[0105]
采用van genuchten模型将非饱和土壤水渗透系数表示成饱和渗透系数和有效饱和度的函数：
[0106][0107]
式中：k为非饱和土壤水渗透系数；ks为饱和渗透系数；se为有效饱和度；m为拟合参数，l为模型的经验参数。将有效饱和度se可表示成土壤体积含水率的函数：
[0108][0109]
式中：θ为土壤体积含水率(cm3/cm3)，θs、θr分别为饱和体积含水率和残余体积含 (cm3/cm3)。
[0110]
建立降雨强度与土壤体积含水率的关系式：
[0111]
θ＝f(k)
[0112]
考虑雨强和阻力对流速的影响，得到流速和雨强、阻力的关系：
[0113]
u～r
1-m
[0114]
根据液态指数m的定义，本文将其简化成地表粗糙度系数n，再结合曼宁公式，得出了在雨强r影响下的地表污染物扩散流速计算公式：
[0115][0116]
4将自适应时间步长算法应用到扩散模型中：
[0117]
本发明根据污染物在不同地表位置的扩散速度设计了一种随污染物质量和坡度衰减的元胞自适应时间步长演变算法，来动态确定一次演变时间的间隔，同时随时间推移动态表示一次演变中的污染先后过程，从而更真实地模拟污染物在地表随时间的扩散过程，该算法的流程如图2所示。
[0118]
主要流程是获取当前所有已经被污染的元胞对象，演变一次，剔除已经被污染的元胞，获得下个时刻的所有污染元胞对象，确定该方向上的演变时间t，设定参考时间t，判断污染元胞扩散时间是否已经大于参考时间，若大于，将元胞放入污染池，并更新它的颜色。若参考时间为真实时间，则扩散时间会以真实世界的时间为基准，若缩短为真实世界的一半，则扩散会以2倍的速率进行。
[0119]
不同邻域方向上的时间步长计算公式为：
[0120][0121]
式中：x、y分别表示邻域元胞的索引；δt
x,y
为污染物在(x,y)方向上扩散的时间步长；l 为元胞边长。某个时刻的推演时间步长为：
[0122]
δt
i,j
＝max(δt
i-1,j-1
,δt
i-1,j
,...,δt
i+1,j+1
)(δt
i,j
＝0)
[0123]
元胞状态转移方程为：
[0124][0125]
式中：s
i,j
为元胞状态；v
i,j
为污染物在i、j方向上的扩散速度；δt
i,j
为时间步长；l为元胞边长。
[0126]
5模拟结果分析与验证：
[0127]
(1)非降雨条件下下垫面为一般性土壤时：
[0128]
污染物因自身的重力作用向地表低的位置流动，扩散距离越远，扩散的速度就越慢。由于污染物扩散速度缓慢，图4a-1和图4a-2是污染物在地表受自身重力作用下流动1天的污染分布图。
[0129]
从图中可以看出污染物的扩散速度十分缓慢，1天的扩散直线距离约为3.56m，并且由于离污染源越远的网格内污染质量越小，流速会越慢，几乎可以认为已经不发生新的扩散；污染物扩散时间随距离增加，扩散0.70m，所需时间约为0.10min；扩散1.41m，所需时间约为 0.80min，扩散2.12m，所需时间约为6.30min，而扩散2.83m，时间已经超过24h。当污染物每扩散单位元胞长度的距离时，其所花的时间都会成倍的增加，在扩散距离超过2.83m左右后，扩散时间几乎成90
°
上升，因此可认为污染物在扩散2.83m后几乎不发生扩散作用，这也与实际情况中污染物扩散距离随质量衰减的情况相符。
[0130]
(2)非降雨条件下下垫面为水泥路面时：
[0131]
由于只改变了下垫面的糙率，地表的坡度未发生改变，因此图4b-1和图4b-2相对图4a-1 和图4a-2，污染物的扩散范围只在左下方的位置上有所增大，但从扩散时间上来看，糙率对扩散速度的影响将更加明显。污染物扩散0.70m，所需时间约为0.04min；扩散1.41m，所需时间约为0.30min，扩散2.12m，所需时间约为2.15min，在扩散距离超过2.12m左右后，扩散时间几乎成直线上升。
[0132]
下垫面对污染扩散有很大的影响，从时间尺度上来看，当曼宁系数从0.035降低到
0.012 时，并略去向下渗入量时，扩散速度将会变为原来的2.7倍左右。
[0133]
(3)降雨时长对污染物扩散的影响：
[0134]
污染物扩散模拟结果如图5a-1到图5a-4所示：当降雨10min时，污染物的扩散距离约为 12.7m；当降雨20min时，扩散距离约为19.8m；当降雨40min时，扩散距离约为31.1m；当降雨60min时，扩散距离约为43.1m。从中可以看出，降雨强度是恒定的情况下，随着降雨时长的增加，污染物的扩散距离并不是时间均匀的，定步长模型不能很好地描述这种情况。
[0135]
(4)降雨强度随时间变化下污染物的扩散模拟：
[0136]
污染物扩散模拟结果如图5b-1到图5b-5所示：当降雨5min时，扩散距离约为7.1m；当降雨10min时，距离约为20.5m；当降雨15min时，距离约为38.2m；当降雨20min时，距离约为51.6m；当降雨25min时，距离约为61.5；当降雨30min时，距离约为70.0m。由于降雨是一个二次抛物曲线，随着时间的变化，污染物扩散距离有一个明显的加速过程，一直到15min左右时，扩散速率最快，之后慢慢变缓。自适应时间步长算法将不同时间下的降雨强度考虑在内，推演出适应的时间步长，用以描述这种外界条件复杂变化的情况。
[0137]
(5)不同地点污染物浓度随时间的变化：
[0138]
污染物扩散模拟结果如图6所示：元胞自动机推演30次时，污染源的附近颜色最深，以污染源为中心向四周的颜色越来越浅；推演50次时，发现污染源附近浓度大于24g/l的元胞越来越少，而处于其他浓度区间的数量越来越多；推演80次时，发现浓度大于24g/l的元胞已经消失；推演100次时，发现所有污染区域的浓度都处于18g/l以下；推演300次后，发现最后所有污染区域的浓度都小于6g/l，此时的扩散趋于平稳。
[0139]
(6)一次演变中污染物不同方向上的扩散差异：
[0140]
为了加快模拟速度，以上述实验中的降雨条件为例，假定已经形成了稳定的坡面水流，水流深度约为0.01m，那么此时在坡度为0.01m的方向上，污染扩散速率约为0.4m/s，本文根据算法设计，设定参考时间为真实时间，进行实时模拟，观察一次演变过程中不同邻域方向上的扩散时间差异，选取了邻域方向之间坡度差异较为明显的第1次、第6次、第10次演变过程，对污染物扩散到不同邻域方向的扩散时刻进行了展示，如图7所示；图中绿色的元胞位置代表该次演变过程中即将发生扩散的位置，但由于一次演变还未完成，需等待其他邻域方向上的到达指定的扩散时间，才会更新这次演变后的所有污染物浓度，并重新计算每个元胞对应的渐变颜色，最后在计算机中绘制。图7展示了污染物在一次演变过程中不同邻域方向上的污染扩散时间差异，这种可变时间步长算法将更加真实地模拟污染物随时间的变化过程，更加符合真实扩散情况。
[0141]
(7)污染扩散模拟结果与验证：
[0142]
将实验区地表高度划分成5个区间，使用渐变色生成了实验区的坡度图，如图8-1到图 8-3所示。从图中可看出实验区的坡向整体自西北向东北和正南，而在东南方向也有一个较缓的坡度；将降雨条件下两个不同时刻污染物的扩散范围与坡度图叠加，从图中可以看出，污染物主要的扩散方向与坡度的方向基本保持一致，主要向地势更低的东北、西南方扩散，速度较快，而在地势差异较小的东南方向的扩散则较为缓慢。