基于FLAC3D不同溶洞情况对风车稳定性分析

Analysis of the Impact of Different Karst Cave Conditions on Wind Turbine Stability Using FLAC3D

摘要: 本研究基于岩土工程数值分析软件FLAC3D，系统探讨了风车地基下方存在不同形态溶洞时，其对风车、地基整体稳定性及变形特性的影响。溶洞作为一种典型的地质灾害，严重威胁上部结构的安全。为量化分析其风险，本研究建立了考虑风车–土体–溶洞相互作用的三维数值模型，并采用莫尔–库仑本构模型来模拟土体的弹塑性行为。研究重点分析了不同溶洞的洞径、形状、空间位置等关键因素。系统获取了地基的位移场、应力场及位移随深度变化的关系与发展规律。模拟结果表明：溶洞距离基础1 m到2 m时变化较大，2 m之后溶洞位置降低对基础的影响就没有明显变化，及在某一个范围内对风车影响较大，过了这个范围影响就有所减小，同时还可以发现溶洞所处深度相同时溶洞体积越小对风车的Z轴影响越小。

Abstract: This study systematically investigates the effects of cavities with different morphologies beneath wind turbine foundations on the overall stability and deformation characteristics of both the turbine and foundation using the geotechnical numerical analysis software FLAC3D. As a typical geological hazard, cavities pose a serious threat to the safety of superstructures. To quantify these risks, a three-dimensional numerical model accounting for wind turbine-soil-karst interaction was established, employing the Mohr-Coulomb constitutive model to simulate soil elastoplastic behavior. Key factors such as cavity diameter, shape, and spatial location were critically analyzed. The study systematically obtained the foundation’s displacement field, stress field, and the relationship and development patterns of displacement with depth. Simulation results indicate: 1) Caves within 1~2 meters of the foundation exhibit significant influence. 2) Beyond 2 meters, reducing cave location has negligible effect on the foundation. 3) Caves within a specific range significantly impact the wind turbine, with influence diminishing beyond this range. 4) At equivalent depths, smaller cave volumes exert lesser Z-axis effects on the turbine.

1. 引言

岩溶地貌在全球范围内分布广泛，尤其在我国西南地区极为普遍。随着该地区基础设施建设的飞速发展，大量公路、铁路及建筑设施不可避免地需要穿越复杂的岩溶地质区域。岩溶地基中存在的溶洞、土洞等不良地质体，对其上部结构的安全构成了严重威胁，这些隐伏于地下的溶洞在自然或工程荷载作用下，可能发生顶板塌落、洞体失稳，进而引发地基不均匀沉降、地面塌陷等灾害，直接导致工程结构损坏，甚至造成生命财产的巨大损失[1] [2]。

溶洞对地基稳定性的影响是一个复杂的岩土力学问题，其稳定性主要受溶洞自身的几何形态(如洞径、截面形状)、空间分布位置以及与基础荷载的相对关系等多种因素的控制[3]。传统的评价方法，如理论公式计算和工程类比法，虽然简便，但难以精确考虑岩土体的非线性本构关系、溶洞的复杂空间形态及其与地基的相互作用，因而在准确性和普适性上存在局限。随着计算机技术的发展大量学者开始利用数值模拟方法来进行分析，Zhen Liu等，开发了一个结合黎曼映射的复杂函数模型，以推导任意洞穴几何形状的解析应力解，量化致密喀斯特洞穴群中的应力相互作用，解决单洞分析的局限性。Sung B K等[4]，应用有限元极限分析及修正的Hoek-Brown，对基础荷载下风化岩石的稳定性进行分析计算，得到水平和异质性扰动被认为是随深度变化的常数或线性变化因素。Li Z等[5]，提出了一种能够有效分析三维岩溶地基承载力的上界法。利用该方法研究了不同桩径、顶板厚度与桩径之比、溶洞宽度与桩径之比对岩溶地基稳定性的影响，结果表明当同时受到拉、压两种水平应力作用时，若最大水平应力超过岩体的抗拉强度，则顶板岩体可能发生破坏。Sheng M等[6]，采用物理模型试验和有限元分析相结合的方法，对岩溶桩基承载特性进行了研究，研究成果为类似岩溶地区嵌岩桩基础设计提供参考。

目前关于地下溶洞对地基稳定性分析主要以有限元方式展开。其中，FLAC3D作为一种基于显式有限差分法的三维数值分析程序，在处理大变形、材料非线性及岩土材料屈服破坏等方面具有独特优势，能够很好地模拟地基与溶洞系统的应力分布、塑性区发展及变形失稳的全过程。同时大量学者利用FLAC3D开展稳定性研究，学者Guorui W等[7]，采用了三维连续体快速拉格朗日分析(FLAC 3D)模拟地下煤矿开采，分析其位移变化等情况。Yongshuai S等[8]，采用FLAC3D三维有限差分软件模拟结合位移监测数据验证方法，对基坑支护过程中，应力、位移等参数进行监测。Haoran L等[9]，基于FLAC3D引入带张力截止的莫尔–库仑屈服准则并推导其有限差分格式，建立了表征失效的横观各向同性模型(称为AN-MC模型)。但目前利用FLAC3D分析不同情况溶洞对地面建构筑物的影响较少。

基于此，本研究旨在利用FLAC3D数值模拟软件，建立一系列地基下伏溶洞的精细化三维模型及地面以风车建筑为例开展数值研究。通过系统改变溶洞的洞径、空间位置及溶洞埋深等参数，深入分析地基的沉降规律、应力重分布特征对地面风车建筑的影响情况。本研究期望通过定量化的分析，揭示不同形态溶洞影响地基及其地面建筑物的稳定性的内在机理，以期为岩溶地区的工程选址、地基稳定性评价及灾害防治提供科学依据和理论支持。

2. 计算原理及模型建立

FLAC3D (三维快速拉格朗日分析)的核心计算原理是采用显式有限差分法和动态松弛技术，通过求解节点运动方程，逐步迭代追踪材料从受力到达到静力平衡的全过程。与基于刚度矩阵和整体平衡方程的传统有限元法不同，这种方法特别擅长模拟岩土材料的大变形、塑性流动和物理不稳定过程。其中显示有限差分法，是对空间导数进行有限差分近似，将控制方程在离散的网格上转化为代数方程；动态松弛技术，引入虚拟的质量和阻尼，将静态问题转化为动态问题，用显式时间步进法求解节点运动方程。

一个案例的计算流程为：由节点运动方程(已知力求运动)驱动网格变形 → 由变形求单元应变 → 由应变增量通过本构模型求应力增量 → 由应力求节点不平衡力 → 再回到运动方程，如此循环直至系统平衡[10]。

1) 求解的基本控制方程是牛顿第二定律，表示为节点(或网格点)的运动方程：

$\frac{\partial {\sigma }_{ij}}{\partial x_j}+\rho g_i=\rho \frac{\text{d}{v}_i}{\text{d}t}$ (1)

其中${\sigma }_{ij}$ 为柯西应力张量分量；$\rho$ 是材料密度；g_i是重力加速度分量；v_i是速度分量；d_v是质点加速度。

动量方程离散：

$F^{\left(l\right)}={\displaystyle \sum _{elements}{T}^{\left(element\right)}}\approx {\displaystyle \sum \left(-\frac{1}{3}\sigma \cdot n^{\left(face\right)}{S}^{\left(face\right)}\right)}$ (2)

T^(element)是单个单元对节点l的贡献；$\sigma$ 是单位应力张量；n^(face)是单元的法向量；S^(face)是对应面的面积.

2) 时间离散：动态松弛法(显式时间递推)，得到节点不平衡力F后，FLAC3D通过引入虚拟质量和阻尼，将静力问题转化为动力问题求解。节点的运动由以下显式中心差分格式更新：

速度及位移与坐标更新：

$v_i^{\left(t+\Delta t/2\right)}=v_i^{\left(t-\Delta t/2\right)}+\frac{F_i^{\left(t\right)}}{m_i}\Delta t$ (3)

$\begin{array}{l}{u}_i^{\left(t+\Delta t\right)}=u_i^{\left(t\right)}+v_i^{\left(t+\Delta t/2\right)}\Delta t\\ x_i^{\left(t+\Delta t\right)}=x_i^{\left(t\right)}+u_i^{\left(t+\Delta t\right)}\end{array}$ (4)

3) 单元应变与应力更新(本构关系)

根据当前循环中所有节点的速度场，计算每个四面体(或子单元)的应变率张量，应力增量则是由应变增量通过本构模型决定，其关系如下：

${\dot{ϵ}}_{ij}=\frac{1}{2}\left(\frac{\partial v_i}{\partial x_j}+\frac{\partial v_j}{\partial x_i}\right)$ (5)

$\Delta {\sigma }_{ij}=C_{ijkl}\left({\sigma }_{ij},\Delta {ϵ}_{kl},\text{Material Parameters}\right)\cdot \Delta {ϵ}_{kl}$ (6)

2.1. 模型建立

利用CAD软件绘制二维风车线框从小到大包含混凝土垫层、风车基础、填土、基岩4个部分，并分别绘制不同埋深、形状及大小的溶洞，其中溶洞最大直径不超过4 m，如下图所示。将二维线框导入Midas GTX软件，将各部分绘制成平面，接着利用旋转、扩展功能分别绘制各部分三维模型如下图1、图2所示。

Figure 1. Schematic of the model wireframe

图1. 模型线框示意

(a) 垫层 (b) 风车基础

(c) 填土 (d) 基岩

Figure 2. Schematic diagram of the 3D model

图2. 三维模型示意图

建立材料属性，划分垫层、风车基础、填土、基岩网格尺寸分别为0.5、0.8、1、3并赋予不同材料属性，得到网格划分图如下图3所示。

(a) 垫层 (b) 风车基础

(c) 填土 (d) 整体网格

Figure 3. Schematic diagram of grid partitioning

图3. 网格划分示意图

2.2. 参数设置及模拟方案

1) 模拟方案

溶洞位于基岩内部，分别研究球体溶洞相同大小不同深度(3种)、长方体溶洞相同大小不同水平位置(3种)、圆柱体溶洞相同深度不同大小及相同大小不同深度(3种)，共9种工况下溶洞对对风车基础的影响情况。

2) 模拟方案及参数设置

本研究基于FLAC3D平台建立数值模拟模型，其主要流程为：首先，对导入的复杂几何模型进行网格组重命名与重组，以清晰界定地基、基础等不同结构区域。随后，利用自编的get_para函数(FISH语言)为各组单元批量赋予相应的材料参数。通过model gravity 0，0，−10命令施加重力场，以模拟自然重力作用；随后是接触面的设置，通过interface property命令在基础与地基等关键接触界面创建接触单元，采用Mohr-Coulomb接触本构模型来模拟界面可能发生的摩擦滑移与脱开行为，并命令流定义其法向/切向刚度、摩擦角、黏聚力等关键参数。

考虑到风车在风力作用下的动荷载效应及风车自重，同样在计算模型中应设置的荷载包括自重与风荷载。风荷载简化为静态的竖向与水平荷载组合，通过zone face apply stress命令并结合range position限定作用区域进行施加。其中，竖向荷载为4000 kN (模拟结构自重与压重)，水平荷载为1500 kN (模拟风推力)。基于此基本荷载状态，后续将开展不同溶洞工况下的稳定性对比分析，得到不同溶洞影响下风车稳定性情况，以评估结构的安全性。

不同结构主要参数如下：垫层混凝土采用弹性模型参数使用经验数值，弹性模量5 × 10⁹ Pa，泊松比0.23，杨氏模量6 × 10⁸ Pa，重度24.5 KN。风车钢筋基础参数分别为，弹性模量5 × 10⁹ Pa，泊松比0.25，杨氏模量2 × 10⁸ Pa，重度25.5 KN。填土为摩尔库伦模型参数分别为，弹性模量2.5 × 10⁸ Pa，泊松比0.2，重度25 KN，粘聚力5 × 10³，内摩擦角30˚，抗拉强度3300 Pa。基岩为中分化灰岩，弹性模量3.1 × 10⁸ Pa，泊松比0.27，重度25 KN，粘聚力5 × 10³，内摩擦角30°，抗拉强度3000 Pa。

3. 不同工况下模拟分析

在前面对风车及其基础、地基进行建模，并设置各自结构间的接触关系，本章将对在风荷载作用下不同溶洞现状、埋深、大小对风车稳定性的影响情况及地基稳定性分析。

3.1. 无地下溶洞工况结果分析

模型X、Z方向的最大负位移为−0.0025238 m、−0.020827 m，最小负位移为−0.0005、−0.0005 mm。最大正位移为0.0024734 m、0 m，最小正位移为0.0005 m，0 m。X轴位移分布在基岩与垫层交接的两侧，填土表面受到与下部的位移量相反。其中正位移为受拉变化产生，负位移为受压作用产生。Z位移从风车基础以辐射状往外围减少，如图4所示。

(a) X位移云图 (b) Z轴位移

Figure 4. Displacement contour map for conditions without underground cavities

图4. 无地下溶洞工况位移云图

模型XX、ZZ最大拉应力为−6.7353 × 10⁵、−2.2848 × 10⁶，最小拉应力为−2 × 10⁵、−2 × 10⁵。最大压应力为4.5624 × 10⁵、1.9867 × 10⁵。其中正为拉应力压为负应力。XX压应力主要集中在风车基础与垫层接触的正中间范围内，拉应力集中在风车基础与风车主体连接的位置处。ZZ方向受到的拉应力主要分布在风车基础与填土连接的位置，如图5所示。

(a) XX应力 (b) ZZ应力

Figure 5. Stress contour map for conditions without cavities

图5. 无溶洞工况下应力云图

3.2. 不同地下溶洞形态工况结果分析

1) 球体(半径2 m)溶洞位于垫层底下1 m处

模型X方向的最大负位移为−0.0040473 m，最小负位移为−0.0005 m。最大正位移为0.0040098 m，最小正位移为0.0005 m。分布在基岩与垫层交接的两侧，填土表面受到与下部的位移量相反。其中正位移为受拉变化产生，负位移为受压作用产生。模型Z方向的最大负位移为−0.13563 m，最小负位移为−0.01 m。没有正位移。负位移从风车主体最高处往下位移逐渐减小。其中正位移为受拉变化产生，负位移为受压作用产生，如下图6所示。

(a) X轴位移云图 (b) X轴随水平方向位移变化曲线(地表)

(c) Z轴位移云图 (d) Z轴位移随深度变化曲线(地下)

Figure 6. Displacement contour map of the spherical cavern (radius 2 m) located 1 m below the cushion layer

图6. 球体溶洞(半径2 m)位于垫层底下1 m处位移云图

模型X面X轴方向的最大拉应力为−1.6296 × 10⁶，最小拉应力为−2.5 × 10⁵。最大压应力为1.2159 × 10⁶，最小压应力为1 × 10⁶。其中正为拉应力压为负应力。压应力主要集中在风车基础与垫层接触的正中间范围内，拉应力集中在风车基础与风车主体连接的位置处。模型Z面Z轴方向的最大拉应力为−5.3541 × 10⁶，最小拉应力为−5 × 10⁵。压应力为3.2795 × 10⁵。其中正为拉应力压为负应力。压应力主要集中在基岩填土表面及风车主体最高位置处，拉应力最大值位置在主体与填土接触的位置附近，如下图7所示。

2) 球体溶洞位于垫层底下2 m

模型X方向的最大负位移为−0.0038208 m，最小负位移为−0.0005 m。最大正位移为0.0037373 m，最小正位移为0.0005 m。分布在基岩与垫层交接的两侧填土表面受到与下部的位移量相反。与第1种情况相比第2种情况下的最大位移有所减小。模型Z方向的最大负位移为−0.13551 m，最小负位移为−0.01 m。没有正位移。负位移从风车主体最高处往下位移逐渐减小。与情况1对比情况2下Z轴的位移量有所减少，减少量为0.00012 m，2工况模拟结果如图8所示。

(a) XX应力 (b) XX应力曲线

(c) ZZ应力 (d) ZZ应力曲线图

Figure 7. Stress variation contour map

图7. 应力变化云图

(a) X轴位移 (b) X轴随水平方向位移变化曲线(地表)

(c) Z轴位移 (d) Z轴随垂直方向位移变化

Figure 8. Displacement variation contour map

图8. 位移变化云图

模型X面X轴方向的最大拉应力为−1.6597 × 10⁶，最小拉应力为−2.5 × 10⁵。最大压应力为1.0323 × 10⁶，最小压应力为2.5 × 10⁵。压应力主要集中在风车基础与垫层接触的正中间范围内，拉应力集中在风车基础与风车主体连接的位置处。与情况1的XX应力对比受到的压应力有所减少，拉应力略微有所增大。模型Z面Z轴方向的最大拉应力为−5.3512 × 10⁶，最小拉应力为−5 × 10⁵。压应力为3.3057 × 10⁵。其中正为拉应力压为负应力。压应力主要集中在基岩填土表面及风车主体最高位置处，拉应力最大值位置在主体与填土接触的位置附近。与第1种情况相比ZZ受到的拉应力值降低，受到的压应力略微升高，结果如图9所示。

(a) XX应力 (b) XX应力曲线

(c) ZZ应力 (d) ZZ应力曲线

Figure 9. Stress contour map

图9. 应力云图

3) 长方体(4 × 2 × 4)溶洞，垫层下1 m正中间位置

模型X方向的最大负位移为−0.0045469 m，最小负位移为−0.0005 m。最大正位移为0.0041256 m，最小正位移为0.0005 m。分布在基岩与垫层交接的两侧填土表面受到与下部的位移量相反。与工况1相比工况3的溶洞体积更大，从X位移云图上发现位移增大了0.0004996 m。可知在相同埋深条件下，溶洞的体积越大风车受到位移的变化越明显。模型Z方向的最大负位移为−0.13574 m，最小负位移为−0.01 m。没有正位移。负位移从风车主体最高处往下位移逐渐减小。与工况1相比工况3的Z轴位移增加了0.00011 m，与X位移变化相同都是相对工况1变化明显，结果如下图10所示。

模型X面X轴方向的最大拉应力为−1.6210 × 10⁶，最小拉应力为−2.5 × 10⁵。最大压应力为1.0153 × 10⁶，最小压应力为2.5 × 10⁵。压应力主要集中在风车基础与垫层接触的正中间范围内，与工况1分布情况相似。对比工况1及工况3的XX受到的应力情况云图发现，相比于工况1工况3受到的应力都有所减小。模型Z面Z轴方向的最大拉应力为−5.3778 × 10⁶，最小拉应力为−5 × 10⁵。压应力为3.3497 × 10⁵。分布位置与工况1相似。与工况1相比ZZ受到的拉应力及压应力都有所增加，如图11所示。

(a) X位移 (b) X轴随水平方向位移变化曲线(地表)

(c) Z轴位移 (d) Z轴位移随深度变化曲线(地下)

Figure 10. Displacement contour map

图10. 位移云图

(a) XX应力 (b) XX应力曲线

(c) ZZ方向应力 (d) ZZ方向应力曲线

Figure 11. Stress contour map

图11. 应力云图

4) 长方体溶洞，垫层下1 m中间向X轴正方向4 m处

模型X方向的最大负位移为−0.0037995 m，最小负位移为−0.0005 m。最大正位移为0.0042100 m，最小正位移为0.0005 m。分布在基岩与垫层交接的两侧填土表面受到与下部的位移量相反。与工况3相比工况4的溶洞体积相同埋深相同，但工况4的溶洞位于X轴正方向4 m处(处于压位移区域内)。通过两个工况的X轴位移云图可以发现X轴的负位移减低，正位移升高。模型Z方向的最大负位移为−0.13557 m，最小负位移为−0.01 m。没有正位移。负位移从风车主体最高处往下位移逐渐减小。与工况3相比工况4的Z轴位移减小了0.00017 m，如下图12所示。

(a) X轴位移 (b) X轴随水平方向位移变化曲线(地表)

(c) Z轴位移 (d) Z轴位移随深度变化曲线(地下)

Figure 12. Displacement contour map

图12. 位移云图

模型X面X轴方向的最大拉应力为−1.6955 × 10⁶，最小拉应力为−2.5 × 10⁵。最大压应力为1.1632 × 10⁶，最小压应力为2.5 × 10⁵。压应力主要集中在风车基础与垫层接触的正中间范围内，与工况3分布情况相似。对比工况3及工况4的XX受到的应力情况云图发现，相比于工况3工况4受到的应力都有所增加。模型Z面Z轴方向的最大拉应力为−5.3630 × 10⁶，最小拉应力为−5 × 10⁵。压应力为3.4379 × 10⁵。分布位置与工况1相似。与工况3相比ZZ受到的拉应力有所降低，但压应力有所增加，如图13所示。

5) 圆柱体(半径1 m，高2 m)溶洞，位于垫层正下1 m处

模型X方向的最大负位移为−0.0040525 m，最小负位移为−0.0005 m。最大正位移为0.0039168 m，最小正位移为0.0005 m。分布在基岩与垫层交接的两侧填土表面受到与下部的位移量相反。模型Z方向的最大负位移为−0.13531 m，最小负位移为−0.01 m。没有正位移。负位移从风车主体最高处往下位移逐渐减小，如图14所示。

(a) XX应力 (b) XX应力曲线

(c) ZZ应力 (d) ZZ应力曲线

Figure 13. Stress diagram

图13. 应力示意图

(a) X轴位移 (b) X位移曲线

(c) Z轴位移 (d) Z轴位移曲线

Figure 14. Displacement contour map

图14. 位移云图

模型X面X轴方向的最大拉应力为−1.6302 × 10⁶，最小拉应力为−2.5 × 10⁵。最大压应力为1.1600 × 10⁶，最小压应力为2.5 × 10⁵。压应力主要集中在风车基础与垫层接触的正中间范围内。模型Z面Z轴方向的最大拉应力为−5.3565 × 10⁶，最小拉应力为−5 × 10⁵。压应力为3.1988 × 10⁵。分布位置与工况1相似，结果如下图15所示。

(a) XX应力 (b) XX应力曲线

(c) ZZ应力 (d) ZZ应力曲线

Figure 15. Stress contour map

图15. 应力云图

6) 圆柱体(半径1 m，高1 m)溶洞，位于垫层正下1 m处

模型X方向的最大负位移为−0.0038111 m，最小负位移为−0.0005 m。最大正位移为0.0039432 m，最小正位移为0.0005 m。分布在基岩与垫层交接的两侧填土表面受到与下部的位移量相反，与工况5相似。对比工况5发现工况6的X方向负位移减小0.0002314 m，正位移增加0.0000264 m。模型Z方向的最大负位移为−0.13517 m，最小负位移为−0.01 m。没有正位移。负位移从风车主体最高处往下位移逐渐减小。与工况5相比工况6的Z轴位移有所降低，减少量为0.00014 m，如下图16所示。

(a) X轴位移 (b) X位移曲线

(c) Z轴位移 (d) Z轴位移曲线

Figure 16. Displacement contour map

图16. 位移云图

模型X面X轴方向的最大拉应力为−1.6190 × 10⁶，最小拉应力为−2.5 × 10⁵。最大压应力为1.1844 × 10⁶，最小压应力为2.5 × 10⁵。压应力主要集中在风车基础与垫层接触的正中间范围内。对比工况5工况6的拉应力降低，压应力增加。模型Z面Z轴方向的最大拉应力为−5.3557 × 10⁶，最小拉应力为−5 × 10⁵。压应力为3.2413 × 10⁵。分布位置与工况5相似。与工况5相比拉应力减小，压应力增大见图17所示。

(a) XX应力 (b) XX应力曲线

(c) ZZ应力 (d) ZZ应力曲线

Figure 17. Stress contour map

图17. 应力云图

4. 结论及建议

工况3、4分别为相同体积、相同埋深、不同水平位置。对比工况3、4发现工况4的X轴的负位移减低，正位移升高，Z轴的位移减小，XX受到的应力都增加，ZZ受到的拉应力有所降低，但压应力有所增加。XX受到的拉应力减小，压应力增加，ZZ受到的应力都有所减小。以上结果可以得出水平力在X面X轴负方向时，溶洞位于左侧的影响比在右侧的影响要明显；工况5，6为相同体积不同埋深，进行对比分析位移及应力变化情况。可以发现，相同深度下溶洞体积减小X方向负位移减小，正位移增加，Z轴位移有所降低，XX拉应力降低，压应力增加，ZZ拉应力减小，压应力增大。相同体积不同埋深下X轴负位移，Z轴位移升高，XX拉应力增加，压应力减小，ZZ应力减小。

通过FLAC3D分析得到6种不同工况下风车的位移及应力云图及曲线图，通过控制相同体积溶洞不同深度埋深可以发现，溶洞距离基础1 m到2 m时变化较大，2 m之后溶洞位置降低对基础的影响就没有明显变化，及在某一个范围内对风车影响较大，过了这个范围影响有所减小。通过控制相同溶洞体积及埋深，改变溶洞的水平位置，可以发现溶洞在受拉区域内对风车的影响要大于溶洞在受压区域内。通过改变不同溶洞体积，控制相同埋深，可以发现溶洞体积越小对风车的Z轴影响越小。同时通过6种不同工况的位移及应力云图发现，位移增加相对应的应力就减少，反之就增加。对比没有溶洞工况，有溶洞时对风车影响十分明显，所以遇见溶洞时要进行填充处理。

致 谢

感谢华北水利水电大学地球科学与工程学院提供试验及计算平台。