Constraint on low-velocity layer using higher mode Rayleigh waves in the shallow structure research
-
摘要: 在实际中常遇到一类含有低速层的介质，其频散曲线在感兴趣的频率范围并无视觉上的交叉现象，也无明显的低速层指示。针对此类介质，基于地震反射资料中的面波信息，综合实例和数值模拟进行分析，结果显示：如果观测的基阶模式频散曲线不包含对低速层深度敏感的频段，单独基于基阶模式的频散曲线可能无法恢复模型的低速特征；但如果同时利用相同频段的多阶模式信息，即使观测的频散曲线没有明显的低速特征指示，也可以通过多模式频散信息反演重建模型的低速特征。Abstract: Due to high sensitivity to S-wave velocity, Rayleigh-wave dispersion curves of the fundamental and higher modes are usually used to invert near-surface S-wave velocities in engineering geophysical exploration. For the model containing a low-velocity layer, the dispersion curves of the fundamental and higher modes show two typical characteristics. One typical characteristic is that the crossover would be observed between different modes, and the fundamental mode shows obvious indication of low-velocity characteristics in interested frequency ranges. For the other kind of model with low-velocity layers, the dispersion curves have no visual crossing phenomenon in the frequency range of interest, and the low-velocity characteristics may not be observed in the measured dispersion curves. For the latter model containing a low-velocity layer, which is often encountered in practice, investigations on the inversion of multi-mode Rayleigh waves are conducted in this paper based on seismic reflection data. The studies show that if the observed fundamental-mode dispersion curve does not include the frequency band sensitive to the depth of the low-velocity layer, the inversion based on the fundamental-mode alone may not be able to recover the low-velocity characteristics of the model. But the low-velocity layer can be reconstructed accurately by inversion considering both the fundamental and higher mode Rayleigh waves even the observed fundamental mode dispersion curve has no obvious indication of low-velocity characteristics.
-
引言
干热岩（hot dry rock，缩写为HDR）一般是指温度为150—650℃、内部不存在流体或仅有少量流体的高温致密岩体，其埋藏深度一般为3—10 km （许天福等，2016）。我国地热资源丰富（蔺文静等，2012，2013；汪集旸等，2012），合理开发利用干热岩不仅可以提供清洁无污染能源，还可以避免传统化石燃料燃烧所带来的环境污染问题。因此，针对干热岩能源的开发在我国具有重要意义和广阔的发展前景，研究如何合理开发利用干热岩资源至关重要。
目前干热岩开发一般通过增强型地热系统（enhanced geothermal systems，缩写为EGS）实现，即由注入井、压裂改造后的热储层、生产井及发电厂形成闭合环形结构，通过流体流动进行热量的传递和提取（图1）。由于干热岩处于地壳深层部位，且温度较高，在注水采热开发过程中储层处于热-流-固（thermo-hydro-mechanical，缩写为THM）耦合之中，各项物理因素相互影响，进而影响采热效率。研究THM耦合下EGS的采热性能，对我国EGS系统的生产建设与优化具有重要的指导意义，并有助于揭示干热岩开发利用过程中地质、生产因素的影响机理。
干热岩处于高温高压地质环境，其EGS地热开采不仅与渗流（H）、力学（M）有关，还与温度（T）和化学反应（C）等复杂的多物理场耦合过程相关（Zang et al，2014；Gaucher et al，2015；Kumari，Ranjith，2019；Rathnaweera et al，2020）。目前EGS地热开采数值模拟主要集中于热-流（TH）耦合（Shaik et al，2011；胡剑等，2014；丁军锋，王世民，2019；Li et al，2019；Han et al，2020；刘汉青等，2023）或热-流-固耦合（Jacquey et al，2016；Pandey et al，2017；Rutqvist，2017；Yao et al，2018；高坤，2020）下热储层各物理场随时空的演化规律，以及采热效率影响因素等方面。TH耦合重点研究温度场和渗流场对采热的影响，而THM耦合则研究温度场、渗流场及应力场之间的相互作用及其对地热开采的影响。尽管THM耦合相对于TH耦合只是增加了力学场的引入，但是从耦合原理来说，TH耦合只有两项耦合项，而THM耦合则是增加了四项耦合项。因此，相对于TH耦合，THM耦合更加完善且模拟精度更高。多数有关TH耦合或THM耦合下的EGS地热开采模拟研究将基岩、裂隙或者流体各项参数设为常数，而在真实地质情况下这些性质是动态变化的。因此，对于目前EGS地热开采模拟的多物理场耦合并不够完善，不能更好地体现地热开采过程中物理场的两两耦合与相互影响。
影响EGS地热开采效率的因素复杂，Ma等（2020）基于多井注入增强型地热系统（EGS）模型，探究了注采参数和裂缝参数对采热性能的影响；崔翰博等（2020a，b）以青海共和盆地地热田为对象研究了注采参数、渗透率和流体损失等因素对地热开采及渗流场演变规律的影响；段云星和杨浩（2020）以云南腾冲热海热田为地质背景，分析了井间距、注入流量、注入温度、储层渗透率等因素对地热开采的影响。刘汉青等（2023）的研究显示压裂改造的储层的渗透率对地热开采效率起着关键作用。而天然裂缝的存在会使水力裂缝复杂化，对开采模式有重要影响（Guo et al，2019），尤其是对于构造复杂的青海共和地区而言（张森琦等，2018），这些隐伏天然断层对地热开采效率的影响不能忽视。Zhang等（2019b）采用TH耦合模型，提出了EGS断层或者裂缝控制模型的最佳生产方案。Zhang等（2023）采用THM模型讨论了断层和储层的非均匀性对EGS采热的影响。以上研究大多认为断层具有较好的渗透性，但是也有部分断层相对储层具有阻水作用（Lockner et al，2000；Mizoguchi et al，2008）。因此，讨论隐伏断层渗透类型及储层渗透率等地质因素对地热提取的影响很有必要。
本文建立了全耦合的热-流-固耦合模型用以描述地热资源（干热岩）开发过程中的注采过程。在该模型中考虑了流体压力变化、温度变化所引起的储层变形，受储层变形影响的储层孔隙率、渗透率和断层裂隙开度，以及受温度控制的流体性质（密度、黏度等）和流体流动引起的热传导。之后以青海共和干热岩GR1井为例，建立三维数值模型，研究了定流量和定压力两种典型生产模式下，储层渗透率和断层渗透结构等多种因素对地热提取的影响。此外，还讨论了引起定压和定流生产两种模式引起的地热提取的差异，两种生产模式可能引起的生产安全问题并给出了共和地区地热开采综合建议。
1. 热-流-固全耦合模型
干热岩地热开发受到温度场、渗流场和应力场三场之间耦合作用的控制，具体体现在：① 温度场-渗流场耦合控制作用，岩石的温度场因低温流体的干扰而发生变化，而流体温度的变化会改变其自身物理性质进而影响流体流动；② 温度场-应力场耦合控制作用，岩石温度变化将会导致岩石变形，而岩石变形过程中所释放的热量又会反过来影响温度场；③ 渗流场-应力场耦合控制作用，高温岩石中的应力分布受流体压力的干扰，岩石产生变形，而岩石的孔隙度和渗透率受岩石变形的影响。
1.1 储层变形方程
注水过程中，储层的变形控制方程为（Izadi，Elsworth，2015）
$$ G\nabla {\boldsymbol{u}}＋\frac{G}{1－2 \nu }\nabla （ \nabla · {\boldsymbol{u}} ） ＝f＋{\alpha }_{{\mathrm{r}}}{p}_{{\mathrm{r}}}＋3{\beta }_{{\mathrm{T}}}K\nabla {T}_{{\mathrm{r}}}\text{，} $$ (1) 式中，$ \boldsymbol{u} $为位移，$\nu $为泊松比，f为体应力，pr为储层流体压力，αr为储层比奥系数，K为体积模量，G为剪切模量，βT为线性热膨胀系数，Tr为储层温度，$\nabla $为哈密顿算子。
1.2 储层渗流方程
采用储水模型描述多孔的岩石储层流体渗流过程（Fox et al，2013）：
$$ {\rho }_{{\mathrm{r}}}S\frac{{\text{∂}} {p}_{{\mathrm{r}}}}{{\text{∂}} t}＋\nabla · （ \rho {\boldsymbol{v}}_{{\mathrm{r}}} ） ＝{Q}_{{\mathrm{m}}}\text{，} \,\,\,\, {\boldsymbol{v}}_{{\mathrm{r}}}\text＝\text－\frac{{k}_{{\mathrm{r}}}}{\mu } （ \nabla {p}_{{\mathrm{r}}}－{\rho }_{{\mathrm{r}}}{g} ） {\text{.}}$$ (2) 采用裂隙流描述断层中流体渗流过程：
$$ {d}_{{\mathrm{f}}}\frac{{\text{∂}} （ {\phi }_{{\mathrm{f}}}{\rho }_{{\mathrm{f}}} ） }{{\text{∂}} t}＋\nabla · （ {d}_{{\mathrm{f}}}{\rho }_{{\mathrm{f}}}{\boldsymbol{v}}_{{\mathrm{f}}} ） ＝{d}_{{\mathrm{f}}}{Q}_{{\mathrm{m}}}\text{，} \,\,\,\, {\boldsymbol{v}}_{{\mathrm{f}}}\text＝\text－\frac{{k}_{{\mathrm{f}}}}{\mu } （ \nabla {p}_{{\mathrm{f}}}－{\rho }_{{\mathrm{f}}}{g} ） \text{，} \,\,\,\,{k}_{{\mathrm{f}}}＝\frac{{d}_{{\mathrm{f}}}^{2}}{32}\text{，} $$ (3) 式中：下标r和f分别表示多孔的干热岩储层和断层；v为达西流速；p为孔隙压力；S，g，ρw，k，μ，df分别为储水系数、重力加速度、流体密度、岩石渗透率、流体黏度和裂隙开度；Qm表示裂缝与岩石基质之间流体的质量交换。
1.3 热传导方程
1.3.1 储层热传导方程
对于不可变形固体，热传导方程为
$$ \rho_{\mathrm{rs}}C_{\mathrm{prs}}\frac{\text{∂}T_{\mathrm{rs}}}{\text{∂}t}＋\nabla · \boldsymbol{q}_{\mathrm{rs}}＝Q_{\mathrm{s}}\text{，} $$ (4) 式中，Cp代表定压比热容，ρ代表密度，T代表温度，q代表热通量，下标r代表储层，下标s代表储层中的固体颗粒。
忽略流体压力做功和流体黏性的耗散能量，流体热传导方程可以表示为（Aliyu，Archer，2021）
$$ {\rho }_{{\mathrm{w}}}{C}_{{\mathrm{pw}}}\frac{{\text{∂}} {T}_{{\mathrm{w}}}}{{\text{∂}} t}＋{\rho }_{{\mathrm{w}}}{{C}_{{\mathrm{pw}}}}{\boldsymbol{v}}_{{\mathrm{r}}}\nabla {T}_{{\mathrm{w}}}＋\nabla · {\boldsymbol{q}}_{{\mathrm{w}}}＝{Q}_{{\mathrm{w}}}\text{，} $$ (5) 式中，v代表从渗流场得到的速度场，下标w代表储层中的流体。
采用混合法则，将固体传热方程乘以固体体积分数θs，流体传热方程乘以孔隙度φr，二者相加（通常情况下有$ {\theta }_{\mathrm{s}}＋{\varphi }_{\mathrm{r}} {\text{≈}} 1 $）。利用局部热平衡理论，假设在某一点处流体温度和固体颗粒温度相等，即$ {T}_{{{\mathrm{w}}}}＝{T}_{{{\mathrm{s}}}}＝T $，则储层多孔介质的传热方程可以写为：
$$ { （ {\rho }_{{\mathrm{r}}}{C}_{{\mathrm{pr}}} ） }_{{\mathrm{eff}}}\frac{{\text{∂}} T}{{\text{∂}} t}＋{\rho }_{{\mathrm{w}}}{{C}_{{\mathrm{pw}}}}{\boldsymbol{v}}_{{\mathrm{r}}}\nabla T＋\nabla · {\boldsymbol{q}}_{{\mathrm{r}}}＝{Q}_{{\mathrm{T}}}\text{，} $$ (6) 式中：ρw为流体密度；Cpw为定压条件下流体比热容；（ρrCpr）eff为恒压下的有效容积热容，其表达式为
$$ { （ {\rho }_{{\mathrm{r}}}{C}_{{\mathrm{pr}}} ） }_{{\mathrm{eff}}}＝{\theta }_{{\mathrm{s}}}{\rho }_{{\mathrm{rs}}}{C}_{{\mathrm{prs}}}＋{\varphi }_{{\mathrm{r}}}{\rho }_{{\mathrm{w}}}{C}_{{\mathrm{pw}}}\text{，} $$ (7) 式中：φr为孔隙度；ρrs为固体颗粒密度；qr为热传导通量，其表达式为
$$ {\boldsymbol{q}}_{{{\mathrm{r}}}}\text＝\text－{{\textit{κ}} }_{{\mathrm{eff}}}\nabla T\text{，} $$ (8) 式中κeff为有效热传导系数。如果导热在固体和流体中并行发生，则有效导热系数为导热系数的加权算术平均值（Liu et al，2021）：
$$ {{\textit{κ}} }_{{\mathrm{eff}}}＝{\theta }_{{\mathrm{s}}}{{\textit{κ}} }_{{\mathrm{s}}}＋{\varphi }_{{\mathrm{r}}}{{\textit{κ}} }_{{\mathrm{w}}}\text{；} $$ (9) 如果导热是串联进行的，所有的热流都通过固体和流体，那么有效导热系数为导热系数的加权调和平均值：
$$ \frac{1}{{{\textit{κ}} }_{\mathrm{e}\mathrm{f}\mathrm{f}}}＝\frac{{\theta }_{{\mathrm{s}}}}{{{\textit{κ}} }_{{\mathrm{s}}}}＋\frac{{\varphi }_{{\mathrm{r}}}}{{{\textit{κ}} }_{{\mathrm{w}}}}{\text{.}} $$ (10) 1.3.2 断层热传导方程
断层热传导方向主要沿断层的切线方向，因此梯度算子可以简化为
$$ \nabla T＝{\nabla }_{{\mathrm{t}}}T\text{，} $$ (11) 式中下标t代表切向方向。
该假设通常适用于较薄的孔隙结构，与相邻区域相比是良好的热导体，热量可以沿切线方向快速对流。基于该假设和能量混合原理，热传导方程可以写为
$$ {d}_{{\mathrm{fs}}}{ （ {\rho }_{{\mathrm{fs}}}{C}_{{\mathrm{pfs}}} ） }_{{\mathrm{eff}}}\frac{{\text{∂}} T}{{\text{∂}} t}＋{d}_{{\mathrm{fs}}}{\rho }_{{\mathrm{w}}}{{C}_{{\mathrm{pw}}}}{\boldsymbol{v}}_{{\mathrm{f}}}{\nabla }_{{\mathrm{t}}}T＋{\nabla }_{{\mathrm{t}}} · {d}_{{\mathrm{fs}}}{\boldsymbol{q}}_{{\mathrm{f}}}＝{d}_{{\mathrm{fs}}}{Q}_{{\mathrm{T}}}＋{q}_{0}\text{，} $$ (12) 式中dfs为断层厚度。在传热模型内相邻域之间的边界处求解时，两个实体之间会进行一定的热量交换，根据$ {q}_{0}＝\boldsymbol{n} · \boldsymbol{q} $，式（12）可以写为（Liu et al，2021）
$$ {d}_{{\mathrm{fs}}}{ （ {\rho }_{{\mathrm{fs}}}{C}_{{\mathrm{pfs}}} ） }_{{\mathrm{eff}}}\frac{{\text{∂}} T}{{\text{∂}} t}＋{d}_{{\mathrm{fs}}}{\rho }_{{\mathrm{w}}}{{C}_{{\mathrm{pw}}}}{\boldsymbol{v}}_{{\mathrm{f}}}{\nabla }_{{\mathrm{t}}}T＋{\nabla }_{{\mathrm{t}}} · {d}_{{\mathrm{fs}}}{\boldsymbol{q}}_{{\mathrm{f}}}＝{d}_{{\mathrm{fs}}}{Q}_{{\mathrm{T}}}＋{\boldsymbol{n}} · {\boldsymbol{q}} $$ (13) 式中：（ρfsCpfs）eff为恒压下的有效容积热容，其表达式与储层热传导控制方程类似；κeff为有效热传导系数，其表达式与上述断层中的表达式相同。
1.4 渗透率控制方程
1.4.1 储层渗透率控制方程
孔隙度计算公式如下：
$$ \varphi_{\mathrm{r}} ＝{\varphi }_{0} ［ 1＋{c}_{{\mathrm{p}}} （ p-{p}_{{\mathrm{ref}}} ） ＋3{\beta }_{{\mathrm{T}}} （ T-{T}_{{\mathrm{ref}}} ） ］ \text{，} $$ (14) 式中，cp为孔隙压缩系数，pref和Tref分别为参考孔隙压力和参考温度。
利用Kozeny-Carman方程可以通过孔隙度变化计算储层的渗透率变化（Xu，Yu，2008）：
$$ {k}_{{\mathrm{r}}}＝{k}_{{\mathrm{r}}\,0}{\left(\frac{1-{\varphi }_{0}}{1-\varphi_{\mathrm{r}} }\right)}^{3}{\left(\frac{\varphi_{\mathrm{r}} }{{\varphi }_{0}}\right)}^{3}\text{，} $$ (15) 式中，kr和kr0分别为储层的渗透率和初始渗透率，φr和φ0分别为储层的孔隙度和初始孔隙度。
1.4.2 断层渗透率控制方程
断层中裂缝开闭与有效正应力的关系为（Faoro et al，2009）：
$$ b＝{b}_{{\mathrm{max}}}-\frac{A{\sigma { \text{′}}_{ {\mathrm{n}}}}}{1＋B{\sigma { \text{′}}_{ {\mathrm{n}}}}} \text{，} $$ (16) 式中，bmax为有效应力为0时的最大裂隙开度，${\sigma { \text{′}}_{ {\mathrm{n}}}}$为作用在裂隙壁上的有效正应力，A和B为裂隙岩石的力学性质参数。
1.5 基本属性变化
1.5.1 流体属性变化
流体属性的变化也是影响耦合作用的一个重要因素。在高温高压作用下，水的密度ρw 不再是常数，一般可以表示为温度的函数（Zhao et al，2015）：
$$ {\rho }_{{\mathrm{w}}}＝\left\{\begin{split}&6.31 {\text{×}} 1{0}^{－5} {T}^{3}－0.06 {T}^{2}＋18.92 T－950.7 \qquad \text{273.15} \text{ {\text{＜}} }{ T } \text{ {\text{＜}} } \text{293.15}\text{，} \\& 1 {\text{×}} 1{0}^{－5} {T}^{3}－0.01 {T}^{2}＋4.97 T＋432.26 \qquad \quad \, \, \text{293.15} \text{ {\text{＜}} }{ T }\text{ {\text{＜}} } \text{553.75} \text{，} \end{split}\right. $$ (17) 动态黏度系数为
$$ {\mu }_{{\mathrm{w}}}＝\left\{\begin{split} &1.38－0.02 {T}＋1.36 {\text{×}}1{0}^{－4} {T}^{2}－4.65 {\text{×}} 1{0}^{－7} {T}^{3}＋ 8.90 {\text{×}}\\& \qquad \,\,\,\,\,\, 1{0}^{－10} {T}^{4}－9.08 {\text{×}}1{0}^{－13} {T}^{5}＋3.85 {\text{×}} 1{0}^{－16} {T}^{6} \qquad \quad \,\,\, \text{273.15} \text{ {\text{＜}} }{ T }\text{ {\text{＜}} } \text{413.15}\text{，} \\& 0.004－2.11 {\text{×}} 1{0}^{－5} {T}＋3.86 {\text{×}} 1{0}^{ \text{－8}} {T}^{2}－2.401{0}^{\text{－11}} {T}^{3} \qquad \text{413.15} \text{ {\text{＜}} }{ T }\text{ {\text{＜}} } \text{553.75}\text{，} \end{split}\right. $$ (18) 比热容为
$$ {C}_{{\mathrm{pw}}}＝12\,\,010－80 {T}＋0.31 {T}^{2}－5.4 {\text{×}} 1{0}^{－4} {T}^{3}＋3.63 {\text{×}} 1{0}^{－7} {T}^{4} \text{，} $$ (19) 热传导系数为
$$ {{\textit{κ}} }_{{\mathrm{w}}}＝－0.87＋0.009 {T}－1.58 {\text{×}} 1{0}^{－5} {T}^{2}＋7.98 {\text{×}} 1{0}^{{－9}} {T}^{3}{\text{.}} $$ (20) 1.5.2 储层性质变化
在以往理论模型的能量守恒方程中，往往将一些动态模拟参数假设为不变量，例如储层岩石导热系数。设置为不变量的作法将显著扩大计算误差，例如花岗岩从0℃升高至500℃时，导热系数由2.95 W/（m·K）降低至1.59 W/（m·K），降低了约1.36 W/（m·K），降低幅度为46.1% （高平，2015），由此可见，温度对花岗岩导热系数的影响是十分明显的。导热系数与温度的关系可以表示为（高平，2015）：
$$ {{\textit{κ}} }_{\mathrm{s}}＝{\left(0.324\,1＋0.000\,\,449\,\,02T－\frac{29.600\,2}{T}\right)}^{－1}{\text{.}} $$ (21) 1.6 耦合关系
数值模拟技术能够考虑复杂地质条件和不同时空尺度下的力学变形、流体流动、能量传导过程，是分析和预测干热岩系统中多物理场耦合过程最有力的工具（王晓星等，2012；孙致学等，2016；唐志伟等，2016；Kumari，Ranjith，2019）。本文采用有限元模拟软件COMSOL Multiphysics （版本5.4）实现以上三个物理场的耦合过程，其中：采用达西模块中的储水模型计算储层注水、抽水以及储层内流体流动过程；采用多孔介质传热模块实现地热开发过程中的热传导过程；采用固体力学模块实现地热开发过程中的固体变形过程。
三个物理场的耦合过程如图2所示。注采引起的流体压力变化（耦合关系④）、温度变化（耦合关系①）引起储层变形；储层变形通过影响储层孔隙率、渗透率和断层裂隙开度控制着储层和断层流体流动过程（耦合关系⑤）；储层、断层中的流体作为传热介质影响着温度场的分布（耦合关系②），温度改变后，流体性质（密度、黏度）发生相应的改变（耦合关系③）。由于在干热岩处于高温状态时，固体变形引起的生热相对较小，因此本文忽略此影响。
![]()
图 2 温度场、渗流场和应力场耦合关系示意图Figure 2. The coupling relationship of temperature field，seepage field and stress field2. 算例分析
青海省共和盆地干热岩资源丰富（严维德，2015；张超等，2018；张森琦等，2020），其中GR1干热岩勘探井（以下简称“GR1井” ）的深度为
3705 m，井底温度为236℃，是中国目前钻获温度最高的干热岩。GR1井位于共和县恰卜恰镇东南5 km处的早中更新世湖积台地上，地处共和盆地二级构造单元切吉凹陷的东缘。该地区地质构造复杂，周边分布多条隐伏断层（张森琦等，2018）。中国地质调查局和青海省国土资源厅在青海恰卜恰地热区钻探的深部钻孔数据显示，GR1井的干热岩分布深度在2900 —4400 m，GR1井在2614 m深处温度超过150℃，在3404 m深处温度达到180.6℃，花岗岩基底的温度梯度估计为（38.9±0.39）℃/km （Zhang et al，2019a，2023，图3）。本文以GR1井为例，结合该地区地质背景及物性参数建立了热-流-固全耦合的三维有限元模型。![]()
Figure 3. The data collected from the GR1 well in the Qiabuqia HDR in Qinghai （Zhang et al，2019a，2023）2.1 几何模型
根据GR1井的实际数据构建含断层储层的三维模型，模型的长、宽、高分别为
3000 m，1500 m和1500 m，深度范围为2.5—4 km，储层位于地质模型中间高度位置处，储层的厚度为500 m （深度范围2.9—3.4 km），储层水平。在水平位置500 m和1500 m处分别设置注水井和抽水井，流体的注入与采出以垂直井方式进行。注水井和抽水井中间设置隐伏断层一处，其中心点距离注水井约300 m、距离抽水井约700 m，断层走向与y轴的夹角约为45°，其几何尺寸与各地质体关系即坐标系如图4所示。![]()
图 4 含断层储层三维地质几何模型示意图Figure 4. Schematic diagram of a 3-D geological geometric model of a faulted reservoir2.2 初始值和边界条件
在模型计算过程中，假设储层初始渗透率为kr0＝1×10−12 m2，其中储层垂直方向的渗透率约为水平方向的1/10，上、下盖层的渗透率约为2.8×10−18 m2。储层孔隙度、弹性模量和泊松比分别为φr＝0.06，E＝48.4 GPa，$ \nu $＝0.15，其它重要参数设置如表1所示。流体密度和动态黏度均与温度相关，具体可见式（17）—（21）。流-固耦合因子比奥系数设置为一个常数，αr＝0.6。
表 1 模型的材料参数Table 1. Material parameters of the model参数 弹性模量
/GPa泊松比 比奥
系数孔隙压缩
系数/Pa−1初始
渗透率/m2孔隙度 导热系数
/（W·m−1·K−1）比热容
/（J·kg−1·K−1）密度
/（kg·m−3）热膨胀
系数/K−1盖层 84.4 0.15 0.6 1×10−12 2.8×10−18 0.03 2.43 1100 2652 3×10−5 储层 84.4 0.15 0.6 1×10−8 1×10−12 0.06 方程（22） 1100 2652 3×10−5 断层 28.4 0.25 0.6 1×10−8 1×10−11 0.5 2.43 900 2500 3×10−5 基岩 48.4 0.15 0.6 1×10−12 2.8×10−18 0.03 2.43 1100 2652 3×10−5 对于固体力学边界条件，底端采用固定约束边界条件，上边界施加63.6 MPa荷载，模拟上覆荷载，四周采用应力进行约束，x方向的应力荷载约与垂直方向相等，y方向的应力荷载约为垂直方向的1.5倍。
流体（水）通过注水井注入，由抽水井抽出，注入流体温度为330 K。计算区域上表面温度为420 K，温度梯度为0.04 K/m （Zhang et al，2023）。
2.3 网格划分
采用自由四面体对三维模型进行网格剖分，本文进行了三组网格剖分下的计算，用以评价网格的敏感度。① 细化网格：网格共计3万
5936 个，平均网格质量为0.66，最大网格尺寸为0.24 km，最小网格尺寸为0.03 km；② 常规网格：网格共计2万1 514 个，平均网格质量为0.65，最大网格尺寸为0.3 km，最小网格尺寸为0.054 km；③ 粗化网格：网格共计1万1506 个，平均网格质量为0.62，最大网格尺寸为0.45 km，最小网格尺寸为0.084 km。在相同地质及注采条件下（注入压力40 MPa，生产压力20 MPa）对上述三种网格分别进行计算，所得结果如图5所示。可见常规网格与细化网格得出的结果非常接近，而粗化网格有较大的偏差。虽然细化网格计算结果更加精确，但相较于常规网格计算时间太长，综合考虑计算精度和计算速度，本文使用常规网格进行网格剖分。![]()
图 5 不同网格剖分下的计算结果（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 5. Calculation results for different meshing（a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency2.4 采热效果评价
为更加方便地评价EGS的采热效果，本文给出平均生产温度、输出热功率、采热率和累积采热量四个评价指标，其中：生产温度指生产井所产出热液的温度；输出热功率指生产井在某一时刻下产出液的热功率，代表EGS的发电能力；采热率指系统运行一段时间后已经开采出来的热能与储层中总的可开采热能的比值；累积采热率指系统运行一段时间后采出的总热量。
平均生产温度的计算方法为
$$ {T}_{{\mathrm{pro}}}＝\frac{{\displaystyle\int }_{V}T （ t ） {\mathrm{d}}V}{V}\text{，} $$ (22) 式中，T（t）为不同时刻生产井处所对应的温度大小，V为生产井的体积。此式适用于三维模型。实际操作中，采用COMSOL Multiphysics中内置的非局部耦合算子中的平均值算子aveop计算温度的平均值。
输出热功率可以由下式计算得到：
$$ {P}_{\mathrm{w}}＝q{\rho }_{\mathrm{w}}{C}_{\mathrm{p}\mathrm{w}} （ {T}_{\mathrm{p}\mathrm{r}\mathrm{o}}-{T}_{\mathrm{i}\mathrm{n}\mathrm{j}} ） \text{，} $$ (23) 式中：q为生产井的平均体积流量，单位为m3/s；ρw为采热流体的密度；Cpw为采热流体的定压比热容；Tinj为注入井中注入流体的温度。
采热率等于从储热岩层采出的总热量与储热岩石储存的总热量之比，可以表示为
$$ \eta ＝\frac{{\displaystyle\int \int \int }_{ V}\,{\rho }_{\mathrm{s}}{C}_{\mathrm{p}\mathrm{s}} ［ {T}_{\mathrm{i}\mathrm{n}\mathrm{i}}-T （ t ） ］ {\mathrm{d}}V}{{\displaystyle\int \int \int }_{ V}\,{\rho }_{\mathrm{s}}{C}_{\mathrm{p}\mathrm{s}} （ {T}_{\mathrm{i}\mathrm{n}\mathrm{i}}-{T}_{\mathrm{i}\mathrm{n}\mathrm{j}} ） {\mathrm{d}}V}\text{，} $$ (24) 式中${T}_{\mathrm{i}\mathrm{n}\mathrm{i}} $为储层初始温度。
累积采热量是对输出热功率进行时间t积分，即
$$ R＝{\displaystyle\int }_{0}^{t}q{\rho }_{\mathrm{w}}{C}_{\mathrm{p}\mathrm{w}} （ {T}_{\mathrm{p}\mathrm{r}\mathrm{o}}-{T}_{\mathrm{i}\mathrm{n}\mathrm{j}} ） {\mathrm{d}}t{\text{.}} $$ (25) 3. 模拟结果
断层的存在将地热储层划分为注水域和抽水域，本文重点考察在定压和定流量两种生产模式下，储层渗透率和断层渗透类型对地热开采的影响。
3.1 定压生产模式
定压生产指的是在地热生产过程中生产井和注水井的流体压力保持不变。假设储层初始孔隙压力为30 MPa，上下盖层的初始孔隙压力分别根据水的重力产生的压力给出。设置注水井的井口压力为40 MPa，生产井的压力为20 MPa （Zhang et al，2019b）。
3.1.1 储层渗透率
本文分析了储层渗透率（1×10−14—1×10−11 m2）对地热开采效果的影响（图6）。其中最小值1×10−14 m2代表了储层原始渗透率，而最大值1×10−11 m2代表了储层改造的最大极限。当储层未经改造时（渗透率为1×10−14 m2），在十年内生产温度下降很小，输出热功率和累积采热量也极小，采热率也只有10%。而当储层经过初步改造，渗透率达到1×10−13 m2时，储层温度下降幅度、输出热功率、累积采热量和采热率都有了显著的提升；当储层经过大规模改造，渗透率被大幅提升至1×10−12 m2时，储层温度大幅下降，输出热功率、累积采热量和采热率得到进一步提升。特别是当储层渗透率增大3个数量级，达到1×10−11 m2时，开采速度达到一个非常大的数值，一年内其采热率基本可以达到80%，但是这并不利于地热的持续性开采。
![]()
图 6 定压生产模式下储层渗透率k （1×10−14—1×10−11 m2）对地热开采效果的影响（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 6. The impact of reservoir permeability k （1×10−14−1×10−11 m2） on geothermal exploitation effect in constant pressure production mode（a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency3.1.2 断层渗透类型
储层中断层的存在也影响着地热的开采效率，在定压生产模式下，本文分析了断层渗透类型（连通型和阻隔型）对于地热产热效果的影响。按照其水力传导性能的差异，断层可以分为连通型断层和阻隔型断层。为此，定义一个新参数α表征储层渗透率kr与断层渗透率kf的比值α＝kr/kf，其中α＜10代表连通型断层，α＞10代表阻隔型断层。从图7中可以看出，断层的水力传导性能在定压模式下对地热开发有着较大影响：阻隔型断层阻碍了注水区域与抽水区域之间的连通通道，降低了其输出功率、累积采热量和采热率，严重影响了地热的开采效率；而连通型断层对水力传导的影响很小，不同连通型断层渗透率对应的输出功率、累积采热量和采热率变化并不明显，说明在这种条件下，储层本身的渗透性能决定了地热储层的采热效果。
![]()
图 7 定压生产模式下断层渗透类型对地热开采效果的影响α表征储层渗透率kr 与断层渗透率kf的比值，α＝kr/kf，α＜10代表连通型断层，α＞10代表阻隔型断层，下同（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 7. The impact of fault types on geothermal exploitation effect in constant pressure production modeα characterizes the ratio of reservoir permeability kr to fault permeability kf，expressed as α＝kr/kf. α＜10 indicates a connected fault type，whereas α＞10 signifies a sealing fault type，the same below （a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency3.2 定流量生产模式
定流量生产指的是在地热生产过程中注采井的注采速率保持不变，本文主要模拟计算了储层渗透率、断层渗透类型及注采速率对地热开发的影响。
3.2.1 储层渗透率
定流量生产模式下不同储层渗透率（1×10−14—1×10−11 m2）对地热开采效果的影响如图8所示。与定压生产模式不同，在定流量生产模式下，当储层渗透率较小时，生产井的平均温度反而降低得较快。而且四种评价指标中，只有最大渗透率（1×10−11 m2）的情况下有明显差异，其它渗透率的情况下得出的结果基本一致，说明在定流量生产模式下，一定范围内的储层渗透率变化对开采效率的影响较小。
![]()
图 8 定流量生产模式下储层渗透率k （1×10−14—1×10−11 m2）对地热开采效果的影响（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 8. The impact of reservoir permeability k （1×10−14−1×10−11 m2） on geothermal exploitation effect in constant flow production mode（a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency3.2.2 断层渗透类型
定流量生产模式下断层渗透类型（连通型和阻隔型）对于地热开采效果的影响如图9所示。与定压生产模式相反，在定流量生产模式下，阻隔型断层的输出热功率、累积采热量及采热率反而相对较高，但其相应数值与连通型断层相差甚小，说明断层的水力传导性能在定流量模式下对地热的开采影响很小。
![]()
图 9 定流量生产模式下断层渗透类型对地热开采效果的影响（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 9. The impact of fault types on geothermal exploitation effect in constant flow production mode（a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency3.2.3 注水速率
除地质因素外，工程因素对地热的开采也会产生不容忽视的影响，本文分析了注采速率对干热岩地热开发的影响，结果如图10所示，注采速率在0.3—30 kg·m−1·s−1范围内变化。从图10中可以看出，在流量为最小值0.3 kg·m−1·s−1时，储层温度在10年内基本保持不变，10年后才有一个较小的变化，而输出热功率、采热量和采热率均处于较低的水平；增大流量后（3 kg·m−1·s−1，10 kg·m−1·s−1），生产温度下降速度加快，输出热功率、采热量和采热率也显著提升；当流量达到最大30 kg·m−1·s−1时，生产温度迅速降低，且在短时间内达到很大的输出功率和采热率，不利于地热的持续性开采。
![]()
图 10 注采速率q对地热生产的影响（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 10. The impact of injection/production rate on geothermal production（a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency4. 讨论
4.1 不同耦合模式影响
上述结果均是在热-流-固全耦合模型下计算得到的，为了探究耦合关系在计算过程中对采热结果的影响，我们分别去除了耦合关系① （温度场对固体场的影响）和耦合关系③ （温度场对流体场的影响）来对比这两种耦合关系对计算结果的影响，结果如图11所示。从中可以明显地看出，与热-流-固全耦合模型对比，不考虑耦合关系①的曲线与其基本重合，但还是稍有偏差，而不考虑耦合关系③则有较大的偏差。若在耦合关系①和③均不考虑的情况下，其曲线与仅不考虑耦合关系③时稍有偏差。这说明耦合关系①在热-流-固耦合模型中对采热结果虽然影响不大，但在情况允许的条件下最好不要忽略；而耦合关系③对采热结果的影响很明显，这表明在模拟干热岩地热开采过程中考虑流体的密度、黏度等物性参数随温度的变化很有必要。
![]()
图 11 考虑不同耦合关系时各评价指标的计算结果（a） 生产温度；（b） 输出热功率；（c） 累积采热量；（d） 采热率Figure 11. Calculation results of various evaluation indicators considering different coupling relationships（a） Production temperature；（b） Thermal output power；（c） Cumulative geothermal extraction；（d） Geothermal extraction efficiency4.2 定压与定流量生产下取热效率差异性分析
两种生产模式下的生产温度、采热率及累积采热量等干热岩采热性能指标对不同储层渗透率的影响明显不同（图6，图8）。为了分析引起这种差异性的原因，给出了两种生产模式下储层内流体达西流速的变化，如图12所示。可见在定压生产模式下储层内流体达西流速随渗透率的变化很大，而定流生产模式下则差异较小。其中较大的达西流速会加速储层中流体温度的传递，这也可以很好地解释在定压生产模式下采热效率受储层渗透率的影响大，而定流模式下几乎不受影响的现象。结合不同渗透率生产温度变化（图6a，图8a），可以明显地看出，在定压模式下生产温度的变化随渗透率的变化很大，而定流模式则影响很小，这也可以很好地与图12中的达西流速变化相对应。因此，在定压生产模式下，储层渗透率是决定采热效率的一个非常重要的参数。
![]()
图 12 定压力生产模式（a）及定流量生产模式（b）下不同渗透率k储层内达西流速变化Figure 12. Changes of Darcy flow velocity in reservoirs with different permeabilities k in constant pressure production mode （a） and constant flow production mode （b）两种生产模式下采热性能指标对不同的断层渗透类型影响也明显不同（图7，图9）。根据不同断层渗透类型储层内达西流速的变化（图13）可知：在定压模式下，阻隔型断层严重降低了储层内流体的流速，因为其几乎隔断了断层两侧的流体交换，这也是阻隔型断层影响采热效率的主要原因；而在定流模式下阻隔性断层的存在几乎对达西流速没有影响，甚至比存在连通型断层的流速稍大，其主要原因是由于持续的注水会导致孔隙压力的升高，形成更大的压力差从而增大达西流速。
![]()
图 13 定压力生产模式（a）及定流量生产模式（b）下不同断层渗透类型储层内的达西流速变化Figure 13. Changes of Darcy flow velocity in reservoirs with different fault permeability types in constant pressure production mode （a） and constant flow production mode （b）根据上述分析，我们认为定压生产和定流生产模式下采热效率之所以不同，主要与储层的流体达西流速密切相关，因此要想提高采热效率，就要合理地提高流体在储层中的流速。如果采用定压生产模式，就必须在生产前期对储层进行充分的压裂改造，使得储层整体的渗透率（包括阻隔型隐伏断层）相对较高，从而可以保持较高的采热效率。除此之外，过高的渗透率会导致地热在短时间内开采速率过大从而造成大量的资源浪费，不利于地热的持续性开采。因此，提前计算出合适的储层改造的渗透率是十分必要的。而采用定流量生产模式，从采热效率的角度来说，前期储层渗透率的压裂改造对地热开采效率的影响并不明显。
4.3 定压与定流量生产模式下的安全性分析
由于过高的流体压力会导致套管变形及损害，也会存在诱发地震的危险，因此，本文分析了两种生产模式下井口的压力变化（图14）。定压生产模式下，井口压力可控且相对较低，最小的渗透率引起的压力也不会超过注水压力40 MPa；而定流量生产模式下，井口压力随储层的渗透率变化非常明显，当储层渗透率最小（1×10−14 m2）时，井口压力在前期就迅速增大到100 MPa以上，后期甚至接近140 MPa，这对于实际的开采是十分危险的。因此，从工业安全的生产角度来看，定流量生产模式也需要在前期对储层进行合理的改造，使其既具有较高的采热效率，又可避免发生安全事故。
![]()
图 14 定压生产模式（a）及定流量生产模式（b）下不同渗透率k对应的井口压力Figure 14. Wellhead pressure with different permeabilities k in constant pressure production mode （a） and constant flow production mode （b）4.4 共和地区地热开采综合建议
通过上述分析，在定压生产模式下采热效率受地质因素影响很大，过小的储层渗透率以及阻隔性断层都会显著降低采热效率，而过大的渗透率则会导致生产温度下降过快，不利于地热的可持续开采。因此在采用定压生产模式时改造合适的储层渗透率（1×10−13—1×10−12 m2）是非常必要的，并且应当避免注水井与生产井之间存在阻隔型断层。而在定流生产模式下虽然地质因素对采热效率的影响较小，但也要尽量改造储层使其具有较大渗透率（大于1×10−13 m2），并避免阻隔型断层的存在，这是由于定流量生产在较低渗透率时会造成很高的孔隙压力的累积，从而导致一些安全问题的发生。
5. 结论
本文建立了地热资源（干热岩）开发过程中的热-流-固耦合模型，并将上述大型偏微分方程组通过商业软件Comsol Multiphysics 进行求解。以青海共和干热岩GR1井为例，建立三维数值模型，研究了定流量和定压力两种典型生产模式下，储层渗透率、断层渗透类型对地热资源（干热岩）开发过程的影响。主要结论如下：
1） 建立了用于描述地热资源注采过程的热-流-固耦合模型，耦合关系主要包括：流体压力变化、温度变化引起储层变形；储层变形通过影响储层孔隙率和断层裂隙开度控制着储层和断层流体流动过程；储层、断层中的流体作为传热介质影响着温度场分布，流体性质（密度、黏度）随着温度变化发生相应改变；
2） 在定压生产模式下，地热的开采受地质因素（包括储层渗透率和断层渗透类型）的影响很大。其中储层渗透率越大，采热效率越高，过低的渗透率导致流体流速低，从而影响采热效率，但是过高的渗透率会导致采热速率过快，从而造成资源的浪费，不利于持续型开采。阻隔型断层也会显著降低采热效率，因此应避免注水井与生产井之间存在阻隔型断层的存在；
3） 在定流量生产模式下，虽然地热的开采效率受地质因素影响很小，注采速率是决定采热效率的主要因素。但是在相同注入速率下，储层低渗透率或阻隔型断层会产生更大的孔隙压力响应，存在套管变形及损害甚至诱发地震的危险；
4） 综合考虑施工安全因素和开采效率，对于青海共和地区的地热资源开采我们给出以下建议：在定压生产模式下，储层渗透率应改造至1×10−13—1×10−12 m2范围内；在定流生产模式下，储层渗透率需改造至1×10−13 m2以上，此外还应充分探明注水井与生产井之间的地质背景，避免存在阻隔型断层。
-
![]()
图 1 跨天津北断裂的地震勘探测线及钻孔位置示意图（陈宇坤等，2013）
Figure 1. Schematic diagram of seismic survey line and borehole location across the north fault of Tianjin （Chen et al，2013）
![]()
图 2 跨天津北某活动断层浅层地震勘探数据的时间序列
Figure 2. Time series of the shallow seismic exploration data across the north fault of Tianjin
![]()
图 3 图2所示时间序列在频率-波数域和频率-相速度域中的能量分布
（a） 未剔除体波的结果；（b） 剔除体波后的结果；（c） 提取的三个模式的频散点
Figure 3. The energy distribution in the frequency-wavenumber and frequency-phase velocity domains for the time series shown in Fig.2
（a） The result before removing body waves；（b） The result after removing body waves； （c） The extraction of surface-wave dispersion curves for three modes
![]()
图 4 基于多模式瑞雷面波频散曲线的反演（左）及其最终模型理论频散曲线与观测频散曲线的拟合（右）
（a） 仅用基阶模式；（b） 基阶、一阶模式联合；（c） 基阶、一阶和二阶模式联合
Figure 4. The inversion based on multi-mode Rayleigh-wave dispersion curves （left） and the fitting between the theoretical and observed dispersion curves （right）
（a） Using only the fundamental mode；（b） Using the fundamental and the first modes；（c） Using the fundamental，the first and the second modes
![]()
图 5 表2中两个模型的理论频散曲线
红点表示的频率范围与本文实际数据采集得到的观测频散曲线的频带范围相同，如图3c所示（a） 模型1：速度递增模型；（b） 模型2：含有低速层的模型
Figure 5. The theoretical dispersion curves calculated from two models 1 and 2 listed in Table 2
The red dots indicate the theoretical dispersion curves with the same frequency range as the observed dispersion curves shown in Fig. 3c （a） Model 1：A normal layered model with S-wave velocity increasing with depth； （b） Model 2：A model with a low-velocity layer
![]()
图 6 基于模型1多模式瑞雷面波频散曲线的反演（左）及其最终模型理论频散曲线与观测频散曲线的拟合（右）
（a） 仅用基阶模式；（b） 基阶和一阶模式联合；（c） 基阶、一阶和二阶模式联合
Figure 6. The inversion based on multi-mode Rayleigh-wave dispersion curves of model 1 （left） and the fitting between the theoretical and observed dispersion curves （right）
（a） Using only the fundamental mode；（b） Using the fundamental and the first modes；（c） Using the fundamental，the first and the second modes
![]()
![]()
图 8 基于模型2多模式瑞雷面波频散曲线的反演（左）及其最终模型理论频散曲线与观测频散曲线的拟合 （右）
（a） 仅用基阶模式；（b） 基阶和一阶模式联合
Figure 8. The inversion based on multi-mode Rayleigh-wave dispersion curves of model 2 （left） and the fitting between the theoretical and observed dispersion curves （right）
（a） Using only the fundamental mode；（b） Using the fundamental and the first modes
![]()
图 8 基于模型2多模式瑞雷面波频散曲线的反演（左）及其最终模型理论频散曲线与观测频散曲线的拟合 （右）
（c） 基阶、一阶和二阶模式联合；（d） 基阶、一阶和二阶模式联合，其中初始模型为图（c）的最终模型
Figure 8. The inversion based on multimode Rayleigh-wave dispersion curves of model 2 （left） and the fitting between the theoretical and observed dispersion curves （right）
（c） Using the fundamental，the first and the second modes；（d） Using the fundamental， the first and the second modes，and the initial model is the final model of Fig.c；
![]()
图 9 基于模型2多模式瑞雷面波频散曲线的反演 （左） 及其最终模型理论频散曲线与观测频散曲线的拟合 （右）
（a） 基阶、一阶和二阶模式联合（仅使用部分二阶频散数据）；（b） 基阶、一阶和二阶模式联合，其中初始模型为图8b的最终模型（仅使用部分二阶频散数据）
Figure 9. The inversion based on multi-mode Rayleigh-wave dispersion curves of model 2 （left） and the fitting between the theoretical and observed dispersion curves （right）
（a） Using the fundamental，the first higher and the second higher modes （only part of the second higher mode dispersion data are used）；（b） Using the fundamental，the first and the second modes，and the initial model is the final model of Fig.8b （only part of the second mode dispersion data are used）
![]()
![]()
图 11 模型1 （a）和模型2 （b）在不同频率和不同模式下瑞雷波相速度的S波速度深度敏感核
Figure 11. The sensitivity analysis of the multimode Rayleigh-wave phase velocities at different frequencies for models 1 （a） and 2 （b）
![]()
图 12 模型1 （a）和模型2 （b）中采用不同模式基于较窄（左）和较宽（右）频段的频散曲线的反演结果
Figure 12. Inversion results of dispersion curves based on narrower （left） and wider （right） frequency bands using different modes for models 1 （a） and 2 （b）
![]()
图 13 S波速度剖面 （a）及钻孔土层 （b）
Figure 13. S-wave velocity profile （a） and borehole soil layer （b）
表 1 初始模型参数和各模式频散曲线联合反演的S波速度
Table 1 The parameters of the initial model and the S-wave velocities inverted from multi-mode dispersion curves
层序 层厚/m 深度/m ρ/（g·cm−3） vP/（m·s−1） 初始vS/（m·s−1） 反演vS/（m·s−1） 1 1.0 1.0 1.60 300 152 135.4 2 1.0 2.0 1.73 400 156 145.3 3 1.0 3.0 1.77 470 163 200.8 4 1.0 4.0 1.79 550 165 221.5 5 1.0 5.0 1.79 570 170 195.4 6 1.0 6.0 1.79 590 180 169.6 7 1.0 7.0 1.79 600 185 155.3 8 1.0 8.0 1.79 700 190 102.4 9 2.0 10.0 1.79 800 200 152.8 10 3.0 13.0 1.79 900 210 254.6 11 4.0 17.0 1.80 1 100 220 259.3 半空间 $ \infty $ $ \infty $ 1.80 1 200 230 250.1 
下载: 导出CSV
表 2 模型1和模型2的模型参数
Table 2 The parameters of models 1 and 2
层序 层厚/m 深度/m 模型1 模型2 $\mathrm{\rho }/ （ \mathrm{g}\cdot {\mathrm{c}\mathrm{m} }^{-3} ） $ ${v}_{{\rm{P}}}/ （ {\rm{m}}\cdot {{\rm{s}}}^{-1} ） $ ${{v} }_{{{\rm{S}}} }/ （ \mathrm{m}\cdot {\mathrm{s} }^{-1} ） $ $\mathrm{\rho }/ （ \mathrm{g}\cdot {\mathrm{c}\mathrm{m} }^{-3} ） $ ${v}_{{\rm{P}}}/ （ {\rm{m}}\cdot {{\rm{s}}}^{-1} ） $ ${{v} }_{{{\rm{S}}} }/ （ \mathrm{m}\cdot {\mathrm{s} }^{-1} ） $ 1 1.0 1.0 1.581 7 249.1 126.5 1.605 3 266.2 135.4 2 1.0 2.0 1.654 6 399.7 156.0 1.629 7 373.2 145.3 3 1.0 3.0 1.728 3 556.2 192.9 1.724 1 580.1 200.8 4 1.0 4.0 1.735 2 659.8 196.8 1.776 2 738.1 221.5 5 1.0 5.0 1.739 5 646.0 199.2 1.732 7 654.9 195.4 6 1.0 6.0 1.739 5 664.8 199.3 1.683 6 554.8 169.6 7 1.0 7.0 1.739 7 653.0 199.4 1.653 0 505.1 155.3 8 1.0 8.0 1.742 1 739.8 200.8 1.508 2 376.8 102.4 9 2.0 10.0 1.754 8 832.9 208.2 1.647 2 611.3 152.8 10 3.0 13.0 1.767 5 926.3 216.0 1.824 7 1 089.9 254.6 11 4.0 17.0 1.780 3 1 120.5 224.1 1.831 0 1 297.7 259.3 半空间 $ \infty $ $ \infty $ 1.792 9 1 210.7 232.3 1.818 5 1 303.7 250.1
下载: 导出CSV
-
陈宇坤,李振海,邵永新,王志胜,高武平,杨绪连. 2008. 天津地区第四纪年代地层剖面研究[J]. 地震地质,30(2):383–399. doi: 10.3969/j.issn.0253-4967.2008.02.005 Chen Y K,Li Z H,Shao Y X,Wang Z S,Gao W P,Yang X L. 2008. Study on the Quaternary chronostratigraphic section in Tianjin area[J]. Seismology and Geology,30(2):383–399 (in Chinese).
陈宇坤, 赵国敏, 闫成国, 李振海, 杨菲, 杨绪连, 王志胜, 刘芳, 刘红艳, 任峰, 纪静, 张春丽, 杨港生, 李文栋. 2013. 天津市活动断层探测与地震危险性评价[M]. 北京: 科学出版社: 171–175. Chen Y K, Zhao G M, Yan C G, Li Z H, Yang F, Yang X L, Wang Z S, Liu F, Liu H Y, Ren F, Ji J, Zhang C L, Yang G S, Li W D. 2013. Detection on Active Faults and Assessment on Seismic Risk in Tianjin[M]. Beijing: Science Press: 171–175 （in Chinese）.
凡友华,刘雪峰,陈晓非,刘家琦. 2009. 瑞雷波勘探的f-k域能量最大模方法[J]. 哈尔滨工业大学学报,41(1):105–107. doi: 10.3321/j.issn:0367-6234.2009.01.024 Fan Y H,Liu X F,Chen X F,Liu J Q. 2009. Max-mode method of Rayleigh wave prospecting in frequency-wavenumber domain[J]. Journal of Harbin Institute of Technology,41(1):105–107 (in Chinese).
李建平. 2018. 浅层地震反射资料的多阶振型面波反演[J]. 地震学报,40(1):24–31. doi: 10.11939/jass.20170116 Li J P. 2018. Inversion of multi-mode surface waves extracted from the shallow seismic reflection data[J]. Acta Seismologica Sinica,40(1):24–31 (in Chinese).
鲁来玉,张碧星,汪承灏. 2006. 基于瑞利波高阶模式反演的实验研究[J]. 地球物理学报,49(4):1082–1091. doi: 10.3321/j.issn:0001-5733.2006.04.021 Lu L Y,Zhang B X,Wang C H. 2006. Experiment and inversion studies on Rayleigh wave considering higher modes[J]. Chinese Journal of Geophysics,49(4):1082–1091 (in Chinese).
鲁来玉,丁志峰,何正勤. 2011. 浅层有限频率面波成像中的3D灵敏度核分析[J]. 地球物理学报,54(1):55–66. doi: 10.3969/j.issn.0001-5733.2011.01.007 Lu L Y,Ding Z F,He Z Q. 2011. Analysis of 3D sensitivity kernels of the finite frequency surface wave tomography in shallow subsurface[J]. Chinese Journal of Geophysics,54(1):55–66 (in Chinese). doi: 10.1002/cjg2.1586
罗银河,夏江海,刘江平,刘庆生. 2008. 基阶与高阶瑞利波联合反演研究[J]. 地球物理学报,51(1):242–249. doi: 10.3321/j.issn:0001-5733.2008.01.030 Luo Y H,Xia J H,Liu J P,Liu Q S. 2008. Joint inversion of fundamental and higher mode Rayleigh waves[J]. Chinese Journal of Geophysics,51(1):242–249 (in Chinese).
王辉,丁志峰. 2006. 浅层地震勘探资料处理中的速度分析参数选取[J]. 地震地质,28(4):597–603. doi: 10.3969/j.issn.0253-4967.2006.04.007 Wang H,Ding Z F. 2006. Parameters selection for velocity analysis in shallow seismic data processing[J]. Seismology and Geology,28(4):597–603 (in Chinese).
夏江海. 2015. 高频面波方法[M]. 武汉: 中国地质大学出版社: 85–86. Xia J H. 2015. High-Frequency Surface-Wave Method[M]. Wuhan: China University of Geosciences Press: 85–86 （in Chinese）.
尹晓菲,胥鸿睿,郝晓菡,孙石达,王芃. 2020. 水平层状模型中多模式瑞雷波和拉夫波相速度频散曲线的灵敏度分析[J]. 石油地球物理勘探,55(1):136–146. Yin X F,Xu H R,Hao X H,Sun S D,Wang P. 2020. Sensitivity analysis of multi-mode Rayleigh and Love wave phase-velocity dispersion curves in horizontal layered models[J]. Oil Geophysical Prospecting,55(1):136–146 (in Chinese).
张碧星,鲁来玉. 2002. 层状半空间中导波的传播[J]. 声学学报,27(4):295–304. doi: 10.3321/j.issn:0371-0025.2002.04.002 Zhang B X,Lu L Y. 2002. Propagation of guided waves in stratified half space[J]. Acta Acustica,27(4):295–304 (in Chinese).
Alleyne D,Cawley P. 1991. A two-dimensional Fourier transform method for the measurement of propagating multimode signals[J]. J Acoust Soc Am,89(3):1159–1168. doi: 10.1121/1.400530
Forbriger T. 2003a. Inversion of shallow-seismic wavefields:I. Wavefield transformation[J]. Geophys J Int,153(3):719–734. doi: 10.1046/j.1365-246X.2003.01929.x
Forbriger T. 2003b. Inversion of shallow-seismic wavefields:II. Inferring subsurface properties from wavefield transforms[J]. Geophys J Int,153(3):735–752. doi: 10.1046/j.1365-246X.2003.01985.x
Forchap E A,Schmid G. 1998. Experimental determination of Rayleigh-wave mode velocities using the method of wave number analysis[J]. Soil Dyn Earthq Eng,17(3):177–183.
Foti S. 2000. Multistation Methods for Geotechnical Characterization Using Surface Waves[D]. Turin: Politecnico di Torino: 63–65, 109–127.
Gabriels P,Snieder R,Nolet G. 1987. In situ measurements of shear-wave velocity in sediments with higher-mode Rayleigh waves[J]. Geophys Prospect,35(2):187–196. doi: 10.1111/j.1365-2478.1987.tb00812.x
Gao L L,Xia J H,Pan Y D. 2014. Misidentification caused by leaky surface wave in high-frequency surface wave method[J]. Geophys J Int,199(3):1452–1462. doi: 10.1093/gji/ggu337
Herrmann R B. 2013. Computer programs in seismology:An evolving tool for instruction and research[J]. Seismological Research Letters,84:1081–1088.
Ikeda T,Matsuoka T,Tsuji T,Hayashi K. 2012. Multimode inversion with amplitude response of surface waves in the spatial autocorrelation method[J]. Geophys J Int,190(1):541–552. doi: 10.1111/j.1365-246X.2012.05496.x
Liang Q,Chen C,Zeng C,Luo Y H,Xu Y X. 2008. Inversion stability analysis of multimode Rayleigh-wave dispersion curves using low-velocity-layer models[J]. Near Surf Geophys,6(3):157–165. doi: 10.3997/1873-0604.2007040
Lu L Y,Zhang B X. 2004. Analysis of dispersion curves of Rayleigh waves in the frequency-wavenumber domain[J]. Can Geotech J,41(4):583–598. doi: 10.1139/t04-005
Lu L Y,Zhang B X. 2006. Inversion of Rayleigh waves using a genetic algorithm in the presence of a low-velocity layer[J]. Acoust Phys,52(6):701–712. doi: 10.1134/S106377100606011X
Lu L Y,Wang C H,Zhang B X. 2007. Inversion of multimode Rayleigh waves in the presence of a low-velocity layer:Numerical and laboratory study[J]. Geophys J Int,168(3):1235–1246. doi: 10.1111/j.1365-246X.2006.03258.x
Luo Y H,Xia J H,Liu J P,Liu Q S,Xu S F. 2007. Joint inversion of high-frequency surface waves with fundamental and higher modes[J]. J Appl Geophys,62(4):375–384.
Maraschini M,Ernst F,Foti S,Socco L V. 2010. A new misfit function for multimodal inversion of surface waves[J]. Geophysics,75(4):G31–G43. doi: 10.1190/1.3436539
Maraschini M,Foti S. 2010. A Monte Carlo multimodal inversion of surface waves[J]. Geophys J Int,182(3):1557–1566. doi: 10.1111/j.1365-246X.2010.04703.x
Matsuzawa H,Yoshizawa K. 2019. Array-based analysis of multimode surface waves:Application to phase speed measurements and modal waveform decomposition[J]. Geophys J Int,218(1):295–312. doi: 10.1093/gji/ggz153
Mcmechan G A,Yedlin M J. 1981. Analysis of dispersive waves by wave field transformation[J]. Geophysics,46(6):869–874. doi: 10.1190/1.1441225
Socco L V,Foti S,Boiero D. 2010. Surface-wave analysis for building near-surface velocity models:Established approaches and new perspectives[J]. Geophysics,75(5):75A83–75A102. doi: 10.1190/1.3479491
Tokimatsu K,Tamura S,Kojima H. 1992. Effects of multiple modes on Rayleigh wave dispersion characteristics[J]. J Geotech Eng,118(10):1529–1543. doi: 10.1061/(ASCE)0733-9410(1992)118:10(1529)
Wang J N,Wu G X,Chen X F. 2019. Frequency-Bessel transform method for effective imaging of higher-mode Rayleigh dispersion curves from ambient seismic noise data[J]. J Geophys Res:Solid Earth,124(4):3708–3723. doi: 10.1029/2018JB016595
Xia J H,Miller R D,Park C B. 1999. Estimation of near-surface shear-wave velocity by inversion of Rayleigh waves[J]. Geophysics,64(3):659–992. doi: 10.1190/1.1444578
Xia J H,Miller R D,Park C B,Tian G. 2003. Inversion of high frequency surface waves with fundamental and higher modes[J]. J Appl Geophys,52:45–57. doi: 10.1016/S0926-9851(02)00239-2
Xu Y X,Xia J H,Miller R D. 2006. Quantitative estimation of minimum offset for multichannel surface-wave survey with actively exciting source[J]. J Appl Geophys,59(2):117–125. doi: 10.1016/j.jappgeo.2005.08.002
Zhang S X,Chan L S. 2003. Possible effects of misidentified mode number on Rayleigh wave inversion[J]. J Appl Geophys,53(1):17–29. doi: 10.1016/S0926-9851(03)00014-4
-
期刊类型引用(1)
1. 刘朝阳. 基于光纤微动的泉水敏感区不良地质探查方法研究. 铁道建筑技术. 2024(10): 184-187+203 . 
百度学术
其他类型引用(0)
计量
- 文章访问数: 846
- HTML全文浏览量: 260
- PDF下载量: 185
- 被引次数: 1
