天堂国产午夜亚洲专区-少妇人妻综合久久蜜臀-国产成人户外露出视频在线-国产91传媒一区二区三区

當(dāng)前位置:主頁 > 科技論文 > 水利工程論文 >

三維自由面湍流場數(shù)值模擬及其在水利工程中的應(yīng)用

發(fā)布時間:2016-11-21 15:06

  本文關(guān)鍵詞:三維自由面湍流場數(shù)值模擬及其在水利工程中的應(yīng)用,由筆耕文化傳播整理發(fā)布。



河海大學(xué) 博士學(xué)位論文 三維自由面湍流場數(shù)值模擬及其在水利工程中的應(yīng)用 姓名:王志東 申請學(xué)位級別:博士 專業(yè):水力學(xué)及河流動力學(xué) 指導(dǎo)教師:汪德爟 20040322

河海大學(xué)博士論文

摘要
自然界中的流體流動問題一般都具有三維性和湍動性,對于最常見的流體一水則
同時還擁有自由表面,水利工程中繞閘墩的泄洪水流、船

舶工程中的船體繞流等均具 有上述流動的三種特性,這三種特性是計(jì)算流體力學(xué)(CFD)領(lǐng)域非常重要的研究內(nèi) 容和方向。


本文針對粘流場數(shù)值模擬技術(shù)的研究現(xiàn)狀和發(fā)展動態(tài)作了比較全面的回顧與展 望,重點(diǎn)研究了網(wǎng)格生成技術(shù)、數(shù)值求解技術(shù)、湍流模型技術(shù)和動邊界模擬技術(shù),在 此基礎(chǔ)上建立了模擬自由面湍流場的數(shù)學(xué)模型,并成功地應(yīng)用N-維溢流壩和帶閘墩
的三維溢流壩過壩水流的數(shù)值計(jì)算。

本文完成的具體工作和主要結(jié)論如下: 1)以代數(shù)網(wǎng)格生成方法為基礎(chǔ)提出了一種簡單的、可獨(dú)立于網(wǎng)格生成方法之外的邊 界正交化技術(shù);針對分區(qū)結(jié)構(gòu)網(wǎng)格系統(tǒng)建立了分區(qū)交界面處的數(shù)據(jù)結(jié)構(gòu)與計(jì)算模型: 2)利用有限體積方法在非正交同位網(wǎng)格系統(tǒng)中建立了SIMPLE求解算法,對非正交 網(wǎng)格系統(tǒng)中的邊界條件、延遲修正技術(shù)及計(jì)算節(jié)點(diǎn)的梯度計(jì)算等專題進(jìn)行了深入討 論。通過對傾斜方腔驅(qū)動流、后臺階繞流、圓柱繞流等典型流場的數(shù)值模擬,表明本 文建立的非正交同位網(wǎng)格系統(tǒng)中的SIMPLE算法是合理的; 3)研究比較了各種湍流模型的特點(diǎn)及應(yīng)用范圍,以水翼繞流場為例完成了RNGk一£ 模型與標(biāo)準(zhǔn)k—s模型的計(jì)算比較,結(jié)果表明RNGk—s模型能夠更好地模擬流場的湍
流特性:

4)研究了運(yùn)動邊界模擬技術(shù)中的VOF方法,詳細(xì)建立了一種精度更高的自由面重構(gòu) 模型,通過對典型運(yùn)動界面的數(shù)值模擬表明,該方法能夠更加準(zhǔn)確地確定運(yùn)動界面的
形狀和位置:

5)利用本文建立的非正交網(wǎng)格系統(tǒng)中自由面湍流場數(shù)學(xué)模型,對二維溢流壩過壩水 流進(jìn)行了數(shù)值計(jì)算,分析了不同壩頂水頭下的水面線位置、壩面壓力及流量、流速等 水流特性,通過與典型實(shí)驗(yàn)資料的對比,水面線高度及壩面負(fù)壓極值均與實(shí)驗(yàn)值具有
非常好的一致性。 6)通過建立溢流壩和閘墩的三維整體模型,完成了對三維過壩水流粘流場的數(shù)值計(jì)

算,針對閘墩的型式及在壩面的布置計(jì)算了閘墩對水面線形狀、壩面壓力、流量系數(shù) 等設(shè)計(jì)參數(shù)的影響,并將不同墩型與布置形式的計(jì)算結(jié)果進(jìn)行比較。結(jié)果表明,閘墩

II

河海大學(xué)博士論文

的存在抬高了水面線的位置,其中在閘墩頭部尤其明顯:墩型對流量系數(shù)和壩面壓力 影響較大.3A型閘墩相對于2型閘墩具有更大的流量系數(shù)和更小的壩面壓力:不同
的墩厚閘寬比f/6對泄流能力也將產(chǎn)生顯著的影響.隨著墩厚閘寬比的增加,壩面壓

力降低,而當(dāng)tlb=0.2時溢流壩具有更大的流量系數(shù)。 7)本文建立的三維自由面湍流場模型能夠比較準(zhǔn)確地模擬帶閘墩溢流壩三維過壩水 流湍流場,可以針對不同的壩型、墩型、壩面布置形式以及不同的壩頂水頭完成系列
水力計(jì)算,為泄水工程建筑物的設(shè)計(jì)提供可靠的分析依據(jù)。

關(guān)鍵詞:數(shù)值模擬,網(wǎng)格,湍流模型,流體體積法,溢流壩,閘墩

II

河海大學(xué)博士論文

ABSTRACT
Three dimension and turbulent
are

common characteristic of fluid flow

problem in such flow
as

nature

and

added flee surface for the best familiar fluid-water,

flood discharge flow around frusta of brake in waterpower in ship engineering,and
SO

project,

around hull

on,all

have the three

characteristics which are the very important study content and aspects in
CFD.
A comprehensive review and expectation about researching status in quo

and development trends of viscous flow numerical simulation technique are given in this paper

and grids generation

methods,numerical
job


computation
are

methods,turbulent model and moving boundary simulation methods emphasis study in the paper.On the basic of above

the

mathematical model

of simulation turbulent flow with free surface is established and applied
to

triumphantly

numerical computation

on

2D spillway

and

3D spillway with frusta

of brake. The material

job and main

conclusion

are

as follows:

(1)A simple boundary orthogonalization procedure independence grids
generation method is put forward
on

the base of

algebraic鰣d

generation

method;Data
aiming
at

structure

and

computational model

on

interface are established

blocks of structured grids. procedure is established in nonorthogonal grids by
on
use

(2)SIMPLE

of
as

the finite volume method and in-depth discusses

special topic such

boundary condition in nononhogonal grids,deferred correction method grads compute
on

and

calculational

nodes,and

SO

on.The

SIMPLE procedure in

nonorthogonal grids established in the paper is proved to be reasonable through representative numerical simulation such as incline square cavity driving flow,flow around back sidestep,flow around cylinder.

(3)Characteristic and studied and compared and

application range of different turbulent model are compute compare RNGk-c turbulent model with
case

k-占model is established by the

of viscous flow through hydrofoil
call

and

t11e result show that beRer characteristic turbulent flow
RNGk一占.

be obtained by

(4)VOF method in moving boundary simulmion method
iv

is studied and



河海大學(xué)博士論文

higher precision free surface reconstruction model is established detailedly. The method is proved motion interface interface.
to

well and truly confirm the shape and location of numerical simulation of representative motion

through

(5)Compute of flow
free surface、pressure

over

2D spillway is done by

BSe

of

established and

free

surface turbulent mathematical model in nonorthogonal grids
on

location of
as

dam

and

flow characteristic

on

the such

flux、 the

velocity of flow of different water height of dam peak

are

analyzed and
on

result of height of free surface and negative extremum of pressure surface are well inosculate
to

dam

result of representative experiment.

(6)3D

integer model of spillway and frusta of
over

brake are established and

numerical computation of 3D viscous flow

spillway is

completed.The
shape of free

influence of frusta of brake for design parameters such surface,pressure
on

aS

dam,flux coefficient is computed aim at type of frusta of
on

brake and

disposal

dam

and compute result of different type of frusta and
of of

disposal model is compared.The result show that the position of free surface is

higher

because of frusta of

brake;flux coefficient frusta

brake,especial obvious in the head of frusta and pressure on dam are quite influenced by type

and 3A type frusta of brake relative to 2 type frusta of brake haS bigger flux coefficient and smaller pressure on dam;Obvious influence of ability of
discharge flow for different ratio of thickness of frusta

and

breadth of brake

and

pressure

on

dam reduces with increasing ration of thickness of frusta

and

breadth ofbrake,while bigger flux coefficient is obtained when t/b=0.2.

(7)The
brake

3D free surface turbulent model established in the paper is

comparatively well

and

truly simulates 3D flow

over

spillway with frusta of

and series waterpower computation can be completed aim at different type of dam、type of frusta and disposal on the dam and different water height
on

peal(of dam,which will supplies reliable

analysis

data for design of

waterpower construction.

keyword:numerical

simulation,酣d,turbulent model,Volume

ofFluids,

spillway,frusta of brake



河海大學(xué)博士論文

.1上?



——1一



本論文以水利工程中泄水工程的水力設(shè)計(jì)為背景,研究了自由面湍流場數(shù)值模擬
技術(shù),并將之應(yīng)用于帶閘墩溢流壩三維過壩水流場的數(shù)值計(jì)算。

隨著計(jì)算機(jī)科學(xué)和數(shù)值計(jì)算方法的發(fā)展,數(shù)值模擬在科學(xué)技術(shù)研究中具有更大的 優(yōu)勢和發(fā)展空間,一方面可以確定物理問題整體舶和局部的細(xì)致行為,另一方面在解 決實(shí)際工程問題的過程中便于進(jìn)行方案優(yōu)選,從而降低費(fèi)用,縮短周期,此外還可以 替代一些昂貴的、危險(xiǎn)的甚至難以實(shí)現(xiàn)的實(shí)驗(yàn),因此數(shù)值模擬在整個科學(xué)技術(shù)進(jìn)步中
起到非常重要的作用。


然而,利用數(shù)值模擬技術(shù)解決實(shí)際工程問題必須要求數(shù)值模擬過程中的每一個環(huán) 節(jié).如:數(shù)學(xué)模型、控制條件、數(shù)值算法等都是高效的j先進(jìn)的,并經(jīng)過充分驗(yàn)證和 確認(rèn)的,它們決定了數(shù)值模擬成果最終的可靠性和實(shí)用性。 在前人研究工作的基礎(chǔ)上,本論文對數(shù)值模擬的關(guān)鍵技術(shù):網(wǎng)格生成技術(shù)、數(shù)值 求解技術(shù)、湍流模型和動邊界模擬等進(jìn)行了較深入的研究,建立了非正交同位網(wǎng)格系 統(tǒng)中自由面湍流場數(shù)學(xué)模型,并成功地應(yīng)用于二維溢流壩和帶閘墩的三維溢流壩過壩
水流數(shù)值模擬。


通過本論文的工作,取得了以下研究成果: 1)提出了一種簡單的、可獨(dú)立于網(wǎng)格生成方法之外的網(wǎng)格邊界正交化方法; 2)利用分區(qū)結(jié)構(gòu)網(wǎng)格技術(shù)和非正交同位網(wǎng)格系統(tǒng)中SIMPLE算法建立了求解復(fù)雜流 場的數(shù)學(xué)模型,提高了計(jì)算效率和計(jì)算精度; 3)詳細(xì)建立了一種精度更高的自由面重構(gòu)方法,能夠更加準(zhǔn)確地確定運(yùn)動界面的形
狀和位置:

4)利用以上模型首次計(jì)算了帶閘墩溢流壩三維過壩水流粘流場。針對不同的壩頂水 頭、閘墩的型式以及在壩面上的布置完成了系列水力計(jì)算,為泄水工程建筑物的水力 設(shè)計(jì)提供可靠的分析依據(jù)。

學(xué)位論文獨(dú)創(chuàng)性聲明:

本人所呈交的學(xué)位論文是我個人在導(dǎo)師指導(dǎo)下進(jìn)行的研究工作 及取得的研究成果。盡我所知,除了文中特別加以標(biāo)注和致謝的地方
外,論文中不包含其他人已經(jīng)發(fā)表或撰寫過的研究成果。與我一同工

作的同事對本研究所做的任何貢獻(xiàn)均己在論文中做了明確的說明并 表示了謝意。如不實(shí),本人負(fù)全部責(zé)任。

論文作者(簽名):
(注:手寫親筆簽名)

£色:壘
墜&:堅(jiān)

旦壘年





衛(wèi)2日

學(xué)位論文使用授權(quán)說明: 河海大學(xué)、中國科學(xué)技術(shù)信息研究所、國家圖書館、中國學(xué)術(shù)
期刊(光盤版)電子雜志社有權(quán)保留本人所送交學(xué)位論文的復(fù)印件或

電子文檔,可以采用影印、縮印或其他復(fù)制手段保存論文。本人電子
文檔的內(nèi)容和紙質(zhì)論文的內(nèi)容相一致。除在保密期內(nèi)的保密論文外, 允許論文被查閱和借閱。論文全部或部分內(nèi)容的公布(包括刊登)授

權(quán)河海大學(xué)研究生院辦理。

論文作者(簽名):
(注:手寫親筆簽名)

互盤壘

竺年

3月2工日

河海大學(xué)博士論文

符號表
松弛因子,控制網(wǎng)格邊界正交化的程度
P,Q


阻尼參數(shù),控制網(wǎng)格正交性的衰減程度
流體密度 流體運(yùn)動的速度分量

U,(f=1,2,3)

P,
-“

粘性流體中的法向應(yīng)力
粘性流體中的切向應(yīng)力
流體的動力粘度



Re=UL/y U


雷諾數(shù)

特征速度
特征長度 流體的運(yùn)動粘度





渦粘系數(shù) 湍動能 湍動能耗散率
壁面切應(yīng)力

F Hs

流體體積分?jǐn)?shù)
溢流壩高 溢流壩頂設(shè)計(jì)水頭 溢流壩頂水頭
溢流壩流量

H1

日。

,竹

溢流壩流量系數(shù)
溢流壩面壓力水頭
閘墩厚度
閘孔寬度



,



IX

河海大學(xué)博士論文

第一章緒論
1.1數(shù)值模擬技術(shù)的工程研究意義
數(shù)值模擬、理論分析和模型試驗(yàn)是研究流體力學(xué)(水力學(xué))的三種主要手段, ~般來說,這三種研究手段和方法必須相互配合,相互補(bǔ)充,相互促進(jìn),才能共同 推進(jìn)流體力學(xué)學(xué)科的發(fā)展并解決各種工程實(shí)際問題。

理論分析方法是在研究流體運(yùn)動規(guī)律的基礎(chǔ)上提出各種簡化流動模型,建立各
種控制方程,在一定的假使和條件下,經(jīng)過一系列的解析推導(dǎo)和運(yùn)算得到問題的解 析解,理論分析的許多方法仍是目前解決實(shí)際問題,主要是在初步設(shè)計(jì)階段中常常 采用的方法,但理論分析法常常無法用于研究復(fù)雜的、以非線性為主的流動現(xiàn)象。

模型試驗(yàn)方法是研究流動機(jī)理、分析流動現(xiàn)象、探討并獲得流動新概念,推動
流體力學(xué)發(fā)展的主要手段,其不足是完成一個完整的試驗(yàn)過程需要解決一系列復(fù)雜 的技術(shù)問題,所需周期長、費(fèi)用高。 數(shù)值模擬方法的特點(diǎn)是在用比模型試驗(yàn)所需花費(fèi)少得多的情況下給出流場內(nèi)部 細(xì)節(jié)的詳細(xì)描述。并且只要數(shù)值模擬的數(shù)學(xué)提法(控制方程、邊界條件、數(shù)值算法 等)是正確的,則可以在較廣泛的流動參數(shù)(如雷諾數(shù),空化數(shù)等)和物理設(shè)計(jì)參

數(shù)范圍內(nèi)較快地給出流場的定量結(jié)果,不受試驗(yàn)中固有約束條件的影響。因此隨著 計(jì)算流體力學(xué)和計(jì)算機(jī)科學(xué)的迅速發(fā)展,數(shù)值模擬在科學(xué)技術(shù)研究中將具有更大的 優(yōu)勢和發(fā)展空間,主要體現(xiàn)在:①數(shù)值模擬在某種意義下比理論和實(shí)驗(yàn)研究對流體
運(yùn)動過程的認(rèn)識更加深刻、細(xì)致,不僅可以得到運(yùn)動的結(jié)果,而且可以了解整體的

和局部的細(xì)致行為;②在解決實(shí)際工程問題的過程中,利用數(shù)值模擬可以選擇最優(yōu) 的設(shè)計(jì)方案或?qū)嶒?yàn)方案,從而減少實(shí)驗(yàn)次數(shù),周期及費(fèi)用:③數(shù)值模擬可以從理論 上探索流體運(yùn)動新的現(xiàn)象和規(guī)律.而且可以替代一些昂貴的、危險(xiǎn)的甚至難以實(shí)現(xiàn)
的實(shí)驗(yàn),比如水壩坍塌造成的災(zāi)害預(yù)報(bào)、水下爆炸與核爆炸的過程和效應(yīng)、大氣運(yùn) 動現(xiàn)象等。因此數(shù)值模擬在流體力學(xué)研究乃至整個科學(xué)技術(shù)進(jìn)步中起到非常重要的
作用。

雖然數(shù)值模擬在流體力學(xué)研究中占有如此重要的地位,然而要使其一方面能夠?yàn)?br />工程設(shè)計(jì)提供高質(zhì)量、短周期、可靠的分析設(shè)計(jì)依據(jù),另一方面為探索流動新概念、 新機(jī)理提供新途徑,則要求數(shù)值模擬過程中的每~個環(huán)節(jié),如:數(shù)學(xué)模型、控制條件、 數(shù)值算法等都是高效的、先進(jìn)的,并經(jīng)過充分驗(yàn)證和確認(rèn)的,它們決定了數(shù)值模擬成 果最終的可靠性和實(shí)用性。

河海大學(xué)博士論文

1.2粘流場數(shù)值模擬關(guān)鍵技術(shù)研究現(xiàn)狀 1.2.1網(wǎng)格生成技術(shù)
網(wǎng)格生成技術(shù)是計(jì)算流體力學(xué)(CFD)發(fā)展的一個重要分支,在CFD高度發(fā)展的

美國,網(wǎng)格生成所需的人力時間占計(jì)算任務(wù)全部人力時間的60%,因此將CFD技術(shù)進(jìn)
一步推向工程應(yīng)用,首先必須完善網(wǎng)格生成技術(shù)“。。對于復(fù)雜程度一般的計(jì)算區(qū)域邊

界,采用結(jié)構(gòu)網(wǎng)格系統(tǒng)完全能夠滿足工程要求的精度,由于算法相對簡單,因而目前 仍被廣泛采用。常用的三維結(jié)構(gòu)網(wǎng)格生成方法大致可分為代數(shù)生成法“。、橢圓微分方 程生成法”’”和雙曲微分方程生成法。”等三類。 代數(shù)生成法是利用已知的邊界值進(jìn)行中間插值來生成網(wǎng)格,其核心是確定插值函
數(shù),晟常用的是“超限插值”(Transfinite interpolation),內(nèi)部網(wǎng)格節(jié)點(diǎn)的分布 可以由邊界上的伸縮函數(shù)來控制,通過一些專門的、不是十分復(fù)雜的技術(shù)能夠?qū)崿F(xiàn)網(wǎng) 格在邊界的正交或近似正交“““,因此,對于幾何形狀不是十分復(fù)雜的流場,代數(shù)生

成方法不失為一種簡單有效的方法。 橢圓型網(wǎng)格生成法是通過求解橢圓型微分方程來實(shí)現(xiàn)網(wǎng)格的生成,該方法的最
早代表為著名的TTM方法”’,即通過求解Laplace方程來生成區(qū)域網(wǎng)格,該方法的 優(yōu)點(diǎn)是所生成的網(wǎng)格是光滑的,能夠處理復(fù)雜的幾何邊界,同時若在方程右端增加

源項(xiàng)交成泊松方程則可以控制網(wǎng)格分布的疏密及傾斜度,其缺點(diǎn)是這樣的源項(xiàng)很難 確定,且較難實(shí)現(xiàn)內(nèi)部節(jié)點(diǎn)的控制。目前最常用的源項(xiàng)控制方法是根據(jù)邊界正交性 和網(wǎng)格間距的要求壹接推導(dǎo)出右端項(xiàng)的表達(dá)式”’…,橢圓型網(wǎng)格生成方法是最常用
的一種單域貼體網(wǎng)格生成方法。 雙曲型網(wǎng)格生成法是通過求解雙益型微分方程來實(shí)現(xiàn)的,該方法的優(yōu)點(diǎn)是能夠保 證網(wǎng)格在邊界處的正交性,容易控制兩格的尺寸,且計(jì)算時間比橢圓方法少卜2個量

級,缺點(diǎn)是邊界上的任何不連續(xù)都會傳播到內(nèi)部網(wǎng)格節(jié)點(diǎn),外邊界位置不能事先給定。 但隨著對一些敏感問題,如邊界條件、物面不連續(xù)和角點(diǎn)的處理等所開展的討論和改
進(jìn)。提高了該方法的適應(yīng)能力,在近年來得到了比較廣泛的應(yīng)用”’”‘“。。

隨著邊界復(fù)雜程度的提高,生成單域貼體計(jì)算網(wǎng)格更加困難,同時網(wǎng)格質(zhì)量不 能得到保證。影響最終的數(shù)值計(jì)算結(jié)果,為此人f|、]提出了分區(qū)結(jié)構(gòu)網(wǎng)格、非結(jié)構(gòu)網(wǎng) 格和結(jié)構(gòu)/非結(jié)構(gòu)混合網(wǎng)格技術(shù)。 分區(qū)結(jié)構(gòu)網(wǎng)格系統(tǒng)是目前比較成熟的、并經(jīng)常采用的構(gòu)造復(fù)雜區(qū)域網(wǎng)格的方 法,分區(qū)結(jié)構(gòu)網(wǎng)格方法的基本思想是根據(jù)邊界的特點(diǎn)將整體流場劃分為若干個子
區(qū),對每一個子區(qū)分別建立網(wǎng)格,并在其中對流動控制方程求解。但在區(qū)域劃分對

需遵循以下原則:盡量使每個子域的邊界簡單以方便結(jié)構(gòu)網(wǎng)格生成;各子區(qū)的大小
也盡量相近以實(shí)現(xiàn)計(jì)算負(fù)荷的平衡,這一點(diǎn)對并行計(jì)算尤為重要。常用的分區(qū)網(wǎng)格



河海大學(xué)博士論文

方法有:對接網(wǎng)格方法和重疊網(wǎng)格方法。
國外分區(qū)網(wǎng)格生成技術(shù)于80年代興起并逐漸形成實(shí)用軟件,具有代表性的有

EAGLE…。和GRIDGEN““。前者是以橢圓型網(wǎng)格生成理論為基礎(chǔ),兼有代數(shù)生成方法 (超限插值)和表面網(wǎng)格生成技術(shù):后者是以代數(shù)網(wǎng)格生成方法為基礎(chǔ),功能不如 前者完善。到了90年代,二者又補(bǔ)充了許多新的技術(shù),如基于弧長的超限插值函 數(shù),提高正交性的Hermit三次型函數(shù),良好的用戶界面(6UI)和加強(qiáng)的與CAD接 口能力等,從而發(fā)展成為功能很強(qiáng)的網(wǎng)格生成系統(tǒng)。其基本思想是依靠超限插值方 法生成表面和空間網(wǎng)格,然后由橢圓方法光順,通過調(diào)整代數(shù)方法中的Hermite超

限插值或橢圓方法中的源項(xiàng)控制函數(shù)來使網(wǎng)格線在交界面處的斜率具有連續(xù)性,以
保證交界面兩側(cè)網(wǎng)格的正交性。 非結(jié)構(gòu)網(wǎng)格技術(shù)于80年代末90年代初得到迅速發(fā)展“““。,其基本思想是:由 于四面體是三維空間最簡單的形狀,任何空間區(qū)域都可以被四面體單元所填滿,即

任何空間區(qū)域都可以被以四面體為單元的網(wǎng)格所劃分。顯然,非結(jié)構(gòu)網(wǎng)格方法舍去
了網(wǎng)格節(jié)點(diǎn)的結(jié)構(gòu)性限制,節(jié)點(diǎn)和單元的分布是任意的,能夠較好地處理復(fù)雜邊界。 非結(jié)構(gòu)網(wǎng)格方法在其生成過程中都采用一定準(zhǔn)則進(jìn)行優(yōu)化判定,因而能生成高質(zhì)量 的網(wǎng)格,且很容易控制網(wǎng)格的大小和節(jié)點(diǎn)的疏密度,一旦在邊界上指定網(wǎng)格分布, 在兩個邊界之間即可自動生成網(wǎng)格,因此非結(jié)構(gòu)網(wǎng)格系統(tǒng)在近年來受到了很大的重

視,有了較大的發(fā)展,如文獻(xiàn)[16]針對二維粘流場構(gòu)造出適應(yīng)于非結(jié)構(gòu)網(wǎng)格的有限 體積時間推進(jìn)格式;文獻(xiàn)[17]利用點(diǎn)自動生成的Delaunay三角化方法構(gòu)造非結(jié)構(gòu)
網(wǎng)格,同時處理成自適應(yīng)網(wǎng)格形式,能夠計(jì)算帶移動邊界的瞬時流動問題等。非結(jié) 構(gòu)網(wǎng)格的主要缺點(diǎn)是所需計(jì)算內(nèi)存較大,計(jì)算效率較低。 結(jié)構(gòu)/非結(jié)構(gòu)混合網(wǎng)格系統(tǒng)則兼顧了二者的優(yōu)點(diǎn),對于幾何形狀復(fù)雜的三維問 題可以由近邊界處的非結(jié)構(gòu)網(wǎng)格和遠(yuǎn)離邊界的結(jié)構(gòu)網(wǎng)格構(gòu)成,具有靈活模擬復(fù)雜邊
界和計(jì)算效率高的特點(diǎn)。

目前分區(qū)結(jié)構(gòu)網(wǎng)格的生成方法和并行計(jì)算方法、結(jié)構(gòu)/非結(jié)構(gòu)混合網(wǎng)格系統(tǒng)是計(jì) 算具有復(fù)雜幾何形狀流場的主要方法,其進(jìn)一步的發(fā)展方向是減少工作量、實(shí)現(xiàn)自動

生成和自適應(yīng)加密,具有良好的人機(jī)界面及可視化,強(qiáng)調(diào)更有效的數(shù)據(jù)結(jié)構(gòu)等,最終
的目的是能夠快速生成適用于粘流場計(jì)算的高質(zhì)量計(jì)算網(wǎng)格。

今后,網(wǎng)格生成技術(shù)的研究重點(diǎn)應(yīng)當(dāng)是網(wǎng)格的專項(xiàng)技術(shù)研究,如:邊界正交化技
術(shù)、分區(qū)結(jié)構(gòu)網(wǎng)格及并行計(jì)算技術(shù)、網(wǎng)格的自適應(yīng)加密技術(shù)以及動態(tài)網(wǎng)格技術(shù)等。

河海大學(xué)博士論文

1.2.2數(shù)值求解技術(shù) 數(shù)值求解技術(shù)實(shí)際上是一種離散近似的計(jì)算方法,各種各樣的流體力學(xué)問題必須
從給定的微分方程或基本定律出發(fā),建立物理上合理、數(shù)學(xué)上適定、能夠在計(jì)算機(jī)上 進(jìn)行計(jì)算的離散化數(shù)學(xué)模型。大多數(shù)數(shù)值求解技術(shù)的基本思想可歸納為:把原來在時 間和空間上連續(xù)的物理量的場(如速度場、溫度場等),用有限個離散點(diǎn)上的值來代 替,然后按一定的方式建立關(guān)于這些值的代數(shù)方程組并求解,從而獲得物理量場的近 似解。目前流場計(jì)算的主要數(shù)值求解方法有:有限差分法(FDM)、有限元法(FEM)、

有限分析法(刪肘)、邊界元法(BEM)、有限體積法(FVM)等,其中前四種數(shù)
值計(jì)算方法的特點(diǎn)及適應(yīng)范圍在眾多的文獻(xiàn)中均可查閱[18,19,20,2l,22,23],下
面重點(diǎn)介紹有限體積法。

有限體積法由McDonald在1971年首次用于求解二維歐拉方程,1972年被 Patanker等人用于SIMPLE算法計(jì)算恒定不可壓流。有限體積法的基本思想是:將計(jì) 算區(qū)域劃分成一系列規(guī)則或不規(guī)則的單元或控制體積,并使每一個控制體積包含一個 網(wǎng)格節(jié)點(diǎn),將待解的微分方程對每一個控制體積分得到一組離散方程,其中的未知數(shù) 是網(wǎng)格節(jié)點(diǎn)上因變量的值,首先計(jì)算出通過每個控制體邊界沿法向輸入(輸出)流量 的動量通量,然后對每一個控制體分別進(jìn)行流量及動量平衡計(jì)算,便可得到計(jì)算時刻

各個控制體內(nèi)流動參數(shù)的平均值(平均流速、平均壓強(qiáng)等)。因此,有限體積法實(shí)際 上是對推導(dǎo)原始微分方程所用控制體方法的回歸,與FDM和FEM的數(shù)值逼近相比 其物理意義更加直接明確,因而在流場計(jì)算中獲得了廣泛的應(yīng)用“。。’“。“。。
若邊界通量的計(jì)算只使用前一個時刻的數(shù)值,為顯式FVM,反之,當(dāng)涉及到下 一個時刻的數(shù)值時,則為隱式FVM。此外還有半隱式FVM“““。。

1.2.2.1數(shù)值方法對離散格式的要求
數(shù)值求解首先要對微分方程進(jìn)行離散,然而不同的離散格式將帶來完全不同的 結(jié)果,因?yàn)橥晃⒎址匠滩煌碾x散格式并不等價(jià),也不一定能保持原來微分方程 的性質(zhì),對漸變光滑的情況其結(jié)果可能相差不大,對迅變以致間斷的情況其解則可

能相差很大。因此利用FVM在建立和使用離散格式時需要滿足以下條件才能保證 數(shù)值求解的精度和可靠性。①相容性:即當(dāng)時空步長趨于零時截?cái)嗾`差(離散方程 與原始方程之差)為零;②守恒性:是指在每一個空間有限的控制體積內(nèi)保持質(zhì)量 和動量平衡,具有守恒性的數(shù)值格式在計(jì)算過程中不會產(chǎn)生守恒誤差,并且計(jì)算精
度和穩(wěn)定性較好,同時只有守恒格式才能處理計(jì)算空間中可能產(chǎn)生的間斷面,即流

場中可能產(chǎn)生的激波:③逆風(fēng)性:即空間導(dǎo)數(shù)離散時根據(jù)流速方向選擇差分方向,
從而使粘性的耗散是合理的:④無虛假振蕩:是指在計(jì)算間斷(或解梯度大)時要

河海大學(xué)博士論文

求格式具有足夠耗散,不出現(xiàn)虛假振蕩,當(dāng)然,無假振要求與計(jì)算精度的要求是相

互矛盾的,目前大多采用冊格式計(jì)算有大梯度和間斷的解;⑤高分辨率地捕捉間
斷:是指要求人工粘性和格式粘性的強(qiáng)度與解的光滑性相適應(yīng),在光滑區(qū)粘性自動

減少以提高精度,在間斷點(diǎn)附近或梯度大的區(qū)域又自動加強(qiáng),以防假振并達(dá)高分辨 率,目前一些優(yōu)良格式的間斷分辯率可達(dá)二到三個網(wǎng)格:⑥計(jì)算穩(wěn)定性:是指當(dāng)數(shù)
值解在受到微小的擾動后不發(fā)生突變現(xiàn)象,目前很難通過數(shù)學(xué)推導(dǎo)給出實(shí)用的穩(wěn)定 性條件,應(yīng)用中幾乎總是通過數(shù)值試驗(yàn)來確定時間步長,對時間步長起決定作用的

因素有時間離散格式、網(wǎng)格的規(guī)則性及邊界條件的處理等;⑦解的收斂性:數(shù)值格
式的準(zhǔn)確解與原微分方程的真解之差成為離散化誤差,當(dāng)時空步長趨于零時,離散 誤差亦趨于零,則稱格式具有收斂性,顯然相容性只保證數(shù)值格式收斂于原微分問

題,更重要的是數(shù)值解收斂于真解;⑨精度:數(shù)值解的精度取決于算法、網(wǎng)格、邊
界處理等諸多因素,除網(wǎng)格形狀要求規(guī)則外,還要求網(wǎng)格光滑,步長漸變,否則要

按照步長比進(jìn)行校正:⑨健全性及通用性:健全性是指能保證數(shù)值解不園所使用的 網(wǎng)格而異,例如可選擇一個簡單的明渠均勻流問題,所用格式在不同網(wǎng)格上的計(jì)算
結(jié)果必須一致,此外格式中不含專門為計(jì)算而引入的、需要人工加以調(diào)整的參數(shù)(如 人工粘性系數(shù)),目的是保證結(jié)果的客觀性,有文獻(xiàn)預(yù)言,隨著格式其它性能的完 善,未來十年內(nèi)對格式的主要要求為健全性;通用性則要求格式適用于各種流態(tài)(恒 定流/非恒定流,低速流/高速流,連續(xù)/間斷)、任意幾何形狀,多種網(wǎng)格布置、多 種邊界條件等盡可能多的情況。 能夠滿足上述大多數(shù)要求的格式可以稱作為高性能計(jì)算格式”…。 為了達(dá)到上述要求,首先必須了解各種離散格式的優(yōu)缺點(diǎn)及其改進(jìn)方法。下面 對微分方程目前主要的空間和時間離散方法進(jìn)行綜述。

(1)空間離散化方法: ①中心格式:中心格式的原理是泰勒級數(shù)展開式,其最大優(yōu)點(diǎn)是計(jì)算簡單、精度
高,有算例表明,在加粘后其計(jì)算量仍只是二階TFD格式的]/5,因而成為計(jì)算流 體動力學(xué)的代表方法之一,常用來計(jì)算光滑流動;缺點(diǎn)是計(jì)算解在計(jì)算空間中會產(chǎn)

生假性振蕩,同時穩(wěn)定性較差,目前實(shí)際計(jì)算中最有代表性的中心型數(shù)值方法是
Jameson““等首先提出的。 ②逆風(fēng)格式:逆風(fēng)格式是根據(jù)流速方向來選擇差分方向,優(yōu)點(diǎn)是穩(wěn)定性較好,如 二階單側(cè)逆風(fēng)格式允許的f為同樣點(diǎn)數(shù)的對稱格式的兩倍,對急流和出流邊界必須 或經(jīng)常使用逆風(fēng)格式,入流邊界使用單側(cè)逆風(fēng)格式是不可取的,但是,一階逆風(fēng)的 耗散往往太強(qiáng),精度低,二階線性逆風(fēng)格式在計(jì)算間斷解時常產(chǎn)生和中心格式相似

的振蕩,這就要求在經(jīng)典逆風(fēng)格式的基礎(chǔ)上研究非線性的單調(diào)性保持或TFD格式。

河海大學(xué)博士論文

文獻(xiàn)[32,33]以中心格式和逆風(fēng)格式為基礎(chǔ)建立了通量修正格式,該格式既防止假
性振蕩又可達(dá)到高階精度,且格式較為健全。

③TVD格式
TVD(Total Variation

Diminishing)的含義是總變差減小,對于微分方程

的數(shù)值解“(一)來說,其總變差對應(yīng)網(wǎng)格點(diǎn)的變差為:TV(u(?))=∑lu(x,)--U(X。)I
i=l

如果對于不同時間步長的近似解有:TV(u”、)≤TV(u“) 則稱總變差是減小的,即具有TVD性質(zhì)。

TVD格式的建立是以方程罷+罷:0的柯西問題為基本對象,在間斷形成前, 拼盤
解的全變差Zy保持常數(shù),以后r礦下降,同時在解的遞推過程中不產(chǎn)生新的極值,
原有極值也不增加,從而初值的單調(diào)性得以保持,據(jù)此可建立嚴(yán)格意義下的TVD格
式。

建立二階TVD格式的途徑有:修改通量(在粘性過大的~階逆風(fēng)格式中加入

受限制的反擴(kuò)散通量):通量限制(通過引入限制因子對反擴(kuò)散項(xiàng)加以限制以提高
精度的階數(shù)):變量插值(在變量插值時對坡度進(jìn)行限制)[18,30,34J。

④NND格式 NND格式即為無波動、無自由參數(shù)的高分辨率格式,由張涵信所構(gòu)造““,其 基本出發(fā)點(diǎn)是在Euler方程和N.s方程的右端附加一個三階導(dǎo)數(shù)項(xiàng),并使該項(xiàng)的系
數(shù)在激波上游為正,在激波下游為負(fù),以此來抑制差分解在激波上下游的非物理波 動,格式中不包含任何需要人工調(diào)整的參數(shù),并具有二階精度的TIeD性,對應(yīng)不同 的時間推進(jìn)形式,可以得到不同的顯式全離散格式NND—l到NND一5,并可構(gòu)造

出隱式時間推進(jìn)無條件穩(wěn)定的NND格式,其中NND一2到NND一5格式對時間均 為二階精度。因此NND格式以其準(zhǔn)確可靠,程序簡單,計(jì)算量小等特點(diǎn)在科學(xué)研 究及工程計(jì)算領(lǐng)域得到了廣泛的應(yīng)用并取得了良好的效果。 需要說明的是中心格式和逆風(fēng)格式可以相互轉(zhuǎn)化,中心格式添加人工粘性后可 解釋為逆風(fēng)格式,中心格式也可改寫為逆風(fēng)格式加上反擴(kuò)散項(xiàng),對后者加以限制即 構(gòu)成TVD格式。反之,逆風(fēng)格式也寫作中心格式加耗散項(xiàng),對后者加以限制也構(gòu) 成TVD格式。因此兩類格式的主要區(qū)別在于如何處理耗散,中心格式的人工粘性
和逆風(fēng)或TVD格式的粘性是相通的,故可以相互借鑒。

(1)時間離散化格式
時間坐標(biāo)軸與空間坐標(biāo)軸不同,它是單向的,只能是過去影響將來,時閶積分

河海大學(xué)博士論文

必須是采用相對r。的單側(cè)格式,且往往與空間離散化分別對待。時間積分主要分為 顯式和隱式兩類。顯式是指可以由已知的f。時的初值直接計(jì)算出t。時刻的未知解,

不必求解方程組。其優(yōu)點(diǎn)是每個時間步的計(jì)算工作量和存儲量均較小,且程序簡單, 缺點(diǎn)是時間步長受到穩(wěn)定性要求的限制,當(dāng)局部網(wǎng)格很密時,時間步長必須非常小
而使整個計(jì)算次數(shù)很大。隱式是指可同時由,。和,。時的解逼近空間導(dǎo)數(shù)和源項(xiàng),

故通常需要求解非線性微分方程組。其優(yōu)點(diǎn)是在線性分析時往往是無條件穩(wěn)定的,
計(jì)算時間步長可以取較大的值從而提高總體效率。缺點(diǎn)是需要對每一個時間步長求 解方程組,因而計(jì)算工作量和內(nèi)存占用量均很大,程序設(shè)計(jì)也比較復(fù)雜,主要用于 恒定流或解的時間變率小、流動時間歷程長、高雷諾數(shù)流等情況。

應(yīng)當(dāng)注意,時間離散格式的選擇對址的確定和計(jì)算量至關(guān)重要,址應(yīng)在保證 穩(wěn)定性的前提下由精度、單調(diào)性等綜合考慮確定,通常需要通過數(shù)值試驗(yàn)來最終確
定。

1.2.2.2線性方程組的解法 線性方程組的解法有兩類:直接法和迭代法。直接法可通過有限步得到準(zhǔn)確解, 求解精度高,但占有內(nèi)存大,僅適用于小型的方程組,而迭代法是從給定初值出發(fā),
通過某種算法反復(fù)加以校正,直至收斂。迭代法主要用于求解大型多維問題的不對

稱稀疏線性方程組,缺點(diǎn)是不一定能保證迭代收斂,且終止迭代的準(zhǔn)則不明確。
常用的迭代法是Jacobi法和Gauss.Seidel(GS)法,在GS方法中為加速收斂一般

均引入松弛因子田,稱為SOR法(逐次超松弛法),若m選取適當(dāng),可大大加速收 斂,∞的取值范圍為[0,2],若∞>l,則稱超松弛法,SOR法通常性能更好,因而
在計(jì)算流體力學(xué)中應(yīng)用廣泛。SOR法的一種變型為SSOR(對稱逐次超松弛法),

其收斂速度常遠(yuǎn)遠(yuǎn)超過SOR,SSOR的每步迭代由兩小步組成,一是向前SOR,一 是向后SOR,實(shí)際計(jì)算中需通過數(shù)值試驗(yàn)來選取甜。 1.213湍流模型技術(shù)
工程實(shí)際中存在的流動大部分都是湍流流動,湍流問題一直是計(jì)算流體力學(xué)中

的一個重要難題,湍流的特點(diǎn)是具有很強(qiáng)的隨機(jī)性、非恒定性,并且無一例外地都 是三維問題。由于湍流現(xiàn)象的復(fù)雜性,湍流運(yùn)動以及與之相聯(lián)系的熱和物質(zhì)的輸運(yùn)

現(xiàn)象都很難描述,從而也極難進(jìn)行理論預(yù)測。如果要通過計(jì)算來解決這一問題,就
不可避免地對湍流輸運(yùn)過程提出各種假設(shè),采用一些經(jīng)驗(yàn)性的結(jié)果和假定,把湍流 輸運(yùn)過程中的各種物理量與時均流場聯(lián)系起來,這就是湍流模型研究的內(nèi)容。所謂

湍流模型就是以雷諾時均運(yùn)動方程與脈動運(yùn)動方程為基礎(chǔ)。根據(jù)現(xiàn)有的理論和經(jīng)驗(yàn) 引進(jìn)一系列假設(shè)而建立起來的一組描述湍流平均量的封閉方程組。湍流模型的選取

河海大學(xué)博士論文

不僅要考慮其模擬復(fù)雜流動的可靠度和精確度,而且還需要確定計(jì)算所需的費(fèi)用和 處理問題所需的時間,后者對于工程應(yīng)用尤為重要。通常評價(jià)湍流模型優(yōu)劣的標(biāo)準(zhǔn) 是:能夠適用于較多種類型的水流現(xiàn)象;具有足夠的模擬精度;人力和計(jì)算機(jī)費(fèi)用 適度;復(fù)雜程度適當(dāng)。顯然,那些適用范圍大,精度高的湍流模型必然復(fù)雜,消耗

費(fèi)用較大,而好的湍流模型正是在這些相互矛盾的要求之間達(dá)到某種理想的平衡。 目前的湍流模型主要有雷諾時均模型、大渦模擬模型(Large Simulation-LESl和直接數(shù)值模擬湍流模型(Direct
Turbulence—DNSn。
Numerical
Eddy Simulation

雷諾時均模型:該模型在總體上可劃分為基于Boussinesq假定的紊動粘性模型

和基于雷諾應(yīng)力輸運(yùn)方程的二階矩封閉模型,前者又可分為混合長度模型、一方程
模型和二方程(如k—e,k—u,k-r等)模型,后者則有代數(shù)應(yīng)力模型(ASM)和輸
運(yùn)方程模型(DSM)兩大類。

Boussinesq假定的基本思想是:認(rèn)為雷諾應(yīng)力各向異性張量與平均場的應(yīng)變率

張量成正比,即湍流應(yīng)力的作用類似于層流應(yīng)力,但其比例系數(shù)為一流動變量一渦
粘性系數(shù),正是對渦粘性系數(shù)不同的計(jì)算方法導(dǎo)出了零方程、一方程、二方程等湍 流模型。不難想象,在渦粘性模型中,雷諾剪切應(yīng)力的方向與平均場的應(yīng)變率相同 這一特點(diǎn)可以說是渦粘性模型的主要缺陷之一,即雷諾應(yīng)力與平均場應(yīng)變率之間的

關(guān)系為線性關(guān)系,使得k-e模型無法再現(xiàn)流動中(如在方管中)出現(xiàn)的二次流“…, 這是由于復(fù)雜流動中,湍流場的變化總是滯后于平均場的變化。 在二階矩封閉模型中,雷諾應(yīng)力與平均速度之間的關(guān)系隱含在一組聯(lián)立的偏微
分方程中,它們之間的相互關(guān)系包括一些復(fù)雜的物理機(jī)制,同時在雷諾應(yīng)力方程中

的壓力一應(yīng)變率相關(guān)項(xiàng)為雷諾應(yīng)力輸運(yùn)過程中提供了快慢兩種不同的變化機(jī)理,在 一定程度上反映了湍流的非線性作用,然而通常形式的二階矩封閉模型對雷諾應(yīng)力 來說基本上是線性的,所包含的非線性作用非常小,還不能令人滿意地用來解決復(fù)
雜的工程問題。

上述各類雷諾時均湍流模型的特點(diǎn)及應(yīng)用可參見文獻(xiàn)[18,33,37,38,39,40],這 里就不再詳細(xì)展開。

大渦模擬模型(LES):通過不斷研究人們認(rèn)識到,在湍流流動中除存在許多隨機(jī)
性很強(qiáng)的小尺度渦運(yùn)動外,還有一些組織很好的大尺度渦結(jié)構(gòu),小尺度渦運(yùn)動受流
動邊界條件的影響很小,且近似是各相同性的,流動中大部分質(zhì)量、動量或能量的

輸運(yùn)主要來自于大渦運(yùn)動。因此,如果將比網(wǎng)格尺度大的大渦運(yùn)動通過數(shù)值求解 N.S方程來計(jì)算,而將比網(wǎng)格尺度小的小渦運(yùn)動利用湍流模型來模擬,使得通過模 型提供的計(jì)算在整個求解中只占很小的比重,降低最終結(jié)果對模型的依賴性,從而 達(dá)到提高計(jì)算精度的目的,這就是大渦模擬的基本思想。目前大渦模擬方面的工作

河海大學(xué)博士論文

開展不多,主要集中在氣象學(xué)的流動問題或基礎(chǔ)性研究及簡單的流動問題,如文獻(xiàn) [4l,42】利用大渦模型分別對方柱繞流和明渠流動進(jìn)行了數(shù)值模擬。關(guān)于LES更詳細(xì)
的介紹可參見文獻(xiàn)[38]

直接數(shù)值模擬模型(DNsT):模型化總是存在缺陷驅(qū)使人們尋求更好的解決
端流問題的途徑,由于N.S方程本身就是封閉的,不需要建立模型,由此提出了不 引入任何湍流模型,直接求解完整的三維非定常N.S方程,對湍流的瞬時運(yùn)動進(jìn)行 直接數(shù)值模擬,各種感興趣的物理量都可通過對瞬時量的平均運(yùn)算來取得,這就是 直接數(shù)值模擬湍流。該方法顯然有許多優(yōu)點(diǎn),如精確、可提供流場的全部信息,能

實(shí)現(xiàn)實(shí)驗(yàn)室中很難做到的流動條件控制且比實(shí)驗(yàn)更為經(jīng)濟(jì),但是直接數(shù)值模擬的計(jì) 算工作量卻非常巨大,由于湍流運(yùn)動所包含的單元(即小渦)比運(yùn)動區(qū)域的尺度要
小得多,典型的數(shù)量級是流動區(qū)域尺度的lO。倍,而數(shù)值計(jì)算所需要的網(wǎng)格必須比

湍動單元的尺度更小,在三維情況下,為了覆蓋流動區(qū)域至少需要lO。個網(wǎng)格點(diǎn),

如此多的網(wǎng)格點(diǎn)和相應(yīng)的計(jì)算次數(shù)是現(xiàn)代計(jì)算機(jī)的儲存能力和運(yùn)算速度所不能達(dá) 到的。目前國際上進(jìn)行的湍流直接數(shù)值模擬還僅限于較低雷諾數(shù)(Re=200)和簡單
幾何邊界條件的問題,主要用于湍流的基礎(chǔ)研究.如發(fā)現(xiàn)新結(jié)構(gòu),揭示新機(jī)理,提 供新概念,檢驗(yàn)和改進(jìn)湍流模型等。 顯然,無論是大渦模擬還是直接數(shù)值模擬,主要的困難不僅在于計(jì)算枧的限制,

還有方法本身的問題,限于目前計(jì)算機(jī)的發(fā)展水平,短期內(nèi)還不能夠利用上述兩種方 法對工程問題中的復(fù)雜流動進(jìn)行數(shù)值模擬,因此對現(xiàn)有的湍流模型加以改進(jìn),提出更
加合理的新的湍流模型依然是解決工程實(shí)踐中湍流問題的主要途徑。 重整化群(Renormalization Group。RNG)是一種用于構(gòu)筑許多物理現(xiàn)象模型的 通用方法,該方法起源于量子力學(xué)和高能物理中對基本粒子場的研究,自70年代后

期開始把RNG方法引入到湍流研究領(lǐng)域,Yakhot和0rszag…。于1986年應(yīng)用RNG 方法建立了第一個湍流模型,簡稱為RNGk-e。通過數(shù)值模擬表明它較之傳統(tǒng)湍流 建模方法具有顯著的優(yōu)越性和發(fā)展?jié)摿。自從基于RNG方法的湍流模型提出以后?br />許多學(xué)者應(yīng)用該模型進(jìn)行了數(shù)值模擬,Speziale G.G.”“和Yakhot““利用該模型求

解了后臺階繞流中的湍流分離現(xiàn)象,蔣莉““將RNG k—e湍流模型與壁面率相結(jié)合, 應(yīng)用于90。彎曲槽道湍流流動的數(shù)值模擬,結(jié)果表明,該湍流模式可以有效地模擬有 曲率影響的湍流流動,李玲”“利用該模型對鈍體繞流的尾流場進(jìn)行了數(shù)值模擬,結(jié) 果優(yōu)于標(biāo)準(zhǔn)的k—e湍流模型。由于利用RNG方法建立的湍流模型中不包含任何經(jīng)驗(yàn) 常數(shù)和可調(diào)節(jié)的參數(shù),其模型參數(shù)是利用RNG理論精確推導(dǎo)出來的,因而是通用的:
該模型適用于各種雷諾數(shù)范圍,包括層流、過渡流以及充分發(fā)展的湍流。此外,RNG

模型能夠較好地反映流動的各向異性,對于與時間相關(guān)的大尺度運(yùn)動也能給出真實(shí)的
模擬。因此RNG湍流模型在湍流計(jì)算中具有進(jìn)一步的推廣和應(yīng)用的價(jià)值。

河海大學(xué)博士論文

1.2.4自由面模擬技術(shù) 自由液面問題是計(jì)算流體力學(xué)(CFD)中非常重要的研究內(nèi)容,自由液面模擬的
合理性直接關(guān)系到含自由液面問題數(shù)值模擬成果的可靠性及精度。由于自由液面形狀 的復(fù)雜性和整體位置的不斷變化,特別是當(dāng)自由液面發(fā)生破碎麗使得該闖題成為強(qiáng)非 線性問題,所以精確計(jì)算和模擬自由液面的難度非常大,可以說,自由液面問題已經(jīng)

成為計(jì)算流體力學(xué)領(lǐng)域的重要分支和前沿課題。目前在該領(lǐng)域的研究處于領(lǐng)先地位的 是日本的東京大學(xué)和日本船研所(SRI),此外還有美國的LAWO大學(xué)和韓國的SRI,,
我國的702研究所、上海交通大學(xué)、大連理工大學(xué)等單位對自由面模擬技術(shù)開展了一 些卓有成效的研究工作,但由于問題的復(fù)雜性,上述結(jié)果還遠(yuǎn)遠(yuǎn)沒有達(dá)到令人滿意的 程度c48-49’5引。

下面對處理自由液面的一般方法及其應(yīng)用作簡單描述,重點(diǎn)介紹VOF方法。 剛蓋假定方法是將自由面看作是一個可移動的固體壁面,直接采用固體壁面的 “不可穿入”條件,僅適用大水體的宏觀運(yùn)動; 靜壓假定法”“則認(rèn)為流體中的壓強(qiáng)沿水深的分布服從靜壓規(guī)律,目前幾乎所有
的大體積水流數(shù)值模擬都以此假定為基礎(chǔ);

高度函數(shù)法”“是將自由面處的水深h定義為坐標(biāo)(x,y,z)和時間f的函數(shù),通過求 解高度函數(shù)h滿足的微分方程來確定自由液面的位置,該方法適用于非恒定自由表面
問題.但無法處理標(biāo)高函數(shù)為多值函數(shù)的情況,如射流、波浪破碎等問題。

流線迭代法”“是直接從Navier-Stokes方程出發(fā),迭代求解壓力分布及流線位置。 文獻(xiàn)[231采用此方法計(jì)算了溢流壩反弧曲線段的自由水面。 除上述幾種方法外,目前處理自由液面問題的典型方法是歐拉模型中的MAC法
和VOF法。MAC(The Marker and Cell Method)是設(shè)想在整個流場中均勻分布沒有

體積和質(zhì)量的微小顆粒一標(biāo)記點(diǎn),這些標(biāo)記點(diǎn)以當(dāng)?shù)亓黧w質(zhì)點(diǎn)的速度隨流體一起運(yùn) 動,顯然只要計(jì)算出標(biāo)記點(diǎn)在空間的位置即可確定自由表面的形狀和位置”‘…。該方 法的主要缺點(diǎn)是需要很大的計(jì)算機(jī)內(nèi)存和計(jì)算工作量,并且對于非均勻流場在網(wǎng)格中 會出現(xiàn)虛假的密度很高或很低的標(biāo)記點(diǎn),造成自由液面形狀的失真,而VOF方法則 克服了上述不足,是目前研究自由液面問題的理想方法。
VOF(Volume

ofFluid)”“方法是定義這樣一個體積分?jǐn)?shù)F,在有流體的點(diǎn)取值

為l,在無流體的點(diǎn)取值為0,這樣,在一個網(wǎng)格單元內(nèi),F的平均值則代表單元內(nèi) 流體所占的份額,若F的均值為1.則網(wǎng)格充滿流體,若F的均值為0,則為空網(wǎng)格,

若F的均值在O和l之閻,則網(wǎng)格包含自由液面,為自由面網(wǎng)格。由于VOF法追蹤
的是網(wǎng)格中的流體體積,而不是追蹤流體質(zhì)點(diǎn)的運(yùn)動,因而具有容易實(shí)現(xiàn)、計(jì)算量小

和精度高等優(yōu)點(diǎn),并且可以處理自由面折疊、自由面入水等強(qiáng)非線性悶題,成為目前


河海大學(xué)博士論文

自由面運(yùn)動模擬的主流方法。

在VOF方法中,關(guān)鍵技術(shù)是建立運(yùn)動界面的輸運(yùn)與重構(gòu)模型,如文獻(xiàn)[571研究了 運(yùn)動界面追蹤技術(shù)中的界面重構(gòu)算法,并討論了由此引起的數(shù)值誤差問題,文獻(xiàn)[581 以VOF理論為基礎(chǔ)建立了三維空間自由面輸運(yùn)及重構(gòu)的一階模型,文獻(xiàn)[591將有限體 積法和動網(wǎng)格技術(shù)應(yīng)用到運(yùn)動界面追蹤,自由面的位置和形狀利用運(yùn)動邊界條件和空
間守恒率確定等,總之70F方法中有許多問題有待進(jìn)一步深入研究。

I'3本論文的工作意義 本論文的主要工作是研究流場數(shù)值模擬中的關(guān)鍵技術(shù),并將之應(yīng)用于水利工程中 帶閘墩溢流壩過壩水流三維流場的水力計(jì)算。 21世紀(jì)人類面臨的主要問題是生態(tài)保護(hù)型可再生能源的安全供電、用于抵抗饑 荒、貧窮和疾病的優(yōu)質(zhì)而充足的供水以及可滿足人類基本需求的能源和食物,這些問 題的解決在很大程度上依靠設(shè)計(jì)開發(fā)新的水電工程來實(shí)現(xiàn),水電工程對世界經(jīng)濟(jì)繁榮 和可持續(xù)發(fā)展起著至關(guān)重要的作用。在我國隨著水電技術(shù)的迅速發(fā)展,高壩的建設(shè)日 益增多,其特點(diǎn)是水流落差大,流量大,功率大,多位于河谷狹窄地帶,帶來的一個 突出問題是泄洪消能、防沖及安全.可以說我國在該領(lǐng)域遇到問題的技術(shù)難度居世界
第一,因此需要更加重視大壩安全問題。

在泄水建筑物一溢流壩的水力設(shè)計(jì)中需要解決的問題不僅包括壩面下游反弧段 上的摻氣、振動、空蝕防蝕及霧化等問題,同時需要解決壩面上游壓力分布、水面線 的位置及形狀、泄洪流量、閘墩對流量及水面線形狀的影響等問題。目前主要的研究 手段是物理模擬(水工模型試驗(yàn))和數(shù)值模擬.物理模擬有著公認(rèn)的實(shí)際意義和科學(xué) 價(jià)值.它可以預(yù)演工程中非常復(fù)雜的三維水流現(xiàn)象,檢驗(yàn)設(shè)計(jì)方案的可靠性和準(zhǔn)確性, 但是它具有費(fèi)用高、周期長、無法準(zhǔn)確測取工程的某些湍流特征量和不便于進(jìn)行多方 案優(yōu)選分析等弱點(diǎn),因而現(xiàn)在物理模擬中已有相當(dāng)部分的工作為數(shù)值模擬所替代。數(shù) 值模擬便于進(jìn)行方案比較和優(yōu)選,同時研究高速水流、空化和空蝕現(xiàn)象已在相當(dāng)程度 上可借助于湍流的數(shù)值模擬來實(shí)現(xiàn)。 在水利工程中,基于水深方向的二維模型在一定程度上能反映出壩面的基本流動 規(guī)律并解決了水利工程中的一些實(shí)際問題,如溢流壩反弧曲線的設(shè)計(jì),因此在80—90 年代得到了廣泛的研究和應(yīng)用。許協(xié)慶利用有限元方法計(jì)算了泄洪洞進(jìn)口水流和溢流 壩反弧水流…。;顧兆勛等利用N—s方程及相應(yīng)的湍流模型模擬了二維非定常帶自由 水面的湍流流動!;王奇峰采用極坐標(biāo)系下的水流控制方程和k-i湍流模型模擬了 反弧水流的水力特性”1;魏文禮利用正交曲線坐標(biāo)系和曲率修正的k—e模型模擬了 二維溢流壩流動,并提出了一種抗蝕性能較優(yōu)的反弧曲線形式“…。但是這些成果還 不能滿足水利工程實(shí)踐的需要,因?yàn)樾顾こ讨械乃骶哂忻黠@的空間三維特性,二

河海大學(xué)博士論文

維模型無法揭示復(fù)雜流場內(nèi)部的空間特性,必須借助于三維數(shù)學(xué)模型。目前國外在水
利工程中河口航道的三維水流數(shù)值模擬方面開展了較多的研究工作,提出了許多新模

型和新算法.國內(nèi)學(xué)者對河口三維水動力學(xué)的研究也給予了很大的關(guān)注……“,開發(fā) 出一些有特色的三維潮流場計(jì)算軟件。與此相比,在泄水工程的水動力學(xué)研究方面, 利用三維湍流模型全面揭示整個流場水力特性的研究工作開展不多,其中馬福喜利用 修正的k-e湍流模型和流體體積分?jǐn)?shù)法模擬自由表面計(jì)算了重力壩及其下游消力池 中的水流流場和一個三維水躍流場…。,金忠青利用k—e模型數(shù)值模擬了帶有自由水 面的三維強(qiáng)紊動流場“…。對溢流壩壩型和平面布置(如:頂部溢流曲線、閘墩的形
狀、位置及尺寸等)的綜合三維數(shù)值模擬研究工作開展不足,而這些因素將直接關(guān)系

到工程的泄流能力、抗蝕體型的優(yōu)化乃至整個大壩的安全。 正是由于泄水工程的水力設(shè)計(jì)對水利大壩工程的泄洪效率及安全營運(yùn)起著至關(guān)
重要的作用,因此必須進(jìn)一步加強(qiáng)對泄水工程復(fù)雜三維流場特性的研究,其中包括建 立具有足夠精度的水動力學(xué)數(shù)學(xué)模型和提供高效、可靠的數(shù)值計(jì)算方法,在此基礎(chǔ)上 則可以在較多的設(shè)計(jì)方案和大壩體型范圍內(nèi)快速而準(zhǔn)確地給出整個流場的定量結(jié)果,

從而把大型水利工程的建設(shè)建立在更扎實(shí)的科學(xué)研究的基礎(chǔ)之上。 1.4本論文的研究內(nèi)容 1)對網(wǎng)格生成技術(shù)進(jìn)行深入研究,詳細(xì)給出了網(wǎng)格的代數(shù)生成方法:研究了邊界正 交化技術(shù),給出了一種簡單的網(wǎng)格邊界正交化方法;建立了非對接分區(qū)結(jié)構(gòu)網(wǎng)格生
成方法,包括分區(qū)交界面上的數(shù)據(jù)結(jié)構(gòu)形式和數(shù)值計(jì)算方法:

2)利用有限體積方法建立了非正交同位網(wǎng)格系統(tǒng)中的SIMPLE算法,詳細(xì)給出了運(yùn) 動方程的離散形式,數(shù)值計(jì)算采用二階精度的通量修正格式,并完成了傾斜方腔驅(qū) 動流、后臺階繞流、圓柱繞流等流場的數(shù)值模擬: 3)研究比較了各種湍流模型的特點(diǎn)及應(yīng)用范圍,建立了RNGk一占湍流模型,并以水 翼粘流場為例完成了與標(biāo)準(zhǔn)k一占模型的計(jì)算比較:

4)對運(yùn)動界面模擬技術(shù)進(jìn)行了比較深入的研究,針對VOF方法中關(guān)鍵技術(shù)一自由面
重構(gòu)方法詳細(xì)給出了利用傾斜直線段模擬運(yùn)動界面的方法,并對旋轉(zhuǎn)流場中的馬蹄
鐵形界面運(yùn)動、平移流場中圓形界面的運(yùn)動、剪切流場中圓形界面的運(yùn)動進(jìn)行模擬:

5)利用建立的數(shù)學(xué)模型對WES型溢流壩過壩水流二維流場進(jìn)行數(shù)值計(jì)算,給出了三 種不同壩頂永頭下的自由水面線位置以及五種壩頂水頭下的壩面壓力分布曲線,并 與實(shí)驗(yàn)結(jié)果進(jìn)行了比較。同時計(jì)算了流量系數(shù)及壩面流速分布。
6)研究了帶閘墩溢流壩三維過壩水流粘流場,針對閘墩的型式及壩面的布置計(jì)算了

閘墩對水面線形狀、壩面壓力、流量系數(shù)等設(shè)計(jì)參數(shù)的影響,并將不同墩型與布置
形式的計(jì)算結(jié)果進(jìn)行比較。



河海大學(xué)博士論文

第二章網(wǎng)格生成技術(shù)研究
在對流動問題進(jìn)行數(shù)值計(jì)算之前,首先要解決如何將流動區(qū)域離散化的問題,目 前已經(jīng)發(fā)展出許多種將物理區(qū)域進(jìn)行離散以生成計(jì)算網(wǎng)格的方法,這些方法統(tǒng)稱為網(wǎng)
格生成技術(shù)(grid
generation

techniques)。必須注意,網(wǎng)格生成技術(shù)在流場數(shù)值模擬過

程中占有非常重要的地位,網(wǎng)格生成質(zhì)量的好壞直接關(guān)系到計(jì)算結(jié)果的精度,因此需
要引起足夠的重視。

由于眾多的網(wǎng)格生成方法各有特點(diǎn),在第一章中己經(jīng)介紹,關(guān)于一般的網(wǎng)格生成
方法可參閱文獻(xiàn)[67,68,69,70,71,721,因此本章主要針對網(wǎng)格生成方法中的專項(xiàng)技術(shù)進(jìn)

行討論,這些討論有利于在數(shù)值計(jì)算之前幫助建立合理的計(jì)算網(wǎng)格系統(tǒng)。 2.1有限體積方法中處理不規(guī)則區(qū)域的網(wǎng)格生成方法 工程中大多數(shù)的流動問題都具有一個不規(guī)則的復(fù)雜幾何邊界,如水利工程中的溢
流壩、泄洪道、水輪機(jī)蝸殼等過水建筑物,船舶工程中的船體及附體,天然河道中的 水流等問題都包含了非常復(fù)雜的幾何邊界,為此需要研究針對不規(guī)則區(qū)域的網(wǎng)格生成
方法,常用的方法有:

2.1.1采用階梯形邊界逼近真實(shí)邊界 當(dāng)流場邊界為曲線或斜直線時,一種最簡單的邊界逼 近方法是采用階梯形直線段逼近,如圖2.1,這種方法的
特點(diǎn)是同一網(wǎng)格線上的網(wǎng)格點(diǎn)(或控制體)的個數(shù)是不等

,_

真實(shí)邊界



的,同時計(jì)算邊界為相互垂直的鋸齒狀粗糙表面,當(dāng)網(wǎng)格
較稀疏時會帶來較大的計(jì)算誤差,隨著網(wǎng)格在曲線邊界處
的不斷細(xì)化,這種誤差會逐漸減小。
/ /


,

由于這種網(wǎng)格構(gòu)造簡單,可以適應(yīng)任何形狀的邊界, 若同時配合網(wǎng)格的局部加密技術(shù),該方法仍然不失為一種
處理不規(guī)則邊界的方法。

,

t t

計(jì)算邊界
圖2.1階梯形邊界逼近曲線邊界

需要注意的是,當(dāng)采用這種方法逼近真實(shí)邊界時,如果把求解區(qū)域限制在計(jì)算邊
界之內(nèi),則對每一個具體問題都需要單獨(dú)開發(fā)計(jì)算程序,如果把計(jì)算區(qū)域擴(kuò)充為一個

規(guī)則區(qū)域,則通過特殊的處理可以針對不同的問題使用同一個規(guī)則區(qū)域內(nèi)的計(jì)算程 序,這是處理形狀不是特別復(fù)雜的計(jì)算區(qū)域的有效方法,稱之為區(qū)域擴(kuò)充法(Domain
extension

method),具體細(xì)節(jié)可參見文獻(xiàn)[73]。

河海大學(xué)博士論文

2.1.2采用曲線坐標(biāo)變換 在水利工程中經(jīng)常采用的一種處理不規(guī)則邊界的方法是引入廣義曲線坐標(biāo)進(jìn)行 曲線坐標(biāo)變換‘68t23,74,'75,76】,即將原來的物理空間(x,Y,z)變換到廣義曲線坐標(biāo)空間 (孝,刁,f),同時保證物理空間的幾何邊界與曲線坐標(biāo)系中的坐標(biāo)軸重合或平行,從而 將(x,Y,z)空間中不規(guī)則的物理區(qū)域變換到({,77,f)空間中形狀規(guī)則的計(jì)算區(qū)域,如圖 2.2所示,坐標(biāo)變換公式為
孝=毒(x,Y,z)’
r/=r/(x,Y,:),

f=f(x,Y,2)

因此如果在變換后的計(jì)算平面上進(jìn)行數(shù)值計(jì)算,則網(wǎng)格的生成就十分簡便。

物理空間

計(jì)算空間

圖2.2曲線坐標(biāo)變換

但是,由于數(shù)值計(jì)算是在變換后的廣義曲線坐標(biāo)系中進(jìn)行的,因此描述流體運(yùn)動 的控制方程也需要進(jìn)行相應(yīng)的變換,從而在方程中會出現(xiàn)一些附加項(xiàng),使原本十分復(fù)
雜的控制方程變得更加復(fù)雜,同時這些附加項(xiàng)的物理意義不明確,而且在對附加項(xiàng)的 導(dǎo)數(shù)進(jìn)行離散時增加了離散誤差,這是曲線坐標(biāo)變換方法的主要不足之處。

2.1.3邊界貼體非正交網(wǎng)格 在計(jì)算復(fù)雜幾何區(qū)域的流體運(yùn)動時,經(jīng)常使用邊界貼體非正交網(wǎng)格系統(tǒng)【3”,其中
網(wǎng)格形式可以是結(jié)構(gòu)網(wǎng)格、分區(qū)結(jié)構(gòu)網(wǎng)格,也可以是非結(jié)構(gòu)網(wǎng)格,這類網(wǎng)格的優(yōu)點(diǎn)是

可以適應(yīng)任意的幾何形狀,同時比正交曲線坐標(biāo)變換更容易實(shí)現(xiàn)。由于非正交網(wǎng)格的 網(wǎng)格線是沿邊界布置的,因而相對于階梯形逼近曲線邊界的情形更容易應(yīng)用邊界條 件,特別是當(dāng)使用分區(qū)結(jié)構(gòu)網(wǎng)格或非結(jié)構(gòu)網(wǎng)格時,還可以沿流動方向布置網(wǎng)格以及在 物理量變化劇烈的區(qū)域進(jìn)行網(wǎng)格加密,因此目前在眾多的商用軟件中均采用邊界貼體
非正交網(wǎng)格系統(tǒng)。

當(dāng)然,非正交網(wǎng)格也存在一些不足,如方程中由于網(wǎng)格非正交而產(chǎn)生的修正項(xiàng)增 加了方程求解的難度;非正交網(wǎng)格可能導(dǎo)致出現(xiàn)非物理解;此外變量在非正交網(wǎng)格系 統(tǒng)中的布置方式對求解精度和算法的有效性都將產(chǎn)生影響。

14

河海大學(xué)博士論文

2.1.4分區(qū)結(jié)構(gòu)網(wǎng)格
結(jié)構(gòu)化網(wǎng)格是指網(wǎng)格系統(tǒng)中節(jié)點(diǎn)排列有序,相鄰節(jié)點(diǎn)的關(guān)系明確;分區(qū)結(jié)構(gòu)網(wǎng)格 是指把一個復(fù)雜的計(jì)算區(qū)域劃分成若干個塊,每~塊內(nèi)均采用結(jié)構(gòu)化網(wǎng)格,但不同塊

內(nèi)的網(wǎng)格系統(tǒng)可以不同。分區(qū)結(jié)構(gòu)網(wǎng)格的優(yōu)點(diǎn)是能夠用于幾何形狀復(fù)雜的流場計(jì)算, 同時可以使用結(jié)構(gòu)網(wǎng)格各個系統(tǒng)中成熟的求解技術(shù),適合于并行計(jì)算。便于提高計(jì)算
效率[77,78 791。

2.1.5非結(jié)構(gòu)兩格 與結(jié)構(gòu)化網(wǎng)格不同,在非結(jié)構(gòu)化網(wǎng)格中節(jié)點(diǎn)的位置無法用一個固定的法則進(jìn)行有 序地定義,但是它可以適應(yīng)任何形狀的幾何邊界,若在非結(jié)構(gòu)網(wǎng)格系統(tǒng)中應(yīng)用有限差 分方法和有限體積方法,則可以使這兩種數(shù)值計(jì)算方法對不規(guī)則邊界的適應(yīng)性增強(qiáng)到 與有限元方法等同的程度,這是非結(jié)構(gòu)網(wǎng)格的主要優(yōu)點(diǎn)1171。非結(jié)構(gòu)網(wǎng)格的缺點(diǎn)是網(wǎng)格
生成的工作量大,離散方程的求解速度較慢。
2.2

網(wǎng)格生成的代數(shù)方法一雙邊界方法 在眾多的網(wǎng)格生成方法中,代數(shù)方法以其簡單、通用性強(qiáng)而被廣泛使用,若同時

配合網(wǎng)格的區(qū)域加密技術(shù)、邊界正交化技術(shù)和分區(qū)網(wǎng)格技術(shù),其適用的范圍和精度將

被大大提高,因此首先介紹一種具有通用意義的代數(shù)網(wǎng)格生成方法一雙邊界方法It,9,舯1。
設(shè)在物理平面上有一不規(guī)則區(qū)域ABCB,其中AB,CD為兩條不直接相聯(lián)的邊界, 首先選定這兩條邊界上的刁值,分別為rib,礬,于是這兩條邊界上的x,Y僅隨f變化,


k=x^(善).y。=Y.(舌)
X,=z,(考), Y。=Y,(孝)

其中:6,t分別表示底邊AB和頂邊CD。 雙邊界方法提供了~種在兩個邊界之問內(nèi)插的方法,從而確定內(nèi)部節(jié)點(diǎn)的位置, 具體公式如下:

x(善,J7)=工r(f)^(77)+xz(亭)也(J7)+!i!墮j;;三塑2h,(ri)+‘!叢j尹^4(叩)
彤,野)=Yl(嬲(椰+Yzg)^:(7)+魚生i導(dǎo)盟縞(叩)+魚塵粵業(yè)‰(彳)
dr/

d17

其中:

河海大學(xué)博士論文

(1)XI

Y-,x2,Y2為邊界1和邊界2上的離散點(diǎn)坐標(biāo):

(2)hi=2礦一3叩2+1,h2=一2r/3+31/2,也=r/3—2r/2+r/,h4=礦一礦

若將^(f=1,2,3,4)表達(dá)式中的叩(2了乞擊)用一維伸縮函數(shù)



剛叩+(I-P。币、th(Q(1r-r/))]

(3)掣---K1皓)粵瞎),掣:KI瞎)掣
cI亡

來代替,則可實(shí)現(xiàn)網(wǎng)格密度的控制,上式中的P,Q即為控制節(jié)點(diǎn)分布的參數(shù)。

dC

d,7

d‘

萼掣“:皓)掣,唑掣一礎(chǔ))訾
oC 0‘
OC

O‘

系數(shù)Kl,K2用于控制正交網(wǎng)格深入?yún)^(qū)域內(nèi)部的程度。
說明:

(1)插值公式中的后兩項(xiàng)導(dǎo)數(shù)項(xiàng)是為了保證區(qū)域邊界處的正交性。

(2)當(dāng)h.=2叩3—3叩2+1,h2=-2r/3+37/2時,可同時實(shí)現(xiàn)網(wǎng)格在兩個邊界處的加密
而當(dāng)h。=1一,7,h:=r/(線性情況)時,生成的網(wǎng)格比較均勻。 (3)當(dāng)P<1時,網(wǎng)格線向J=l的邊界處加密 當(dāng)P>I時,網(wǎng)格線向J=JM的邊界處加密 (4)Q用于控制善與刁的關(guān)系偏離線性的程度,對網(wǎng)格影響不大。

(5)若將紅表達(dá)式中的叩(2孟矗)用下式替代,也能夠?qū)崿F(xiàn)網(wǎng)格的加密

其中:成為大于l的常數(shù),若玎接近于1,則可使網(wǎng)格在叩=0和玎=l(即J=l
和J=JM)處加密。

2.3網(wǎng)格邊界的正交化技術(shù) 網(wǎng)格的正交性,特別是在物理區(qū)域邊界附近的正交性,對于整個流場的計(jì)算精度 是至關(guān)重要的。實(shí)現(xiàn)網(wǎng)格正交化的方法有許多種,比如:在微分方程網(wǎng)格生成方法中,
改善正交性的途徑是選取合適的控制函數(shù),根據(jù)網(wǎng)格的正交條件g.,=0給出了控制

函數(shù)及控制方程,并實(shí)現(xiàn)了網(wǎng)格在邊界上的近似正交”4”。:在網(wǎng)格生成過程中引入
Hermite插值,通過求解一個泊松方程和一個拉氏方程來實(shí)現(xiàn)網(wǎng)格的正交化”。。 在上述方法中,網(wǎng)格的正交性算法隱含在整個網(wǎng)格生成過程當(dāng)中,并且與網(wǎng)格的

河海大學(xué)博士論文

生成是同時實(shí)現(xiàn)的,此外由于控制函數(shù)是按正交性的要求而設(shè)定的,故無法同時實(shí)現(xiàn)
網(wǎng)格的加密。還有一種正交化的方法是通過移動邊界點(diǎn)來實(shí)現(xiàn)的,但是當(dāng)邊界曲線的 曲率較大時容易引起網(wǎng)格的變形“…。 文獻(xiàn)【6]介紹了一種網(wǎng)格邊界正交化的代數(shù)方法,本節(jié)以該方法為基礎(chǔ),對其進(jìn)

行了改進(jìn)和完善,給出了更具有一般性,容易實(shí)現(xiàn)且獨(dú)立于網(wǎng)格生成過程之外的代數(shù) 網(wǎng)格邊界正交化方法,通過控制參數(shù)的設(shè)定很容易控制網(wǎng)格正交化的程度和范圍。 2.3.1正交化算法
.、j

本算法的基本思想是給定邊界點(diǎn)的分
一——

布,將區(qū)域內(nèi)部網(wǎng)格節(jié)點(diǎn)向?qū)?yīng)邊界點(diǎn)的法 矢量方向移動,從而實(shí)現(xiàn)網(wǎng)格在邊界區(qū)域附
近的正交化,下面以J=1邊界為例說明該正 交性算法。如圖2.3,利用超限插值方法生
v、m

、、、j.^
。,.

一一。J-
(卜I,1)

,

、.蕊、~
,i

成的區(qū)域內(nèi)初始網(wǎng)格點(diǎn)為x:,Y:,首先求
出邊界曲線J:l在(i,1)點(diǎn)處的法向矢量 月,則可以得到初始網(wǎng)格點(diǎn)在n上的投形點(diǎn)
圖2.3網(wǎng)格在邊界處的正交化

《,Y:。將初始網(wǎng)格點(diǎn)向投影點(diǎn)移動,則
可以實(shí)現(xiàn)網(wǎng)格在邊界處的正交性,移動的程度可通過一個權(quán)重?cái)?shù)∞來描述。

x。=(1一%)x:+%x; YⅣ:(1一曲Ⅳ)y:+D口Y;
顯然,權(quán)重出,,應(yīng)滿足0≤∞,,≤l 甜。=0表示網(wǎng)格點(diǎn)沒有移動 ∞。=1表示網(wǎng)格點(diǎn)無阻尼地移動到對應(yīng)邊界點(diǎn)的法線上。

在算法的具體實(shí)現(xiàn)中需要解決兩個問題: 第一個問題是網(wǎng)格點(diǎn)在對應(yīng)邊界點(diǎn)法線上投影點(diǎn)的計(jì)算。通過解析幾何的方法推得投 影點(diǎn)的計(jì)算分式為:

x:=隊(duì),砘)+吉h+吲他+%,) y:=一÷x:+y。+古x。
其中I.為邊界線在點(diǎn)(i,1)處的切線斜率

第二個問題是權(quán)重甜。的給定,本文采用文獻(xiàn)[6]推薦的公式,并進(jìn)行改進(jìn),給出了更

河海丈學(xué)博士論文

一般形式的權(quán)重計(jì)算公式。

oa,,=exp{枷卜老)-2_l】卜,烈‰酬_。一p】。
其中:d,j=√(xg—z叫)2十(j,Ⅱ一y,{)2
P,Q稱為阻尼參數(shù),P值控制網(wǎng)格正交性在J方向上衰減程度,P值越大,則沿
J方向向區(qū)域內(nèi)部傳播的正交性衰減越快,Q值控制網(wǎng)格正交性沿I方向的衰減程度, Q值越大,則表明只有在邊界的中心區(qū)域網(wǎng)格是正交的。 在解決了上述兩個問題之后,即可實(shí)現(xiàn)網(wǎng)格的正交化。

2.3.2正交網(wǎng)格生成算例
為檢驗(yàn)正交算法的可靠性,給出了正交網(wǎng)格生成的幾個算例…,為便于比較,在

圖2.4中給出了利用雙邊界法,但未采用正交化技術(shù)生成的初始網(wǎng)格。從圖2.4、2。5 中可以看出,使用正交化技術(shù)能夠?qū)崿F(xiàn)網(wǎng)格在邊界處的完全正交,區(qū)域內(nèi)部正交化的 程度可根據(jù)情況進(jìn)行調(diào)整,周時網(wǎng)搐的加密與正交化可同時進(jìn)行。

鐾蓁 。。冀熏
量垂肇薹蘩蕊瀏
。一一 0 20 0.40

三三三至至量耋善摹?’膏誓◆:曩

審器至毒害耋薹輦j尊毫享三量 !鞭睫啡浪颍橹`:摹霉 蔓芝三三主三三妄;}■÷如i、‘‘:一.j

’∞一一…~——一一一一一…

~——

?一O¨一-一
1∞0.00

O瑚0.80

4∞0.40

~.一一
O∞



0∞

1∞

(1)

(2)

圖2.4網(wǎng)格生成算例1一未使用正交化技術(shù),2--fl!hqi交化技術(shù)

《三~
圖2.5正交網(wǎng)格生成算例


I∞

河海大學(xué)博士論文

2.4分區(qū)結(jié)構(gòu)網(wǎng)格技術(shù) 2.4.1分區(qū)結(jié)構(gòu)網(wǎng)格的一般特點(diǎn)
實(shí)際流動問題往往具有非常復(fù)雜的幾何物理區(qū)域,傳統(tǒng)的處理方法是通過曲線坐 標(biāo)變換將復(fù)雜的物理區(qū)域轉(zhuǎn)化為規(guī)則的計(jì)算區(qū)域,由此導(dǎo)致變換后的方程組更加復(fù)

雜,方程中各項(xiàng)的物理含義也模糊不清,計(jì)算速度和精度大大降低,目前流行的處理 方法是采用非結(jié)構(gòu)網(wǎng)格技術(shù)和分區(qū)結(jié)構(gòu)網(wǎng)格技術(shù),由于非結(jié)構(gòu)網(wǎng)搐的節(jié)點(diǎn)和單元分布 是任意的,因此具有非常優(yōu)越的幾何靈活性,可應(yīng)用于任意形狀流場的模擬,其主要 缺點(diǎn)是需要更多的計(jì)算機(jī)內(nèi)存和CPU時間,結(jié)構(gòu)網(wǎng)格中已經(jīng)成熟的算法在非結(jié)構(gòu)網(wǎng) 格中無法應(yīng)用,同時在復(fù)雜粘性流場中仍存在一些尚需解決的問題。分區(qū)結(jié)構(gòu)網(wǎng)格的 基本思想則是把幾何外形復(fù)雜的計(jì)算區(qū)域劃分為若干個盡可能規(guī)則的子區(qū)域,在每個 子區(qū)域中獨(dú)立生成網(wǎng)格并進(jìn)行流場計(jì)算,其最大的優(yōu)點(diǎn)是可以使用結(jié)構(gòu)網(wǎng)格中非常有 效、可靠的求解技術(shù),適合于并行計(jì)算,計(jì)算效率高,因此分區(qū)網(wǎng)格技術(shù)越來越受到
人們的重視。 分區(qū)網(wǎng)格系統(tǒng)一般分為相鄰子區(qū)有重疊部分的重疊網(wǎng)格(overlapping grids)18
2'馴

和相鄰子區(qū)無重疊網(wǎng)格的拼接網(wǎng)格【32】(patched grids),分區(qū)拼按網(wǎng)格系統(tǒng)又分為相鄰 予區(qū)網(wǎng)格線在交界面上對接和非對接兩種情況,非對接情況比對接情況要復(fù)雜,但其 應(yīng)用的空間和優(yōu)勢也更大,因?yàn)榭梢詫?shí)行逐塊加密而實(shí)現(xiàn)每一區(qū)的網(wǎng)格最優(yōu),同時允 許沿交界面移動網(wǎng)格,因此交界面處網(wǎng)格線非對接的分區(qū)結(jié)構(gòu)網(wǎng)格擁有非結(jié)構(gòu)網(wǎng)格的 幾何靈活性和局部加密能力。同時可以利用可靠的結(jié)構(gòu)網(wǎng)格求解技術(shù),使計(jì)算過程更
加簡單、高效.因而對于大多數(shù)流動問題都是非常有吸引力的。

2.4。2分區(qū)結(jié)構(gòu)網(wǎng)格及其數(shù)據(jù)結(jié)構(gòu) 分區(qū)結(jié)構(gòu)網(wǎng)格的優(yōu)點(diǎn)從圖2.6中可以看出,如果采用單區(qū)網(wǎng)格技術(shù),則不論是利 用代數(shù)生成方法還是解橢圓方程生成方法,都將會產(chǎn)生網(wǎng)格線在物面處非正交、網(wǎng)格 扭曲、尖點(diǎn)多、網(wǎng)格走向與流動方向不一致等缺點(diǎn),從而嚴(yán)重影響計(jì)算精度及效果,
若將流場按照其幾何外形和流動特點(diǎn)劃分為圖中所示的I,lI。IIl區(qū),則問題迎刃而
解。

為了保證在對分區(qū)結(jié)構(gòu)網(wǎng)格數(shù)值計(jì)算時能夠使用 單區(qū)結(jié)構(gòu)網(wǎng)格求解器,在分區(qū)結(jié)構(gòu)網(wǎng)格中需要采用兩 種數(shù)據(jù)結(jié)構(gòu),即在每一區(qū)中的規(guī)則數(shù)據(jù)結(jié)構(gòu)(f√)和
分區(qū)交界面的數(shù)據(jù)結(jié)構(gòu),典型的分區(qū)交界面如圖2.7。


在規(guī)則數(shù)據(jù)結(jié)構(gòu)中,整個流場的網(wǎng)格節(jié)點(diǎn)和計(jì)算節(jié)點(diǎn)
按逐塊和逐列的方法利用一維數(shù)組來標(biāo)記和存放。




河海大學(xué)博士論文

,=Ⅳ。k。‘+(f-1)N;+J

k=1,.,NB,i=l,..,Ⅳ?,_,=1,.,Nj

其中:NB~總的分區(qū)數(shù);Ⅳ?,川k一第k?yún)^(qū)網(wǎng)格在I?J方向的網(wǎng)格節(jié)點(diǎn)數(shù)
Ⅳ:=N?×N:一第k?yún)^(qū)網(wǎng)格的節(jié)點(diǎn)總數(shù)。
此外,在交界面數(shù)據(jù)結(jié)構(gòu)中,需要對交界 面左右兩側(cè)控制體中心以及交界面上所有網(wǎng)格 節(jié)點(diǎn)進(jìn)行序號標(biāo)記。由于交界面上兩側(cè)子區(qū)的 網(wǎng)格線是非對稱的,故同一個控制體積在交界
面上可能擁有幾個控制面,如圖2.7中典型的
r,tr冉 工 月 月 衛(wèi)

舡卉J

分區(qū)交界面,控制體積L在交界面上有3個控 制面.因此需要確認(rèn)交界面上每一個網(wǎng)格面所 對應(yīng)的左右兩側(cè)的控制體。關(guān)于分區(qū)結(jié)構(gòu)網(wǎng)格

的數(shù)值計(jì)算方法在第三章中單獨(dú)給出。
2.5小結(jié)

圖2.7分區(qū)網(wǎng)格交界面

本章以通用形式的代數(shù)網(wǎng)格生成方法為基礎(chǔ)研究了網(wǎng)格的邊界正交化技術(shù),提出 了一種簡單的、能夠獨(dú)立于網(wǎng)格生成技術(shù)之外的邊界正交化方法,給出了詳細(xì)的數(shù)學(xué)

模型并編制了計(jì)算程序:通過網(wǎng)格生成算例表明,該方法能夠大大改善網(wǎng)格在邊界處 的正交性。此外對分區(qū)結(jié)構(gòu)網(wǎng)格技術(shù)進(jìn)行深入的討論,并建立了非對接分區(qū)結(jié)構(gòu)網(wǎng)格
的數(shù)據(jù)結(jié)構(gòu)形式。

河海大學(xué)博士論文

第三章非正交同位網(wǎng)格系統(tǒng)中的SIMPLE算法
3.1交錯網(wǎng)格系統(tǒng)與同位網(wǎng)格系統(tǒng)
數(shù)值計(jì)算的基本思想就是將求解區(qū)域劃分成許多小的區(qū)域,即網(wǎng)格,在每一個網(wǎng) 格中定義求解變量并進(jìn)行計(jì)算,網(wǎng)格中計(jì)算變量所在的點(diǎn)稱之為計(jì)算點(diǎn),計(jì)算點(diǎn)的布 置方式有多種,可以在網(wǎng)格節(jié)點(diǎn),可以在網(wǎng)格中心點(diǎn),也可以在網(wǎng)格線中心點(diǎn),在圖 3.1中給出了三種情況計(jì)算節(jié)點(diǎn)對應(yīng)的控制體。

——一~一一…一1一。一一_’一
一一~~一+÷—¥^’'}~一一一J-2L2』‘“‘一—— ∽j臻:麓


j-:E, 7∥,■

_■一j—j.

1『一再,7

~M

7,≯謄n,≯ 。;:A j夠

幽3.I計(jì)算節(jié)點(diǎn)及對應(yīng)的控制體

(a)交錯網(wǎng)格系統(tǒng)

(b)同位網(wǎng)格系統(tǒng)

圖3.2非正交網(wǎng)格系統(tǒng)中變量的布置形式

如果求解變量定義在不同的計(jì)算節(jié)點(diǎn)且擁有不同的控制體,則稱為交錯布置,對 應(yīng)的網(wǎng)格系統(tǒng)稱為交錯網(wǎng)格(staggered grid),采用交錯網(wǎng)格系統(tǒng)的主要目的是解 決N—s方程中壓力梯度項(xiàng)及連續(xù)方程的離散困難,同時避免出現(xiàn)波動的壓力場和速度 場14”,但是在進(jìn)行編程計(jì)算時,需進(jìn)行各物理變量的定位及復(fù)雜的插值運(yùn)算。 如果求解域中所有的物理變量均定義在同一個計(jì)算點(diǎn)上且擁有相同的控制體,則 稱為同位布置,對應(yīng)的網(wǎng)格系統(tǒng)稱為同位網(wǎng)格(colocated grid)。典型的交錯網(wǎng)格
系統(tǒng)與同位網(wǎng)格系統(tǒng)見圖3.2。

顯然,只要能解決數(shù)值解的波動問題,同位網(wǎng)格在編程的簡易性和通用性方面更 具有優(yōu)勢,因?yàn)樗械淖兞繐碛型粋控制體,它們對應(yīng)離散方程中的許多項(xiàng)是~樣 的,計(jì)算量及存儲量均減少,但是由于網(wǎng)格面上沒有布置速度變量,在計(jì)算通過控制 面的通量時其速度需要由相鄰兩個計(jì)算節(jié)點(diǎn)上的速度插值得到,從而在程序中增加了

河海大學(xué)博士論文

更多的插值運(yùn)算,盡管如此,同位網(wǎng)格系統(tǒng)在應(yīng)用于復(fù)雜計(jì)算區(qū)域的非j_匕父嗍格和非

結(jié)構(gòu)網(wǎng)格中更具有交錯網(wǎng)格無法比擬的優(yōu)判2 51。
3.2

Navier-Stokes方程及其在非正交同位網(wǎng)格系統(tǒng)中的離散方法

3.2.1積分形式的Navier-Stokes方程
應(yīng)力形式的N.S方程的一般形式為【33,84l:

判+旦照趔:反一堡+魚
Ot ax.

…蘇.

(3.1)

ax.

為便于理解,寫出分量形式

掣ot塒?嘞=六一警+C等+等] 掣塒?刪=六一警+(魯+魯]


d.


@z,

d,J

@s,

掣塒?咖=正一警+[等+等]
其中:P,=一P+2∥警一姜腳v礦為粘性流體中的法向應(yīng)力; 戚


慨a,

。=∥(等+斟“一燁一Ⅲ粘性流體中的切向魷溈切向應(yīng)
力所在面的法線方向,J表示應(yīng)力本身的方向。

代入(3.1)式,并注意到研究對象是不可壓縮的牛頓流體:

掣+掣=∥一若+毒卜考]
西=一grad(pgh)。
其中:h=z一‰,zo為參考點(diǎn)的垂向坐標(biāo)。

c。.s,

其中:,為重力加速度g在三個坐標(biāo)軸方向上的分量。顯然,若p=constant,則

則P+觸=蘆
從而可以將質(zhì)量力項(xiàng)合并到壓力項(xiàng)中處理,為方便起見,在以下各式中,將聲仍 記作P,從而可以得到微分形式的N?S方程:

河海大學(xué)博士論文

——————卟一; 善十毒卜善j 融,I’良,j
8t

a(∥,),ab.“,)
瓠i

(3.6)

蘇。

將上述方程對控制體Q積分,并利用Gauss定理即可得到積分形式的N.S方程:

曇£艘地+肛,礦.蕊=一fPi,?贏+p}?蕊

(。。,)

其中:i(f-1,2,3)為三個坐標(biāo)軸方向上的單位矢量,等式左端第一項(xiàng)為隨時間變化的
非定常項(xiàng),第二項(xiàng)為對流項(xiàng),等式右端第一項(xiàng)為包含質(zhì)量力的壓力項(xiàng),第二項(xiàng)為粘性 力項(xiàng),又稱擴(kuò)散項(xiàng)。 積分形式的連續(xù)方程直接寫出:

杪?蕊=o
利用有限體積法求解流場的原型N-S方程即為(3,7)、(3.8)式。

(3.8)

需要說明的是,關(guān)于有限體積方法的特點(diǎn)在文獻(xiàn)f30,33,55】中有非常詳細(xì)的介紹,在 這里就不詳細(xì)展開,直接利用有限體積方法推導(dǎo)非正交同位網(wǎng)格系統(tǒng)中的離散模型。
3.2.2

Navier-Stokes方程的離散方法

3.2.2.1網(wǎng)格點(diǎn)的布置
本文采用非『E交同位網(wǎng)格系統(tǒng),計(jì)算點(diǎn)布置在網(wǎng)格中心,對應(yīng)的控制體即為網(wǎng)格
本身。

3.2.2.2動量方程的離散 以2D情況為例,推導(dǎo)非正交同位網(wǎng)格系統(tǒng)中動量方程離散形式,離散過程中涉
及到的計(jì)算節(jié)點(diǎn)、控制體、控制面等信息參見圖3.3
z.


~~

f~

●E



,i



一~~—鑫&蘭
?

;




圖3.3計(jì)算節(jié)點(diǎn)與控制體

河海大學(xué)博士論文

1)對流通量Ipu,V?rids
僅以控制面中的e面為例說明,其它控制面的形式完全相同,速度分量“,11--I殳 符號≯替代。同時為保證離散形式具有二階精度,并避免產(chǎn)生波動解,采用中心格式
與上風(fēng)格式聯(lián)合的方法,具體方法在3.4中介紹。

①顯式對流通量(中心差分格式)

k。)=…=L卻曠?ids=L矽?hds。丸zme'以
其中:

(3?9)

卅,=£p礦-贏=b礦?二k s。=pp。?以+s7丸)為通過e面的質(zhì)量通量;
s。,Sy為控制面積S。沿X,Y方向的投影面積: 以,≯,為速度在x,Y方向的分量:

以為控制面e中心處e點(diǎn)的速度,由線性插值確定,如圖3.4所示。

吮噸+(覿o。--。)+㈣(y。M)

锨埔㈣)+㈣。"(新一丑忙x。,) +[(詈]。丸+(考],o一九)]。。一
值的方法確定。

(3.10)

其中,相鄰計(jì)算節(jié)點(diǎn)連線與e面的交點(diǎn)e處的坐標(biāo)、遞廈及遠(yuǎn)廈梯度i電遼線住艷

x。.=xE,,l。+斗(1-2c) y,=y,五+¨(1一五。)

丸咆+罷∽--。p等帆誡.)
砘”辦(1--五e)+罷(Xe--Xe,)+詈(壙虬)

(罷]。=(罷]。以+(罷],‘(1吐)

c。?Ⅲ

河海大學(xué)博士論文

其中:

丸=I÷蘭i為線性插值矢量;計(jì)算節(jié)點(diǎn)上的速度梯度計(jì)算在后面單獨(dú)給出。

圖3.4控制面上的變量插值

②隱式對流通量(上風(fēng)差分格式)


k。尸=上廁蕊:卜辦…
Im!汀

礦(y?扦)。>0

礦艫.’nl<0
(3.12)

=min【m。,0 I毋E+max(m。,o)≯P
\ /

則總的對流通量采用下式計(jì)算

c。:k。尸+Ic)一一k‘嚴(yán)r

(3.13)

其中【】中的項(xiàng)稱為延遲修正項(xiàng)(defferedcorrection),放在方程右端的源項(xiàng)之中, 即為離散方程中的Q‘項(xiàng),并利用上次的迭代值計(jì)算。采用延遲修正的原因是一方面 可使迭代數(shù)值解具有二階精度,另一方面可以避免因單純上風(fēng)格式造成系數(shù)矩陣的非 對稱性,從而導(dǎo)致迭代算法的發(fā)散。

2)擴(kuò)散通量[∥挈i淼
。掀

仍以控制面e為例說明,動力粘性系數(shù)∥及速度”,分別用一股變量r和礦替代,
于是擴(kuò)散通量的一般形式為:

【Fgradq)?rids

①顯式擴(kuò)散通量

河海大學(xué)博士論文

忙“廣。=f,Fgrn彬-礎(chǔ)=(Fgradq6?;I‘s。

=L伊+-。妫妫蟆

(3 14)

上式中需要計(jì)算物理量在網(wǎng)格面上的梯度值(—蘭)。,一般近似采用中心差分格式,即 佩
直接由控制體中心處的導(dǎo)數(shù)插值,只是需要對網(wǎng)格的非正交性影響而引入修正。

昏*磬=鳊re re—c刪烈隔吐,
出, 廢,。
l—1 %一,尸|

其中當(dāng)網(wǎng)格在e面處正交時,上式第二項(xiàng)等于0。

②隱式擴(kuò)散通量

k。尸=£(咧螂肛L(豢)。遼

2等¨¨鄧。(刪∥(禹以) 一C¨so(a彘#axm+霧知”卜蹦刪刪“a隔rE-rp噸)(3
其中:f,2radd,)“有網(wǎng)格中一I)處的梯度線性插值,

15)

禹噸。與掣
則擴(kuò)散通量

(,。。一_y。)i+(】‰一石。)_7

扛乏面=

F:

(F—mra+k一廣’一k“嚴(yán)r
correction)

(3.16)

上式[】中的項(xiàng)稱為延遲修正(deffered

放在離散方程的源項(xiàng)之中,即

為∥項(xiàng),并利用上一次的迭代值計(jì)算。
3)源項(xiàng)計(jì)算
首先給出壓力源項(xiàng)的計(jì)算

Q’=一p?淼=一.fagradP?№

(3.17)

河海大學(xué)博士論文

對于“勰酣=一£駟=一(新△Q
具甲:Af2為控制體體積,任意形狀控制體體積的計(jì)算方法見3.4.4 對于v方程,情形類似。
非定常項(xiàng)的計(jì)算方法如下:

(曇脅堿=害@?+1--4U~一)
=4;“i;1一Q=.

煦恥等卅,=害㈣?叫)
此外,在對流項(xiàng)與擴(kuò)散項(xiàng)的計(jì)算中,分別出現(xiàn)了對流項(xiàng)修正和擴(kuò)散項(xiàng)修正,也作
為離散方程的源項(xiàng)進(jìn)行處理,即

Q‘=k‘xp 7一k。尸

g“=k哼…一k“嚴(yán)
以“方程中的e面為例給出具體的表達(dá)式:

餅=川。【丸無+辦(1一t)+X。--Xe,OX 一min(m。,O)屯+max(m。,O)紼

)+警(咒一以)】 刪


Q?2

L(罷s‘+警s’)~等e警缸。+芳匈。)。

最后給出由于網(wǎng)格非正交而引起的修正項(xiàng)

Qg=呱刪H囂高詞】[1蠅】
3.2.3離散形式的Navier-Stokes方程
將(g.13)、(3.16)式代入(3.7)式,得動量方程的離散形式:

4辦+∑鉑=g


(3.18)

其中:I=E,W,N,S為與計(jì)算點(diǎn)P相鄰的計(jì)算點(diǎn)

河海大學(xué)博士論文

。





鵬一k

一 鵬一‰ 一 叢‰







髓一k

AP=√;一Af—Aw~AⅣ一AJ

Q一:鐘+Q;+彰+Q;+Q8.該式右端分別為壓力項(xiàng)、對流項(xiàng)修正項(xiàng)、擴(kuò)散
項(xiàng)修正項(xiàng)、非定常項(xiàng)以及非正交網(wǎng)格修正項(xiàng)對源項(xiàng)的貢獻(xiàn)。 3.3非正交同位網(wǎng)格系統(tǒng)中的壓力校正方程及其離散
由于隱式離散方法對時間步長沒有嚴(yán)格的限制。故正如上節(jié)中介紹的那樣,動量

方程采用隱式方法進(jìn)行離散,為推導(dǎo)壓力校正方程,我們將動量方程中的壓力項(xiàng)從源
項(xiàng)中分離出來,以u方程為例進(jìn)行說明。

群砷+∑群“,=(Q一餅)+餅
,

(3.19)

若動量方程的第m次迭代值為“,,則

彳;“;’+∑鐘“?+=(Q。一餅)+鐘


(3.20)

”f=_—一+÷餅
(包一餅)一∑掣礦,



群一

(3.21)




記右端第一項(xiàng)為印同時咖一11駟--(孰垃
蠆亥JP
△Q,印、
(3.22)

河海大學(xué)博士論文

顯然,速度“:’滿足動量方程,但不滿足連續(xù)方程,為此需要對該速度進(jìn)行修正,
以滿足連續(xù)方程(3.8),在SIMPLE算法中,是通過修iE壓,J7場來實(shí)現(xiàn)的。

咖夥一等嗜Opm)
其中:

(3 23)

“;=”;++“;為滿足連續(xù)方程的速度,“j為校正速度。
P”=P+P’為相應(yīng)的壓力場,P’為校正壓強(qiáng)。

若將(3.23)式代入連續(xù)方程塑竺衛(wèi):o,可得壓力泊松方程

砉【等(iOpm m=‘TO(pHT")】,

(3.2a)

通過直接求解壓力泊松方程(3.24),即可得到滿足連續(xù)方程的速度場Ⅳ?,但此 時的速度場及壓力場并不滿足動量方程,因此需要將迭代進(jìn)行下去,上述的迭代過程,

即通過壓力場校正來滿足速度場,稱為外部迭代,而內(nèi)部迭代是指線性代數(shù)方程組的
內(nèi)部求解迭代,這種直接求解壓力泊松方程來確定壓強(qiáng)的方法稱為Projection
method。

在實(shí)際計(jì)算中最常用的一種方法是壓力校正法,(3.23)式減去(3.22)式可得
到校正壓力P’與校正速度”’之間的關(guān)系

啦即籌(參
值之后控制面e上的修正速度與修正壓強(qiáng)梯度依然保持上述關(guān)系,即

(3.25)

為了應(yīng)用連續(xù)方程(3.8),需要知道控制面上的速度修正和壓強(qiáng)修正,可以通過 插值的方法來確定。由于在修正速度的插值計(jì)算中包含了修正壓強(qiáng)梯度的插值,故插

咖∥一等《)。
∑A,ul
中,該項(xiàng)忽略不計(jì)。


(3 26)

式中q=一々,但由于在求解壓力校正方程之前,珥未知,故在SIMPLE算法
虬一萬‘玄’e 吩一等c拳
(3.27)

河海大學(xué)博士論文

同理

吩一等(答

(3 28)

將(3.27)、(3.28)式代入離散型連續(xù)方程m。+m。+m,+m。=0,并注意到由動量方 程第m次迭代解出的速度場并不一定滿足連續(xù)方程,因此

捕:+確:+rh:+廓:+△確:=0

其中△祝為迭代速度場對應(yīng)的質(zhì)量通量誤差。

即(psu’)。一(psu’)。+(psu’)。-(psu’),+Ath:=0
在非正交網(wǎng)格系統(tǒng)中,修正質(zhì)量通量計(jì)算采用下面的公式:

(3.29)

咖(∥卜一(pAc2s)。(刁1。(拳

=一器咖喊啊)-(倒p從¨)】 鼽(拳z南№h卜(gradp’H瑤刊】
(gradpl)。=(gradp‘)E丑+(gradp’)P(1~20)

(3 30)

(gra印’)。.(瑤一斥)=【(譬L)。丑。+(芒L),(1一九)】o。一x,)

“(嬰)。丑。+(孕),(1一五。)慨一yp) 們 卻
應(yīng)當(dāng)注意,式(3.30)的第二項(xiàng)是由于網(wǎng)格的非正交性而引入的,當(dāng)網(wǎng)格非正交 性較弱時,該項(xiàng)影響不大。但當(dāng)網(wǎng)格是強(qiáng)非正交的,則必須計(jì)入該項(xiàng)的影響,將另外

三個控制面上的質(zhì)量通量代入(3.29)式,并只進(jìn)行一次修正,可得非正交網(wǎng)格中的 壓力校正方程為:

爿;p;+∑彳,|D;=一△廊:+Q:

(3.31)

舯班一(等)。而1
30

河海大學(xué)博士論文

儼一(等)。瓦1麗
爿:=_(、pA雒DS)。百殺

班一c等,,赤
爿f=一∑鐘 a;=∑(gradp。);(不一‘),k=e,w,s,n,K=E,∥,S,N
△痢:=一(嘶:+嘲:+,fl:+嘲:)”1為前一次迭代速度場對應(yīng)的控制體質(zhì)量通量。
至此,給出了非正交同位網(wǎng)格系統(tǒng)中SIMPLE算法中的動量方程和壓力校正方程 的離散形式,下面將SIMPLE算法的計(jì)算步驟簡述如下:

(1)將n時刻的計(jì)算值“?,p”作為初始值計(jì)算n+l時刻的∥、p∥;
(2)求解動量方程得到“?‘(該速度場不滿足連續(xù)方程,要通過修正壓力場來修正速
度場)

(3)求解壓力校正方程得到P’

(4)修正速度及壓力,得到滿足連續(xù)方程的速度場和壓力場 “●=“?++“’,P”=P+P。 (5)重復(fù)(2).利用妒和P”作為對∥”,P”’的修正,計(jì)算直到所有的修正量足夠
小。

(6)推進(jìn)到下一個時間步長。
(2)~(5)稱為外部迭代。

3.4非正交同位網(wǎng)格系統(tǒng)中的專題討論 3.4.1邊界條件 3.4.1.1壁面邊界條件 壁面邊界條件通常采用非滑移邊界條件-(Dirichlet條件),即流體速度等于壁面速 度,在有限體積方法中一般利用法向應(yīng)力等于零的Neumann邊界條件,以二維情況
的南邊界為例進(jìn)行說明,見圖3.5

河海大學(xué)博士論文



—7一—77_,一,—!,——≥—_———7—7—7——7———,——.r /////7/,/一。,/,//7/ /////
圖3.5壁面邊界條件

“=2∥(罷)。=0 砂
擴(kuò)散通量的貢獻(xiàn)為

(3.32)

由于壁面邊界是不可穿透的,故對對流通量項(xiàng)沒有貢獻(xiàn),只有對擴(kuò)散通量的貢獻(xiàn)。對

吩I,T、ads=j磚㈣s,(孰
Yv—ys

=以S,盟 :盟0,一“。):4,一“,sS。t


(3.33)

其中A,6,S。6分別為邊界條件對離散方程中系數(shù)及源項(xiàng)的貢獻(xiàn)。 對于3D情況,與2D情況相類似,下面以底面(Bottom)為例.

對“方程:由于在底畫上r。=娑+_Ov:0,所以 洲o'x

‘4=fp,+k這 =I‘r,瓠

=∥。£(筆+警產(chǎn)。
Mj,摯,

M&警

河海大學(xué)博士論文

2警曠-以L‘--“s
對于v方程,類似可得:

(3 34)

f.。=警V,一警Vs
其中,底面對W方程中的擴(kuò)散通量沒有貢獻(xiàn)

(3.35)

3.4.1.2動量方程中的入13及出13邊界條件
由于入口邊界處的速度是已知的,故對流通量可直接計(jì)算

m,=pVns,=p吣,+"vs,J
如果采用隱式離散格式,則入口邊界對對流通量的貢獻(xiàn)為

(3.36)

F。=min(m,,oh。.+max(m,,o-,

(3.37)

其中:第一項(xiàng)為入口邊界對u方程源項(xiàng)的貢獻(xiàn),max■,,o)為入13邊界對離散方程中系
數(shù)一.的貢獻(xiàn)。入口邊界對擴(kuò)散通量的貢獻(xiàn)可類似于一般情形擴(kuò)散通量的計(jì)算,只是 速度梯度采用單向插值計(jì)算:
Fd

fFgradu‘贏=r(黔,=號¨訓(xùn)

㈣。s,

出口邊界條件對對流及擴(kuò)散通量的影響,其處理方法類似于入口邊界,只是出口

邊界要放在距離入口較遠(yuǎn)的下游,同時出口邊界上的速度采用一階上風(fēng)近似,即 廬f=辦,≯=“,V,‘w 3.4.1.3壓力校正方程中的入口及出口邊界條件
壓力校正方程中的邊界條件同樣需要特別注意。若通過某邊界的質(zhì)量通量是給定 的,則壓力校正方程中的質(zhì)量通量修正等于0,即采用零梯度的Neumann條件.

在出口,若入口質(zhì)量通量是給定的,對于定常運(yùn)動,則出口邊界上的速度可采用 外推法而獲得(即“,=",),但是需要對由外推法得到的出口速度進(jìn)行修正,以保證
各個流場內(nèi)的質(zhì)量守恒(即出口質(zhì)量通量與入口質(zhì)量通量相等),出口速度修正方法

見圖3.6。利用修正后的速度進(jìn)行下一次的外部迭代,同時出口邊界上的質(zhì)量通量修
正值為零。

為保證解的唯一性,需要選取一個固定點(diǎn)的壓強(qiáng)作為參考壓強(qiáng),所有其它點(diǎn)的壓 強(qiáng)修正計(jì)算都要減去同一個參考點(diǎn)的壓強(qiáng),但是需要注意,在分區(qū)結(jié)構(gòu)網(wǎng)格中,參考

點(diǎn)只能取一個,不能每區(qū)各取一個。

河海大學(xué)博士論文

圖3.6出口邊界的速度修正

3.4.2延遲修正技術(shù)(deffered correction)
在有限體積方法中.網(wǎng)格控制面上的變量值都需要通過插值來確定,若利用高階 插值的方法會導(dǎo)致計(jì)算量迅速增大,如對二維情況,利用四階辛普生(Simpson)公式, 則每個通量需要計(jì)算15個節(jié)點(diǎn)上的變量值,并且每個控制體對應(yīng)的代數(shù)方程包含25 個節(jié)點(diǎn)上的變量值,對于非線性方程組來說,其求解工作量是非常大的。

我們知道,利用前一步的迭代值顯式計(jì)算高階通量是簡單的,而隱式的低階通量 僅需要與計(jì)算節(jié)點(diǎn)最鄰近的節(jié)點(diǎn)上的變量值,計(jì)算工作量小。若能將二者結(jié)合起來,
則將同時兼顧計(jì)算的精度與簡單性,這種思想由Khosla
出【8卯,具體形式如下:
and

Rubin于1974年首先提

F二=F?9+tF,i—F?smpl

lold

t3.39)

其中:F州為某種低階隱式離散格式; r州為某種高階顯式離散格式; 式中第二項(xiàng)由前一步的迭代值計(jì)算,與第一項(xiàng)的隱式部分相比是很小的,因而顯 式處理不會影響計(jì)算的收斂性。這里將高階與低階離散格式聯(lián)合的目、的有兩個:一是
避免出現(xiàn)物理上不合理的波動解:二是改善某一迭代求解算法的收斂性[32,331。

河海大學(xué)博士論文

3.4.3控制體中心的梯度計(jì)算
在對流通量,擴(kuò)散通量及壓力源項(xiàng)計(jì)算中 計(jì)算方法如下,以擴(kuò)散通量為例:
Fed

都需要知道控制體中心P處的梯度

t莩(黔。J

(3.40)

其中交界面e處的梯度可以通過控制體積中心點(diǎn)上的梯度插值得到

(考]。=[考],屯+[等],o以,
控制體積中心處梯度可利用Gauss定理來確定

c。.t?,

E等m肼ha,s
如缸. 壚‘一‘
△Q△Q

z莩箬2面1融+丸S。'-#iw吖一∥)
3.4.4三維任意形狀網(wǎng)格控制體體積計(jì)算

(3.t2)

任意形狀的控制體積計(jì)算的原理是Gauss定理, 即將體積分轉(zhuǎn)換成面積分
由于

擊v(xi)=1,故控制體體積

△Q=肛=pVG彳妞=Ix7ds z∑tS。x

(3.43)

其中c為控制體擁有的控制面,疋。為控制體的面積矢量在x方向上的分量,即

控制面積S。在yz坐標(biāo)面上的投影面積

s。=s。。i+s。yj+墨。i
下面以空間六面體為例,c=l、2、3、4、5、6,工?扇。憧刂

面上四個頂點(diǎn)X坐標(biāo)的平均值計(jì)算,如圖3.7,S。1為第c個控
制面在yz坐標(biāo)面上的投影面積
圖3.7控制面投影

s。J=圭I舅一i)×E—iI=吉[(y.一y,)(Z2--Z4)一(2'I--Z3)(y:一y。)】
應(yīng)為負(fù)號。

(3.44)

將每個控制面上的Xc和S!(jì)算出來之后,利用(3.43)式即可計(jì)算出該控制體體積。 注意:若控制面外法線方向與X坐標(biāo)軸方向的夾角大于90。,則(3.43)式中該項(xiàng)前

河海大學(xué)博士論文

3.5分區(qū)結(jié)構(gòu)網(wǎng)格中SIMPLE算法的實(shí)現(xiàn)

動量方程的離散形式為:Ap辦+∑Ak≯k=Q尸,女=E,∥,N,S,在單區(qū)網(wǎng)格內(nèi)部,


每一個控制體擁有4個控制面,而交界面處的控制體則可能擁有多于4個的控制面, 故實(shí)際計(jì)算時每個網(wǎng)格面上的對流通量和擴(kuò)散通量需單獨(dú)處理,然后補(bǔ)充到交界面兩 側(cè)控制體的離散方程系數(shù)中去,如計(jì)算節(jié)點(diǎn)上對應(yīng)的e面為分界面,則取A,=0,在 離散方程的系數(shù)A,中補(bǔ)充交界面的貢獻(xiàn)一R,A。,,A馬,同理,若控制體R:的∥面

為分界面,取A,=O,在對應(yīng)的系數(shù)A,中補(bǔ)充交界面的貢獻(xiàn)A。。 在計(jì)算對流通量時,需要知道物理量在控制面上的值,對于正交網(wǎng)格可直接有網(wǎng)
格面兩側(cè)的計(jì)算點(diǎn)插值得到。但在一般情況下,交界面上的網(wǎng)格線是非正交的,如圖 3.8所示,即,與,’不重合,因此交界面上每一個網(wǎng)格面的變量不能由兩側(cè)計(jì)算節(jié)點(diǎn)直 接插值得到,而需要進(jìn)行修正,即

辦叫+(》“,一)
bJok^

?N




?

其中:

,

一一’一7L——’-

力=≯.-+≯。(1一五,-)

“?≮~
圖3.8

一…i!;i

學(xué)),=(缸”嘗M,1)
以上兩式均由網(wǎng)格面兩側(cè)計(jì)算節(jié)點(diǎn)線 性插值確定。
非對接分區(qū)交界面

A,=;}—魯為線性插值矢量,


h一■l

計(jì)算節(jié)點(diǎn)的梯度值采用下式計(jì)算

。耠挲



其中S:包括該控制體在分區(qū)交界面上的所有網(wǎng)格面。

上述格式可以保證交界面上的通量計(jì)算具有二階精度。 在計(jì)算交界面上擴(kuò)散通量的貢獻(xiàn)時,需要知道交界面上的梯度矢量,若采用對流 通量的處理格式,則將出現(xiàn)高階偏導(dǎo)數(shù)從而增加了計(jì)算工作量。常用的方法是直接采

河海大學(xué)博士論文

用中心差分格式,即由控制體中心的導(dǎo)數(shù)線性插值,同時對交界面處網(wǎng)格的非正交性 進(jìn)行修正。

c等。憧迹唬察o焉一c等,產(chǎn)c音!,,
其中當(dāng)網(wǎng)格在交界面處正交時右端的第二項(xiàng)等于0。

非正交網(wǎng)格中的壓力校正方程的離散形式為:爿;p:+∑彳,p,,i=一△;+鮮
其中需要注意的是對交界面上各個網(wǎng)格面上的質(zhì)量通量進(jìn)行修正,以,控制面為 例:

妒(棚=_(腳)『(扣挈,

一器(扣小礎(chǔ)如dp'M¨)]
其中:(!蘭),z=_。剑椋郏ǎ穑阂唬穑海┮唬ǎ纾颍幔洌穑,(乓一五)】 靠(rR一丘妒
(gradp’),=(gradp’)R^+(graclp’)£(1—4)

(gradp‘),.(不一五)=【冬)。丑+(冬)。(1一刪‰一吒) 盤 劣 +【冬)。^+(嬰)。(1-丑)艦一yL) 卯
cry

應(yīng)當(dāng)注意,上式的第二項(xiàng)是由于網(wǎng)格的非正交性而引入的,當(dāng)網(wǎng)格非正交性較弱時, 該項(xiàng)影響不大,但當(dāng)網(wǎng)格是強(qiáng)非正交的,則必須計(jì)入該項(xiàng)的影響。 3.6非正交同位網(wǎng)格系統(tǒng)中SIMPLE算法求解算例 3.6.1傾斜方腔頂板驅(qū)動流 頂板移動的方腔內(nèi)流場計(jì)算常被用來檢驗(yàn)數(shù)值算法的穩(wěn)定性和精度,為了便于比 較,本文首先計(jì)算了正交網(wǎng)格下方腔的內(nèi)流場,網(wǎng)格數(shù)為80×80,迭代過程中速度 及壓強(qiáng)的欠松弛矢量為0,8,0.8,0.3,給出了速度分布曲線,并與文獻(xiàn)[33】的結(jié)果進(jìn) 行比較,非常一致,見圖3.9。為了檢驗(yàn)非正交網(wǎng)格系統(tǒng)下的SIMPLE算法,對方腔 壁傾斜角為45。時的傾斜方腔內(nèi)流場進(jìn)行了數(shù)值模擬,網(wǎng)格數(shù)仍為80×80,見圖

河海大學(xué)博士論文

3.10,但迭代過程中速度及壓強(qiáng)的欠松弛矢量為0.7。0.7,0.2。在圖3.1l中給出了雷

諾數(shù)Re=100時的斜方腔中縱剖面及中橫剖面上的速度分布曲線,其中Re:墮,“為
,,

頂板移動速度,,為方腔長度,y為流體的運(yùn)動粘性系數(shù),圖3.12,3.13中分別給出 了Re=100和Re=1000時對應(yīng)的流線矢量圖和壓力等值線分布圖,可以看出,在斜方

腔右下方產(chǎn)生的次渦中心隨著雷諾數(shù)的增加向上移動。
’∞]
:—4

y州
o∞一

,.一一一


…jt忡mI





『『-
/.

。。礦\
4∞_/
、'、\

¨、j



。。j/
D∞ a∞
a∞


1∞


0∞ 口棚o∞
tM

圖3.9

二i一一一~.一7型”一…一一.一一Ⅲ


1∞

8--90。,Re=1000,方腔中縱斷面及中橫斷面上的速度分布曲線: 圖中離散點(diǎn)為文獻(xiàn)【25】中的計(jì)算結(jié)果

圖3.10傾斜方腔中的非正交網(wǎng)格

河海大學(xué)博士論文

一一㈠
圖3
,



l 1

B=450,Re=100,方腔中縱斷面及中橫斷面上的速度分布曲線

一…一三!型.
0帥 0∞


20

●∞

¨ ∞ ¨ ∞

o 0 1 1 5 2

(a)

Re=100

,

¨ % ¨ 艙

o 0 2

(b)

Re=1000

圖3.12 13=45。傾斜方腔內(nèi)的流線矢量曲線

河海大學(xué)博士論文

圖3 13

B-=45。傾斜方腔內(nèi)的流線矢量曲線

3.6.2后臺階繞流
后臺階繞流是典型的流動問題,如文獻(xiàn)[86,87]分別對三維后臺階的層流流動及后

臺階繞流的出流邊界條件進(jìn)行討論。為檢驗(yàn)分區(qū)結(jié)構(gòu)網(wǎng)格算法的可靠性,對計(jì)算流體 力學(xué)中的典型算例一后臺階繞流進(jìn)行分區(qū)計(jì)算,臺階的尺寸為入口高度5.2mm,臺階 高度4.9ram,臺階上游長度200ram,臺階下游長度500ram,計(jì)算結(jié)果給出了不同Re下

40

河海大學(xué)博士論文

的回流區(qū)長度的變化曲線,其中Re中的特征速度和特征長度分別為臺階上游的平均 速度及臺階高度,同時給出了Armaly
et

a1.的試驗(yàn)曲線【78】和文獻(xiàn)【79】的計(jì)算結(jié)果見圖

3.14,可以看出,在Re<400時,計(jì)算結(jié)果的一致性很好,而當(dāng)Re>400時與試驗(yàn)曲 線偏離較大,其主要原因是試驗(yàn)數(shù)據(jù)是在三維流動的情況下獲得的。

0——?-,-----r,----———————-———,——————---—-———-—-————————一,——








Re









圖3.14二二維后臺階回流區(qū)長度計(jì)算結(jié)果比較

3.6.3非對稱平板間圓柱繞流 當(dāng)研究位于兩個平板間的圓柱繞流時【E”,如圖3.15,重點(diǎn)關(guān)注的是柱面壓力與尾 渦,因而需要對柱面附近和尾部區(qū)域使用較密的網(wǎng)格,此時使用單區(qū)網(wǎng)格是不合適的,

要么整體利用非結(jié)構(gòu)網(wǎng)格,或者使用分區(qū)結(jié)構(gòu)網(wǎng)格.本文采用六分區(qū)網(wǎng)格對其進(jìn)行了
數(shù)值計(jì)算,計(jì)算條件為:圓柱直徑為D,圓柱中心至上下平板的距離分別為2.08D和 2,OD,圓柱上游距離為8D,下游距離為20D,來流速度為u,雷諾數(shù)Re=UD/,,, 這里僅給出Re=100和Re=1000對應(yīng)的壓力等值線,速度等值線和流線圖,見圖3.16。 從圖中可以清楚地看到圓柱后尾渦的形成和脫落過程,隨著Re的增加,流線的波動 逐漸加劇。

≈≯ 誓1_㈡-] 一 一 |;纛一
,■≈





圖3 15平板間圓柱繞流

4I

河海大學(xué)博士論文

圈3 16

Re=100(a)和Re=1000(b)對應(yīng)的壓力,速度等值線與流線圖

3.7小結(jié) 本章在非正交同位網(wǎng)格系統(tǒng)中推導(dǎo)出Navier-Stokes方程的離散模型、壓力校正方 程及其離散模型,對非正交同位網(wǎng)格系統(tǒng)中的邊界條件、延遲修正技術(shù)、控制體中心 的梯度以及任意形狀控制體體積計(jì)算等進(jìn)行專項(xiàng)討論,建立了非正交同位網(wǎng)格系統(tǒng)中 的SIMPLE算法,并以此為基礎(chǔ)建立了非對接分區(qū)結(jié)構(gòu)網(wǎng)格數(shù)值求解模型。利用上述 算法編制計(jì)算程序,對傾斜方腔驅(qū)動流、后臺階繞流、非對稱平板間圓柱繞流進(jìn)行計(jì) 算,通過與實(shí)驗(yàn)結(jié)果的比較表明,本文建立的數(shù)學(xué)模型是可靠的和正確的。

河海大學(xué)博士論文

第四章
4.1湍流研究的工程背景和意義

湍流模型

大多數(shù)工程中的流體力學(xué)問題都是湍流問題,應(yīng)當(dāng)用湍流理論進(jìn)行處理,但是由

于實(shí)際工程中的流動甚為復(fù)雜,目前的湍流理論和計(jì)算技術(shù)都未達(dá)到可以大量解決實(shí)
際工程中流體力學(xué)問題的水平,人們?yōu)榱藵M足工程應(yīng)用與設(shè)計(jì)方面的要求,往往顧及

不上湍流場的分辨率,而僅對湍流影響作簡單化處理。有些涉及到流體內(nèi)部結(jié)構(gòu)的工 程流體力學(xué)問題,不用湍流理論不可能較好地得到解決,如環(huán)境工程、水利工程中的
高速水流和挾沙水流、空氣動力學(xué)及船舶動力學(xué)中的繞流問題等都是如此,因而近

20年來,國內(nèi)外都在研究如何用湍流理論解決工程中的湍流問題。其主要內(nèi)容包括: ①利用目前比較成熟的湍流理論解釋一些特定的流動現(xiàn)象, ②應(yīng)用不同的湍流模型預(yù)測工程中的流體力學(xué)問題。, ③以湍流理論為依據(jù),提出一些解決工程問題的新理論和新方法。 上述研究內(nèi)容緊密聯(lián)系,相互促進(jìn),一方面湍流理論已被廣泛地用于解決工程問
題,同時在工程應(yīng)用中又促進(jìn)了湍流理論本身的發(fā)展,如直接服務(wù)于工程的湍流模式
理論等。

根據(jù)目前湍流研究的趨勢,除對湍流基本結(jié)構(gòu)大力進(jìn)行研究外,結(jié)合工程中的問 題進(jìn)行工程湍流研究,己成為國際湍流研究的重要方面,如湍流模式、湍流控制、湍
流測量、湍流與傳質(zhì)傳熱等。

因此,以工程為背景開展湍流研究,對于國民經(jīng)濟(jì)建設(shè)和科學(xué)技術(shù)的發(fā)展都具有 十分重要的意義,為此,必須以發(fā)展和利用已有的湍流知識解決實(shí)際技術(shù)問題為目標(biāo), 在工程中廣泛應(yīng)用湍流理論,為解決工程湍流問題提供實(shí)踐基礎(chǔ),同時積極促進(jìn)湍流
理論及其相關(guān)技術(shù)的發(fā)展,正如周培元教授指出的那樣需要“積極開展流體湍流運(yùn)動 的實(shí)驗(yàn)、理論和應(yīng)用研究”13 71


4.2湍流模式理論 由于湍流的復(fù)雜性,目前要想獲得復(fù)雜的計(jì)算結(jié)果幾乎是不太可能的,需要引入 各種假定,這些假定是對實(shí)際湍流流動的某種近似模擬。所謂湍流模式就是指將真實(shí) 的湍流運(yùn)動;癁槿藶樵O(shè)計(jì)的模型,這種模型的平均行為應(yīng)與實(shí)際湍流統(tǒng)計(jì)平均行為
基本一致。


雷諾時均模型在總體上可劃分為基于Boussincsq假定的湍動粘性模型和基于 雷諾應(yīng)力輸運(yùn)方程的二階矩封閉模型,前者又可分為混合長度模型、一方程模型和二 方程(如k.£,k.m,k-r等)模型,后者則有代數(shù)應(yīng)力模型(ASM)和輸運(yùn)方程模型(DSM)
兩大類。

河海大學(xué)博士論文

Bottssinesq假設(shè)的基本思想在于,認(rèn)為湍流應(yīng)力的作用類似于層流應(yīng)力

f=∥罷,但其比例系數(shù)為~流動變量——渦粘系數(shù)以,即f;∥.—d_u。其中^是
uy

8y

與,完全不同的量,它不是流體的一種性質(zhì),而是決定于湍流的流動狀態(tài),因而幾在 不同流動狀態(tài),甚至在同一流場中的不同點(diǎn),其值都很不相同。根據(jù)對渦粘系數(shù)一的 近似計(jì)算方法的不同,即可構(gòu)造出不同的湍流模型,如零方程,一方程、二方程等湍 流模式。需要注意,在渦粘性模式中,雷諾剪切應(yīng)力的方向與平均場的應(yīng)變率相同, 即渦粘系數(shù)是各向同性的標(biāo)量,這一特點(diǎn)是Boussinesq假設(shè)模式的主要缺陷之一。 與Boussinesq假設(shè)模式不同,在二階矩封閉模式中,雷諾應(yīng)力與平均速度之間的 關(guān)系隱含在一組聯(lián)立的偏微分方程中,在一定程度上反映了湍流的非線性作用。然而 通常形式的二階矩封閉模式對雷諾應(yīng)力來說基本上是線性的,所包含的非線性作用非 常小,還不能令人滿意地用來解決復(fù)雜工程問題,可以說非線性湍流模式在理論和應(yīng)
用上目前都已成為工程湍流計(jì)算研究中的核,Ii,問題之一【89】。

4.3湍流運(yùn)動的基本方程 在工程計(jì)算中,往往是流動參數(shù)的平均量對工程設(shè)計(jì)更具有重要作用。因此,流 體的運(yùn)動方程都建立在雷諾平均的概念下,即對N.S方程取某種形式的平均,如時間
平均,空間平均、系統(tǒng)綜合平均(指對相同初始條件和邊界條件下大量實(shí)驗(yàn)結(jié)果求平
均)。

對不可壓流體來說,平均場的動量方程(即雷諾方程RANS=Reynolds
Navier-Stokes

Averaged

Equations)的形式如下:
【q.1 J (t.1)
?。

魯五籌一羆+毒(,謄dx一瓦)壇 a。教k礎(chǔ):融|、.! —二+“I一一十一ly一一“.”}+P.

其中:“∽為雷諾應(yīng)力,是方程中的未知量,U.和JP分別為平均速度分量和壓力, 島為體積力,,,為流體的動力粘度。當(dāng)平均方式為時間平均時,-拋'-;7。:0。如果流動 是非定常的,則不能取是平均。
不可壓縮流體運(yùn)動的連續(xù)方程為:

墮:0
敏.

(4.2)

由于在方程中出現(xiàn)了雷諾應(yīng)力,因而方程(4.I)、(4.2)是不封閉的,即未知量

河海大學(xué)博士論文

的個數(shù)大于方程的個數(shù),因此求解控制方程(4.1)、(4.2)的首要問題是確定雷諾應(yīng)

力(一“,“,),在湍流模式理論中是通過引進(jìn)低階的關(guān)系或平均量來近似表示這些湍動 量的,由這些微分形式或代數(shù)關(guān)系表示的湍動量模擬關(guān)系式與控制方程(4.1)、(4.2)
共同組成封閉方程組。

需要注意的是,雷諾應(yīng)力與粘性應(yīng)力是不同的,粘性應(yīng)力是流體分子運(yùn)動與分子 相互作用的效應(yīng),而雷諾應(yīng)力則是流體宏觀的湍動效應(yīng),不是流體的物質(zhì)本性,是流

體運(yùn)動到一定階段后的產(chǎn)物,從力的觀點(diǎn)看,雷諾應(yīng)力是遷移慣性力(非線性項(xiàng))因 取時均而出現(xiàn)的,因此在實(shí)質(zhì)上是屬于慣性力的范疇,只有在統(tǒng)計(jì)意義下雷諾應(yīng)力才
會出現(xiàn)。 在雷諾平均框架下的流動控制方程中,湍流的全部信息可以被認(rèn)為是包含在二階 速度相關(guān)項(xiàng)“。“,之中,若“,“,可以被忽略不計(jì)時,則方程(1)回歸到層流的控制方 程。由此可見,層流可被看作是湍流的一個極限特例,在大部分工程湍流問題中,甜.“, 不僅不能忽視,而且它對流動混合,擴(kuò)散等起決定性作用。 在建立“弘的封閉模式之前,有必要討論雷諾應(yīng)力輸運(yùn)的控制方程,該方程的具 體推導(dǎo)過程可參見文獻(xiàn)[901,這里僅給出表達(dá)式:

警:導(dǎo)(,學(xué)一而一叢小絲西)
一—-_ 一—-_

_—_.

_2v當(dāng)當(dāng) 面 生% 麗 墮奶 卜 ‰p 贏一% 巫魄
+ + 0冀.O%.

(4.3)

該方程等式左邊為應(yīng)力的對流輸運(yùn)項(xiàng)(q),右邊的第一個括號中的項(xiàng)(吒)代表
湍流的擴(kuò)散特性,控制湍流的空間分布,由三個不同的物理機(jī)理組成:粘性擴(kuò)散、三 階速度相關(guān)和壓力與速度的相互作用。第二個括號表示湍流的生成機(jī)制(只),即湍
流在通常情況下由雷諾應(yīng)力和速度梯度的相互作用所產(chǎn)生。第三個括號為湍動壓力與

應(yīng)變率的相關(guān)項(xiàng),最后一項(xiàng)為應(yīng)力的耗散率(毛)。為了方便起見,(4.3)可寫作

q=吒+弓+m,一毛
通過對方程中的毛,中。和61i作適當(dāng)模化即可得到不同的湍流模式。
4.4

(4?4)

t一£兩方程模型。
對方程(4.4)中的f√進(jìn)行縮并,根據(jù)縮并原理及連續(xù)條件有中。=0,可以得到

湍動能≈=三瓦的控制方程:

河海人學(xué)博士論文

嬖:§y婺一一k'ut一土麗卜瓦婺一y挈晏 面2面咿面一 一石掣p“1,嘶葛吖菇瓦
項(xiàng)為湍動擴(kuò)散,右端第二項(xiàng)為湍動能發(fā)生項(xiàng),最后一項(xiàng)為湍動能耗散項(xiàng)。

(4.5)

其中:左端為動能變化率,右端第一項(xiàng)為擴(kuò)散項(xiàng),【】中內(nèi)第一項(xiàng)為分子擴(kuò)散,后兩

在k一£模型中,為了封閉湍流運(yùn)動方程,還需要增加耗散率占的控制方程:

坐dt=旦&l悟L罟磊一萬]_
2y墮(墮t墮蘭)一2掣笪型一 zH等矗一
‘融|、敏j a]c? axl瓠t。 。a)ct瓠i瓠l
口f

:隔2
(4.6)

式中宰為耗散能變化率,右端第一項(xiàng)為擴(kuò)散(包括紊流擴(kuò)散和分子擴(kuò)散)項(xiàng),第二
及第四項(xiàng)為產(chǎn)生能項(xiàng),第三及第五項(xiàng)為破毀能項(xiàng)

方程(4.1)一(4.6)共12個方程,若取Bi,M:“:,k,s共12個量作為未知量?
則構(gòu)成一個近似封閉的方程組,為了能夠求解,還需要將方程組中其余的未知量與 12個選定的量建立聯(lián)系,因此還需要對方程(4.4)、(4.5)、(4.6)進(jìn)行摸化,給出

季:水譬]掣+,期扣¨。㈢阿卻
‘(巴專妒),
cs個槲
(1個方程)
(4.7)

筆;=毒[c。(譬]善+y考]一占

一“



面i

警=言[c。晤)毒+y爿一ci(妻]雨≤‘:(守c,個方程,

韓只=_(而蕁+麗習(xí)~再等。

河海大學(xué)博士論文

Ck=0.225,Cl=2.2,C2=o.4,e=o.13,巳=1.45,e,=1.92。
這些常數(shù)都是由實(shí)驗(yàn)或特例的計(jì)算給定。 方程(4.1)、(4.2)、(4.7)構(gòu)成了可以求解的封閉方程組,需要指出的是,方

程(4.7)中出現(xiàn)的六個常數(shù)G,cl,c:,c。,t.,c。是在模化過程中采用近似的結(jié)果, 顯然在;胁捎貌煌募俣〞玫讲煌某(shù)。
在標(biāo)準(zhǔn)k一占模型中,雷諾應(yīng)力的計(jì)算并不采用它的摸化方程,而是直接采用 Boussinesq假設(shè),即雷諾應(yīng)力與平均場應(yīng)變率應(yīng)有如下關(guān)系:

一瓦=u譬+等)一;毛七
其中:渦粘系數(shù)_在k—s模式中取

(4.s)

v,=C。生

而湍流動能(膏=寺蠔“,)和湍能的耗散率s采用模化方程:

面Dk=毒((¨嘗)毒)+pt—s 嘗=言ccV+≯毒,+妻cc以一t:s,
其中仇為由于平均速度梯度引起的湍動能生成項(xiàng),既:一p麗挈

(4.9)

ca.?。,

方程(4.9)(4.10)與方程(4.1)(4.2)即構(gòu)成了標(biāo)準(zhǔn)k—s湍流模型,在該模型中 的常數(shù)取值如下:




c, 。.。9

C引 1.44

C;2 l,92

盯★

盯‘

l,O

1.3

k一£兩方程模型是目前應(yīng)用比較廣泛的一種湍流模型。在推演過程中,采 用了以下幾項(xiàng)基本處理:(1)用湍能七反映特征速度;(2)用湍能耗散率占反映特征 長度尺度;(3)引進(jìn)了■=C。k2/F的關(guān)系;(4)利用了Boussinesq假定進(jìn)行簡化。 因此,k一占模型具有以下優(yōu)點(diǎn):(1)通過求偏微分方程考慮湍流物理量的輸運(yùn)過程, 即通過求解偏微分方程確定湍動特征速度與平均場速度梯度的關(guān)系,而不是直接將兩 者聯(lián)系起來:(2)特征長度不是由經(jīng)驗(yàn)確定,而是以耗散尺度作為特征長度,并由求 解相應(yīng)的偏微分方程得到。由于脈動特征速度和特征長度是通過解相應(yīng)的微分方程求 得,因而k一占模型在一定程度上考慮了流動場中各點(diǎn)的湍能傳遞和流動的歷史作用。
47

河海大學(xué)博士論文

計(jì)算結(jié)果表明,它能比較好地用于某些復(fù)雜的流動,例如環(huán)流、渠道流、邊壁射流和 自由湍流,甚至某些復(fù)雜的三維流等。但是它也有局限性,主要有以下幾點(diǎn):(1)采 用Boussinesq假定,即采用了梯度型和u各向同性的概念,因而使k一占模型難以準(zhǔn)

確地模擬剪切層中平均流動方向的改變對湍流場的影響;(2)引入了一系列的經(jīng)驗(yàn)系 數(shù),而這些系數(shù)都是在一定試驗(yàn)條件下得來的,因而也限制了模型的使用范圍。 4.5RNG女一占模型 RNGk—s模型由Yakhot和Orszag于1986年應(yīng)用重整化群(Renormalization
Group)的方法導(dǎo)出,對于高雷諾數(shù)湍流,RNGk一占模型具有和標(biāo)準(zhǔn)k一£模型相同

的形式,只是在s方程中出現(xiàn)了一個附加生成項(xiàng),當(dāng)流動快速畸變時,這一項(xiàng)將顯著 增大。RNGk一占模型與標(biāo)準(zhǔn)k一占的主要不同之處在于它們的模型系數(shù)不同,標(biāo)準(zhǔn) k一£模型系數(shù)由經(jīng)驗(yàn)給定,而RNG^一£模型系數(shù)均由理論精確計(jì)算得出。此外
Yakhot和Orszag認(rèn)為RNG k一占模型不需要采用壁面率,可直接積分到壁面,而文

獻(xiàn)【46】作者在應(yīng)用中發(fā)現(xiàn),若采用RNGk—E模型直接積分到壁面,必須在壁面附近 劃分足夠細(xì)的網(wǎng)格,這必然增加計(jì)算工作量,尤其在三維工程湍流計(jì)算中巨大的網(wǎng)格 數(shù)量往往是難以接受的,因此類似于文獻(xiàn)[46]的作法,本論文采用考慮壁面率的
RNG女一占模型,并將RNGk—s湍流模型推廣應(yīng)用于三維湍流模擬。

RNGk—s模型與標(biāo)準(zhǔn)的k—s模型具有相似的形式,該模型中的兩個標(biāo)量方程
k方程與s方程形式如下:

iOk M婺:魯((y+互)婺)+A一占 ol
出I 盤I ok呶‘

(4.11)

?0動6'+u,毒=云((,+與O"。罷O'Xi)一矗+{(e.n—e:占)
程中的R項(xiàng)為RNGk一£模型與標(biāo)準(zhǔn)k~£的主要區(qū)別,R=

(4.12)

其中:^方程中P.為由于平均速度梯度引起的湍動能生成項(xiàng),Pk:一p麗挈,占方
。O'X.

—百礦T’
巴礦(1一qlq。)£2

平均應(yīng)變鄹㈣珊7=厝孝船均湖擱瞄湍流時間尺趔E,

河海大學(xué)博士論文

湍流粘性系數(shù)∥。=p?G。七二/5,叩。=4.38是印在剪切流中的典型值, 口=0.015。此外.與標(biāo)準(zhǔn)女一s模型的不同之處還有方程中的常數(shù)不是用經(jīng)驗(yàn)方法 確定的,而是采用RNG理論推導(dǎo)出的精確值,具體數(shù)值如下:


CF:O.845,口k=0.717 9,O-。=0?717 9,Ccl-1.42,c。2

1?68

與標(biāo)準(zhǔn)女一s模型相比,RNGk一£模型具有以下特點(diǎn):

①RN6k一£模型的占方程中有一附加項(xiàng)R,該項(xiàng)代表平均應(yīng)變率對£的影響,可以 改善迅交流中的模擬精度: ②在RNG≈一s模型中包含了旋渦對湍流的影響,提高了對旋渦流的模擬精度; ③方程中的常數(shù)不是用經(jīng)驗(yàn)方法確定的,而是利用RNG理論推導(dǎo)出的精確值; ④k一£模型為高雷諾數(shù)模型,而RNGk一£模型則同時適應(yīng)于低雷諾數(shù)的情形,
此時需要在近壁面處作特別處理

綜上所述,RNG七一f模型比標(biāo)準(zhǔn)七一占模型具有更高的精度和更廣的應(yīng)用領(lǐng)域。
4.6

k—s湍流模型中的邊界條件
對t一£方程求解需要給出合理的邊界條件,這些邊界條件的應(yīng)用對所有的方程

基本上是相類似的,只是對于固體壁面。問題可能變得比較困難。在壁面上取七=0,
j6

但£≠0,由娑=0替代,其中”為壁面的法向。


對于高見流動問題,由于粘性底層很薄,無法在壁面處采用非常細(xì)密的網(wǎng)格模 擬,解決的方法是采用壁面函數(shù)(walI Functions),該函數(shù)與速度分布的對數(shù)有關(guān)。湍
流邊界層的速度分布在圖4.1中給出,從圖中可以看出,在對數(shù)區(qū)域的速度分布為:
,,



“+:二土:二lnn++B
“.

(4。13)



∞礦俗


,

o 2



10

20

50

100礦

圖4.1湍流邊界層速度分布

河海大學(xué)博士論文

斯∞與壁面平行嗍溉“。黼應(yīng)力戤或稱為摩擦斌鏟層,
f。為壁面切應(yīng)力,K稱為卡門常數(shù)(Kaman constant)(r:OAl),B表示粘性底層厚
度的實(shí)驗(yàn)常數(shù)(平板邊界層中B 化。


5.2),”為網(wǎng)格點(diǎn)到壁面的距離,n+為n的無因次

n+:型


(4.14)

由于通常假定流體是局部平衡的,即湍動能的產(chǎn)生與耗散是近似相等的,因此,

“,=c:4√.i}

(4.15)

由(4。14)(4.15)兩式可以推導(dǎo)出壁面上第一個網(wǎng)格點(diǎn)的速度與壁面切應(yīng)力的關(guān)系

鏟倒;=斛4r壓赤
建立更加精確的方程。

(4.16)

其中:E=Po,由于緊靠壁面的控制體的一個控制面在壁面上,因此在動量方程中

對于平行于壁面的控制體,其壁面上的切應(yīng)力值是需要的,因此利用該邊界條件可以

當(dāng)采用壁面率邊界條件時,女在壁面上的擴(kuò)散通量取為o,即婺:o,關(guān)于耗
散率占的邊界條件可以根據(jù)局部平衡的假定,即壁面處產(chǎn)生與耗散近似相等推導(dǎo)出 來。在壁面區(qū),湍動能的產(chǎn)生項(xiàng)可由下式計(jì)算:

只一。魯
aH

(4.17)

同時該區(qū)內(nèi)切應(yīng)力近似為常數(shù),為了計(jì)算出控制體中心處的耗散率,通過對(1)式
的對數(shù)速度分布求導(dǎo),可得:

舀O,:÷:掣


n)8i他~'’jo J ’P槲, 肼p

該式與(4.16)式聯(lián)立即可計(jì)算出耗散率。
當(dāng)采用上述邊界條件時,近壁面處的控制體中是不需要求解占方程的,而是由下

河海大學(xué)博士論文

式直接計(jì)算

。P

,:竺豎
砌口 必須注意,上述邊界條件是使用只有當(dāng)靠近壁面的第一個網(wǎng)格點(diǎn)處于對數(shù)區(qū),即

月:>30時才有效。 此外,在入口邊界處,k,6通常未知,一般給膏賦一個小量值,如10一云2,而£的

值一般由式占。譬確定。
4.7算例:水翼流場模擬 利用非正交網(wǎng)格同位網(wǎng)格系統(tǒng)中的SIMPLE算法對NACA0015型二維水翼粘流 場進(jìn)行數(shù)值計(jì)算,湍流模型分別采用標(biāo)準(zhǔn)k一£模型和RNGk—e,雷諾數(shù)Re=1.97e~, 在圖4.2中給出了不同攻角下升力系數(shù)曲線,并與文獻(xiàn)[91]的實(shí)驗(yàn)結(jié)果進(jìn)行比較,在 小攻角下兩種湍流模型的計(jì)算結(jié)果與實(shí)驗(yàn)值均比較一致,RNGk-£模型略優(yōu)于標(biāo)準(zhǔn)
k一£模型。
1.4 1.2 甚1 0.8 0.6
0.4

0.2 0










1,14 g(dea)

1(1

16

18

20

22

圖4.2
4.8小結(jié)

NACA0015水翼升力系數(shù)一攻角變化曲線

本章討論了湍流模式理論研究的工程背景和意義,詳細(xì)給出了標(biāo)準(zhǔn)k一占模型和 RaNGk一£模型,比較了兩種模型的特點(diǎn)及應(yīng)用范圍,為更好地處理固壁邊界,采用 與壁面率相結(jié)合的RNGk一占模型模擬湍流流動,通過對NACA0015型水翼湍流場的 模擬表明,RNGk—s模型能夠更好地模擬出流場的湍流特性。

河海大學(xué)博士論文

第五章運(yùn)動界面數(shù)值模擬技術(shù)
5.1運(yùn)動邊界模擬技術(shù)簡介 自由面模擬技術(shù)在水利工程、船舶工程、海洋工程、機(jī)械工程等領(lǐng)域有著深刻的 現(xiàn)實(shí)意義和應(yīng)用價(jià)值,隨著Hirt等人提出用VOF方法模擬自由面追蹤技術(shù),目前已 經(jīng)發(fā)展和形成了許多實(shí)用軟件,如VOF.3D,NASA.VOF,SOLA.VOF等,同時推廣 到一般介質(zhì)面的運(yùn)動界面模擬,該理論已經(jīng)成功應(yīng)用于半導(dǎo)體工業(yè)、噴墨技術(shù)等現(xiàn)代
高新技術(shù)產(chǎn)業(yè)。據(jù)不完全統(tǒng)計(jì),1999年度的J.Comput.Phys.雜志涉及運(yùn)動界面追蹤問

題的文章不少于30篇,可見運(yùn)動邊界數(shù)值模擬的研究和分析方法受到越來越多的關(guān) 注,但是由于運(yùn)動界面的形狀和整體位置不斷變化,所以精確模擬的難度非常大。 運(yùn)動界面的模擬模型一般分為兩大類:歐拉模型和拉格朗日模型。歐拉模型的特 征是坐標(biāo)系或者固定或者按給定的方式移動以適應(yīng)不斷變化的計(jì)算域形狀,其主要優(yōu) 點(diǎn)是可以適應(yīng)內(nèi)部界面的大變形和大幅度移動而能保持一定的計(jì)算精度。利用歐拉模 型需要解決的主要問題是用什么形式的控制方程最方便,如何選取合適的數(shù)值離散方
法和最佳的內(nèi)界面跟蹤技術(shù)。拉格朗日模型的特征是坐標(biāo)隨流體雨運(yùn)動,但每個網(wǎng)格 單元內(nèi)包含的流體質(zhì)點(diǎn)保持不變。該模型的優(yōu)點(diǎn)是可以精確跟蹤內(nèi)界面并準(zhǔn)確描繪出 運(yùn)動界面的形狀,便于使用界面的運(yùn)動學(xué)和動力學(xué)邊界條件;拉格朗日模型還可以在

運(yùn)動界面附近進(jìn)行網(wǎng)格加密,以精確模擬運(yùn)動界面附近的流場內(nèi)部結(jié)構(gòu),同時精細(xì)網(wǎng) 格區(qū)可隨界面移動,大大提高了運(yùn)動界面的模擬精度。拉格朗日模型的主要不足是隨 著計(jì)算的進(jìn)行,網(wǎng)格會出現(xiàn)重疊、纏繞等現(xiàn)象,隨著網(wǎng)格越來越不規(guī)則,計(jì)算精度也
越來越低。 在歐拉模型中,模擬運(yùn)動邊界的典型方法有三種:即MAC法,VOF法,Level
Set

方法。MAC法的基本原理是在計(jì)算區(qū)域中包含流體的網(wǎng)格內(nèi)布設(shè)“標(biāo)記點(diǎn)”(Mal(er),

標(biāo)記點(diǎn)沒有質(zhì)量,不參與流場計(jì)算,只是用來標(biāo)記流體質(zhì)點(diǎn)的運(yùn)動軌跡,隨著流場計(jì) 算的進(jìn)行,可獲得每個標(biāo)記點(diǎn)的速度,并由速度確定出標(biāo)記點(diǎn)的最新位置,根據(jù)標(biāo)記 點(diǎn)的最新位置即可確定哪些網(wǎng)格單元被流體充滿,哪些網(wǎng)格單元被流空,從而確定運(yùn) 動界面或自由液面在下一時刻的位置;VOF法是Hirt和Nichols[921于1981年首先提

出,通過對主場的計(jì)算采用了Donor--Acceptor(施主一受主)的概念、Upwind特征
思想和補(bǔ)償效應(yīng),成功地模擬出潰壩和涌浪(brokendam,breakingbore)自由面,避 開了對于自由面問題通常采用的高度函數(shù)法和計(jì)算量十分浩大的Marker點(diǎn)的方案,
引入了一個流體體積函數(shù)方程,與流體運(yùn)動基本方程聯(lián)合求解,開辟了自由面計(jì)算的

新途徑。盡管當(dāng)初的格式、方法是比較粗糙的,但其基本思想至今一直受到重視,特 別是對運(yùn)動界面的重構(gòu)方法開創(chuàng)了一個新的研究方向:1988年,Osher[931等人提出了

河海大學(xué)博士論文

一種零等值面方法(Level Set方法),在許多復(fù)雜的界面追蹤問題中獲得成功,該方

法的主要思想是將運(yùn)動界面定義為一個函數(shù)的零等值面(線),即妒(x,,)=0,然后始 終保持它是零等值面,并且在界面附近保持單調(diào)。
5.2

VOF方法 VOF方法的基本思想是定義一個函數(shù)F,在有流體的點(diǎn)取值為1,在無流體的點(diǎn)

取值為0,這樣在一個單元內(nèi),F的平均值就可以代表單元內(nèi)流體所占的份額,如果 F的均值為1,則單元內(nèi)充滿流體,如果單元內(nèi)無流體,則F的值為零。而F值在0

和1之間則表示該單元必然包含自由液面,所以,VOF法提供了與MAC法類似的有 關(guān)內(nèi)部界面的信息。由于VOF方法只需對每一個單元提供一個內(nèi)存就可以了,不像 MAC法那樣需要對標(biāo)志點(diǎn)提供內(nèi)存,從而節(jié)省的計(jì)算空間,提高了計(jì)算效率。 5.2.1控制方程
流體體積分?jǐn)?shù)F隨時間的變化過程由下面的方程控制:

~OF+“笪+v堡:0
af



(5.1)



對于不可壓縮流體,控制方程的守恒形式為:

堡+O(F/2)+型:0




(5.2)



自由面模擬問題的控制方程一般由上面的流體體積函數(shù)方程與流體運(yùn)動基本方
程構(gòu)成。

5.2.2流體體積分?jǐn)?shù)F的輸運(yùn)模型 以Hirt和Nichols建立的施主一受主模型(donor—acceptor)為例說明F輸運(yùn)的基本思
想。

考慮網(wǎng)格的右側(cè)邊界,國時段內(nèi),通過該網(wǎng)格面的流體和空氣的體積通量 V=“田,其中”為網(wǎng)格面上的法向速度,其符號決定了網(wǎng)格是“施主網(wǎng)格”還是“受 主網(wǎng)格”,若“的符號為正,則網(wǎng)格面上游是施主網(wǎng)格,下游是受主網(wǎng)格。在一個時 間步長內(nèi)通過網(wǎng)格面的F的通量等于6F與網(wǎng)格面積的乘積。8F由下式計(jì)算:

b'F=min眈DM+CF,R衄D}

(5.3)

其中第一項(xiàng)表示輸運(yùn)的流體的量,第二項(xiàng)表示施主網(wǎng)格所擁有的流體的量。

CF=max{(1.0一只D)lyI一(1.o—FD)缸D,0.0}

(5.4)

河海大學(xué)博士論文

式中,下標(biāo)A代表受主網(wǎng)格,下標(biāo)D代表施主網(wǎng)格,下標(biāo)AD則表示可以是施

主網(wǎng)格,也可以是受主網(wǎng)格,根據(jù)邊界面上流動的方向確定,比如,同一個網(wǎng)格相對
于不同網(wǎng)格邊界的情況。(5.3)式中的Mill是為了防止通過界面的通量超過施主網(wǎng)格

所能提供的F的量,(5.4)式中的max是為了保證當(dāng)空氣的通量超過可能提供的量時,
需要補(bǔ)充流體通量,這個補(bǔ)充的量就是CF。

在圖5.1中給出了施主一受主模型通量輸運(yùn)的基本形式,當(dāng)流動通過一個垂直的單 元邊界面時,在(a)中給出了施主網(wǎng)格與受主網(wǎng)格的定義:當(dāng)下游網(wǎng)格為受主網(wǎng)格

即AD=D時,則凹=0,F的通量為F:Folvl,其中V=。叮,F值用于確定網(wǎng)格
邊界面通過流體的那部分?jǐn)嗝婷娣e,當(dāng)然要求川<Ax以保證施主網(wǎng)格內(nèi)的流體不會
流空,見(b):當(dāng)AD=A?xí)r,則施主網(wǎng)格內(nèi)的F值用來確定網(wǎng)格邊界面上通過流體的
斷面面積(c),此時施主網(wǎng)格內(nèi)的所有流體都將流失,即所有存儲在虛線和邊界面之 間的流體和空氣均進(jìn)入受主網(wǎng)格;在(d1中表明,當(dāng)AD=A,但通過網(wǎng)格面的流體通

量超過YoIVI時,則需要補(bǔ)充額外的流體以滿足需要的通量,其中虛線和邊界面之間
的額外流體就等于CF。

(b)

(c) 圖5.1
(a)

(d)

施主.受主模型中F通量的計(jì)算方法舉例

施主?受主模型網(wǎng)格的配置,(b)、(c)、(d)中雙陰影線區(qū)域表示實(shí)際傳輸?shù)模频牧?br />
事實(shí)上,使用施主網(wǎng)格中的F值還是受主網(wǎng)格中的F值來確定邊界面上通過流體 的斷面面積完全取決于自由液面的平均定向,這可分為兩種情況①當(dāng)自由液西基本上

河海大學(xué)博士論文

沿其法向移動時,AD=A,②當(dāng)自由液面基本上沿其切向移動時.AD=D。此外如
果受主網(wǎng)格是空的,或者施主網(wǎng)格的上游是空的。則一律有AD=A,也就是說,施 主網(wǎng)格必須首先補(bǔ)充流體,然后才能向下游的空單元輸運(yùn)流體。 用通量值5F乘以流體通過邊界的斷面面積就可以得到輸運(yùn)流體的量,從施主網(wǎng)

格中減去該通量值添加到受主網(wǎng)格,然后對所有的網(wǎng)格邊界重復(fù)上述步驟即可完成F
的輸運(yùn)計(jì)算。 在上述計(jì)算過程中可能會出現(xiàn)F值小于零和大于l的情況,因此必須對F的 計(jì)算值歸整到0和1之間,同時將由此引起的誤差進(jìn)行累計(jì),以保證總的流體體積守

恒。此外由于在數(shù)值計(jì)算中產(chǎn)生的舍入誤差使得F的值不可能正好是1或0,為此通 常設(shè)定一個很小的數(shù)如占,=10一,當(dāng)F小于卸時,網(wǎng)格是空的,當(dāng)F大于1一£時,
網(wǎng)格是滿的。

5。2。3運(yùn)動界面的重構(gòu)方法 Hirt和Nichols的貢獻(xiàn)主要體現(xiàn)在:~是作者采用的施主一受主(Donor-Acceptor) 方法給出了界面通量輸運(yùn)計(jì)算的初始格式,同時該方法不僅考慮本網(wǎng)格單元上的流體
體積分?jǐn)?shù),而且用到相鄰網(wǎng)格上的流體體積分?jǐn)?shù),從而提高了計(jì)算精度:作者的另一 個重大貢獻(xiàn)是提出了一種自由面的重構(gòu)方法,因此掌握運(yùn)動界面的重構(gòu)方法必須從

H—N模型出發(fā)。

一赫
.÷一1I歷黝嬲勉
圖5.2施主.受主模型中的自由面重構(gòu)方法

該模型重構(gòu)運(yùn)動界面的基本思想是:將流體自由面看作是局部的單值函數(shù) r(x)或X(y),計(jì)算中需要用到9個網(wǎng)格的信息,如圖5.2。首先計(jì)算f一1,f,f+l網(wǎng)格

列的r值和J一1,,,J+1網(wǎng)格列的z,值,估算出每個網(wǎng)格上自由面的斜率值譬,譬, 靠聊
然后根據(jù)流體體積分?jǐn)?shù)和斜率的大小確定網(wǎng)格(f,J)上的自由面位置和方向,即

河海大學(xué)博士論文

K=∑F,kSyI,

?_i一1,i,i+1

(5 5)

(氅,:o暨斗
、dr“

(5 6) 、…7

蘇…+2融,+蘇H

x,=∑r,k4v^,

,=產(chǎn)1,J,j+l

(5 7)

(箏,:善婆蘭業(yè) I一}
、咖。 ¥H+2。妫А海

理論上講,以上兩個導(dǎo)數(shù)應(yīng)該互為倒數(shù),但數(shù)值計(jì)算結(jié)果則不然,作者的意圖是比較
它們絕對值的大小,來決定自由面在該網(wǎng)格內(nèi)是“水平”的還是“垂直”的。即如果

=一

(5。,

l:,.6, 、 。

引爿
則定義自由面是水平的;否則自由面是垂直的,如圖5.2所示。

(5.9)

盡管H.N運(yùn)動界面重構(gòu)方法比較簡單,但它卻起到了拋磚引玉的開創(chuàng)性作用, 在這一思想的啟迪下,自由面的重構(gòu)技術(shù)得到了迅速的發(fā)展,比如:Chofinl94】引入

l,4和3/4單元以連接水平和垂直的自由面單元;Puckett[951用離散的線段來表示自由
面,通過相鄰網(wǎng)格的體積分?jǐn)?shù)預(yù)測自由面的法向;Hongfangwenl96JNN斜線段來代替 自由面,其中斜線段方程中的系數(shù)由自由面網(wǎng)格周圍8個網(wǎng)格的體積分?jǐn)?shù)確定:

Ashgriz和Poo的FLAIR方法p。】則利用每兩個水平相鄰網(wǎng)格的體積分?jǐn)?shù)求出斜線段來 模擬自由液面,但只適用于均勻網(wǎng)格情況。Rudmarll9”對上述兩類自由面重構(gòu)算法進(jìn)
行了比較計(jì)算指出,利用斜線段來逼近自由液面比其它方法精度更高,這些方法的基

本思想都是首先確定網(wǎng)格內(nèi)運(yùn)動界面的具體位置及其運(yùn)動方向,然后確定運(yùn)動界面在 下一個時刻的新的位置,同時更加注意運(yùn)動界面與網(wǎng)格線的交點(diǎn)的位置,因此得到的 運(yùn)動界面也更加精確,普遍采用的格式是用網(wǎng)格內(nèi)具有法向量的斜直線段Y=ax+b 來逼近,也有用二次曲線Y=O.X2+bx+c來構(gòu)造運(yùn)動界面1981,只是考慮得更加精細(xì)而
已。

本論文利用FLAIR方法重構(gòu)自由面,并對該方法進(jìn)行了推廣。
FLAIR自由面重構(gòu)方法的基本思想是在兩個相鄰的網(wǎng)格內(nèi)構(gòu)造一條傾斜的直線

段來近似該網(wǎng)格邊界的界面,然后計(jì)算透過該網(wǎng)格邊界的流體逶量,作為修正流體體
積分?jǐn)?shù)的數(shù)值通量,由于界面構(gòu)造僅涉及到兩個網(wǎng)格,故需要根據(jù)網(wǎng)格內(nèi)的流體體積

分?jǐn)?shù)劃分為許多種情況分別計(jì)算,各個網(wǎng)格內(nèi)的流體體積分?jǐn)?shù)可能有以下四種情況,
見圖5.3。

河海大學(xué)博士論文

(I)施主網(wǎng)格是滿網(wǎng)格,即瓦=1 (2)施主網(wǎng)格是空網(wǎng)格,即屹=0 (3)施主網(wǎng)格是半網(wǎng)格,受主網(wǎng)格也是半網(wǎng)格,即0<FD<1,0<只<1 (4)施主網(wǎng)格是半網(wǎng)格,受主網(wǎng)格滿網(wǎng)格或半網(wǎng)格,即0<%<1,只=0,只=0

(c)

i-,?T T二



jji彩
(e) (f)

∥%/

圖5.3相鄰網(wǎng)格間的關(guān)系

情況(1)、(2)比較簡單,情況(3)可以分為16種,經(jīng)過轉(zhuǎn)換可以歸結(jié)為4種,

見圖5.4,即£≥厶,且目標(biāo)流體在運(yùn)動界面以下。此時只要確定了直線的斜率及位
置即可計(jì)算出通過網(wǎng)格邊界的流體輸運(yùn)量,如對情況(a),設(shè)左邊是施主網(wǎng)格,則舀 時間內(nèi)輸運(yùn)到相鄰網(wǎng)格的體積量V=”.西,即圖中的斜線陰影區(qū)域的面積;如果右邊 是施主網(wǎng)格,則應(yīng)取右邊網(wǎng)格內(nèi)靠近網(wǎng)格邊界對應(yīng)區(qū)域面積,施主網(wǎng)格和受主網(wǎng)格由 網(wǎng)格邊界上的速度方向確定。

河海大學(xué)博士論文



]赫勿
(c, (d)

圖5.4

FLAIR模型運(yùn)動界面重構(gòu)方法

情況(4)最為復(fù)雜,以圖5.5(a)為例,施主網(wǎng)格(f√)是半網(wǎng)格,而受主網(wǎng)格(i+l,,)

是空網(wǎng)格,可以用網(wǎng)格(f,『)垂直方向的兩個相鄰網(wǎng)格分別按照第三種情況構(gòu)造界面,
然后求通過網(wǎng)格邊界流入受主網(wǎng)格的流體體積輸運(yùn)量,同樣可以歸結(jié)為圖5.5(b)所示 的四種情況:若受主網(wǎng)格是空網(wǎng)格,則取靠右邊的寬度為“.西的區(qū)域面積作為輸運(yùn)

量,如果受主網(wǎng)格是滿網(wǎng)格,則取左邊的寬度為It.毋的區(qū)域面積作為輸運(yùn)量,即保
證流體分布靠近滿網(wǎng)格一邊。

E瓢 簍貔 遴~絲擻
(a)

(b)

圖5.5運(yùn)動界面的四種重構(gòu)方法

下面以情況(3)為例給出FLAIR模型構(gòu)造運(yùn)動界面的方法,見圖5.4,并將該方法

推廣到非均勻網(wǎng)格,丘,兀分別表示左、右網(wǎng)格內(nèi)的體積分?jǐn)?shù),見圖5.6。

河海大學(xué)博士論文

Case 1

Case 3

Case 4

圖5.6
Casel

運(yùn)動界面結(jié)構(gòu)形式

不妨假設(shè)兩個相鄰網(wǎng)格中的體積分?jǐn)?shù)滿足:L≥以,自由面認(rèn)為是與網(wǎng)格邊界相交 的直線段,故自由面可由Y=刪+6來描述,其中系數(shù)口,b由已知的體積分?jǐn)?shù)Z,五確
定。因此可得到下面的方程組:

厶?血,緲=裊,(ax+6協(xié)=i1口(x}x2I)+b(xi-x,_I) 兀-Axi+l緲=璧+1(ax+b)dx=去n(毒l一#)+6(鼉+J—xi)
聯(lián)立求解,得到

。:!亟二型心
b=厶.緲+立二生L(厶一廠6).緲
^i+I—j卜I
case



類似于casel可列出如下方程

厶-Ax。Ay=去Ⅱ(z}x三1)+6(xf—Xi-I)
/j?缸,+l?Ay=:I口(x;一x;)+6(x6一Xi)
由于出現(xiàn)了新的未知量%,tt需奉t.充一個方程

河海大學(xué)博士論文

0=aXb+b

聯(lián)立上述三個方程,求解可得

口:—2(fn.a—y-b)
Zt+Xi一1

其中,。:Xi,,:壘』型
kxi

b=2l帆+如)一



缸?

case3

類似于case 2可列出方程組
Yj。IT-Xa+b

厶缸:緲=x。一xH)Ay+罷b; x。2)+60。一工。)

厶止Ⅲ每=曇cc磊,一x;)+b(x“
聯(lián)立求解可得:

口:—2[Aa—y-b]
b=,_【(2s+2t一1)一2吮一2(s—p)A

—z石i面刁‘面習(xí)F麗i西可了了】.緲
其中5:二L

血.

”—■疊一
。一kxf+l
LU.

x,+z』_1)

,:塾!!二蘭二塾=!
Xf+Xi一1

P:

!型二。。剑
2xJ+I—xf一3xf—

case4

Yj=ax日+b
Y,一1=c“6+b

xn+x6=2(fa?血i+兀kxf+1)+2xj

河海大學(xué)博士論文

^緲.缸M=曇僻一x;)+6(Xb--Xi)
聯(lián)立求解得:

口:坐
xd一工6

e~川一i△xi.xb
{。b 2——五『_一
} B+√B2—4AC

lXa=2(/:-Axf+.,j?血J+1)+2xf一1一x6

其中:爿=2yj-I+等
B=2(_‘?△xf十fb?△x“l+xl“+xt)Y』一1+Gf+2fb△xf+1)?卻

c=2眈?缸,+^-缸i+1+Xi+1)(XiYj-I+厶血州?緲)+等緲
容易驗(yàn)證,對于尺度為h的方形網(wǎng)格,文獻(xiàn)[6】中的結(jié)果即為本方法的特例。 顯然在確定了n,b之后,利用Y=似+b即可完成自由液面的重構(gòu)。

本文研究的運(yùn)動界面是自由液面,自由液面單元的壓強(qiáng)P。需要利用內(nèi)插(或外

p“2(1一_)p~+秘s
這里,77=d。/d,d。是自由液面單元和鄰近


的內(nèi)部單元中心間的距離,d是自由液面至內(nèi)部
鄰接單元中心的距離,用于內(nèi)插的鄰接單元的選 取標(biāo)準(zhǔn)是,聯(lián)接自由液面單元中心和所選取內(nèi)部 單元中心的聯(lián)線最靠近自由液面的法線,該單元 又稱內(nèi)插鄰接單元。

、平
/p¨

—、~

王i塵
p_

圖5 7自由面上的壓強(qiáng)邊界條件

河海大學(xué)博士論文

5.2.5運(yùn)動界面的寬度與等值線的給定方法

由于運(yùn)動界面的初始位置是由流體體積分?jǐn)?shù)定義的㈣,隨著流場的變化和時間的
推移,計(jì)算所得到的自由面輪廓是從流體體積分?jǐn)?shù)為0的高度變到1的高度,即計(jì)算 流場中得到的自由面為F∈[O,11的模糊界面,其寬度可以跨越幾個阿格。如圖5.8

圖5.8運(yùn)動界面的寬度

所謂等值線,就是根據(jù)不同問題的實(shí)際情況,流體體積分?jǐn)?shù)取[o,l】區(qū)間內(nèi)的某個 值,如令F=O 5的等值線作為自由面的位置,當(dāng)然該值不是固定的,有可能在其它 情況下。疲剑埃陡线m。這是采用單一定值等值線確定自由面的方法。此外,還可 以利用F值的百分比寬度來給出邊界面的形狀,如96%寬度的界面;顯然,如果有 某種方法可以進(jìn)一步將具有一定寬度的界面重構(gòu)為單值輪廓線,則得到的自由西位置 將具有更高的精度。 5.2.6數(shù)值穩(wěn)定性條件 (1)對流穩(wěn)定條件 為了正確的模擬自由液面或內(nèi)界面的變化,我們只考慮相鄰單元之間的通量,因 此,時間步長必須滿足下面的不等式,

趴幽協(xié)尚}
上式要求在計(jì)算網(wǎng)格確定之后,時I'uq步長西的選取必須使流體在該時段內(nèi)流過的 距離不能超也網(wǎng)格的尺寸,由于rain是對所有的單元取的,故在實(shí)際計(jì)算中的毋一 般取上式右端的l/4~1/3。

河;饘W(xué)博士論文

(2)擴(kuò)散穩(wěn)定性 對于粘性流體,要求在每一個時間步內(nèi),動量的擴(kuò)散同樣不超過一個網(wǎng)格單元。
線性穩(wěn)定分析給出如下的限制條件:

V藤2篇Ax



.2+△y.

5.3運(yùn)動邊界模擬算例 算例1:Zalesak模型 為檢驗(yàn)本文所建立的的自由面重構(gòu)算法的可靠性,對經(jīng)典的Zalesak界面模型在 平面旋轉(zhuǎn)流場中的運(yùn)動進(jìn)行了計(jì)算模擬,計(jì)算狀態(tài)為:計(jì)算區(qū)域【0,1]×【0,l】,網(wǎng)格
節(jié)點(diǎn)數(shù)200x 200.旋轉(zhuǎn)速度場為:u(x,Y)=—石(y—O.5)'v(x,力=x(x—O.5),初始位置

為一個下方有缺口的圓周;計(jì)算時間步長為母=O.0005。在圖5.9中分別給出了經(jīng)過 四個時間步長f=0.5s,f=ls,f=1.5s,f=2s對應(yīng)的二維界面圖形,同時給出了文獻(xiàn)[1001
中采用Hirt和Nichols模型的計(jì)算結(jié)果?梢钥闯,Hirt和Nichols方法給出的界面

比較粗糙,無法準(zhǔn)確地描述運(yùn)動后的界面形狀,而本文所給出的方法則能夠獲得非常
好的結(jié)果。

,=0.5s

f=1.0s

f=1.5s

t=2.0s

(a)本方法

t=0.5s

f=1.0s

t=1.5J

f=2.0s

(b)Hirt&Nichols方法 圖5.9

Zalesak模型的計(jì)算結(jié)果

算例2:均勻流場中圓形界面模型
計(jì)算區(qū)域及網(wǎng)格設(shè)置與算例l相同,圓的初始位置為:圓,11,(0.3,0‘3),半徑r=0.2

河海大學(xué)博士論文

均勻速度場為U=0.2,v=0.2,經(jīng)過四個時間步長圓的形狀及位置如圖5.10。相對于 初始圓周,面積的變化量分別為as=O.000154,0.000154,0.000153,0.000158,其中每一 步的面積計(jì)算是根據(jù)相應(yīng)的面積分?jǐn)?shù)相加而得到的。

f=0.5s

t=1.05

,=1.5s

f=2.Os

圖5.10均勻流場中圓的移動

算例3:剪切流場中圓形界面的運(yùn)動 平移和旋轉(zhuǎn)流不會使運(yùn)動界面的形狀發(fā)生變化,而在真實(shí)的流體問題中會發(fā)生界

面不斷變形的復(fù)雜物理現(xiàn)象,為了考察本文建立的算法對各種現(xiàn)實(shí)流場的模擬效果,
計(jì)算了剪切流場中圓形界面的變形運(yùn)動,剪切速度場為
u(x,y)=7rcos(x(x—Xo))sin(石(y—Yo)) V(x,y)=一,rsin(x(x—Xo))cos(x(y一),0))

計(jì)算區(qū)域?yàn)椤荆,1]×[O,l】,初始界面為圓心在(O.3,O.3)、半徑為0.2的圓周,見 圖5.11(a),計(jì)算結(jié)果給出了r_O 5s,忙ls,t=1.55,,=2s圓形界面的變化情況(圖 5.12),同時以f=2s的結(jié)果作為計(jì)算初值,將速度反號進(jìn)行反剪切2秒,結(jié)果在圖5.11
(b)中給出。





£丫∞,。3一
、、\ j?f

\j

(a)

(b)

圖5.1l

界面初始位置(a)及反剪切后的恢復(fù)位置(b)

河海大學(xué)博士論文

r=0.5s

r=1.Os

f=1.5s

f=2.0s

圖5.12剪切流場中界面的位置和形狀

5.4小結(jié)

本章討論了運(yùn)動界面的模擬技術(shù),重點(diǎn)研究了VOF方法中運(yùn)動界面的高精度重
構(gòu)算法,以FLAIR方法為基礎(chǔ),建立了可適應(yīng)于非均勻網(wǎng)格的自由面重構(gòu)模型,通 過對旋轉(zhuǎn)流場中的馬踢鐵形界面運(yùn)動、平移流場中圓形界面的運(yùn)動、剪切流場中圓形 界面的運(yùn)動進(jìn)行模擬表明,本章建立的運(yùn)動界面重構(gòu)模型相對于其它方法能夠更加準(zhǔn)
確地確定運(yùn)動界面的形狀和位置。

河海大學(xué)博士論文

第六章溢流壩過壩水流二維流場數(shù)值模擬
6.1溢流壩過壩水流的研究方法
溢流壩的水力設(shè)計(jì)對于整個水力工程建設(shè)具有重要意義,其中溢流壩的流量系 數(shù)、自由水面線的形狀和位置、壩面壓力及壩面流速等都是水力設(shè)計(jì)中的重要參數(shù), 也是選擇合理溢流壩型的重要依據(jù)。這些數(shù)據(jù)多年來主要通過水力學(xué)模型試驗(yàn)或由經(jīng) 驗(yàn)公式而獲得,一方面要耗費(fèi)大量的人力、物力和財(cái)力,同時受到測量條件的限制不

可能給出整個流場區(qū)域的全部信息。隨著計(jì)算技術(shù)的迅速發(fā)展,利用數(shù)值模擬方法求
解溢流壩過壩水流場的各種計(jì)算模型相繼產(chǎn)生,如有限差分方法、邊界元方法[60,107】、 有限元方法[103,104,105]、邊界擬臺曲線坐標(biāo)變換方法【‘吲等。
有限差分法簡易明了,易于理解,其主要不足是不能很好地處理不規(guī)則邊界,因

此在過壩水流數(shù)值模擬中的應(yīng)用受到了限制。

有限元法比有限差分法能更好地?cái)M合復(fù)雜邊界,因而在過壩水流的數(shù)值模擬中得
到廣泛而成功的應(yīng)用。其中最早將有限元方法應(yīng)用到過壩水流求解的是日本的池川昌 弘lI“J,并對自由面的迭代提出了兩種方法。該方法的基本思想是首先假設(shè)一個流量, 分別用緩流和急流的迭代方法進(jìn)行自由面迭代,如果上游緩流邊界能很好地與下游急 流邊界交匯并光滑地過渡到下游,則該流量即為正確流量:如果不能很好地交匯,則

需要重新假設(shè)流量進(jìn)行試算。然而通過計(jì)算,他認(rèn)為兩種流態(tài)邊界不能很好地匯合,
因而未給出最終成果。文獻(xiàn)[109]從可變區(qū)域的概念出發(fā),把決定自由水面形狀的坐 標(biāo)和各節(jié)點(diǎn)的流函數(shù)作為未知量形成方程組同時求解。文獻(xiàn)『1 lo]用有限元法計(jì)算了

過壩水流,提出了表面流帶的概念,將一維水流的臨界流量概念和水面線計(jì)算方法應(yīng)
用于二維,找出了過壩流量和自由表面之間的聯(lián)系,克服了尋找臨界點(diǎn)的困難,文獻(xiàn)

[111】分別給出了用有限元泛函極值加權(quán)余量法離散格式,并用此格式求解過壩水流問
題的方法,得到了與實(shí)驗(yàn)較一致的結(jié)果。文獻(xiàn)f112]提出了一種求解具有未知能頭的
過壩水流迭代方法,根據(jù)自由表面能量極值條件,建立了確定未知水頭的數(shù)學(xué)關(guān)系,

從而得以實(shí)現(xiàn)水頭與自由面同步迭代求解。這一方法既適用于高壩,也適用于低堰。

由此可見,用有限元法對過壩水流進(jìn)行數(shù)值計(jì)算已是碩果累累,占有很重要的位 置。用該方法對流量函數(shù)和自由表面位置的計(jì)算也是相當(dāng)成功的,所得結(jié)果能夠基本
滿足水利工程的設(shè)計(jì)要求。 邊界元法起步較晚,在近20年來才得到了較快的發(fā)展,最早應(yīng)用于勢流計(jì)算, 同時也是應(yīng)用最為成功的領(lǐng)域,但把邊界元法直接應(yīng)用于過壩水流計(jì)算目前還不多

見,文獻(xiàn)[113]求解了溢流壩頂自由出流問題,文獻(xiàn)[114]求解了水力學(xué)中的不定邊界
流場問題。

河海大學(xué)博士論文

邊界擬合曲線坐標(biāo)變換法在近年來被更多地應(yīng)用到過壩水流的數(shù)值模擬,該方法 是通過邊界擬合坐標(biāo)變換,將不規(guī)則的物理流動區(qū)域變換成規(guī)則的矩形計(jì)算區(qū)域,同 時將控制方程作相應(yīng)變換,然后采用比較簡單的有限差分法進(jìn)行求解。文獻(xiàn)[1061干0 用該方法計(jì)算了過壩水流問題,計(jì)算結(jié)果與實(shí)驗(yàn)吻合良好,文獻(xiàn)[23]N用曲線坐標(biāo)變 換對溢流壩反弧段水流進(jìn)行了模擬。 需要注意的是,關(guān)于過壩水流的計(jì)算主要還是以勢流理論為基礎(chǔ),本章利用前面 建立的非正交網(wǎng)格系統(tǒng)中的SIMPLE算法對溢流壩過壩水流二維粘流場進(jìn)行計(jì)算,采 用分區(qū)網(wǎng)格技術(shù)對計(jì)算區(qū)域進(jìn)行網(wǎng)格劃分,湍流模型采用RNGk一占模型,自由水面 線的形狀及位置采用VOF方法確定,分別對高水頭溢流壩的水面線形狀、壩面壓力、 流量系數(shù)及壩面流速分布等進(jìn)行計(jì)算。 6.2計(jì)算模型 目前世界上廣泛使用的壩面溢流曲線是由美國陸軍工程兵團(tuán)水道實(shí)驗(yàn)站(WES) 提供的WES曲線‘“”,該型曲線具有工程量小、無有害的負(fù)壓及泄流能力大的特點(diǎn), 一般表達(dá)式為:

日Y-c(≯
其中:爿一不包括行近流速水頭的溢流壩上的水頭高度 c,6一與上游壩面傾斜坡度有關(guān)的系數(shù),
x,Y一以溢流壩最高點(diǎn)為原點(diǎn)的坐標(biāo)值。 在坐標(biāo)原點(diǎn)上游,壩面曲線一般由兩個或三個圓弧線組成,其中三圓弧曲線的性能優(yōu) 于二圓弧曲線。

本文選用的壩面曲線為三圓。祝牛忧X‘”=2.0川”Y,見圖6.1,6.2,其中
上游壩面坐標(biāo)在表6.1中給出,溢流壩高H。=76.2米,設(shè)計(jì)水頭H。=9.14米,計(jì)算 區(qū)域與網(wǎng)格劃分見圖6.3。

67

河海大學(xué)博士論文

圖6.1

三圓弧溢流壩上游壩面

圖6.2三圓弧溢流壩壩面曲線 表6.1
x l Ha
.0.O

三圓弧WES壩面上游坐標(biāo)(坐標(biāo)原點(diǎn)取在壩頂中心線上)
.0.05 .0.10 .O.15 .O.175 .0.20 .O.22

Y|H d

0.O

0.0025

O.010l

0.023

0.0316

0.043

0.0553

x f Hd

.O.24

.O.26

.0 276

.0.278

.O.28

.0.2818

Y|H d

0 0714

0 0926

0.11 53

0.1190

0.124l

0.1360

河海大學(xué)博士論文

翻6.3計(jì)算區(qū)域與網(wǎng)格劃分

從圖6-3中可以看出,整個流場分為三個區(qū)域進(jìn)行網(wǎng)格劃分,網(wǎng)格類型均為結(jié)構(gòu) 化網(wǎng)格,其中I區(qū)網(wǎng)格數(shù)為7500,II區(qū)網(wǎng)格數(shù)為4000,為保證壩面壓力的計(jì)算精度, 對III區(qū)網(wǎng)格進(jìn)行加密,網(wǎng)格數(shù)為8000。 入13邊界速度按設(shè)定流量計(jì)算給出,出口邊界條件取物理量沿水流方向梯度為 零:所有的氣體邊界都定義為壓力邊界,即壓力為0,湍動動能k和耗散率e按經(jīng)驗(yàn) 公式確定,固體邊界速度按無滑移條件給定,并引入壁面函數(shù)技術(shù),自由面的位置利 用VOF方法確定:時間步長△f=O.002,迭代次數(shù)20000。 6-3自由水面線計(jì)算 水面線位置的計(jì)算按表6.2進(jìn)行,水面線形狀在圖6.4中給出,與試驗(yàn)結(jié)果的比 較在圖6.j中給出,可以看出,計(jì)算出的水面線位置與試驗(yàn)值非常一致。
表5.2水面線計(jì)算
H。/HJ
(初值)

0.50

1.O

1.33

1.5

Q{Qd =(H。/HJ)1



0.329877

1.O

1.578202

1.913137

(計(jì)算值)
H。/Hd
L計(jì)算值)

80,85

85.38

88.22

89.67

0.5087

1.0044

1.313

1.478

河海大學(xué)博士論文

圖6.4不同壩頂水頭下的水面線(紅色為水)

80

75

70

60 lO 0 10 20 30

X(111) (a)He/HF0.5

河海大學(xué)博士論文

90

85

80




75

70

65

60
0 10

20

30

X(m)

(b)He/Ha=i.0
90

85

80

占75


70

65

60 0 10 20

X(m)

圖6.5

水面曲線的形狀

6.4水位一流量關(guān)系
文獻(xiàn)[97】根據(jù)實(shí)驗(yàn)結(jié)果,利用圖解法給出溢流壩的水位一流量關(guān)系曲線計(jì)算公式

若=(魯)l

6,該實(shí)驗(yàn)曲線可以作為數(shù)值計(jì)算結(jié)果的驗(yàn)證比較。

首先計(jì)算設(shè)計(jì)流量Q=CLH。3”米3,秒,其中C=o.5524瓦?m為設(shè)計(jì)流量系數(shù),
其值根據(jù)實(shí)驗(yàn)回歸曲線確定。L為有效壩寬,二維情況下取值為1。 為了確定溢流壩的水位與流量的關(guān)系,給定8種狀態(tài)下的不同流量,計(jì)算出對應(yīng)



河海大學(xué)博士論文

的壩頂水頭高度,繪制出水位一流量關(guān)系曲線,并與實(shí)驗(yàn)結(jié)果進(jìn)行比較,分別列于表
6^3和圖6 6中,可以看出,本文的計(jì)算方法能夠比較準(zhǔn)確地給出溢流壩的水位與流

量關(guān)系,從而一方面可以在溢流壩設(shè)計(jì)時準(zhǔn)確計(jì)算出設(shè)計(jì)流量,另一方面,能夠在全
部溢流水頭范圍內(nèi),對理論結(jié)果或模型實(shí)驗(yàn)值進(jìn)行校核。
表6.3水位一流量關(guān)系 Q|Qd
H。/Hd
(計(jì)算)
O,05 O.1 0.2 0.4 0.6 O.8 1.O 1.33

o 153173

o.229759

0.36105

0 579869

0.71116

0.853392

1.004376



181619

H:/H




0 153765

0 237137

0.365716

0 564011

0.726682

0 869824

1.0



195108

(實(shí)驗(yàn)) l相對誤差
0 000592 0 007378 0.004666

0.01586

0.015522

0 015432

0.00438

0 013489















l'@a















圖6.6水位.流量關(guān)系曲線

6.5流量系數(shù)
流量系數(shù)是決定溢流壩泄洪流量的重要參數(shù),對于同一壩面曲線,不同的壩頂水 頭將具有不同的流量系數(shù),為了便于比較,本文在給定流量Q/Qd=O.33,1.o,1.29, 1.58,1.91時計(jì)算出對應(yīng)的實(shí)際壩頂水頭分別為

He/Hd=0.5087,1.0044,1.166,1.313,1.478,根據(jù)公式m2了南即可得到對應(yīng)的流量

河海大學(xué)博士論文

系數(shù),計(jì)算和實(shí)驗(yàn)結(jié)果分別在表6.4和6.5中給出,同時將結(jié)果繪于圖6.7中,二者
基本一致。
表6.4流量系數(shù)(計(jì)算)

Q}Qd
H:}H



O.33

1.O

1.29

1.58

1.91

0.5087

1.0044

1.166

1.313

1.478

0.455647

0.497958

O.511585

0.525832

0.533604

表6 5流量系數(shù)(實(shí)驗(yàn))

Q舊d
Ⅳ,/HJ

0.33

1.O

1.29

1.58

1.91

O.50

1.0

1.17

1.33

1.5

門,

0.455

0.50l

O.513

0.526

0.536

1.6
















0.8

0.6

0.4

0.2
0.5 m 0.55

0.6

圖6.7流量系數(shù)

需要說明的是,關(guān)于壩頂水頭高度的計(jì)算是在水庫上游面按流體體積分?jǐn)?shù)F=O.5 時對應(yīng)的水面線高度進(jìn)行的。
6.6壩面壓力 壩面壓力的計(jì)算是溢流壩水力設(shè)計(jì)的重要內(nèi)容,為此,本文分別計(jì)算了 日。/H。=O.5,1.O,1.17,1.33,1.5共5種情況下的壩面壓力,繪制成曲線并與文獻(xiàn)【17】中的

實(shí)驗(yàn)點(diǎn)進(jìn)行比較,見圖6.8,壓力分布的基本趨勢是一致的,其中工程中關(guān)心的負(fù)壓

河海大學(xué)博士論文

極值,其出現(xiàn)的位置與大小均與實(shí)驗(yàn)符合良好。
0 5
o 4

o 3

。弋
;\ , I
f L ,j

、

0 2 0 l

\. j
√.
,


"▲U

—I --J
、一 ——-一



0 l






3 4 5 6 T

I 1
\。


,,一




f一
。、

—,一J

≯’

\_.
l l

<?


、 、 /
\/
-0


“’“‘叮。r

叼∞加吣加加叼叼

8 -O.了




U.73

U b





0 9

27Hd

圖6.8壩面壓力分布曲線(圓點(diǎn)為實(shí)驗(yàn)值)

6.7壩面水流速度分布
為考察溢流壩過壩水流的流速變化規(guī)律,在圖6.9中分別給出了斷面 x=0.0,2.576,5,10,20,30,壩頂水頭分別為H。/Hd=O.5,1.0,1.5時速度分布曲線圖

90
??+??X=2 576

85
、‘

———o—一X=5
80
75
———●●-—-Y=1 n








,

—tr!唬兀剑玻 一


70 一 占65
)、

一】
60 55

50
● 45

40
O O








惦弘∽

河海大學(xué)博士論文

90

90

85

85

80

75

。蔫_|||_II.。一

80

佰 陽 ∞ ∞ 踮

70

65

60


576 10 15

55

50

45


0 5



Xl02 ?-?--x=

,

蠣 蚰 %

-…一?x=2.5
—-——-a_--—-X=5 ———啼÷——一X=10

76

40

———■r—一x。20

35 20 2 5 30





10

l 5

20

25

30

V(m/s)

V(m/s)

(b)

(c)

圖5.9壩面水流速度分布(a)H。/Hd=O 5(b)H。/Hd=1.0(c)H。/Hd=1.5

為進(jìn)一步考察水流沿壩面的速度變化,給出了個斷面上最大速度的沿程變化曲線 見圖6.10。從圖5.9,6.10可以看出,水流速度沿壩面逐漸增加,且不同水頭下的最
大速度的差值沿程減小。為考察流線的整體變化情況,在圖6.1l中給出了 H。/H。=1.0時整個流場內(nèi)的流線圖。

圖6’IJ

6.8本章小結(jié)

鎰流壩流線矢量圖

河海大學(xué)博士論文

第七章帶閘墩溢流壩三維過壩水流數(shù)值模擬
帶閘墩溢流壩三維過壩水流數(shù)值模擬的工作開展不多,主要原因計(jì)算區(qū)域增大, 閘墩區(qū)域形狀復(fù)雜,并且包含自由液面問題,這些都使得數(shù)值計(jì)算較為困難。因此過 壩水流的研究工作更多地集中在二維流場計(jì)算。盡管二維過壩流場計(jì)算具有計(jì)算工作 量小、計(jì)算成本低,并能夠在一定程度上指導(dǎo)工程設(shè)計(jì),比如文獻(xiàn)『1161通過計(jì)算過 壩水流與閘墩繞流兩個二維模型的方法近似模擬了三維帶閘墩的溢流壩流場,可以作 為解決三維過壩水流問題的一條途徑,但是將二維數(shù)值模型應(yīng)用于實(shí)際的帶閘墩的溢 流壩過壩水流計(jì)算還有一定的差距。主要表現(xiàn)在閘墩的存在對水面線位置、壩面壓力、 溢流流量等均產(chǎn)生顯著的影響。文獻(xiàn)i117]利用有限體積方法對溢流壩三維湍流場進(jìn) 行數(shù)值模擬,但沒有考慮閘墩的影響,本文直接建立溢流壩面與閘墩的三維整體模型, 對三維湍流場進(jìn)行數(shù)值求解,其中湍流模型仍然采用RNGk—s兩方程模型,自由水 面線的位置采用VOF方法確定,針對不同的墩型和布置形式進(jìn)行計(jì)算,比較準(zhǔn)確地
給出了過壩水流的水面線位置及壩面壓力。 7.1閘墩型式

溢流壩壩面曲線選用三圓。祝牛忧f¨8】z‘”=2.0H≯Y,見圖6,l,6.2。 溢流壩設(shè)計(jì)的重要參數(shù)是流量計(jì)算,對流量影響最大的是閘墩,因此必須對包含 閘墩的溢流壩三維流場進(jìn)行數(shù)值模擬。閘墩的作用是將溢流壩前緣分隔為若干個孔
口,并支承閘門、啟閉機(jī)和橋梁等傳來的載荷,閘墩的斷面形狀應(yīng)盡量減小孔口水流 的側(cè)收縮,使水流能夠平順通過,因此關(guān)鍵在于墩頭和墩尾形狀的選擇,頭部主要影 響側(cè)收縮,從而影響泄洪流量,尾部主要影響下游流態(tài),如可能產(chǎn)生水冠或沖擊波,
影響下游的消能。閘墩墩頭的形狀及閘墩的主要尺寸在圖7.1中給出。




乇 州 o



半徑0 133Hj

、

O 2f
● —、‘

7只、


壩項(xiàng)軸線L


j§ 習(xí) 磊 薄麗彰 悟
0.2t

7Ha.J
7I


f.02I 『片。。I壩頂軸線}.02f
’i

7H.,、I


—』乙?

C、/、





’_r



i¨

2型
圈7.1閘墩的型式

3A型4型

河海大學(xué)博士論文

7.2計(jì)算區(qū)域與網(wǎng)格劃分 含閘墩溢流壩三維流場計(jì)算區(qū)域在圖7.2中給出,可以看出在閘墩區(qū)計(jì)算區(qū)域的 幾何形狀最為復(fù)雜,具有雙向曲率變化,這也正是研究該類問題的難點(diǎn)所在。為了保 證計(jì)算精度,對該區(qū)域進(jìn)行網(wǎng)格劃分需要特別注意,本文采用分區(qū)網(wǎng)格系統(tǒng),將整個 計(jì)算區(qū)域劃分為3個子區(qū),見圖7t3。由于閘墩區(qū)域邊界變形較大,采用結(jié)構(gòu)網(wǎng)格無 法生成可用于計(jì)算的網(wǎng)格系統(tǒng),因此在該區(qū)內(nèi)采用非結(jié)構(gòu)網(wǎng)格,同時在壁面處網(wǎng)格進(jìn) 行加密:在I區(qū)和II區(qū)采用結(jié)構(gòu)網(wǎng)格,由于II區(qū)處于自由表面區(qū),為提高自由面的模 擬精度,同樣進(jìn)行了網(wǎng)格加密。此外。為了保證I區(qū)、II區(qū)與閘墩區(qū)網(wǎng)格的平順過渡, 在閘墩區(qū)的上游和左下方增加一部分區(qū)域歸并到閘墩區(qū)構(gòu)成IU區(qū)進(jìn)行網(wǎng)格劃分,三個 區(qū)的網(wǎng)格數(shù)目約為27000,14000,70000,共計(jì)約1l萬個網(wǎng)格節(jié)點(diǎn)。在圖7.4中給出 了壩面及閘墩區(qū)的網(wǎng)格圖。

.100

圖7.2含閘墩溢漣壩流場計(jì)算區(qū)域

河海大學(xué)博士論文

圖7.3

含閘墩溢流壩流場網(wǎng)格的分區(qū)劃分

圖7.4壩面及閘墩區(qū)的網(wǎng)格劃分

河海大學(xué)博士論文

7.3計(jì)算初始條件與邊界條件

壩高H=76.2米,設(shè)計(jì)水頭H。=9.14米,設(shè)計(jì)流量傷=CL蟛3”,C為設(shè)計(jì)流
量系數(shù),上為壩的凈寬=壩寬.閘墩厚度,分別對2型墩和3A型墩進(jìn)行計(jì)算,其中墩 長均為1.77Hd,墩厚分別為0.267Hd和0.173Hd,墩厚與閘寬之比均為0.2。假定壩 項(xiàng)水頭高度及流量,分別對三種不同水頭高度Ⅳ。/Hd=0.5,1.O,1-33下的自由水面線、 壩面壓力、流量系數(shù)進(jìn)行計(jì)算,計(jì)算初值見下表
表7.1水頭高度及流量初值 H!/Hd
(初值)

O.50

1.O

1.33

Q『Q。 =(H。/Hd)‘6

0.329877

1.O

1.578202

入口邊界的速度己知,k,s值按公式々:0.005“。2,£:燮給出,出口處的速度
)P

及k和£值為未知.計(jì)算根據(jù)與出口斷面相鄰的上游兩個斷面的變量值,線性外延求 得出口斷面之值,作為下一次迭代的邊界條件.計(jì)算過程中,每次迭代都需根據(jù)上一 次迭代計(jì)算結(jié)果修正出口邊界條件,固體壁面上流速按非滑移條件給定,固壁附近的 k和占采用用壁函數(shù)技術(shù);自由表面上的壓力取為大氣壓,而k和£采用零梯度條件, 即令自由表面網(wǎng)格中的k和s的值等于相鄰的流體內(nèi)部網(wǎng)格的k和s的值。計(jì)算開始 流動為非定常流,經(jīng)過一段時間的計(jì)算,流場各水力要素均趨不變,流動趨于定常, 計(jì)算即告結(jié)束。 7.4自由水面線 按照上節(jié)的計(jì)算條件分別對2型墩和3A型墩水面線形狀及位置進(jìn)行計(jì)算,
7,4.1

2型墩水面線計(jì)算 圖7.5中給出了沿2型墩閘孔中心線水面線位置計(jì)算與實(shí)驗(yàn)結(jié)果的比較,三條曲

線分別對應(yīng)日/H。=0.5,1.O,1.33三種情況。二者是比較一致的,在壩頂水頭較低時, 誤差較大,隨著壩頂水頭的增加,計(jì)算與實(shí)驗(yàn)結(jié)果的一致性更好。與沿閘孔中心線的 計(jì)算相比,沿閘墩表面水面線位置計(jì)算結(jié)果與實(shí)驗(yàn)值更加接近(見圖7.6),說明本計(jì) 算模型是準(zhǔn)確、可靠的。

80

河海大學(xué)博士論文

90

85

80

— 5

h 70

75

65

60 10 —5 0 5 10 15 20 25 30

x(m)

圖7.5沿2型墩閘孔中心線水面線位置比較
(從上至下分別對應(yīng)He/Hd=l
90
33.1 0.0.5)

85

80

.、

縣75
>、

70

65

60

10

—5





10

15

20

25

30

x(Ⅲ)

圖7.6沿2型墩表面水面線位置比較
(從上至下分別對應(yīng)He/Hd=I
33,1.0,0 5)

為了比較閘墩對水面線位置的影響,在圖7.7和圖7.8中分別給出了無閘墩、沿 閘孔中心線和沿閘墩表面的自由水面線位置的計(jì)算曲線,顯然閘墩的存在抬高了自由 水面線的位置,且隨著水位的增加,抬高的幅度逐漸增加。 而沿閘墩表面的變化則

是閘墩頭部水面抬升,在閘墩的平行段水面低于閘孔中心線處的水面。 在圖7.9中比較直觀地給出了三種水頭情況下通過2型墩閘孔的過壩水流。圖

河海大學(xué)博士論文

7.10中給出了不同高度處閘墩附近的流線矢量圖。

90

85

80




75

65

60 lO 一5 0 5

x(m30

15

20

25

30

圖7.7無閘墩與2型墩沿閘孔中心線的自由水面線
(從上至下分別對應(yīng)He/Hd=l
33,1

0,0.5)

85一

占75,


70

60!
一10 0 5 10 15 20 25

x(面 圖7.8沿2型墩閘孔中心線與沿閘墩表面的自由水面線
(從上至下分別對應(yīng)He/Hd=I
33,1

0,0.5)

河海大學(xué)博士論文

(a)

(b) 圖7.9 H/Hd=0.5,1.O,1.33時2型墩過壩水流 (紅色為水,黃色為閘墩和壩面)

(c)

圖7.10日/H。=1.0時2型墩過壩水流不同高度處的流線矢量圖 (a)h=77m,(b)h=78m,(c)h=79m,(d)h=80m
7.4.2

3A型墩水面線計(jì)算 圖7.1l,7.12中給出了設(shè)計(jì)水頭H/Hd=O.5,1.0,1.33三種情況下沿3A型墩閘孔

中心線和沿閘墩的水面線位置計(jì)算結(jié)果,并與無閘墩情況進(jìn)行了比較。與圖7.7,7.8 比較可以看出,兩種墩型對孔中心處水面線位置影響不大,為了比較3A型墩與2型 墩對水面線位置的影響,在圖7.13中給出了沿兩種墩型表面的水面線位置計(jì)算結(jié)果, 2型墩使墩前水面產(chǎn)生更高的抬升。

河海大學(xué)博士論文

90

85

80




75

70

65

60 5 0 5


,10





20

25

30

LmJ

圖7.11無閘墩與3A型墩沿閘孔中心線的自由水面線
(從上至下分別對應(yīng)He/Hd=1.33,1
90 0,0 5)

85

80


)、

75

70

65

60 10 —5 0 5

x(拶

15

20

25

30

圖7.12沿3A型墩閘孔中心線與沿閘墩表面的自由水面線
(從上至下分別對應(yīng)He/Hd=l
33,1

0.0.5)

河海大學(xué)博士論文

90

85

80

占75


70

65

60 10 —5 0 5 10 15 20 25 30

X(m)

酗7.13沿2型墩及3A型墩表面的自由水面線
(從上至下分別對應(yīng)He,/Hd=I.33,1 0,0.5)

7.5壩面壓力 壩面壓力計(jì)算,特別是壩面最小壓力的數(shù)值及出現(xiàn)位置是溢流壩設(shè)計(jì)的重要參 數(shù),在圖7.14和圖7.15中給出了三種不同水頭下2型閘墩沿閘孔中心線及沿閘墩的 壩面壓力分布曲線,在圖7.16和圖7.17給出了三種不同水頭下3A型閘墩沿閘孔中

心線及沿閘墩的壩面壓力分布曲線,并與實(shí)驗(yàn)結(jié)果進(jìn)行對比,可以看出,沿閘孔中心

線壩面壓力的計(jì)算結(jié)果與實(shí)驗(yàn)值非常一致,只是在以/z‘=1.0時在壩頂后段計(jì)算值
偏低,而壩面最小壓力值的大小及出現(xiàn)位置均與實(shí)驗(yàn)吻合:對于沿閘墩面的壓力,總

的趨勢是計(jì)算結(jié)果略低于實(shí)驗(yàn)值,并且隨著水頭的增加其偏差量增大。在 H。/H。,=1.O、1.33兩種情況下沿閘墩的壩面上的最小負(fù)壓出現(xiàn)位置與實(shí)驗(yàn)一致,均在
壩頂前0.1H。處。

將2型墩與3A型墩的壩面壓力計(jì)算結(jié)果進(jìn)行比較可以發(fā)現(xiàn),在相同的墩厚閘寬

比情況下,2型墩將使壩面出現(xiàn)更低的負(fù)壓,當(dāng)水頭H。/H。=1.33時負(fù)壓極值為7.3 米水柱,超過了工程中建議的不低于一6.1米水柱以保證壩面不出現(xiàn)空蝕現(xiàn)象的要求,

而3A型墩在相同水頭下對應(yīng)的最低負(fù)壓僅為一4.6水拄,不會出現(xiàn)壩面空蝕現(xiàn)象。

河海大學(xué)博士論文

0.6 0,4 0,2



鼉0

-0,2
—0,4

—0,6 —0,8 -i -o,3一O,1 O,l

O,3

O,5

0,7

0,9

1,1

1,3

X/Hd
圖7.14沿2型墩閘孔中心線上的壩面壓力分布




0.4






.C






一0

—0 4
—0 6

—0 8 1 —0.3—0.1 0.1 0.3 o.5 0.7 0.9 I.1 I.3

X/Hd 圖7.15沿2型閘墩表面的壩面壓力分布
0.6 O.4 0.2 弓

姜0
工 一0.2 —0.4 —0.6 —0.8
—1

降喪釜≤。。銎琴 ■':~:i。辏。楹
0.1
0.3

-0 3—0.1

0.5

0.7

0.9

1.1

1.3

X,,I-Id

圖7.16沿3A型墩閘孔中心線上的壩面壓力分布比較
(從上至下分別對應(yīng)He.,/l-Id=0 5,I.0,l 33)

河海大學(xué)博士論文

圖7.17沿3A型閘墩表面的壩面壓力分布比較
(從上至下分別對應(yīng)HeJHd=0.5,1.0,l 33)

7.6流量系數(shù) 流量系數(shù)關(guān)系到工程的泄洪規(guī)模和具體布置,但流量系數(shù)隨工程的體型而變化, 由于水流均為三維,其流態(tài)和邊界條件都比較復(fù)雜,故在一般的手冊中難以找到準(zhǔn)確 的數(shù)據(jù)。在圖7.18中比較了墩厚閘寬比t/b=0.2時2型墩和3A型墩在三種不同水頭 下的流量系數(shù),顯然,3A型墩具有更大的流量系數(shù)。
m 6



呲lIJ
n 5

¨5
m 4 O,4 0.6 0.8


1.2

1.4

He/lid

圖7.18有閘墩影響的流量系數(shù)(r/b=0.2)

7.7墩厚閘寬比對壩面壓力及流量系數(shù)的影響
除壩型和閘墩型式對泄流流量、壩面壓力、水面線位置產(chǎn)生較大影響外,閘墩在 壩面上的布置同樣至關(guān)重要,特別是對壩面壓力影響更為顯著,因此為了考察閘墩布

置對壩面壓力等的影響,針對3A型閘墩,給定初值HePr-Id=1.0時情況,計(jì)算了墩厚

f與閘寬b之比t/b卸.14,0.2,O.5三種情況下壩面壓力分布,對應(yīng)的實(shí)際水頭高度

河海大學(xué)博士論文

分別為0.968271,1.007659,1.017505,分別給出了沿iYq孑t,ee心線和沿閘墩的壩面壓 力分布曲線(圖7,19,7.20),并給出了對應(yīng)的流量系數(shù)曲線(圖7.21)。
0,6
0,4

0,2
-。

暮0


一0,2
—0,4

-0,6
—0,8 —1

—0.3—0,1

0,l

0.3

0,5

0,7

0,9

i,1

1,3

X/Hd 圖7.19

He/Hd=1.0時不同墩厚閘寬比沿閘孔中心線的壓力分布

0,6

0,4 0,2


罩0
上 一0,2

.撥
一¨_V

r._?。翕殻 垃 -1lt-1一一:

—-?一t/b=0 5 ---II---t/b=0.2 ——?廣_t/b=0 14
r—一二


—0。4 —0,6

—0,8
。 . :

一1
一0.3—0.1



;

0,1

0,3

0,5

0,7

0.9

1,1

l,3

X/Hd

圖7.20

He/Hd=1.0時不同墩厚閘寬比沿閘墩面的壓力分布

88

河海大學(xué)博士論文

0,55

0。53

0,51


0 49

0,47

0.45 0.4 0,6 0,8

t/b 圖7.21

He/Hd=1.0時不同墩厚閘寬比對應(yīng)的流量系數(shù)

7.8本章小結(jié) 本章在二維過壩水流流場計(jì)算的基礎(chǔ)上,直接建立溢流壩和閘墩的三維整體模

型,對三維過壩水流粘流場進(jìn)行數(shù)值計(jì)算,其中采用了分區(qū)網(wǎng)格系統(tǒng)來處理復(fù)雜的幾
何邊界,針對2型墩和3A型墩計(jì)算了閘墩型式對水面線形狀、壩面壓力、流量系數(shù)

等的影響,通過與實(shí)驗(yàn)結(jié)果的對比表明,本文建立的模型能夠比較準(zhǔn)確地給出過壩水 流的水面線位置及壩面壓力,計(jì)算指出,閘墩的存在抬高了水面線的位置,在閘墩頭
部尤其明顯,3A型閘墩相對于2型墩具有更大的流量系數(shù)和更小的壩面壓力。此外, 不同的墩厚閘寬比對泄流能力也將產(chǎn)生顯著的影響,結(jié)果表明,隨著墩厚閘寬比的增 加,壩面壓力降低,而當(dāng)t/b=0.2時溢流壩具有更大的流量系數(shù)。

河海大學(xué)博士論文

第八章結(jié)論與展望
本文對粘流場數(shù)值模擬技術(shù)的研究現(xiàn)狀和發(fā)展動態(tài)作了比較全面的回顧和展望, 重點(diǎn)研究了網(wǎng)格生成技術(shù)、數(shù)值求解技術(shù)、湍流模型技術(shù)和動邊界模擬技術(shù),在此基
礎(chǔ)上建立了模擬自由面湍流場的數(shù)學(xué)模型,并成功地應(yīng)用到二維溢流壩和帶閘墩的三

維溢流壩過壩水流場數(shù)值計(jì)算。

8.1本文完成的主要工作與結(jié)論 1)對網(wǎng)格生成技術(shù)進(jìn)行了較深入的研究,建立了一種簡單的、能夠獨(dú)立于網(wǎng)格生成 技術(shù)之外的邊界正交化方法,給出了詳細(xì)的數(shù)學(xué)模型并編制了計(jì)算程序;通過網(wǎng)格生 成算例表明,該方法能夠大大改善網(wǎng)格在邊界處的正交性;同時針對分區(qū)結(jié)構(gòu)網(wǎng)格系 統(tǒng)建立了分區(qū)交界面處的數(shù)據(jù)結(jié)構(gòu)和計(jì)算模型: 2)利用有限體積方法在非正交同位網(wǎng)格系統(tǒng)中建立了SIMPLE求解算法,并對非正 交同位網(wǎng)格系統(tǒng)中的邊界條件、延遲修正技術(shù)及計(jì)算節(jié)點(diǎn)的梯度計(jì)算等專題內(nèi)容進(jìn)行
了深入討論。通過對傾斜方腔驅(qū)動流、后臺階繞流、非對稱平板間的圓柱繞流等典型

流場的數(shù)值模擬,表明本文建立的非正交同位網(wǎng)格系統(tǒng)中的SIMPLE算法是合理的,
具有更加廣泛的適用范圍:

3)研究比較了兩種湍流模型的特點(diǎn)及應(yīng)用范圍,以水翼粘流場為例完成了RNGk一£ 湍流模型與標(biāo)準(zhǔn)k一£模型的計(jì)算比較,表明RNGk一占模型的模擬結(jié)果更優(yōu); 4)研究了運(yùn)動邊界模擬技術(shù)中的VOF方法,詳細(xì)建立了一種利用斜直線構(gòu)造自由面
的高精度方法,通過對典型運(yùn)動界面的數(shù)值模擬表明,該方法能夠比其它自由面重構(gòu)
方法更加準(zhǔn)確地確定運(yùn)動界面的形狀和位置:

5)完成了對二維溢流壩過壩水流進(jìn)行計(jì)算,分析了不同壩頂水頭下的水面線位置、

壩面壓力及流量、流速等水力特性,通過與典型實(shí)驗(yàn)資料的對比,計(jì)算的水面線高度
及壩面負(fù)壓極值均與實(shí)驗(yàn)值具有非常好的一致性。

6)直接建立溢流壩和閘墩的三維整體模型,對三維過壩水流粘流場進(jìn)行數(shù)值計(jì)算, 針對閘墩的型式及在壩面的布置計(jì)算了閘墩對水面線形狀、壩面壓力、流量系數(shù)等設(shè) 計(jì)參數(shù)的影響,并將不同墩型與布置形式的計(jì)算結(jié)果進(jìn)行比較,結(jié)果表明,閘墩的存
在抬高了水面線的位置,其中在閘墩頭部尤其明顯;墩型對流量系數(shù)和壩面壓力影響 較大,3A型閘墩相對于2型墩具有更大的流量系數(shù)和更小的壩面壓力;不同的墩厚 閘寬比對泄流能力也將產(chǎn)生顯著的影響,隨著墩厚閘寬比的增加,壩面壓力降低,而

河海大學(xué)博士論文

當(dāng)f/b=0.2時溢流壩具有更大的流量系數(shù)。

7)本文建立的三維自由面湍流場模型能夠比較準(zhǔn)確地模擬帶閘墩溢流壩三維過壩水 流湍流場,可以針對不同的壩型、墩型、壩面布置形式以及不同的壩頂水頭完成系列
水力計(jì)算,為泄水工程建筑物的設(shè)計(jì)提供可靠的分析依據(jù)。

本論文的主要創(chuàng)新點(diǎn): 1)提出并建立了網(wǎng)格邊界正交化算法模型:
2)建立了非正交同位網(wǎng)格系統(tǒng)中的SIMPLE算法,并以此為基礎(chǔ)建立了非對接分區(qū) 結(jié)構(gòu)網(wǎng)格數(shù)值求解模型; 3)提出了一種精度更高的自由面重構(gòu)算法模型,通過典型運(yùn)動界面的模擬驗(yàn)證了本 算法是精確的:

4)建立了帶閘墩溢流壩三維整體計(jì)算模型,首次完成了帶閘墩溢流壩三維過壩水流
粘流場數(shù)值模擬,得出了一些有意義的結(jié)論。

8.2展望 雖然本文對帶自由面湍流場數(shù)值模擬技術(shù)開展了一定的研究工作,并成功地應(yīng) 用到泄水工程中三維自由面湍流場的數(shù)值模擬,得到了一些有意義的結(jié)論。但是由 于三維自由面湍流是一個十分復(fù)雜的流體運(yùn)動現(xiàn)象,其復(fù)雜的運(yùn)動規(guī)律和內(nèi)在機(jī)理
需要作進(jìn)一步的分析和研究,因此,本論文只是在幾個方面做了一些探討性的研究, 今后仍有大量的工作要做,具體內(nèi)容包括: 1)進(jìn)一步完善非正交同位網(wǎng)格系統(tǒng)中的SIMPLE算法,使之適用于網(wǎng)格強(qiáng)非正交的
情形;

2)對基于VOF方法的運(yùn)動界面追蹤技術(shù)中的輸運(yùn)模型和重構(gòu)方法仍需作更深入的研

究,進(jìn)一步提高計(jì)算精度,以便更好地模擬水流的破碎、空化等現(xiàn)象; 3)開展帶閘墩溢流壩上下游流場整體模型數(shù)值計(jì)算,探討泄洪流量、空化空蝕等復(fù) 雜現(xiàn)象的內(nèi)在機(jī)理。

9l

河海大學(xué)博士論文





作者在攻讀博士論文的過程中得到了許多的指導(dǎo)、支持和幫助,在此表示衷心的
感謝。

首先感謝我的導(dǎo)師汪德燧教授,正是在導(dǎo)師一步一步的悉心指導(dǎo)和督促下,才得
以順利完成畢業(yè)論文,導(dǎo)師深厚的學(xué)術(shù)造詣、嚴(yán)謹(jǐn)?shù)闹螌W(xué)態(tài)度、踏實(shí)的工作作風(fēng)給我 留下了深刻的印象并將使我受益終生,在論文即將完成之際,特向?qū)煴硎咀钫\摯的 謝意和最美好的祝愿。

感謝華東船舶工業(yè)學(xué)院副校長王自力教授和科技處長蔣志勇教授,他們在我的
工作中給予了非常大的幫助,使我克服了工作中的許多困難,為順利完成論文提供
了良好的工作和學(xué)習(xí)條件。

感謝楊松林教授、朱仁慶副教授、姚震球副教授在論文工作中提供的無私幫助
和建議。

感謝師弟賴錫軍,在與他日常的討論交流中作者獲得了許多有益的幫助和啟

迪,他的聰穎與豁達(dá)也給我留下了深刻的印象。感謝同學(xué)童朝峰、陳雄波,在與他 們的交往過程中獲得了許多幫助。
在論文撰寫過程中,得到了研究生周林慧的幫助,在此表示感謝。 還要感謝我的家人,他們幫助我分擔(dān)了許多家務(wù),提供了大力的支持,正是在

他們的幫助和鼓勵下,我才得以將更多的時間、精力和飽滿的熱情投入到論文工作
中:感謝母親對我的關(guān)心和牽掛,感謝并深深懷念我的父親。 最后再一次向關(guān)心、支持和幫助過我的所有同仁表示最誠摯的謝意。

河海大學(xué)博士論文

參考文獻(xiàn)
1.Thompson J.F.,Weatherill N.P,Aspects of numerical grid genertion:current and art.

A1AA(paper),93—3539
2.Shih T.I—Eand Bailey R.T,Algebraic grid generation for complex geometries,

Int.J.Num


  本文關(guān)鍵詞:三維自由面湍流場數(shù)值模擬及其在水利工程中的應(yīng)用,由筆耕文化傳播整理發(fā)布。



本文編號:184731

資料下載
論文發(fā)表

本文鏈接:http://sikaile.net/kejilunwen/shuiwenshuili/184731.html


Copyright(c)文論論文網(wǎng)All Rights Reserved | 網(wǎng)站地圖 |

版權(quán)申明:資料由用戶1fdd5***提供,本站僅收錄摘要或目錄,作者需要刪除請E-mail郵箱bigeng88@qq.com