三角形布置圓柱體群繞流特性與流動機理研究

涂佳黃，呂海宇，胡 剛，2，王貫宇

（1.湘潭大學土木工程學院，湖南湘潭 411105；2.重慶市城市建設發展有限公司，重慶 400014）

0 引 言

在實際工程中，相較于單、雙柱體形式結構，多柱體結構群形式運用更廣泛，如海洋平臺立柱、輸電線路和橋橋墩柱等。柱體結構群的繞流問題會更復雜，既要考慮結構與流體的作用，還需考慮結構群之間復雜的互擾效應。以往的研究大多集中于串、并列雙柱體結構形式[1-3]。Sumner[1]和Alam[2]研究串列雙圓柱結構繞流時，發現了三種流型：拓展體流、再附流和共同脫渦流，三種流型中，下游柱對上游柱的影響在間距比L/D<5.0 時較顯著[3]。在對18≤Re≤100 范圍內的并列雙圓柱結構研究中，Liu 等[4]將其近尾流區分為穩定與不穩定模式。

部分學者對低雷諾數下多圓柱群繞流問題進行了研究[5-10]。Bao 等[5]采用數值計算方法研究了間距比和來流角度對等邊三角布置圓柱群繞流特性的影響；Zheng等[6]模擬三圓柱群的繞流，發現小間距比工況下，上游布置并列圓柱比下游布置并列圓柱時的流場流動模式更加穩定；同時，Kang 和Su[7]也分析了兩種布置方式下三圓柱群渦脫落頻率和流體力系數的變化，也發現上游布置并列圓柱時的尾流呈周期性變化流場更加穩定；楊梟梟等[8]詳細研究間距比L/D（1.0≤L/D≤6.0）、Re=100 下的等邊三角形布置圓柱的尾流模式，發現6 種模式，按間距比分單體（1.0≤L/D≤1.4）、偏斜（1.4

另外一些學者用實驗的方法去研究高雷諾數工況下不同參數對三圓柱群繞流特性的影響[11-15]。Pouryoussefifi[11]測量了Re=6.08×104和1.26×104工況下等邊三角布置三圓柱表面壓力分布；Igarashi 和Suzuki[12]實驗研究了串列三圓柱的繞流特性，將上游圓柱剪切層對下游圓柱影響方式劃為六種；Lam和Cheung[13]在2.1×103≤Re≤3.5×103范圍內對三圓柱群繞流開展了實驗研究，在1.27

綜上所述，柱體繞流的研究使用了不同的數值方法或者實驗方法，相關的參數涉及了不同的Re、間距比、布置方式等，探究了其尾流模式、流體力特性、渦脫落特征，取得了諸多成果。本文運用浸沒邊界的多松弛格子Boltzmann（IB-MRT-LBM）數值方法研究Re=100 下等邊三角形布置三圓柱的流場特性，該方法可以滿足計算高效率、高準確度和較好適用性的要求，同時數值模擬方法可以彌補實驗方法對低雷諾數的模擬能力的缺陷，便于更加精細刻畫流場的動態過程；主要分析不同來流角度和間距比兩個參數（其中來流角度α選取0°、15°、30°、45°和60°，間距比Kd范圍為1.2≤Kd≤6.0）對等邊三角形布置三圓柱體群尾流模式和流體力系數的影響，闡明流場的轉變過程。

1 數值方法和驗證

1.1 數值方法

1.1.1 多松弛格子Boltzmann方程

格子Boltzmann 方程在速度、時間、空間上表現出特殊的離散形式，可將微觀粒子的速度ξ離散成有限維的速度空間{e0,e1,...,eN}，N表示粒子編號，分布函數f也可進行離散{f0,f1,...,fN}，其中fα=fα(r,eα,t),α= 0,1,...,N。在粒子碰撞過程中，可在不同方向上的分布函數中采用不同的松弛時間進行求解，即多松弛的格子Boltzmann模型（MRT-LBM）：

式中：m是mk組成的向量，mk=φk f,k=1,2,...,b，b等于DrQs模型中速度模型的向量個數s；φk是粒子速度的線性無關多項式；M是φk確定的變換矩陣；meq=Mfeq,m=Mf；feq為平衡態的密度分布函數；S為松弛矩陣，引入D2Q9 模型后，松弛矩陣轉換為S=MΩM-1=dig(0,se,sε,0,sq,0,sq,sυ,sυ)-1，其中Ω是碰撞矩陣，由無量綱松弛時間τ=τ0/δt組成。

MRT-LBM算法中，不同于實際有量綱的物理系統參數，密度ρ、u和p均與分布函數有關，即

式中，cs是格子聲速，定義為cs=c/ 3，在D2Q9模型中，有cs=1/ 3。運動粘度的定義也與實際的物理系統定義不同：

式中，δt=1/Nt，Nt為時間步數，sυ為松弛矩陣中關于運動粘度的一項。

實際流場中存在固體物時，流固邊界的流體質點受到邊界力作用使得流場密度分布函數發生變化，郭照立等[16]將額外作用力施加到流場中，結合不含外力項的MRT-LBM方程，得到新的方程：

式中，δt為時間步長，F?為邊界作用力分布函數，I是單位矩陣。

1.1.2 浸沒邊界-格子Boltzmann方法對外力項的速度進行修正，可以改善傳統IB-LBM 方法在流線穿透方面的不足，從而嚴格滿足無滑移邊界的要求，使流場演化過程中分布函數精度得到提高[17]。本文將浸沒邊界與格子Boltzmann 方法相結合，將節點額外作用力加入格子Boltzmann方程，即公式（4）。按照文獻[18]的計算思路，使用二階拉格朗日插值法對修正速度進行求解，提高作用力f的精度。如圖1所示，曲線代表實際的固體邊界，固體采用拉格朗日坐標系描述，流體則采用歐拉坐標系描述，P為拉格朗日節點，和歐拉節點A之間的距離為Δ1，和歐拉節點B之間的距離為Δ2，通過A1、B1、P點對A和B速度進行修正。拉格朗日節點P的速度（up，vp）為已知，由三點拉格朗日插值可計算出節點A修正后的速度：

圖1 流體和固體邊界節點示意圖Fig.1 Schematic diagram of fluid and solid boundary nodes

同理，歐拉節點B的速度可以由點A1、B1、P進行修正。豎向歐拉節點C和D同樣需要D1、C1、Q三點插值修正，如果同一個節點被修正了兩次，如圖重合的節點B和節點C，那么最終修正速度取兩次修正的平均值：

通過對修正速度的求解，可以求出作用力f：

式中，uˉ為修正速度，u*為中間過渡速度。

本文的IB-MRT-LBM 計算流程如下：

（a）設定流場的初始值ρ和u，利用公式（1）計算出初始的平衡態密度分布函數meq；

（b）求出流場域各格點的宏觀量ρ和u*，進行邊界處理；

（c）由修正速度法，求出邊界周圍歐拉點的修正速度uˉ，根據公式（7）計算流體與固體邊界發生碰撞后的作用力f；

（d）通過作用力f求出F?，代入公式（4）求得下一時刻新的分布函數；

（e）重復上述步驟（b）與（c），直到結果收斂。

1.2 參數有效性驗證

為了驗證算法參數設置的準確性，選取單圓柱模型做繞流模擬，對模型的邊界處理格式與網格分辨率進行驗證。無量綱參數如阻力系數CD、升力系數CL與斯托羅哈數St分別定義為

式中，Fx與Fy分別是浸入邊界所受到的阻力與升力，fs是圓柱漩渦脫落頻率。選取Re=100，無量綱的入口速度u0=0.05，邊界條件以及流場域大小設置如圖2所示。

圖2 單圓柱繞流計算域以及邊界條件示意圖Fig.2 Schematic diagram of computational domain of flow around a single cylinder and the boundary conditions

對單位D內包含的格子數分別為20、30 和40 三種工況的單圓柱繞流問題進行模擬，并將其流體力系數和St的計算結果與文獻[19-20]進行對比，結果見表1。計算結果之間的誤差主要是由于邊界格式的精度、網格數量與計算域大小等方面共同影響產生的。另外，在其它參數設置不變情況下，流場網格的分辨率會對數值結果產生一些影響，本節對不同格子數工況下的單圓柱繞流問題進行了計算分析。格子數20工況的CD-mean、CL-rms和St值相對于格子數40工況的誤差均明顯大于格子數為30工況的數值，但是格子數為40工況計算所需時間會明顯增加，因而單位直徑D包含的格子數選擇為30，可以較好地解決計算精確性與效率之間的問題。

表1 網格分辨率驗證Tab.1 Grid resolution verification

1.3 問題描述

1.3.1 計算模型和邊界條件

均勻來流作用下等邊三角形布置三圓柱繞流的計算模型以及邊界條件如圖3 所示，三個圓柱分別記為圓柱1（C1）、圓柱2（C2）、圓柱3（C3），其圓心圍成的等邊三角形的重心為坐標原點O，三個圓柱的直徑均為D，任意兩圓柱圓心之間的距離與圓柱直徑的比為Kd=L/D，OO2是C2的圓心O2與流場域坐標原點O的連線，來流角度為來流在入口處與流場域上下邊界形成的角度，大小等于直線OO2與X軸的夾角α。計算區域選取為[-25D，40D]×[-20D，20D]。模型參數選擇如下：Re=100，1.5≤Kd≤6.0，α=0°、15°、30°、45°、60°。為保證計算的精確性，設置無量綱速度ux=u0=0.05，uy=0，選擇非平衡外推邊界格式為流場區域邊界，流場出口的密度ρ=1，對應壓力p==1/3。

圖3 等邊三角形布置三圓柱繞流計算域以及邊界條件示意圖Fig.3 Schematic diagram of the computational domain and boundary conditions for the flow around three cylinders in equilateral triangle arrangement

1.3.2 算例驗證

為驗證IB-MRT-LBM 算法在模擬等邊三角布置三圓柱繞流問題的正確性，本小節選擇Re=100，1.5≤Kd≤6.0，α=0°工況進行對比分析。圖4展示了等邊三角形布置三圓柱繞流計算的圓柱阻力系數平均值CD-mean與升力系數脈動值CL-rms隨間距比的變化規律，與文獻結果進行了比較，發現本文的計算結果與文獻的結果基本吻合，存在一些誤差可能是計算方法不同導致的。通過以上驗證，說明本文算法在求解等邊三角形布置三圓柱群繞流問題方面具有正確性與適應性。

圖4 等邊三角布置三圓柱群流體力系數計算結果對比Fig.4 Comparison of calculation results of fluid force coefficients of three cylinders in equilateral triangle arrangement

2 計算結果與分析

2.1 尾流模式

本節對不同來流角度與間距比下的三圓柱體群的尾流模式進行了總結，見表2，主要模式有七種：單漩渦模式（SV）、單漩渦體向雙漩渦過渡模式（SV-DV）、不規則的雙漩渦模式（IDV）、規則的雙漩渦模式（RDV）、雙漩渦向三漩渦過渡模式（DV-TV）、不規則的三漩渦模式（ITV）和規則的三漩渦模式（RTV）。

表2 各工況下三圓柱群尾流模式分布Tab.2 Wake pattern distribution of three-cylinder group at different cases

小間距比工況下，間隙流發展較弱，間隙內部不存在剪切層，三圓柱群上下側的剪切層以穩定頻率脫落漩渦，因而各來流角度工況下尾流模態均為單漩渦模式（SV）。

當Kd=1.5 時，間隙流流量增多，如圖5（a）所示，除了自由流側有明顯漩渦脫落之外，C2 后方負號剪切層和C3自由流側的剪切層融合后脫落漩渦，呈現單漩渦向雙漩渦過渡模式（SV-DV）。間距比增大到2.0，當α=45°時，C1尾流和C2尾流發生了融合，如圖5（b）所示，C1和C2間隙側沒有渦脫落，間隙流被夾在中間，流動方向不斷變化，C2 和C3 間隙側的尾流被壓縮，呈粗細不均的形態，間隙側剪切層有漩渦脫落但周期性不明顯，為不規則雙漩渦模式（IDV）。然而，其它來流角度下，流場會轉變為規則的雙漩渦模式（RDV），如圖5（c）所示，C1 尾流在間隙內部伸長，影響C2、C3 渦脫落，由于沒有充足的動量提供，C2和C3間隙側尾流產生的漩渦被拉伸成條狀，但整體的下游渦街正負交替、規則有序。

圖5 不同尾流模式下流場瞬時渦量圖Fig.5 Instantaneous vorticity diagram of the flow field for different wake patterns

中間距比工況下，C1 自由流側的剪切層有漩渦脫落趨勢，但間隙側剪切層還未脫落漩渦，如圖5（d）所示，層間互擾作用以及下游漩渦變形融合顯著，C3 下游負號漩渦向右下方移動，在對流過程中正負融合后，正號漩渦消失，整個流場下游的漩渦形態極不規則，尾流模式是雙漩渦向三漩渦脫落的過渡模式（DV-TV）。另外一種尾流模式如圖5（e）所示，C1 尾流受間隙流影響很小，C1 下游可以找到負號的漩渦，表明C1尾流完全具備了脫落漩渦的能力，自由流側的剪切層在脫渦過程中和C2剪切層融合，但此時C2 下游漩渦變形程度顯著，使得圓柱受力不穩定性增加，為不規則三漩渦脫落模式（ITV）。

大間距比工況下，如圖5（f）所示，各圓柱剪切層均發生擺動，C1 和C2 剪切層產生的漩渦發生融合，與C3 后方的漩渦形成兩列正負交替的獨立渦街，彼此不受干擾，流場脫渦模式具有周期性的特點，因此尾流模式為規則的三漩渦脫落模式（RTV）。

2.2 流體力系數

2.2.1 時均特性

來流角度的變化會使得流體流動的方向出現偏差，尾流區漩渦也將受此影響，從而柱體間存在不同的繞流特性，然而不同角度下升阻力時均系數（CD-mean，CL-mean）隨間距比的變化規律類似，如圖6和圖7 所示。當α≤30°時，C1 的CD-mean在Kd=1.5 處取得最大值后逐漸減小后穩定，如圖6（a）所示。當α=45°和60°時，C1的CD-mean分別在Kd=2.5與1.5處取得最小值。之后CD-mean均隨著Kd增大而趨于平穩。需注意的是，在α≤30°工況下C1的CD-mean均會大于α=45°和60°工況，由于當α=15°時C1和C3在橫流向接近并列排布，柱體距離較大，C1周邊流速會增大，如圖8（a）所示，于是C1后方形成了低壓區，流線在后方較遠處才開始匯合，如圖8（b）所示，但是該特征在ɑ=45°時并不明顯。同時，與α=45°工況相比，α=15°時C1的剪切層范圍會縮小，表明C1后方流體只在較小范圍內做旋轉運動，其余動量較大的流體則徑直向下游運動，見圖8（c），導致α=15°工況C1承受了更大的阻力。

圖6 不同來流角度下三圓柱體群阻力系數平均值隨間距比的變化Fig.6 Variation of time-averaged value of the drag coefficient for three-cylinder group with the spacing ratio in different incoming flow angles

圖7 不同來流角度下三圓柱體群升力系數平均值隨間距比的變化Fig.7 Variation of time-averaged value of the lift coefficient for three-cylinder group with the spacing ratio in different incoming flow angles

圖8 ɑ=15o，Kd=4.0工況下三圓柱體群時均流場特性圖Fig.8 Time-averaged flow field characteristic diagram of three-cylinder group for ɑ=15o,Kd=4.0

在小間距比（Kd=1.2、1.5）下，當來流角度較大（α=45°、60°）時，C2前表面偏上方的區域存在一個壓力駐點，CD-mean遠大于其余工況下的結果，如圖6（b）所示。α=30°時，C1 和C2 串列，C2 所受阻力會減小。當2.0≤Kd≤3.0與α=30°時，C2的CD-mean會逐漸減小。隨間距比和來流角度改變，CD-mean會出現突變：當Kd=3.0 時，由0.4（α=30°）增加至0.82（α=45°），另外當α=30°時，由0.4（Kd=3.0）增大至0.82（Kd=3.5），此處的機理是流線的分布改變導致C2 的受力方式也發生變化。進一步對其時均流線和壓力圖進行分析，如圖9（a）所示，C1 和C2 的自由流側流線重合，部分的流線偏轉產生了小范圍的回流，向下流動過程中與具有較高動量的間隙流相遇，受流體間剪切力作用回流流體方向再次偏轉，呈“Z”狀的流線，最后流線匯合，C1 對C2 屏蔽作用最明顯，C2 所受的時均阻力最小。然而，ɑ=45°工況下，C2 后方形成“S”形流線，前方自由流側形成一個明顯的壓力駐點，如圖9（b）所示；另外Kd=3.5 時的間隙流量增大，C2后方流體回流后再附著顯著減弱，如圖9（c）所示，造成了突變的形成。

圖9 不同工況下三圓柱體群時均流線和壓力圖Fig.9 Time-averaged streamline and pressure diagram of three-cylinder group at different cases

隨著來流角度的增加，C3 逐漸從上游圓柱轉變為下游圓柱，受力發生變化。當Kd≤2.0 時，隨著來流角度的增大，C3 的CD-mean會減小，如圖6（c）所示。當Kd≥2.5、α≥30°時CD-mean均有明顯的上升幅度。特別是α=45°工況下，當Kd=3.5 和4.0 時，C3 的CD-mean達到最大值，由于C2 和C3 在橫流向近似并列排布，圓柱體間距離大，來流在經過C2和C3間隙時獲得了充足的動量，因而間隙流速度較大，在C3周圍形成低壓區，如圖10（a）～（b）所示。另外，與ɑ=15°工況相比，C3 的剪切層長度會減小，如圖10（c）所示，因此C3受較大的時均阻力作用。Kd≥4.5后變化逐漸平穩，α=0°和α=15°所受的CD-mean較小。

圖10 α=45°，Kd=4.0工況下三圓柱體群時均流場特性圖Fig.10 Time-averaged flow field characteristic diagram of three-cylinder group for α=45°,Kd=4.0

對于升力系數而言，C1 的升力系數平均值變化趨勢受來流角度的影響較小。隨著間距比的增大，其CL-mean值會逐漸減小并趨近于0。在Kd=1.2 處不同來流角度工況CL-mean值相差較大，如圖7（a）所示。然而，當α=60°時，由于三圓柱體群關于來流中心線對稱，此時CL-mean維持在0。

在小間距比工況下（Kd≤1.5），隨著來流角度的增大，C2 的CL-mean會增加，如圖7（b）所示，自由流側的渦脫落對C2 的影響更加明顯。然而，在小來流角度工況下，C2 會被自由流側的剪切層包裹，C2 的下游渦街遠離其表面使得CL-mean減小。在α=15°和α=30°工況下，由于自由流側脫落的漩渦強度增大，因而Kd=1.5～2.0 范圍內CL-mean會增大。在α=45°工況下，CL-mean出現負值，這可能是由于C2 下游渦街較C3更加雜亂，導致渦強度損失大，從而C2受C3下游渦街影響顯著。但Kd=4.5時，由于C3下游渦街也受干擾，致使其對C2的影響減弱，使C2的CL-mean減小。

當Kd=1.2時，間隙流沒有充分發展，外側自由流對C3的影響占主要地位，分離剪切層的作用使得C3 承受向下的吸力，CL-mean的絕對值較大，如圖7（c）所示。同時，隨著來流角度的減小，CL-mean會增大。隨著間距比的增大，C3間隙側也出現分離的剪切層，CL-mean有所減小。然而，當α=60°時，在2.0≤Kd≤4.0范圍內，C3又受到向上的流體力作用，使得CL-mean為正值。

2.2.2 脈動特性

本節分析三圓柱體所受流體力系數脈動特性的變化情況。圖11 和圖12 展示了不同來流角度工況下三圓柱的流體力系數均方根值（CD-rms，CL-rms）隨間距比的變化特性。隨來流角度的增大，C1 的CD-rms會減小，如圖11（a）所示。另外，隨著間距比的增大，C1 的CD-rms逐漸減小趨近于0，當Kd=1.2 時，C1 的CD-rms會達到最大值，當來流角度較小時，自由流側剪切層受到C2 的干擾較小，較少附著在C2 的表面，導致發展越充分，C1的阻力系數波動范圍會增大。

圖11 不同來流角度下三圓柱體群阻力系數均方根值隨間距比的變化Fig.11 Variation of root mean square value of the drag coefficient for three-cylinder group with the spacing ratio in different incoming flow angles

圖12 不同來流角度下三圓柱體群升力系數均方根值隨間距比的變化Fig.12 Variation of root mean square value of the lift coefficient for three-cylinder group with the spacing ratio in different incoming flow angles

在小間距比（Kd≤1.5）下，α=0°工況流場關于順流向中心線對稱，下游渦街分布規則，如圖13（a）與13（b）所示，因而C2的CD-rms大于其它角度下的數值。然而，當α=45°時，Kd增加，間隙側剪切層形成，漩渦脫落位置遠離C2，如圖13（d）所示，對C2 的影響會減弱，此時CD-rms減小。隨著Kd的進一步增加，CD-rms逐漸增大，再次達到極大值后小幅回落。當Kd≥3.0 時，由于尾流模式的改變，C2 的CD-rms出現突變。在α=30°與Kd=3.5 工況下，尾流呈ITV 模式。另外，在α=60°與Kd=4.0 工況下，C2 的CD-rms發生陡降。間隙流進一步發展，切斷了間隙內動量以及渦量的輸送，如圖13（e）所示，C1 尾流拉伸較長而未卷起，間隙流由于層間的相互擠壓速度加快，打斷了剪切層的相互作用，脫落的漩渦變形程度大，強度損失較大，因而單看C2 或者C3 下游渦街正負交替特征并不明顯，同時C2 和C3 渦脫落相位差約為180°，呈反相模式[9]。然而，當Kd增大至4.5后，C1尾流區有了足夠的空間，如圖13（f）所示，能夠形成完整的渦結構，C2 和C3 渦脫落同相，呈同相模式[9]，尾流模式從反相轉變為同相時，也會導致C2 的CL-rms發生較大變化。

圖13 不同工況下的流場瞬時渦量圖Fig.13 Instantaneous vorticity diagrams of the flow field at different cases

當Kd=1.2時，各來流角度下尾流模式均為SV模式，渦街距C3較近，且漩渦脫落最規則，造成C3的CD-rms最高，Kd增至1.5 后，間隙側剪切層初步發展，尾流模式轉變，漩渦位置距C3 較遠，導致其CD-rms減小，Kd≥2.0后變化幅度均較小。

升力系數方面，C1的CL-rms在Kd=1.2處取得最大值，逐漸減小至0附近后趨于穩定。隨著間距比的增加，尾流模式有所轉變，導致CL-rms發生突變。在α=30°與Kd=3.0 工況下CL-rms曲線的躍升，這是由尾流模式從DV-TV 轉變為ITV所引起，C1、C2剪切層融合產生更大的漩渦，導致CL-rms達到0.289，之后的CL-rms均未出現顯著變化情況。

如圖12（b）所示，α=30°與Kd=3.5工況下C2的CL-rms達到1.05。圖13（g）和（h）從機理上展示了突變的原由，Kd=3.0 增至3.5，C1 剪切層卷起后脫渦，C2 后方漩渦在下游形成了更加明顯的渦街，對柱體受力影響顯著。分析兩個工況C2 的CD與CL相位，如圖14（a）所示，Kd=3.0 時，相位關系紊亂，CD、CL整體變化幅度較小，CL在0 附近變化。然而，當Kd=3.5 時，尾流模式從DV-TV 轉變為ITV，其相位曲線呈規則的“月牙形”，CD、CL的變化幅度顯著增加。類似地，在α=45°與Kd=4.5 工況下的CL-rms達到0.75，比較圖13（i）和圖15（j），C1 剪切層區別很大，Kd=4.0 工況下C1 剪切層一直被拉伸而不脫落，而Kd=4.5 工況下C1 剪切層能夠脫落漩渦，受此影響，C2、C3 的剪切層被扭轉和拉伸，脫落的漩渦變成長條狀，尾流模式也由DV-TV 轉變為ITV。從相位分析，如圖14（b）所示，Kd=4.0時，CD、CL相位變化呈長條狀，變化幅度較小，CL均小于0，表明C2只承受向下的流體力作用。另外，CD增大時CL呈減小趨勢，CD減小時則相反，反映了C2受力機理，即CD增大時C2下游的渦吸力效應占主導，CL增大時C3下游渦吸力效應占主導。當Kd=4.5時，由于受下游渦街的影響顯著，相位曲線呈“半月形”，CL近似關于零點對稱，變化幅度顯著增大。

圖14 不同工況下C2的CD、CL相位圖Fig.14 Phase diagram of CD and CL for C2 at different cases

圖15 不同工況下三圓柱群斯托羅哈數隨間距比的變化情況Fig.15 Variation of the Strouhal number for three-cylinder group with the spacing ratio at different cases

如圖12（c）所示，不同來流角度工況下，C3 的CL-rms呈先減小再增大至穩定的變化規律，其中α=0°和α=60°下C3的CL-rms分別與C1、C2一致。值得關注的是，在α=45°與Kd=4.5工況下，C3的CL-rms突然減小，這是由于隨間距比增加，下游圓柱對上游C1尾流發展的抑制效果會減弱，三圓柱群間尾流互擾效應導致C3的繞流發生變化，渦街的演化受到影響。

綜上所述，剪切層間的互擾效應影響了漩渦的脫落，同時在間隙流的干擾下，形成了不同的漩渦脫落模式，對圓柱體的脈動力系數產生了不同的影響。

2.3 渦脫落頻率

圖15 為三圓柱體斯托羅哈數隨間距比和來流角度的變化情況。當間距比較小時（Kd≤1.5），流場漩渦數量較少，尾流模式為SV 和SV-DV 模式，三圓柱體的渦脫落頻率均會逐漸減小。當α≤15°時，三圓柱體結構的渦脫落頻率在Kd=2.0處達到最小，然而，當α>15°時St最小值出現在Kd=1.5處，隨著Kd的增大，三圓柱體群結構的渦脫落得到加強，St值均會有較大幅度的提升。圖16展示了五個工況的能量譜密度（PSD），在圖16（a）～（c）中，C3主峰對應的St值變化不大（分別為0.19、0.185、0.175），表明C3受C1和C2的影響較小。在α=30°與Kd=3.0工況下，C2出現了兩個較大的PSD 主峰，其St值分別為0.127與0.19（和C1 相同），α=45°時C2 的兩個PSD 主峰接近，并且其中一個St仍和C1 相同，但C1 峰值PSD接近于0，在α=30°與Kd=3.5工況下C1和C2都僅有一個主峰，但C1仍接近于0，α=30°與Kd=3.5工況和α=45°與Kd=3.0 工況下，C1 和C2 間隙內不會發生漩渦脫落。α=45°、Kd=4.5 時C3 的St出現一次突變，發生該現象的原因是C3和C2尾流之間的互擾效應致漩渦無法正常脫落。當α=60°時，Kd=4.0與4.5時C2 和C3 主峰對應的St幾乎一致，如圖16（d）～（e）所示，Kd從4.0 變化至4.5 后，C1 剪切層尾部能夠卷起，以C2、C3 相同的頻率脫落漩渦，同時，三圓柱主峰所對應的能量均提升較大，表明當Kd=4.5 時C1尾流能夠干擾圓柱的漩渦脫落方式。隨Kd繼續增大時，剪切層間與下游渦相互作用逐漸減弱，后方渦街更加規則，呈現RTV 模式，渦脫落頻率最終趨于穩定。α=0°工況下渦街演化更規則有序，圓柱間的互擾效應以及流場下游的漩渦融合變形最弱，因此在α=0°與Kd≥3.5下，C1和C3的St值均較大。除α=60°工況，C2 始終位于順流向的下游，其和上游圓柱的相互作用更強，相較于C1 和C3，C2 的渦脫落更容易受到間距比和來流角度的影響。

3 結 論

本文基于IB-MRT-LBM 研究了Re=100下等邊三角布置三圓柱體結構群的繞流特性，分析了不同來流角度和間距比對三圓柱體群尾流模式及其流體力系數的影響，闡述了其流動機理。主要結論如下：

（1）不同來流角度下尾流模式分為SV、SV-DV、RDV、IDV、DV-TV、RTV、ITV 七種模式。α=30°時下游渦街更加規則，且尾流模式轉變更快，α=45°和α=60°時尾流模式則轉變較慢。各來流角度下尾流模式均存在從不規則向規則轉變的過程，Kd≥5.5后尾流模式均為RTV。

（2）流體力系數時均值受間隙流流速影響較大，除了小間距比，各來流角度下隨間距比變化趨勢基本一致，α=45°時在較大間距比范圍內出現了升力系數均值為負的情況。在Kd≥3.5 后，尾流模式的轉變對C1 的流體力系數均方根值有明顯影響。同時，剪切層的擺動、融合等互擾效應也會改變流場的脈動特性。

（3）不同來流角度工況下，在Kd≥3.5 時C1 和C3 的St受Kd的影響較小，數值變化相對穩定。α=0°工況下，當Kd≥3.5時，C1和C3的St值均較大。除α=60°工況，C2的渦脫落更容易受到間距比和來流角度的影響。