水利工程論文-非正交曲線坐標下二維水流計算的SIMPLEC算法.doc

水利工程論文-非正交曲線坐標下二維水流計算的SIMPLEC算法摘要：本文采用Laplace方程坐標變換方法生成正交曲線網(wǎng)格，并對淺水流動的控制方程進行坐標變換，方程離散時采用B型交錯網(wǎng)格。利用“水位掃描法”結(jié)合壁面函數(shù)法來處理移動邊界，用SIMPLEC算法解非正交曲線坐標下的k-雙方程紊流模型，修正了由網(wǎng)格的非正交性引起的誤差。通過對美國Colorado洲FallRiver的資料進行流場驗證，計算結(jié)果與實測資料基本符合，顯示了本模型在不規(guī)則水域計算中的實用價值。關鍵詞：坐標變換k-紊流模型水位掃描法壁面函數(shù)SIMPLEC算法隨著經(jīng)濟發(fā)展和社會進步，水利工程建設的步伐也在進一步加快，其中港航建設、大壩建設中的泥沙問題以及近來倍受世人關注的水污染問題已經(jīng)成為制約水利發(fā)展的瓶頸問題，弄清河流、湖泊、海洋中水動力因素，是解決以上問題的重要基礎。近年來，數(shù)學模型已逐步取代物理模型實驗成為研究水流的重要手段，而淺水流動模型是處理大區(qū)域流場的一種非常有效的模型。它屬于非線性方程組，在目前只能用數(shù)值方法求解，因此，有必要研究一種簡單、高效的方法來求解淺水流動問題。自Patankar和Spalding1發(fā)展了SIMPLE算法以來，該方法被廣泛應用于不可壓縮流場的數(shù)值模擬，而且該方法還得到了進一步的發(fā)展，主要有SIMPLER算法2、SIMPLEC算法3、SIMPLEX算法4和SIMPLET算法5等。這些模型均成功地應用于速度壓力耦合的流場計算，深度平均的淺水流動模型是在靜壓假定下導出的，一般流體模型中的速度壓力耦合也就轉(zhuǎn)換成淺水流動模型中的速度水深耦合6。天然河流、海灣的邊界曲折、地形復雜，采用坐標變換是解決問題的途徑之一。目前多數(shù)N-S方程的坐標變換中，流程全部采用逆變分量，這樣就增加了方程的復雜程度。于是忽略掉方程中的非正交項，利用正交變換下的方程進行數(shù)值求解7,8。對于具有復雜邊界的海灣及彎曲的河流，坐標變換中很難保證每個點都正交，特別是邊界附近。水位變化是水力計算中難點之一，在目前的紊流數(shù)學模型中，多簡單的利用“凍結(jié)法”，這樣做將失去對邊界出流動模擬的準確性。本文研究中，采用正交曲線坐標變換生成數(shù)值網(wǎng)格，而數(shù)值計算中采用非正交曲線坐標下的k-雙方程紊流模型，這樣可以自動修正網(wǎng)格生成中的非正交項。流速除對流項中采用逆變分量，在其余各項中均采用原始分量，這樣使得方程書寫簡單，有利于將各方程寫成通用形式，編寫的程序變得更規(guī)范。作者受JianYe同位網(wǎng)格9的啟發(fā)，對普通交錯網(wǎng)格做了修改，即采用B型交錯網(wǎng)格，使得u,v,k,的計算布置在一個節(jié)點上，有利于節(jié)省計算程序代碼，使程序書寫更加規(guī)范。引入動邊界掃描技術，結(jié)合紊流模型的壁面函數(shù)法，使壁面隨著真實邊界而變化。數(shù)值求解時，采用控制體積法離散方程，運用SIMPLEC算法，使計算的流場更符合實際流場。1數(shù)值網(wǎng)格本文對計算區(qū)域用Laplace方程實施坐標變換，生成正交的貼體網(wǎng)格，控制方程：(1)方程組(1)可采用SOR方法求解，然后用正交曲線網(wǎng)格邊界正交化處理的三次樣條插值邊界滑動法10，處理其邊界處的網(wǎng)格，以提高其正交性。與其它方法相比，此方法雖在一定程度上提高了邊界處網(wǎng)格的正交性，但難以保證完全正交，如果再把曲線坐標系下N-S方程中的非正交項去掉，由斜交網(wǎng)格引起的計算誤差是不容忽視的。為了修正此項誤差，本文采用非正交曲線坐標下的平面二維k-雙方程紊流模型，模型能自動修正網(wǎng)格生成中的非正交項。2數(shù)學模型笛卡兒坐標下深度平均的k-雙方程紊流模型的通用微分方程：(2)方程(2)轉(zhuǎn)換到曲線坐標(、)下，僅在對流項中使用流速的逆變分量，而在其它項中使用原始變量，這樣既簡化了方程，又使所有方程仍可寫為曲線坐標下的通用方程，模型的微分方程可寫為如下通用形式：(3)式中：為所求問題的因變量；U和V分別為直角坐標下流速u和v的逆變分量，僅在對流項中出現(xiàn)；為擴散系數(shù)；S為源項。當表示某一特定量時，和S對此特定量有特定的意義和表達式，這時方程(3)亦賦予特定的意義，模型的控制方程組如表1所示。式中：u和v分別為笛卡兒坐標下x和y方向的深度平均的流速分量；k和分別為深度平均的紊動動能和紊動動能耗散率；H為水深；zs為水位；zb為河床高程；為分子運動粘性系數(shù)；t為漩渦運動粘性系數(shù)；bx和by分別為x和y方向上的底摩擦力eff為有效粘性系數(shù)；C、C1、C2、k和為經(jīng)驗常數(shù)，其中C為0.09，C1為1.44，C2為1.92，k為1.0，為1.3。3模型離散如圖1所示，本文采用B型交錯網(wǎng)格，即將u,v,k，及各標量布置在網(wǎng)格中心，將水深相關項H、zs、zb項布置在網(wǎng)格的左下方節(jié)點上。這樣布置有以下兩個優(yōu)點：有利于邊界條件的引入；在編寫程序u,v,k,各方程共用一套系數(shù)，可節(jié)省程序代碼，使程序簡單明了。本文用控制容積法，采用混合格式來導出離散方程。其邊界條件為：在進口邊界，所有邊界條件都按本征條件給出，即：u=u0,v=0,k=k0,=0；在出口邊界，給定出口處的水位zs，令v=0，其余變量都取法向?qū)?shù)為零，即：在水陸邊界，在計算河流、河口、海灣的非恒定流場時，水位不是恒定不變的，而且每一時間步長的水位可能都不一樣，這就帶來了水陸邊界的移動，如果每一時間步長都進行坐標變換而生成新的網(wǎng)格是非常不經(jīng)濟的。于是本文在前人“凍結(jié)法”的基礎上，提出了適合k-雙方程紊流模型的“移動邊界的壁面函數(shù)法”。如圖2、圖3，水位變化后，使邊界處的網(wǎng)格干出，原來的岸邊界12、34變成了新的岸邊界PQ、ST，而壁面函數(shù)仍然布置在12、34岸邊界處；水域內(nèi)島嶼由于水位的變化而出露，其島嶼內(nèi)部按“凍結(jié)法”來處理，其邊界ABCD周圍也需用壁面函數(shù)法來處理。本文采用“水位掃描法”來判斷新的水陸邊界：從邊界12向34掃描，如水深不為零，則令其為起始點；從邊界34向12掃描，如水深不為零，則令其為終點；然后掃描水域內(nèi)部，確定島嶼的邊界ABCD。之后對水陸邊界實施壁面函數(shù)法。在壁面邊界，k-雙方程紊流模型是高Reynolds數(shù)模型，只適用于充分發(fā)展的紊流，而在近壁處，粘性效應起主導作用。另外，由于近壁處物理量變化很快，需要非常密的網(wǎng)格，才能真實地模擬實際流動。因此本文采用壁面函數(shù)法11，即用半分析的方法得到的解來近似由壁面到紊流核心區(qū)之間的流速、紊動動能、紊動動能耗散率的分布規(guī)律，將壁面的影響(如壁面應力)附加到差分方程中(應先將有關的邊界系數(shù)置為零)12。4SIMPLEC算法求解模型SIMPLE算法的基本步驟如下：(1)計算坐標變換的相關系數(shù)；(2)根據(jù)進、出口邊界的水位，給定全場水位，計算全場初始水深H*；(3)解動量方程求un,vn；(4)解k-紊流模型求k、；(5)解水深校正方程求H，并修正水深H=H*+hH；(6)修正流速u、v及有效粘性系數(shù)eff，然后重復(2)(5)步直至收斂。模型求解中，為了利于非線性迭代的收斂，計算中采用亞松弛技術。計算中采用ADI技術和TDMA算法，以連續(xù)方程的誤差小于給定的值作為判斷收斂的依據(jù)。在恒定流計算中，以相鄰時間步之間全場各點流速誤差的最大絕對值小于給定的值作為達到定常狀態(tài)的依據(jù)。5模型應用為了使本文數(shù)學模型具有代表意義，本文采用天然連續(xù)彎道河流美國Colo