南海油气资源丰富，已成为我国深海油气资源开发的主战场^[1]，但南海内孤立波活动频繁。1990年,在流花油田就曾发生过因内孤立波而导致缆绳拉断、船体碰撞，甚至发生拉断和挤破漂浮软管的事故，严重影响了钻井平台的安全作业^[2]。合理确定各种海洋环境条件下立柱式钻井平台的载荷是深海浮式结构物设计过程中的关键环节。直立圆柱型结构是立柱式钻井平台等各种深海结构物的主体结构形式，目前海洋工程界一般采用Morsion公式这种工程简化方法来计算直立圆柱体的内孤立波载荷，如：Y.Cheng等^[3]和S.Cai等^{[4, 5]}将Morison公式与KdV理论结合，J.Xie等^[6]将Morison公式与MCC(Miyata-Choi-Camassa)理论^[7]结合，研究了内孤立波作用下直立圆柱体的载荷特性问题。尤云祥等^{[8, 9]}将Morison公式与eKdV(extended KdV)理论^[10]结合，研究了内孤立波作用下张力腿和半潜式平台的载荷与动力响应问题，而宋志军等^[11]则将Morison公式与KdV(Korteweg-de Vries)理论^[10]结合，研究了内孤立波作用下Spar平台的载荷与动力响应问题。然而，在上述文献中，Morison公式中的惯性力和拖曳力系数都是参照表面波的方法选取的，缺乏理论和试验依据。为此，黄文昊等^[12]以系列试验为依据给出了圆柱型结构惯性力和拖曳力系数的选取方法，但上述这些简化方法是否适用于现实大尺度海洋环境条件下的载荷计算仍然不确定。此外，南海实际环境条件中，内孤立波的相速度可达1.6～2.0 m/s，最大振幅可达100 m多^[13]，这是无法在内孤立波试验水槽中实现的，受制于试验条件及分析工具的限制，目前对于内孤立波载荷试验的尺度效应尚未进行深入研究，模型试验得到的内孤立波载荷经换算后，在实际大尺度条件下直接使用是否合理，目前仍然并不清楚。计算流体力学(computational fluid dynamics，简称CFD)方法为深入认识内孤立波尺度效应提供了一条有效的途径，采用CFD方法可以在不同尺度条件下直接获得内孤立波与浮式结构物相互作用过程中浮式结构物的水动力特性及载荷构成，从而为验证大尺度条件下工程简化计算方法的适用性以及分析模型试验的尺度效应提供了可能。

鉴于此，笔者采用Navier-Stokes方程，建立了振幅及其波形可控的内孤立波CFD数值模拟方法。在此基础上，对不同尺度条件下内孤立波与立柱式钻井平台的相互作用进行数值模拟，以验证目前立柱式钻井平台内孤立波载荷简化方法在大尺度条件下的适用性，同时对模型试验的尺度效应进行研究和分析。

考虑两层流体中内孤立波与立柱式钻井平台的相互作用，设两层流体均为不可压流体，上层流体的深度与密度分别为h₁和ρ₁，下层流体的深度与密度分别为h₂和ρ₂。内孤立波为平面前进波，界面位移为ζ，沿Ox轴正方向传播，立柱式钻井平台直径为D，吃水深度为d。建立直角坐标系Oxyz，其中Oxy平面位于流体静止时两层流体的界面上，Oz轴与平台中心轴重合且以竖直向上为正。

笔者采用CFD方法模拟内孤立波诱导流场，进而获得内孤立波对立柱式钻井平台的载荷。在此基础上，研究立柱式钻井平台内孤立波载荷的尺度效应。为此，采用求解Navier-Stokes方程的方法模拟内孤立波诱导流场，其中流场控制方程为：

式中：u_i为速度矢量，m/s；p为动压力矢量，Pa；t为时间，s；f_i为重力矢量，N；ν为运动黏性系数，m²/s；ρ为流体密度，kg/m³，当ζ < z < h₁时，ρ=ρ₁；当－h₂ < z < ζ时，ρ=ρ₂。

将平台直立置于两层流体中，平台壁面取无滑移不可穿透边界条件，计算域顶部及底部要求满足如下壁面条件：

流场计算的控制区域包括内孤立波生成传播区和消波区2个区域(见图1)。采用速度入口方法生成内孤立波，当在生成传播区中形成稳定的内孤立波后，对所生成内孤立波的传播特性进行监测分析，并对立柱式钻井平台的内孤立波载荷进行计算。

图1 数值水槽示意 Fig.1 The sketch of the numerical water flume

图选项

设内孤立波振幅为a，相速度为c，则其诱导上下层流体中的层深度平均水平速度分别为^[2]：

式中：为上层流体层深度的平均水平速度，m/s；为下层流体层深度的平均水平速度，m/s。

在内孤立波生成与传播过程中，两层流体的界面会发生变化，采用VOF(volume of fluid)方法追踪两层流体界面的变化^[14]。利用海绵层消波方法对水槽尾部的内孤立波进行消波处理，该方法在消波区通过在动量方程右端添加源项－μ(x)u_i的方式实现，其中μ(x)为海绵层衰减系数^[15]。

立柱式钻井平台的内孤立波水平力F_x及垂向力F_z由摩擦力和压差力2部分组成，即：

式中：S为平台湿表面积，m²；两式中右端的第一项为平台侧表面和底部的载荷摩擦力，N；第二项为平台侧表面的载荷压差力，N；(n_x,n_y,n_z)为平台表面法向矢量，方向指向平台外部。

下面给出立柱式钻井平台力矩的计算方法，设为立柱式钻井平台中心轴上某点到力矩转动中心的距离，F_x为在该点附近作用在单位长度立柱式钻井平台侧表面上的水平力，则在内孤立波作用下立柱式钻井平台的力矩可表示为：

式中：为力矩转动中心至平台底部距离，m；M_y为立柱式钻井平台的力矩，N·m。顺时针时，M_y为正；逆时针时，M_y为负。

依据文献^[16]确定3类内孤立波(KdV、eKdV和MCC)理论的适用性条件，从而计算得到入口边界处上下层流体层深度的平均水平速度，具体计算方法如下：

记h=h₁+h₂为总水深，ε=|a|/h为非线性参数，μ=(h/λ)²为色散参数，其中λ为内孤立波特征宽度。根据内孤立波的振幅等条件，利用有限深水KdV、eKdV 和MCC的理论解，分别计算相应内孤立波的非线性参数ε和色散参数μ。如果μ＜μ₀且ε≤μ，则选择KdV的理论解作为数值模拟的入口边界条件；如果μ＜μ₀且μ＜ε≤，则选择eKdV的理论解作为数值模拟的入口边界条件；如果μ＜μ₀且ε＞或μ≥μ₀，则选择MCC的理论解作为数值模拟的入口边界条件。

文献^[12]利用大型密度分层水槽，对内孤立波作用下立柱式钻井平台的载荷特性进行了系列试验。模型试验参数：平台直径D为0.15 m，吃水d为0.535 m；水槽长L为30 m，宽度B为0.6 m；水深h为1 m；平台中心轴距离速度入口端9 m；上下层流体密度分别为998和1 025 kg/m³；上下层流体深度比分别为1∶9，2∶8和3∶7。结合深海Spar平台的实际尺寸，相对于模型尺寸，笔者选取150∶1，200∶1和300∶1等3种尺度比的数值模型研究平台内孤立波载荷的尺度效应。表1为上述3种尺度比条件下的数值模型尺寸。

力矩转动中心至平台底部距离为=d+0.15k。由于作用于平台的力矩可以由水平力沿垂向积分得到，其载荷特性与水平力载荷特性相类似，所以主要讨论水平力和垂向力的尺度效应。定义F_x=F_x/(ρ₁gD²d)和=F_z/(ρ₁gD³)为无因次内孤立波水平力和无因次内孤立波垂向力，同时定义特征时间=t/k，其中k为相对于试验模型的尺度比。为分析黏性影响的需要，定义雷诺数Re=U_maxD/ν，其中U_max为内孤立波诱导的最大水平速度。

尺度效应是模型试验必须考虑的关键问题之一。根据文献^[12]的试验结果，采用数值水槽方法，对立柱式钻井平台内孤立波载荷的模型试验尺度效应进行分析。为此，在数值模拟中，选择Case A和Case B两种工况进行比较，其中Case A工况的数值水槽及柱体尺寸与文献^[12]中的模型试验一致，水槽上下层流体深度比h₁∶h₂=3∶7，内孤立波无因次振幅a_d/h=0.101；而Case B工况选择基于试验模型尺度比为200∶1的数值模型，数值水槽及柱体尺寸见表1，而上下层流体深度比及无因次振幅均与Case A工况相同。

图2给出了Case A和Case B两种工况下的内孤立波无因次水平力时历数值模拟结果，并与文献^[12]的试验结果进行了比较。由图2可知：Case A工况下，模型尺度下的立柱式钻井平台内孤立波无因次水平力数值模拟结果与试验结果吻合良好，两者相对误差6.25%，且均在特征时间=61达到峰值；而在Case B工况下，无因次水平力时历曲线同样在特征时间=61到达峰值，与Case A工况相符。比较Case A和Case B两种工况下的无因次水平力时历曲线，可以发现Case B工况下无因次水平力幅值明显小于Case A工况下模型试验结果，与模型水平力无因次幅值的相对误差达到31.25%，这说明在求解大尺度浮式结构物内孤立波水平力载荷幅值时，尺度效应非常显著。如果忽略尺度效应，直接对模型试验得到的载荷直接进行尺度换算，将高估平台内孤立波水平力。

图2 Case A和Case B两种工况下内孤立波无因次水平力时历曲线 Fig.2 Time-variant curve for dimensional horizontal forces of internal solitary waves in Case A and Case B

图选项

图3给出了Case A和Case B两种工况下的内孤立波无因次垂向力时历数值模拟结果，并与文献^[12]中的试验结果进行了比较。由图3可知：Case A工况下，模型尺度下的立柱式钻井平台内孤立波垂向力数值模拟结果与试验结果吻合良好，两者相对误差4.17%，且均在特征时间=61.7达到峰值；而在Case B工况下，无因次垂向力时历曲线也在=61.7到达峰值，与Case A工况相符。比较Case A和Case B两种工况下的无因次垂向力时历曲线，可以发现Case B工况下无因次垂向力幅值与Case A工况模型试验结果相差较小，相对误差为5.56%。这说明对于垂向力，尺度效应较小，在求解大尺度浮式结构物内孤立波垂向力载荷幅值时，可以忽略尺度效应，直接进行尺度换算。

图3 Case A和Case B两种工况下内孤立波无因次垂向力时历曲线 Fig.3 Time-variant curve for dimensional vertical forces of internal solitary waves in Case A and Case B

图选项

由上述分析可知，就尺度效应而言，平台内孤立波水平力和垂向力两种载荷存在很大不同，水平力的尺度效应非常显著，不可忽略，而垂向力的尺度效应较小，可以忽略。因此，需要从载荷构成的角度分析出现这种差异的原因，由式(5)和式(6)可知，水平力及垂向力均由摩擦力和压差力两部分构成。根据黏性的影响，可以将立柱式钻井平台内孤立波的水平和垂向压差力进一步划分为波浪压差力和黏性压差力两部分：波浪压差力被认为与内孤立波诱导的水质点波动有关，可以采用无黏条件下Euler模拟方法得到；黏性压差力归因于流体黏性效应导致的压差力增量，通过N-S模拟所得压差力减去波浪压差力得到。这样水平力或垂向力均可以划分为摩擦力、波浪压差力和黏性压差力3部分。

图4为Case A和Case B两种工况下的构成内孤立波无因次水平力3种载荷的时历曲线。由图4可以看出：无因次摩擦力f_x与无因次波浪压差力F_x^pw或无因次黏性压差力F_x^pv相比均非常小，可以忽略；无论是Case A工况还是Case B工况，波浪压差力和黏性压差力的量级相当，这意味着黏性对立柱式钻井平台的内孤立波水平力影响是非常显著的；比较Case A和Case B两种工况的无因次水平波浪压差力幅值，两者的相对误差为0.88%，这说明水平波浪压差力的尺度效应是非常微弱的；比较两种工况的无因次水平黏性压差力幅值，两者的相对误差为26.60%，这说明水平黏性压差力的尺度效应是非常显著的。综合上述分析，造成立柱式钻井平台水平力尺度效应显著的原因是黏性压差力，即黏性效应的影响。另外，通过对比Case A和Case B两种工况可以看出，雷诺数Re从Case A工况的9.7×10³增大到Case B工况的2.7×10⁷时，Case B工况的无因次水平黏性压差力幅值明显小于Case A工况的，同时黏性压差力在水平力中的比重也随之明显降低。

图4 Case A和Case B两种工况下内孤立波水平力3种无因次载荷成分时历曲线 Fig.4 Time-variant curve for three dimensional loads of horizontal forces of internal solitary waves in Case A and Case B

图选项

图5为Case A和Case B两种工况的构成内孤立波无因次垂向力3种载荷的时历曲线。由图5可以看出：无因次摩擦力f_z与无因次波浪压差力F_z^pw或无因次黏性压差力F_z^pv相比均非常小，可以忽略；无论是Case A工况还是Case B工况，黏性压差力与波浪压差力相比也非常小，这意味着黏性效应较弱，立柱式钻井平台的内孤立波垂向力主要由波浪压差力构成；比较Case A和Case B两种工况的无因次波浪压差力幅值，两者的相对误差为5.63%，这说明波浪压差力的尺度效应是非常微弱的，这说明由于构成垂向力主要成分的波浪压差力的尺度效应较小，因而立柱式钻井平台垂向力的尺度效应也较小。

图5 Case A和Case B两种工况下内孤立波垂向力3种无因次载荷成分时历曲线 Fig.5 Time-variant curve for three dimensional loads of vertical forces of internal solitary waves in Case A and Case B

图选项

由此可见，模型试验水平力、垂向力尺度效应差异显著的根本原因在于黏性效应的不同，受黏性影响较大的水平力的尺度效应显著，而受黏性影响较弱的垂向力的尺度效应则可以忽略。

通过上述分析可知，采用模型试验方法将内孤立波水平力试验结果直接换算到实际尺度是不可行的。为此，进一步研究了一种工程简化方法的尺度效应。其具体做法是分别采用Morison公式和傅汝德-克雷洛夫力公式求解平台内孤立波的水平力和垂向力。这种做法在模型试验尺度雷诺数低(Re≤1×10⁴)的情况下被证明是适用的^[12]，但在实际海洋环境中雷诺数高(Re>1×10⁷)的情况下的适用性尚不清楚。因此，笔者采用CFD模拟方法对大尺度条件下该简化方法的适用性进行分析。

记U₁和W₁分别为当ζ < z < h₁时内孤立波诱导的水质点瞬时水平速度和垂向速度，记U₂和W₂分别为当－h₂ < z < ζ时内孤立波诱导的瞬时水平速度和垂向速度，则可得^[17]：

其中

定义U和W如下：当ζ < z < h₁时，U=U₁，W=W₁；当－h₂ < z < ζ时，U=U₂，W=W₂。那么，利用式(7)和式(8)，可得单位长度立柱式钻井平台侧表面内孤立波载荷的Morison公式为^[18]：

式中：C_m为惯性力系数；C_d为拖曳力系数；V=(U,0,W)为水质点速度，m/s；为水质点加速度；m/s²；V_n为垂直于平台侧表面的速度分量，m/s。

文献^[12]根据系列试验结果，建议Morison公式(式(10))中的惯性力系数和拖曳力系数采用下式计算：

另一方面，由伯努利方程可得内孤立波诱导的流体动压力为^[12]：

根据式(13)对平台底部进行压力积分，可得作用其底部的傅汝德-克雷洛夫力为：

式(12)中，C_m和C_d是根据内孤立波水槽试验数据回归得到的，由于内孤立波水槽主尺度等因素，模型试验中的雷诺数Re最大只能做到10⁴量级，而实际情况的雷诺数最大可达10⁸量级。因此，在利用上述工程简化方法估算立柱式钻井平台内孤立波载荷时，同样存在尺度效应问题。

图6给出了尺度比为200∶1，上下层流体深度比h₁∶h₂=3∶7，无因次振幅为a_d/h=0.101工况下，利用简化方法计算所得立柱式钻井平台内孤立波无因次水平力、垂向力与CFD数值模拟结果。由图6可知，简化方法计算得到的无因次水平力与CFD数值模拟结果趋势大体一致，两者幅值相对误差为10.9%；而无因次垂向力与CFD数值模拟结果相吻合，两者幅值相对误差为4.4%。这说明使用简化方法在大尺度、高雷诺数(该工况Re=2.7×10⁷)条件下，求解立柱式钻井平台内孤立波载荷仍然是可行的。

图6 简化方法计算和CFD数值模拟的内孤立波无因次水平力、垂向力幅值 Fig.6 The amplitudes of dimensional horizontal and vertical forces of internal solitary waves calculated with simplified methods and simulated with CFD numerical values

图选项

1) 在大尺度、高雷诺数条件下，采用Morison和傅汝德-克雷洛夫公式计算其内孤立波水平力和垂向力是可行的。

2) 立柱式钻井平台内孤立波水平力及垂向力的尺度效应因受流体黏性影响的不同而出现差异，受黏性影响较大的水平力的尺度效应显著，而受黏性影响较弱的垂向力的尺度效应则可以忽略。

3) 在实际尺度情况下，由于计算量等因素，CFD数值模拟方法仍受到限制。而在实际尺度雷诺数下通过模型试验所得到Morison公式的拖曳力系数和惯性力系数计算公式仍适用，因此可以利用由模型试验获得的拖曳力系数和惯性力系数的计算公式计算出这2个系数，再利用工程简化方法估算实际尺度下的内孤立波载荷。