CN103455667B - 充气法治理承压含水层海水入侵的数值模拟方法 - Google Patents
充气法治理承压含水层海水入侵的数值模拟方法 Download PDFInfo
- Publication number
- CN103455667B CN103455667B CN201310366837.3A CN201310366837A CN103455667B CN 103455667 B CN103455667 B CN 103455667B CN 201310366837 A CN201310366837 A CN 201310366837A CN 103455667 B CN103455667 B CN 103455667B
- Authority
- CN
- China
- Prior art keywords
- kappa
- beta
- phase
- boundary
- water
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 76
- 239000013535 sea water Substances 0.000 title claims abstract description 69
- 238000005273 aeration Methods 0.000 title claims abstract description 30
- 230000009545 invasion Effects 0.000 title claims abstract description 30
- 238000004088 simulation Methods 0.000 title claims abstract description 12
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims abstract description 94
- 239000013505 freshwater Substances 0.000 claims abstract description 35
- 238000004364 calculation method Methods 0.000 claims abstract description 23
- 230000005514 two-phase flow Effects 0.000 claims abstract description 16
- 230000000694 effects Effects 0.000 claims abstract description 12
- 238000012795 verification Methods 0.000 claims abstract description 6
- 239000012071 phase Substances 0.000 claims description 100
- 239000007791 liquid phase Substances 0.000 claims description 38
- 238000009792 diffusion process Methods 0.000 claims description 25
- 230000008859 change Effects 0.000 claims description 23
- 230000014509 gene expression Effects 0.000 claims description 18
- 238000002347 injection Methods 0.000 claims description 18
- 239000007924 injection Substances 0.000 claims description 18
- 239000011148 porous material Substances 0.000 claims description 16
- 230000003068 static effect Effects 0.000 claims description 16
- 239000012267 brine Substances 0.000 claims description 15
- HPALAKNZSZLMCH-UHFFFAOYSA-M sodium;chloride;hydrate Chemical compound O.[Na+].[Cl-] HPALAKNZSZLMCH-UHFFFAOYSA-M 0.000 claims description 15
- 239000007788 liquid Substances 0.000 claims description 14
- 230000009471 action Effects 0.000 claims description 13
- 230000035699 permeability Effects 0.000 claims description 13
- 239000002689 soil Substances 0.000 claims description 13
- 150000003839 salts Chemical class 0.000 claims description 11
- 230000004907 flux Effects 0.000 claims description 10
- 230000001186 cumulative effect Effects 0.000 claims description 8
- 229920006395 saturated elastomer Polymers 0.000 claims description 8
- 239000011159 matrix material Substances 0.000 claims description 7
- 230000008569 process Effects 0.000 claims description 6
- 230000002829 reductive effect Effects 0.000 claims description 6
- FAPWRFPIFSIZLT-UHFFFAOYSA-M Sodium chloride Chemical compound [Na+].[Cl-] FAPWRFPIFSIZLT-UHFFFAOYSA-M 0.000 claims description 5
- 230000005012 migration Effects 0.000 claims description 5
- 238000013508 migration Methods 0.000 claims description 5
- 239000011780 sodium chloride Substances 0.000 claims description 5
- 230000001133 acceleration Effects 0.000 claims description 4
- 238000012512 characterization method Methods 0.000 claims description 4
- 239000012530 fluid Substances 0.000 claims description 4
- 239000000203 mixture Substances 0.000 claims description 4
- 230000004069 differentiation Effects 0.000 claims description 3
- 239000006185 dispersion Substances 0.000 claims description 3
- 238000012821 model calculation Methods 0.000 claims description 3
- 238000012545 processing Methods 0.000 claims description 3
- 238000003672 processing method Methods 0.000 claims description 3
- 238000004458 analytical method Methods 0.000 abstract description 3
- 239000007789 gas Substances 0.000 description 44
- 238000011160 research Methods 0.000 description 8
- 239000003673 groundwater Substances 0.000 description 6
- 230000003247 decreasing effect Effects 0.000 description 5
- 230000002706 hydrostatic effect Effects 0.000 description 3
- 238000004519 manufacturing process Methods 0.000 description 3
- 238000005086 pumping Methods 0.000 description 3
- IYLGZMTXKJYONK-ACLXAEORSA-N (12s,15r)-15-hydroxy-11,16-dioxo-15,20-dihydrosenecionan-12-yl acetate Chemical compound O1C(=O)[C@](CC)(O)C[C@@H](C)[C@](C)(OC(C)=O)C(=O)OCC2=CCN3[C@H]2[C@H]1CC3 IYLGZMTXKJYONK-ACLXAEORSA-N 0.000 description 2
- 230000004888 barrier function Effects 0.000 description 2
- 230000006866 deterioration Effects 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000000670 limiting effect Effects 0.000 description 2
- 238000013178 mathematical model Methods 0.000 description 2
- IYLGZMTXKJYONK-UHFFFAOYSA-N ruwenine Natural products O1C(=O)C(CC)(O)CC(C)C(C)(OC(C)=O)C(=O)OCC2=CCN3C2C1CC3 IYLGZMTXKJYONK-UHFFFAOYSA-N 0.000 description 2
- VREFGVBLTWBCJP-UHFFFAOYSA-N alprazolam Chemical compound C12=CC(Cl)=CC=C2N2C(C)=NN=C2CN=C1C1=CC=CC=C1 VREFGVBLTWBCJP-UHFFFAOYSA-N 0.000 description 1
- 239000008346 aqueous phase Substances 0.000 description 1
- 230000007797 corrosion Effects 0.000 description 1
- 238000005260 corrosion Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000002950 deficient Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 238000006073 displacement reaction Methods 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 239000007792 gaseous phase Substances 0.000 description 1
- 230000001617 migratory effect Effects 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 239000002002 slurry Substances 0.000 description 1
- 239000000243 solution Substances 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 230000001052 transient effect Effects 0.000 description 1
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Investigation Of Foundation Soil And Reinforcement Of Foundation Soil By Compacting Or Drainage (AREA)
Abstract
本发明涉及防止海水入侵技术领域。为治理承压含水层海水入侵,本发明采用的技术方案是:充气法治理承压含水层海水入侵的数值模拟方法,包括如下步骤:步骤一:建立地下水气-液二相流及溶质运移模型,包括基本控制方程及辅助方程,不考虑温度对***的影响;步骤二:模型求解;步骤三:模型边界条件确定:步骤四:模型验证:运用地下水气-液二相流及溶质运移模型,模拟分析承压含水层受咸、淡水共同作用的静力平衡情况,并与已有文献的计算结果进行对比验证;步骤五:分析充气法治理承压含水层海水入侵的效果。本发明主要应用于防止海水入侵。
Description
技术领域
本发明涉及防止海水入侵技术领域,具体讲,涉及充气法治理承压含水层海水入侵的数值模拟方法。
背景技术
近年来,随着沿海地区经济的快速发展,地下水超采越来越严重,造成地下水位持续下降,海水侵入地下含水层,形成海水入侵现象。我国海岸线长达1.8×104Km,出现海水入侵的城市有十几座,入侵面积己超过10000Km2,由此引发的耕地土壤盐碱化、地下设备腐蚀、水质恶化等问题严重阻碍了沿海地区的经济和社会发展。因此,开展海水入侵深入研究,提出防治和减轻海水入侵的措施,具有重要的科学价值和社会意义。
目前,防治海水入侵的措施主要有:限制地下水开采量、迁移抽水井、人工回灌补给、人工抽取咸水、修建地下帷幕。限制地下水开采量并不能从根本上解决海水入侵,有可能会加剧供需水矛盾。迁移抽水井的成本较高,也可能面临其它问题,如建筑层或含水层的尺寸条件。人工回灌需要大量的淡水资源,对于缺水地区,需兴建大量引水工程,势必增加成本费用。人工抽取咸水能够减少海水入侵的程度,但主要问题在于如何处理抽取的咸水。地下帷幕包括实体帷幕、水力屏障和地下充气墙等,实体帷幕和水力屏障的造价昂贵,可能会引起水质恶化和污染,而修建地下充气墙的方法造价低,无需注水或泥浆即可形成挡水帷幕,也不会引起二次污染,较为理想。然而,尽管充气法防治海水入侵在理论上是可行的,但针对充气法防治效果展开试验的研究费用较为昂贵,且耗时耗力。
发明内容
为克服现有技术的不足,治理承压含水层海水入侵,为达到上述目的,本发明采用的技术方案是:充气法治理承压含水层海水入侵的数值模拟方法,包括如下步骤:
步骤一:建立地下水气一液二相流及溶质运移模型,包括基本控制方程及辅助方程,不考虑温度对***的影响,具体如下:
模型的基本控制方程为:
式中,Mκ表示κ组分(纯水w、参考盐水b和空气a)的累积质量密度,Fκ为κ组分的流量密度,包括平流流量密度和扩散流量密度qκ为组分κ的源汇项;
步骤二:模型求解:以TOUGH2/EOS7为工具,空间上采用积分形式的有限差分方法IFDM进行离散,时间上采用一阶向后差分的全隐式方法进行离散,模型的线性化采用Newton-Raphson迭代方法,最后得到大型稀疏系数矩阵的线性方程组;具体如下:
(1)空间上采用积分有限差分法IFDM进行离散
首先将计算域离散成子单元,其性质由形心点代表,分别对各个单元的质量平衡方程的积分形式进行空间离散。对于任意单元n,单元体积为Vn,边界面为Γn,单元的质量平衡方程的积分形式如下:
式中,n为表面单元dΓn的单位法向量,指向控制单元体内为正。
引入适当的体积平均值,有
式中, 为Mκ,gκ在Vn上的平均值。
Γn上的面积分可近似为其所包含的各个表面Anm的面积分的平均值之和,有
式中,m为与单元n相邻的所有单元,Anm是单元n和m相邻的交界面,是Fκ在面Anm上的内法线方向的平均值;
将式(3)、(4)和(5)代入到式(2)中,得到一组关于时间的一阶微分方程组
(2)时间上采用一阶向后差分方法进行离散
对式(6)的时间微分采用一阶向后差分方法,得到任意单元的全隐式非线性方程组,见式(7):
式中,引入了组分κ=w,b,a的余量Δt为时间步长,上标k和k+l表示两相邻的时间步长指标。其中, 分别表示k、k+1时刻Mκ在单元体积Vn上的平均值,表示k+1时刻Fκ在面Anm上的内法线方向的平均值,表示k+1时刻qκ在单元体积Vn上的平均值;右端的流量项和源汇项均采用新的时间步长值。
(3)Newton-Raphson迭代方法
运用Newton-Raphson迭代方法进行线性化,引入迭代指标p,对式(7)中的余量在迭代步p+1处进行泰勒级数展开,只保留一阶项,得到包含3×NEL即计算域内单元数个方程的线性方程组,并且以两迭代步的增量为未知量;最后得到大型稀疏系数矩阵的线性方程组,如式(8):
步骤三:模型边界条件确定:模型计算中的边界条件包括Dirichlet边界条件和Neumann边界条件两种,其数学处理方法如下:
(1)Dirichlet边界条件
Dirichlet边界条件上,边界条件单元的主要变量在计算过程中保持不变,为此,设定边界条件单元的体积非常大,当边界单元的体积相对于土体单元很大时,与土体单元的流量交换将不会改变边界单元的主要变量值;
①对于空气边界,其边界条件单元的相态为仅有气相状态,主要变量为孔隙气压力pg、参考盐水占相态的质量百分数Xb、空气占气相的质量分数和温度T,对于与大气接触的边界上pg=patm,对于有空气超压作用的边界上,如人工充气墙边界上,pg=patm+Δp,其中patm为大气压力,Δp为充气压力;Xb=0.0;T为气温;
②对于已知水头边界,包括地下淡水水头边界和海水水头边界,其边界条件单元的相态均为仅有液相状态,主要变量为孔隙气压力pg、参考盐水占液相的质量百分数Xb、空气占液相的质量分数和温度T,由于液相饱和状态下毛细压力为零,因此有pg=patm+ρlgh,其中ρl为淡水或海水的密度,h为边界上的水深;地下淡水水头边界上Xb等于零,海水边界上的Xb可根据参考盐水的特性以及海水的盐度计算出来;
(2)Neumann边界条件
Neumann边界条件描述的是***与外界的流量交换情况,边界条件单元的单位流量对应式中的源汇项,流入为正,可以是常量,也可以随时间变化,对于不透水边界,是一类特殊的Neumann边界条件,边界上的流量为零,积分形式的有限差分法的处理很简单,就是不设置与之相邻的边界条件单元,则与不透水面单元的流量交换为零;
步骤四:模型验证:运用地下水气-液二相流及溶质运移模型,模拟分析承压含水层受咸、淡水共同作用的静力平衡情况,并与已有文献的计算结果进行对比验证;
步骤五:根据地下水气-液二相流及溶质运移模型,以咸-淡水静力平衡情况作为初始条件,施加充气作用,分析充气法治理承压含水层海水入侵的效果:
(1)盐度变化:充气结束时,相对盐度等值线与不透水底板的交点向海水方向撤退,海水入侵范围减小;
(2)水、气相压力及流场变化:水相压力和气相压力的变化规律基本一致,分为三个阶段:迅速增大至峰值、由峰值迅速减小、缓慢减小趋近于相对稳定值;另外,由于毛细压力的作用,观察点的水相压力略小于气相压力;在注入区及含水层顶部,水、气相压力增大的较多,形成指向海水侧的水力梯度,从而驱替入侵的咸水退出含水层;
(3)空气损失变化:充气初始时期的空气损失较小,然后迅速上升,之后逐渐趋于平稳,达到相对稳定值。
(1)累积质量密度Mκ等于β相中组分κ的质量之和,其表达式为:
式中,β表示流体相(气相或液相),液相为纯水和参考盐水的混合物,φ表示孔隙率,Sβ为β相的饱和度,ρβ为β相的密度,为κ组分占β相的质量百分数;
(2)平流流量的表达式为:
式中,Fβ为β相的平流流量,符合达西定律,其多相流形式的表达式为:
其中,k为固有渗透系数,krβ为β相的相对渗透系数,μβ为β相的粘滞性系数,pβ为β相的孔隙压力,g为重力加速度矢量;
液相的孔隙压力pl为气相压力pg和毛细压力pc(负值)之和:
pl=pg+pc(12)
关于毛细压力~饱和度的关系,采用Leverett模型来表征:
pc=p0γ[1.417(1-Se)-2.120(1-Se)2+1.263(1-Se)3](13)
式中,p0为进气值,γ为表面张力,有效水饱和度Se的表达式为Se=(Sl-Slr)/(Sls-Slr),Sl为液相饱和度,Slr为液相剩余饱和度,Sls为液相最大饱和度;
关于相对渗透率~饱和度的关系吗,采用FattandKlikoff模型来表征:
krg=(1-Se)3(气相)(15)
(3)扩散流量的表达式为:
式中,为κ组分在β相中的分子扩散系数;τ0τβ(Sβ)为迂曲度,τ0是多孔介质特性参数,τβ=τβ(Sβ)是饱和度的函数;
分子扩散系数随压力p和温度T变化的计算式为:
式中,θ为温度系数,通常取值为1.8;
TOUGH2中提供了三种迂曲度模型:分子扩散系数τ0一般取1.0
relativepermeabilitymodel:τ0τβ=τ0krβ(Sβ)
Millington-QuirkModel:
constantdiffusivity:τ0τβ=τ0Sβ(18)
Relativepermeabilitymodel为相对渗透率模型,Millington-QuirkModel是由Millington-Quirk提出的非饱和土壤导水率模型,constantdiffusivity为常扩散系数模型。
本发明可带来如下效果:
(1)建立了地下水气-液二相流及溶质运移模型,可以用来模拟充气作用下,地下水中气-液二相流动过程及其对咸-淡水交界面运移的影响。
(2)本发明避免了由于开展实验性研究的高费用和长工期问题,模拟结果合理可信,为进一步展开充气法治理海水入侵的试验性研究提供了科学指导和分析依据。
附图说明
图1Henry’s问题的研究范围与边界条件。
图2稳态下相对盐度为0.5的等值线位置及其与其他文献结果的对比图。
图3稳态下海水入侵的盐度分布。
图4稳态下地下水流速分布,水平流速为零的等值线。
图5稳态下的水压力分布。
图6稳态下的孔隙水压力分布。
图7充气作用的网格剖分及充气位置。
图8相对盐度为0.5的等值线与不透水底板的交点随时间的变化情况。
图9相对盐度为0.25、0.5和0.75的等值线与不透水底板的交点随时间的变化情况。
图10含水层中总盐度随时间的变化情况。
图11充气结束时水相和气相压力、水流和气流及气体饱和度的分布情况。
图12含水层中各观察点在充气过程中的水相和气相压力变化情况。
图13充气过程中的空气损失率变化情况。
图14本发明的模拟方法流程图。
具体实施方式
本发明采用的技术方案是:
步骤一:建立地下水气-液二相流及溶质运移模型,包括基本控制方程及辅助方程,不考虑温度对***的影响。具体如下:
模型的基本控制方程为:
式中,Mκ表示κ组分(纯水w、参考盐水b和空气a)的累积质量密度,Fκ为κ组分的流量密度,包括平流流量密度和扩散流量密度qκ为组分κ的源汇项。
(1)累积质量密度Mκ等于β相中组分κ的质量之和,其表达式为:
式中,β表示流体相(气相或液相),液相为纯水和参考盐水的混合物,φ表示孔隙率,Sβ为β相的饱和度,ρβ为β相的密度,为κ组分占β相的质量百分数。
(2)平流流量的表达式为:
式中,Fβ为β相的平流流量,符合达西定律,其多相流形式的表达式为:
其中,k为固有渗透系数,kγβ为β相的相对渗透系数,μβ为β相的粘滞性系数,pβ为β相的孔隙压力,g为重力加速度矢量。
液相的孔隙压力pl为气相压力pg和毛细压力pc(负值)之和:
pl=pg+pc(5)
关于毛细压力~饱和度的关系,采用Leverett模型来表征:
式中,p0为进气值,γ为表面张力,有效水饱和度Se的表达式为Se=(Sl-Slr)/(Sls-Slr),Sl为液相饱和度,Slr为液相剩余饱和度,Sls为液相最大饱和度。
关于相对渗透率~饱和度的关系吗,采用FattandKlikoff模型来表征:
krg=(1-Se)3(气相)(8)
(3)扩散流量的表达式为:
式中,为κ组分在β相中的分子扩散系数;τ0τβ(Sβ)为迂曲度,τ0是多孔介质特性参数,τβ=τβ(Sβ)是饱和度的函数。
分子扩散系数随压力p和温度T变化的计算式为:
式中,θ为温度系数,通常可取值为1.8。
TOUGH2中提供了三种迂曲度模型:分子扩散系数τ0一般取1.0
τ0τβ=τ0krβ(Sβ)(relativepermeabilitymodel))
(Millington-QuirkModel)
τ0τβ=τ0Sβ(constantdiffusivity)(11)
Relativepermeabilitymodel为相对渗透率模型,Millington-QuirkModel是由Millington-Quirk提出的非饱和土壤导水率模型,constantdiffusivity为常扩散系数模型。
步骤二:模型求解:以TOUGH2/EOS7为工具,空间上采用积分形式的有限差分方法(IFDM)进行离散,时间上采用一阶向后差分的全隐式方法进行离散,模型的线性化采用Newton-Raphson迭代方法,最后得到大型稀疏系数矩阵的线性方程组;具体如下:
(1)空间上采用积分有限差分法(IFDM)进行离散
首先将计算域离散成子单元,其性质由形心点代表,分别对各个单元的质量平衡方程的积分形式进行空间离散。对于任意单元n,单元体积为Vn,边界面为Γn,单元的质量平衡方程的积分形式如下:
式中,n为表面单元dΓn的单位法向量,指向控制单元体内为正。
引入适当的体积平均值,有
式中, 为Mκ,qκ在Vn上的平均值。
Γn上的面积分可近似为其所包含的各个表面Anm的面积分的平均值之和,有
式中,m为与单元n相邻的所有单元,Anm是单元n和m相邻的交界面,是Fκ在面Anm上的内法线方向的平均值。
将式(13)、(14)和(15)代入到式(12)中,得到一组关于时间的一阶微分方程组
(2)时间上采用一阶向后差分方法进行离散
对式(16)的时间微分采用一阶向后差分方法,得到任意单元的全隐式非线性方程组,见式(17):
式中,引入了组分κ=w,b,a的余量Δt为时间步长,上标k和k+1表示两相邻的时间步长指标;右端的流量项和源汇项均采用新的时间步长值,这种全隐式方法提高了模型求解的数值稳定性。
(3)Newton-Raphson迭代方法
运用Newton-Raphson迭代方法进行线性化。引入迭代指标p,对式(17)中的余量在迭代步p+1处进行泰勒级数展开,只保留一阶项,得到包含3×NEL(计算域内单元数)个方程的线性方程组,并且以两迭代步的增量为未知量;最后得到大型稀疏系数矩阵的线性方程组,如式(18):
步骤三:模型边界条件确定:模型计算中的边界条件包括Dirichlet边界条件和Neumann边界条件两种,其数学处理方法如下:
(1)Dirichlet边界条件
Dirichlet边界条件上,边界条件单元的主要变量在计算过程中保持不变,为此,设定边界条件单元的体积非常大,如1×1040m3。当边界单元的体积相对于土体单元很大时,与土体单元的流量交换将不会改变边界单元的主要变量值。
①对于空气边界,其边界条件单元的相态为仅有气相状态,主要变量为pg,Xb,T,对于与大气接触的边界上pg=patm,对于有空气超压作用的边界上,如人工充气墙边界上,pg=patm+Δp;Xb=0.0;T为气温。
②对于已知水头边界,包括地下淡水水头边界和海水水头边界,其边界条件单元的相态均为仅有液相状态,主要变量为pg、Xb、T,由于液相饱和状态下毛细压力为零,因此有pg=patm+ρlgh,其中ρl为淡水或海水的密度,h为边界上的水深;地下淡水水头边界上Xb等于零,海水边界上的Xb可根据参考盐水的特性以及海水的盐度计算出来。
(2)Neumann边界条件
Neumann边界条件描述的是***与外界的流量交换情况,边界条件单元的单位流量对应式中的源汇项,流入为正,可以是常量,也可以随时间变化。对于不透水边界,是一类特殊的Neumann边界条件,边界上的流量为零,积分形式的有限差分法的处理很简单,就是不设置与之相邻的边界条件单元,则与不透水面单元的流量交换为零。
步骤四:模型验证:运用地下水气-液二相流及溶质运移模型,模拟分析承压含水层受咸、淡水共同作用的静力平衡情况,并与已有文献的计算结果进行对比验证;
步骤五:根据地下水气-液二相流及溶质运移模型,以咸-淡水静力平衡情况作为初始条件,施加充气作用,分析充气法治理承压含水层海水入侵的效果:
(1)盐度变化:充气结束时,相对盐度等值线与不透水底板的交点向海水方向撤退,海水入侵范围减小;
(2)水、气相压力及流场变化:水相压力和气相压力的变化规律基本一致,可分为三个阶段:迅速增大至峰值、由峰值迅速减小、缓慢减小趋近于相对稳定值。另外,由于毛细压力的作用,观察点的水相压力略小于气相压力。在注入区及含水层顶部,水、气相压力增大的较多,形成指向海水侧的水力梯度,从而驱替入侵的咸水退出含水层;
(3)空气损失变化:充气初始时期的空气损失较小,然后迅速上升,之后逐渐趋于平稳,达到相对稳定值。
下面结合附图,针对充气法治理承压含水层海水入侵的数值模拟方法进行具体说明,其步骤如下:
(1)建立数学模型:该数学模型描述的是承压含水层未受污染的淡水受海水入侵的情况,即经典的Henry’s问题。模型选取200米长、100米深、1米厚的含水层为研究对象,其上、下两层不透水顶板、底板之间的含水层均质、各向同性,内陆侧边界上有恒定的淡水流入量(或者等效的淡水作用下的静水压力),海水侧边界上分布着密度较大的海水作用下的静水压力,是一个理想化的准二维模型。其研究范围及边界条件见图1,模型相关参数取值见表1。
假定含水层***处于恒温状态(T=25℃)。
表1模型相关参数取值
模型的基本控制方程为:
式中,Mκ表示κ组分(纯水w、参考盐水b和空气a)的累积质量密度,Fκ为组分κ的流量密度,包括平流流量密度和扩散流量密度qκ为组分κ的源汇项。
①累积质量密度Mκ等于β相中组分κ的质量之和,其表达式为:
式中,β表示流体相(气相或液相),液相为纯水和参考盐水的混合物,φ表示孔隙率,Sβ为β相的饱和度,ρβ为β相的密度,为κ组分占β相的质量百分数。
②平流流量的表达式为:
式中,Fβ为β相的平流流量,符合达西定律,其多相流形式的表达式为:
其中,k为固有渗透系数,krβ为β相的相对渗透系数,μβ为β相的粘滞性系数,pβ为β相的孔隙压力,g为重力加速度矢量。
液相的孔隙压力pl为气相压力pg和毛细压力pc(负值)之和:
pl=pg+pc(5)
关于毛细压力~饱和度的关系,采用Leverett模型来表征:
式中,p0为进气值,γ为表面张力,有效水饱和度Se的表达式为Se=(Sl-Slr)/(Sls-Slr),Sl为液相饱和度,Slr为液相剩余饱和度,Sls为液相最大饱和度。
关于相对渗透率~饱和度的关系吗,采用FattandKlikoff模型来表征:
krg=(1-Se)3(气相)(8)
③扩散流量的表达式为:
式中,为κ组分在β相中的分子扩散系数;τ0τβ(Sβ)为迂曲度,τ0是多孔介质特性参数,τβ=τβ(Sβ)是饱和度的函数。
分子扩散系数随压力p和温度T变化的计算式为:
式中,θ为温度系数,通常可取值为1.8。
TOUGH2中提供了三种迂曲度模型:分子扩散系数τ0一般取1.0
Relativepermeabilitymodel为相对渗透率模型,Millington-QuirkModel是由Millington-Quirk提出的非饱和土壤导水率模型,constantdiffusivity为常扩散系数模型。
(2)模型求解:对计算域进行网格剖分,含水层被划分为5m×5m×1m的8节点的均匀六面体单元,共800个单元,1722个节点。以TOUGH2/EOS7为工具,空间上采用积分形式的有限差分方法(IFDM)进行离散,时间上采用一阶向后差分的全隐式方法进行离散,见式(12);模型的线性化采用Newton-Raphson迭代方法,最后得到大型稀疏系数矩阵的线性方程组,见式(13):
式中,引入了组分κ=w,g,a的余量 表示Mκ在控制体积Vn上的平均值;表示qκ在控制体积Vn上的平均值;m为与单元n相邻的所有单元;Anm是单元n和m相邻的交界面的面积;是Fκ在面Anm上的内法线方向的平均值;△t为时间步长;上标k和k+1表示两相邻的时间步指标。
引入迭代指标p,对式(12)中余量在迭代步p+1处进行泰勒级数展开,只保留一阶项,得到包含3×NEL(计算域内单元数)个方程的线性方程组,并且以两迭代步的增量为未知量:
(3)模型边界条件确定:假定***处于恒温状态T=25℃,在初始条件下,计算域均为液相饱和状态,毛细压力可以忽略,相对渗透系数取值1.0,主要变量为pg、Xb、T。因此有pg=patm+ρlgh,其中ρl为淡水或海水的密度,h为边界上的水深。对于淡水侧(左侧),作用淡水源汇项,将Qin=Vindρw(d为含水层深度100m,ρw为淡水密度)等于7.639×10-3kg/s,均匀分布在左侧边界上,主要变量pg=patm+ρwg(100-z),Xb=0.0,和T=25℃;对于海水侧(右侧),受海水作用,主要变量pg=patm+ρsg(100-z),Xb=1.0,和T=25℃;其中,ρs为海水密度,patm为大气压力,等于1.013×105Pa。初次计算的初始条件为:所有单元液相饱和,压力均为大气压力与静水(淡水)压力之和即pg=patm+ρwg(100-z),Xb=0.0,
以咸、淡水静力平衡为初始条件,施加充气作用后,***的右侧、顶部和底部的边界条件同稳态状态下,对于左侧边界,pg=ρwg(102.5-z),Xb=0.0,和T=25℃;对于充气井边界,主要变量pg=patm+9.0×105,Xb=0.0,和T=25℃。
(4)模型验证:
在静力平衡情况下,相对(海水)盐度Xb=0.5的等值线位置及其与已有文献计算结果的对比如图2。由图可见,本模型的计算结果与其他文献的结果基本一致,从而验证了计算模型模拟海水入侵问题的有效性。
模拟过程中,可得到静力平衡情况下滨海含水层的相对盐度、地下水流速、水压力水头和孔隙水压力水头的分布情况,详见图3~图6。由图可见,含水层在淡水和咸水的共同作用下达到静力平衡,海水沿临海(右)边界上的下半部分大约0~48m的范围内流入含水层,与地下淡水混合形成楔形的咸水入侵体,并在淡水压力作用下沿临海边界上的上半部分大约49~100m的范围流出,形成一个海水入侵环流。在含水层的右上角处地下水流速最大,约为1.59×10-6m/s。地下水流速为零的等值线反映了海水的入侵范围,其与不透水底板的交点约为x=92.0m,即海水入侵含水层的最远距离约为108.0m。
(5)根据地下水气-液二相流及溶质运移模型,分析充气法治理承压含水层海水入侵的效果:
以咸-淡水静力平衡情况作为初始条件,利用钻孔施加充气,当土体水力传导率的范围在10-8~10-4m/s时,可以采用充气法来阻止地下水侵入含水层,有表1数据显示,计算域满足条件。
根据初始的盐度分布情况(图3),海水入侵淡水层的最远距离约为108.0m,将充气井布置在距离海水边界的水平距离为120m处,直径取1.0m。由于空气的密度小于水的密度,空气进入含水层后将上升,因此将注气位置布置在含水层的下半部,选取x=80.5m,15≤z≤30m,如图7。根据压气法的相关经验,只有当钻孔的注气压力大于注气范围内的最大孔隙水压力时,采用充气法驱替咸水才是完全可能的。有图6可知,注气区域的最大孔隙水压力为84.6mH2O(≈8.3×102kpa),因此模型采用的注气压力为9.0×102kpa。注气区的钻孔位置、网格单元如图7所示,其顶部是不透水边界,不透水区域以上的部分属性与含水层相同,充气时间设置为一年(365天)。从三维模型看,钻孔占地下水气流的空间很小,产生的阻碍可以忽略。计算过程中,液相剩余饱和度Slr取值0.15,表面张力γ在25℃条件下为0.072N/m,进气值
为了研究充气过程中水相和气相压力随时间的变化情况,取5个观察点进行研究,如图7。其中①、②、③三个点的水平位置相同,用于研究各相压力随垂向位置的变化规律,③、④和⑤的垂向坐标相等,用于研究各相压力随水平位置的变化规律。
①盐度变化
图8为相对盐度0.5的等值线与不透水底板的交点在充气过程中随时间的变化过程。图9为相对盐度为0.25、0.5和0.75的等值线与不透水底板的交点在充气过程中的瞬态变化,由图可见,相对盐度为0.25、0.5和0.75的等值线分别向海水方向迁移了53.1m、51.0m、43.1m,其平均速度为0.145m/d、0.140m/d、0.118m/d。含水层的整体盐度变化如图10,在整个时间段内从176.5mg/l减少到了47.2mg/l。结果表明充气作用下,海水入侵的范围逐渐减少。
②水、气相压力及流场变化:
图11a-b为充气结束时,水相和气相压力及水流和气流的分布情况,由图可见,由于空气的注入引起了各相压力的增大,在含水层顶部和空气注入区附近增大较多;由于空气的密度较小,空气进入含水层后先垂直向上流动,在含水层不透水顶板处聚集,并向两侧扩散;水相压力的变化形成了由空气注入区指向海水边界的水力梯度,使得含水层中的地下水向海水边界流动,从而驱替入侵的咸水退出含水层。图11c所示为充气结束时,气体饱和度的分布情况。另外,由于空气的注入,土体由饱和过渡为非饱和产生了毛细压力,造成水相压力增大值均略小于气相压力增大值。
图12为各个观察点在充气过程中水相和气相压力的变化情况。由图可见,水相压力和气相压力的变化规律基本一致,可分为三个阶段:迅速增大至峰值、由峰值迅速减小、缓慢减小趋近于相对稳定值。图12a所示的点①、②、③的水平坐标相同,气相和水相压力增大值的峰值随着垂向位置的升高依次增大,且含水层顶部的压力增大值的相对稳定值较大。图12b所示的点③、④和⑤的纵坐标相同,均紧靠近含水层不透水顶板,随着距离充气井的水平距离的增大,各相压力增大值的峰值和相对稳定值逐渐减小。
③空气损失变化
如图13所示,空气损失速率由初始时的q=0.23m3/min逐渐增大,最终在t=90d时达到稳定为q=0.824m3/min;在注气时间(一年)内,总空气损失量大约为Q=2.61×107m3。然而充气结束后,作为水源地,需从含水层底部抽水,如此势必打破新建立的气体平衡,导致空气损失和注气成本增加。因此,开采井应布置在空气扩散的***。
综上所述,对本实施例总结如下:
(1)针对滨海地区淡水资源比较匮乏的特点,采用数值模拟的方法研究了充气法治理承压含水层海水入侵的效果。首先利用地下水气-液二相流及溶质运移模型模拟分析了经典的Henry问题的咸-淡水静力平衡情况,通过与以前文献结果的对比,验证了所采用模型的有效性。然后基于咸-淡水静力平衡情况,模拟分析了充气作用对含水层盐度、各相压力及空气损失的影响。
(2)对模拟结果的分析表明,由于空气的密度小于地下水的密度,空气注入含水层后向上流动,并聚集在含水层顶板进而向两侧流动,因此水相和气相压力在注入区及含水层顶部增大的较多,并形成了指向海水侧的水力梯度,驱替入侵的咸水向含水层外流动,从而减小了海水入侵的范围。
(3)含水层各相压力的增大值及空气损失均在充气作用一段时间后达到某一相对稳定值,说明充气法防治海水入侵的效果及单位成本是比较稳定的;但是充气法应用于实际工程还需要很多研究工作,包括充气法的影响因素、成本及其可能存在的风险研究等。
Claims (1)
1.一种充气法治理承压含水层海水入侵的数值模拟方法,其特征是,包括如下步骤:
步骤一:建立地下水气-液二相流及溶质运移模型,包括基本控制方程及辅助方程,不考虑温度对***的影响,具体如下:
模型的基本控制方程为:
式中,Mκ表示κ组分(纯水w、参考盐水b和空气a)的累积质量密度,Fκ为κ组分的流量密度,包括平流流量密度和扩散流量密度qκ为组分κ的源汇项;
步骤二:模型求解:以TOUGH2/EOS7为工具,空间上采用积分形式的有限差分方法IFDM进行离散,时间上采用一阶向后差分的全隐式方法进行离散,模型的线性化采用Newton-Raphson迭代方法,最后得到大型稀疏系数矩阵的线性方程组;具体如下:
(1)空间上采用积分有限差分法IFDM进行离散
首先将计算域离散成子单元,其性质由形心点代表,分别对各个单元的质量平衡方程的积分形式进行空间离散,对于任意单元n,单元体积为Vn,边界面为Γn,单元的质量平衡方程的积分形式如下:
式中,n为表面单元dΓn的单位法向量,指向控制单元体内为正;
引入适当的体积平均值,有
式中,为Mκ,qκ在Vn上的平均值;
Γn上的面积分可近似为其所包含的各个表面Anm的面积分的平均值之和,有
式中,m为与单元n相邻的所有单元,Anm是单元n和m相邻的交界面,是Fκ在面Anm上的内法线方向的平均值;
将式(3)、(4)和(5)代入到式(2)中,得到一组关于时间的一阶微分方程组
(2)时间上采用一阶向后差分方法进行离散
对式(6)的时间微分采用一阶向后差分方法,得到任意单元的全隐式非线性方程组,见式(7):
式中,引入了组分κ=w,b,a的余量△t为时间步长,上标k和k+1表示两相邻的时间步长指标;其中,分别表示k、k+1时刻Mκ在单元体积Vn上的平均值,表示k+1时刻Fκ在面Anm上的内法线方向的平均值,表示k+1时刻qκ在单元体积Vn上的平均值;右端的流量项和源汇项均采用新的时间步长值;
(3)Newton-Raphson迭代方法
运用Newton-Raphson迭代方法进行线性化,引入迭代指标p,对式(7)中的余量在迭代步p+1处进行泰勒级数展开,只保留一阶项,得到包含3×NEL即计算域内单元数个方程的线性方程组,并且以两迭代步的增量为未知量;最后得到大型稀疏系数矩阵的线性方程组,如式(8):
步骤三:模型边界条件确定:模型计算中的边界条件包括Dirichlet边界条件和Neumann边界条件两种,其数学处理方法如下:
(1)Dirichlet边界条件
Dirichlet边界条件上,边界条件单元的主要变量在计算过程中保持不变,为此,设定边界条件单元的体积非常大,当边界单元的体积相对于土体单元很大时,与土体单元的流量交换将不会改变边界单元的主要变量值;
①对于空气边界,其边界条件单元的相态为仅有气相状态,主要变量为孔隙气压力pg、参考盐水占相态的质量百分数Xb、空气占气相的质量分数和温度T,对于与大气接触的边界上pg=patm,对于有空气超压作用的边界上,如人工充气墙边界上,pg=patm+△p,其中patm为大气压力,Δp为充气压力;Xb=0.0;T为气温;
②对于已知水头边界,包括地下淡水水头边界和海水水头边界,其边界条件单元的相态均为仅有液相状态,主要变量为孔隙气压力pg、参考盐水占液相的质量百分数Xb、空气占液相的质量分数和温度T,由于液相饱和状态下毛细压力为零,因此有pg=patm+ρlgh,其中ρl为淡水或海水的密度,h为边界上的水深;地下淡水水头边界上Xb等于零,海水边界上的Xb可根据参考盐水的特性以及海水的盐度计算出来;
(2)Neumann边界条件
Neumann边界条件描述的是***与外界的流量交换情况,边界条件单元的单位流量对应式中的源汇项,流入为正,可以是常量,也可以随时间变化,对于不透水边界,是一类特殊的Neumann边界条件,边界上的流量为零,积分形式的有限差分法的处理很简单,就是不设置与之相邻的边界条件单元,则与不透水面单元的流量交换为零;
步骤四:模型验证:运用地下水气-液二相流及溶质运移模型,模拟分析承压含水层受咸、淡水共同作用的静力平衡情况,并与已有文献的计算结果进行对比验证;
步骤五:根据地下水气-液二相流及溶质运移模型,以咸-淡水静力平衡情况作为初始条件,施加充气作用,分析充气法治理承压含水层海水入侵的效果:
(1)盐度变化:充气结束时,相对盐度等值线与不透水底板的交点向海水方向撤退,海水入侵范围减小;
(2)水、气相压力及流场变化:水相压力和气相压力的变化规律基本一致,分为三个阶段:迅速增大至峰值、由峰值迅速减小、缓慢减小趋近于相对稳定值;另外,由于毛细压力的作用,观察点的水相压力略小于气相压力;在注入区及含水层顶部,水、气相压力增大的较多,形成指向海水侧的水力梯度,从而驱替入侵的咸水退出含水层;
(3)空气损失变化:充气初始时期的空气损失较小,然后迅速上升,之后逐渐趋于平稳,达到相对稳定值;
其中,(1)累积质量密度Mκ等于β相中组分κ的质量之和,其表达式为:
式中,β表示流体相,其包括气相或液相,其中液相为纯水和参考盐水的混合物,φ表示孔隙率,Sβ为β相的饱和度,ρβ为β相的密度,为κ组分占β相的质量百分数;
(2)平流流量的表达式为:
式中,Fβ为β相的平流流量,符合达西定律,其多相流形式的表达式为:
其中,k为固有渗透系数,krβ为β相的相对渗透系数,μβ为β相的粘滞性系数,pβ为β相的孔隙压力,g为重力加速度矢量;
液相的孔隙压力pl为气相压力pg和毛细压力pc之和:
pl=pg+pc(12)
关于毛细压力~饱和度的关系,采用Leverett模型来表征:
pc=p0γ[1.417(1-Se)-2.120(1-Se)2+1.263(1-Se)3](13)
式中,p0为进气值,γ为表面张力,有效水饱和度Se的表达式为Se=(Sl-Slr)/(Sls-Slr),Sl为液相饱和度,Slr为液相剩余饱和度,Sls为液相最大饱和度;
关于相对渗透率~饱和度的关系吗,采用FattandKlikoff模型来表征:
液相:
气相:krg=(1-Se)3(15)
(3)扩散流量的表达式为:
式中,为κ组分在β相中的分子扩散系数;τ0τβ(Sβ)为迂曲度,τ0是多孔介质特性参数,τβ=τβ(Sβ)是饱和度的函数;
分子扩散系数随压力p和温度T变化的计算式为:
式中,θ为温度系数,通常取值为1.8;
TOUGH2中提供了三种迂曲度模型:分子扩散系数τ0一般取1.0
relativepermeabilitymodel:τ0τβ=τ0krβ(Sβ)
constantdiffusivity:τ0τβ=τ0Sβ(18)
Relativepermeabilitymodel为相对渗透率模型,Millington-QuirkModel是由Millington-Quirk提出的非饱和土壤导水率模型,constantdiffusivity为常扩散系数模型。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310366837.3A CN103455667B (zh) | 2013-08-20 | 2013-08-20 | 充气法治理承压含水层海水入侵的数值模拟方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310366837.3A CN103455667B (zh) | 2013-08-20 | 2013-08-20 | 充气法治理承压含水层海水入侵的数值模拟方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN103455667A CN103455667A (zh) | 2013-12-18 |
CN103455667B true CN103455667B (zh) | 2016-06-15 |
Family
ID=49738023
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201310366837.3A Expired - Fee Related CN103455667B (zh) | 2013-08-20 | 2013-08-20 | 充气法治理承压含水层海水入侵的数值模拟方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN103455667B (zh) |
Families Citing this family (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104537232A (zh) * | 2014-12-23 | 2015-04-22 | 天津大学 | 考虑Lisse现象的浅层地下水位预测方法 |
CN104778356B (zh) * | 2015-04-08 | 2017-06-16 | 重庆交通大学 | 一种对流‑扩散传质过程的数值模拟方法 |
CN106503463A (zh) * | 2016-10-27 | 2017-03-15 | 天津大学 | 模拟海平面上升情况下海水入侵内陆边界的处理方法 |
CN107038272A (zh) * | 2016-11-11 | 2017-08-11 | 福建工程学院 | 一种重力作用下的盐岩动水溶蚀参数模型的创建方法 |
CN108844881B (zh) * | 2018-08-06 | 2020-08-07 | 湖北工业大学 | 一种基于vg模型预测非饱和土相对渗透系数的方法 |
CN108918390B (zh) * | 2018-08-10 | 2023-11-07 | 河海大学 | 泥浆成膜及测定泥膜固结量、进气量的装置及方法 |
CN111455927A (zh) * | 2020-01-21 | 2020-07-28 | 河海大学 | 一种增加海岛地下淡水储量的海岛外环弱透水层设计方法 |
CN115266521B (zh) * | 2022-07-01 | 2023-05-09 | 中国海洋大学 | 考虑温度影响的海岸带地下水渗流模拟***的工作方法 |
CN115587705B (zh) * | 2022-10-18 | 2023-06-16 | 华中科技大学 | 一种城市气候环境快速评估方法及*** |
CN118114592A (zh) * | 2024-03-07 | 2024-05-31 | 河海大学 | 一种不确定性条件下预测咸水入侵动态的替代建模方法 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN1718946A (zh) * | 2005-06-09 | 2006-01-11 | 大连久鼎特种建筑工程有限公司 | 一种阻止海水入侵的地下截潜与集水方法 |
CN102587451A (zh) * | 2012-03-13 | 2012-07-18 | 中国海洋大学 | 一种防治海水入侵的选择性渗透反应墙技术 |
-
2013
- 2013-08-20 CN CN201310366837.3A patent/CN103455667B/zh not_active Expired - Fee Related
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN1718946A (zh) * | 2005-06-09 | 2006-01-11 | 大连久鼎特种建筑工程有限公司 | 一种阻止海水入侵的地下截潜与集水方法 |
CN102587451A (zh) * | 2012-03-13 | 2012-07-18 | 中国海洋大学 | 一种防治海水入侵的选择性渗透反应墙技术 |
Non-Patent Citations (4)
Title |
---|
充气法解决海水入侵问题的探讨;肖江等;《勘察科学技术》;20000630(第6期);第15-19页 * |
水力帷幕防治海水入侵的数值模拟研究;高学平等;《环境污染与防治》;20061130;第28卷(第11期);第872-875页 * |
非饱和带水–气二相流数值模拟研究;孙冬梅等;《岩土工程学报》;20070430;第29卷(第4期);第560-565页 * |
预测压气法隧道施工中的漏气量的模型研究;孙冬梅等;《岩土工程学报》;20090731;第31卷(第7期);第1030-1036页 * |
Also Published As
Publication number | Publication date |
---|---|
CN103455667A (zh) | 2013-12-18 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN103455667B (zh) | 充气法治理承压含水层海水入侵的数值模拟方法 | |
Salimzadeh et al. | A novel radial jet drilling stimulation technique for enhancing heat recovery from fractured geothermal reservoirs | |
Guo et al. | Experimental study and numerical simulation of hydraulic fracturing tight sandstone reservoirs | |
CN105740563B (zh) | 一种成熟油田二次开发之优势通道识别方法 | |
CN111581854A (zh) | 一种考虑非平衡各向异性相对渗透率的油藏状态预测方法 | |
CN106503463A (zh) | 模拟海平面上升情况下海水入侵内陆边界的处理方法 | |
Li et al. | Ensemble-based relative permeability estimation using B-spline model | |
Possemiers et al. | Application of multiple-point geostatistics to simulate the effect of small-scale aquifer heterogeneity on the efficiency of aquifer thermal energy storage | |
Dapeng et al. | An independent fracturing water-flooding development method for shallow low-permeability thin oil layers in multi-layer sandstone reservoirs | |
Younes et al. | An efficient discontinuous Galerkin-mixed finite element model for variable density flow in fractured porous media | |
Wang et al. | Effects of capillary pressures on two-phase flow of immiscible carbon dioxide enhanced oil recovery in fractured media | |
Chuanzhi et al. | Identification and quantitative description of large pore path in unconsolidated sandstone reservoir during the ultra-high water-cut stage | |
CN116467958A (zh) | 盐湖卤水数值模型构建与补水溶矿效率计算方法 | |
Yuan et al. | Modeling multiphase non-Newtonian polymer flow in IPARS parallel framework | |
Guo et al. | Vertically integrated dual-porosity and dual-permeability models for CO2 sequestration in fractured geological formation | |
US11501043B2 (en) | Graph network fluid flow modeling | |
Benham et al. | Axisymmetric gravity currents in anisotropic porous media | |
Li et al. | Impacts of relative permeability hysteresis on the reservoir performance in CO2 storage in the Ordos Basin | |
Dahle et al. | A model-oriented benchmark problem for CO2 storage | |
Potashev et al. | Numerical modeling of local effects on the petroleum reservoir using fixed streamtubes for typical waterflooding schemes | |
Noble et al. | Fluid Flow through a Fractured Porous Reservoir Using CFD Modeling | |
Awag et al. | Comparison between downdip and updip groundwater flow on early CO2 migration in dipping storage aquifers | |
Bánki et al. | Aiding Reservoir Simulation and Maintenance Using analytical Formulae | |
CN106529199A (zh) | 一种砾岩油藏化学驱井距的确定方法 | |
Olalotiti-Lawal et al. | Streamline-Based Simulation of Carbon Dioxide Sequestration in Saline Aquifers |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20160615 Termination date: 20210820 |