一种耦合地表水污染物完整地球化学反应和运移过程的模拟方法

1.本发明涉及环境模拟技术领域，尤其是一种耦合地表水污染物完整地球化学反应和运移过程的模拟方法。


背景技术：

2.水是生命之源，地球表面约有71％被水覆盖。其中地表水是指陆地表面上动态水和静态水的总称，包括各种液态的和固态的水体，主要有河流、湖泊、沼泽、冰川、冰盖等。它是人类生活，生产用水的重要来源之一。然而随着社会的不断发展进步，自然水资源却不断的遭到破坏。工业废水、生活污水的大量排放，河流、湖泊、水库等等污染问题变得越来越严重，对人类的生产、生活造成巨大的影响，时刻威胁着人类的健康。
3.地表水水质模型是描述污染物质在地表水体中混合、随水流输移和迁移转化的数学方程，是地表水环境污染治理规划决策分析的重要工具和有效手段。当前，已有许多地表水水质模型，如wasp、efdc等。但是目前已有的水质模型大都具有较强的针对性，只能对特定种类元素进行模拟，如氨氮、磷酸盐，使用具有一定的局限性。


技术实现要素：

4.本发明所要解决的技术问题在于，提供一种耦合地表水污染物完整地球化学反应和运移过程的模拟方法，通过耦合地表水动力、物质输运与热力学平衡机制，建立研究地区的地表水水质模型，实现了地表水动力场与地球化学场的高效模拟和多类型污染物过程模拟。
5.为解决上述技术问题，本发明提供一种耦合地表水污染物完整地球化学反应和运移过程的模拟方法，包括如下步骤：
6.步骤1、收集研究区域的地表水水文数据；
7.步骤2、构建研究区域的地表水动力模型；
8.步骤3、构建研究区域的物质输运模型；
9.步骤4、收集研究区域的地表水水质数据；
10.步骤5、基于热力学平衡理论构建研究区域的多类型污染物模型；
11.步骤6、将步骤2中的地表水动力模型、步骤3中的物质输运模型与步骤5中的多类型污染物模型耦合，得到地表水多类型污染物运移模型。
12.优选的，步骤1中，收集研究区域的地表水文数据包括对该区域的地形、水位、盐度、流速、流量、降雨量的监测。
13.优选的，步骤2中，地表水研究区域包括河流、湖泊及近岸海域。
14.优选的，步骤2中，地表水动力模型的控制方程如下所示：
15.16.式中：v为速度，p为自由表面压力，r为水体容重，ρm为水体密度，z为水深，υ为运动粘滞系数。
17.优选的，步骤3中，物质输运控制方程如下所示：
[0018][0019]
式中：h为水深；c为等深的物质浓度；t为时间；u和v分别为x、y方向上的速度；k
xx
、k
xy
、k
yx
以及k
yy
为二维扩散系数张量，s是源汇项；式(2)采用有限体积法进行差分求解，利用算子分裂算法，将上式分为对流项和扩散项两部分进行分步计算，模型采用fortran90语言编写。
[0020]
优选的，步骤4中，收集研究区域的地表水水质数据包括水的ph值、pe值、温度和准备研究元素的浓度。
[0021]
优选的，步骤5中，基于热力学平衡理论构建的多类型污染物模型是指以地球化学包phreeqcrm为基础建立的水质模型。
[0022]
优选的，步骤6中，将地表水动力模型、物质输运模型与多类型污染物模型耦合是指代码层面的全耦合，将使用fortran90编写的地表水动力、物质输运程序与fortran版本的phreeqcrm进行统一编译，在每个时间步长内同步求解水动力模型，物质输运模型与水质模型。
[0023]
优选的，步骤6中，每一时间步长的计算流程如下：先求解当前时刻的水动力场，然后基于当前计算所得水动力场进行物质输运计算，得到当前时刻各物质的浓度分布场，最后，将当前各物质的浓度结果输入水质模型中，进行地球化学反应过程的计算，得到下一时间步长的各物质初始浓度场。
[0024]
本发明的有益效果为：本发明以phreeqcrm为反应引擎，可模拟溶解和沉淀反应、离子交换反应、表面络合反应、动力学反应等多种反应过程，该模型基于fortran90程序语言，具有可移植性和可扩展性强的特点，极大的提高了地表水水质模型的适用性。
附图说明
[0025]
图1为本发明的方法流程示意图。
[0026]
图2为本发明一个实例提供的一种基于热力学平衡理论的地表水多类型污染物运移模型的流程示意图。
[0027]
图3为本发明一个实例中考虑顺序一级降解反应的地表水迁移模型的第一顺序物质模拟结果与解析解的验证图。
[0028]
图4为本发明一个实例中考虑顺序一级降解反应的地表水迁移模型的第二顺序物质模拟结果与解析解的验证图。
[0029]
图5为本发明一个实例中考虑顺序一级降解反应的地表水迁移模型的第三顺序物质模拟结果与解析解的验证图。
[0030]
图6为本发明一个实例中考虑化学反应的醋酸在海水中泄露后引起当地海域ph值变化范围图。
[0031]
图7为本发明一个实例中考虑化学反应的醋酸在海水中泄露后引起当地海域co2分压值变化范围图。
[0032]
图8为本发明一个实例中考虑化学反应的醋酸在海水中泄露后引起当地海域caco3值变化范围图。
具体实施方式
[0033]
如图1所示，一种耦合地表水污染物完整地球化学反应和运移过程的模拟方法，包括如下步骤：
[0034]
步骤1、收集研究区域的地表水水文数据；
[0035]
收集研究区域的地表水文数据包括对该区域的地形、水位、盐度、流速、流量、降雨量等的监测。
[0036]
步骤2、构建研究区域的地表水动力模型；
[0037]
地表水动力模型的控制方程如下所示：
[0038][0039]
式中：v为速度，p为自由表面压力，r为水体容重，ρm为水体密度，z为水深，υ为运动粘滞系数。
[0040]
对所述的地表水动力模型控制方程(1)进行离散，确定研究区域的边界，采用三角网格剖分的方式对研究区域进行剖分，同时进行地形数据的插值与区域边界的定义，模型采用fortran90语言编写。
[0041]
步骤3、构建研究区域的物质输运模型；
[0042]
物质输运模型控制方程如下所示：
[0043][0044]
式中：h为水深；c为等深的物质浓度；t为时间；u和v分别为x、y方向上的速度；k
xx
、k
xy
、k
yx
以及k
yy
为二维扩散系数张量，s是源汇项；扩散系数的推导如下；
[0045]kxx
＝k
l cos2(θ)+k
t
sin2(θ)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0046]kxy
＝k
yx
＝(k
l-k
t
)cos(θ)sin(θ)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
[0047]kyy
＝k
l sin2(θ)+k
t
cos2(θ)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0048]
式中：θ流速方向与x轴的夹角；k
l
、k
t
分别是横向、纵向扩散系数；顺流方向及垂直流动方向的扩散系数，可利用如下elder公式：
[0049][0050][0051]
式中：c为谢才系数；α,β为常数，其理论值分别取5.93和0.23，当水体紊动较强时，可取到13和1.2，横、纵方向上物质扩散特征差异非常明显。
[0052]
进一步的，利用算子分裂算法，将上式(2)分为如下对流项和扩散项方程：
[0053]
对流项：
[0054]
扩散项：
[0055]
采用有限体积法对上面对流项和扩散项进行差分计算，其算子计算形式如下：
[0056]
(hc)
n+1
＝l
adv
l
diff
(hc)nꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(10)
[0057]
(hc)
n+2
＝l
diff
l
adv
(hc)
n+1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)
[0058]
式中：l
adv
表示对流算子，l
diff
表示扩散算子。
[0059]
步骤4、收集研究区域的地表水水质数据；
[0060]
收集研究区域的地表水水质数据包括水的ph值、pe值、温度和准备研究元素的浓度等。
[0061]
步骤5、基于热力学平衡理论构建研究区域的多类型污染物模型；基于地球化学包phreeqcrm为基础建立多类型污染物模型。
[0062]
具体的，构建多类型污染物模型部分包括两大部分，初始化phreeqcrm部分和执行反应部分。初始化phreeqcrm部分包括创建phreeqcrm，初始化phreeqcrm，在phreeqcrm中设置模型的初始与边界条件。执行反应部分主要执行数据的传递与实现溶液的反应。
[0063]
具体的，初始化phreeqcrm部分主要涉及下列函数：需要使用rm_create函数创建phreeqcrm线程，使用rm_setcomponenth2o函数设置该实例中h2o的组分形式，使用rm_setfileprefix函数设置输出文件名，使用rm_setunitssolution函数设置溶液单位，使用rm_runfile函数载入初始反应文件，使用rm_settimestep函数设置反应时间步长，使用rm_settime函数设置当前反应时间，使用rm_initialphreeqc2module函数进行水质模型初始化，使用rm_runcells函数进行初始步运行以获得初始浓度场。
[0064]
具体的，执行反应部分主要涉及下列函数：需要使用rm_setconcentrations函数将物质输运浓度结果传递到水质模型中，使用rm_runcells函数对当前时间步进行化学反应模拟，使用rm_getconcentrations函数将反应计算结果传回物质输运模型中。
[0065]
步骤6、将步骤2中的地表水动力模型、步骤3中的物质输运模型与步骤5中的多类型污染物模型耦合，得到地表水多类型污染物运移模型；
[0066]
具体的，将地表水动力模型、物质输运模型与多类型污染物模型进行代码层面的全耦合，将使用fortran90编写的地表水动力、物质输运程序与fortran版本的phreeqcrm进行统一编译，在每个时间步长内同步求解水动力模型，物质输运模型与水质模型。每一时间步长的计算流程如下：先求解当前时刻的水动力场，然后基于当前计算所得水动力场进行物质输运计算，得到当前时刻各物质的浓度分布场，最后，将当前各物质的浓度结果输入水质模型中，进行地球化学反应过程的计算，得到下一时间步长的各物质初始浓度场。
[0067]
如图2所示，为本发明一个实施例提供的一种基于热力学平衡理论的地表水多类型污染物运移模型的流程图示意图。
[0068]
此外，本发明通过两个案例开展了方法可靠性与实用性验证。分别是考虑顺序一级降解反应的地表水迁移模型，以及一个考虑动力学反应的醋酸在海洋中泄露的模型。结果如图3到图8所示。图3到图5是本技术一个实施例提供的考虑顺序一级降解反应的地表水迁移模型验证结果示意图；图6到图8是本发明一个实施例提供的考虑动力学反应的醋酸在
海洋中泄露的模型模拟结果图。由图可见，进一步说明了本发明具有较好的可靠性与实用性。