一種基于多維黎曼解的任意拉格朗日歐拉方法
【專利摘要】本發(fā)明公開了一種基于二維黎曼解的任意拉格朗日歐拉方法，用以求解二維可壓縮流體力學(xué)方程組及以之為基礎(chǔ)的相關(guān)問題。整個算法的流程包括：網(wǎng)格生成與物理量的初始分布；確定下一時刻網(wǎng)格的策略和方式；確定網(wǎng)格單元邊界流量的二維黎曼解法器算法；網(wǎng)格單元邊界流量的表達形式；獲得下一時刻的計算結(jié)果的數(shù)值方法。新的發(fā)明提供了一種確定網(wǎng)格單元邊界流量的新的二維黎曼解法器算法，能夠克服傳統(tǒng)的一維黎曼解法導(dǎo)致的數(shù)值激波不穩(wěn)定、網(wǎng)格易扭曲等現(xiàn)象，改正了現(xiàn)有二維黎曼解方法過于復(fù)雜而難以實施的缺點。本發(fā)明是一種簡潔、健壯和精確的數(shù)值算法組件，適合多介質(zhì)大變形問題的數(shù)值模擬，適合現(xiàn)在流行的有限差分、有限體積和有限元方法的形式。
【專利說明】一種基于多維黎曼解的任意拉格朗日歐拉方法
【技術(shù)領(lǐng)域】
[0001]本發(fā)明屬于流體力學(xué)數(shù)值模擬【技術(shù)領(lǐng)域】，涉及一種基于多維黎曼解的任意拉格朗日歐拉方法，具體地說，涉及一種利用二維黎曼解模擬多介質(zhì)可壓縮流動的任意拉格朗日歐拉方法。
【背景技術(shù)】
[0002]多介質(zhì)大變形問題的數(shù)值模擬，是流體力學(xué)問題數(shù)值模擬中具有挑戰(zhàn)性的課題。它需要求解涉及多種介質(zhì)相互作用、并發(fā)生劇烈運動的問題(如，介質(zhì)撞擊和作用的彈塑性問題，爆轟問題，高溫高壓的輻射流體力學(xué)問題)。在這類問題中，介質(zhì)常常具有復(fù)雜的物態(tài)性質(zhì)，狀態(tài)量(密度、壓強)會發(fā)生突然的改變(激波、爆轟波出現(xiàn))，作用區(qū)域有可能會發(fā)生幾百倍的變化。這種問題的計算方法包含了復(fù)雜的流程，如網(wǎng)格運動、介質(zhì)界面的追蹤和可壓縮流體力學(xué)方程組的數(shù)值解法。
[0003]在工程中目前應(yīng)用的主要方法有歐拉方法、拉格朗日方法和任意拉格朗日歐拉方法。歐拉方法使用固定網(wǎng)將，實施簡便。主要困難一個是難以精確地計算區(qū)域變動巨大的問題，另一個是難以正確地描述不同介質(zhì)處于同一個網(wǎng)格內(nèi)的狀態(tài)(特另是當(dāng)介質(zhì)的物態(tài)性質(zhì)相差較大、或者存在輻射等多物理作用過程時)。在彈塑性問題和輻射流體力學(xué)問題中，人們較少使用這一方法，更多的是采用拉格朗日方法和任意拉格朗日歐拉方法。拉格朗日方法使用的網(wǎng)格隨著介質(zhì)的運動而運動，可以清晰地分辨介質(zhì)界面，不需要計算介質(zhì)在人為混合后的狀態(tài)。但這種方法非常容易造成網(wǎng)格扭曲，導(dǎo)致計算中斷。任意拉格朗日歐拉方法可以在計算過程中調(diào)節(jié)網(wǎng)格的運動方式，例如僅在介質(zhì)界面處讓網(wǎng)格隨著介質(zhì)運動，而在其它地方則調(diào)整網(wǎng)格，使之保持較好的幾何品質(zhì)。該方法在理論上同時兼具歐拉與拉格朗日方法的優(yōu)點，是一種較理想的解決方案。但在實際的計算過程中，計算結(jié)果還是會出現(xiàn)嚴(yán)重的錯誤，如虛假的渦度誤差和數(shù)值激波不穩(wěn)定現(xiàn)象。這些錯誤扭曲了真實的物理圖像，并常常導(dǎo)致計算過程的中斷。計算出現(xiàn)錯誤的主要原因是流體力學(xué)方程組的解法器出現(xiàn)了問題。
[0004]流體力學(xué)方程組的解法器是一種利用當(dāng)前時刻的物理量，獲得下一時刻的物理量的算法實施組件。其關(guān)鍵之處是設(shè)計當(dāng)流體通過網(wǎng)格邊界時數(shù)值通量(流量)的算法。目前在工程應(yīng)用中廣泛使用的流體力學(xué)方程組解法器，都是沿網(wǎng)格邊界法線方向，求解一維的黎曼問題，并以獲得的結(jié)果作為流量的近似值使用的。這種方法對一般流體力學(xué)問題是適用的，但在多介質(zhì)大變形等極端條件下，由于缺少多維信息，對多維效應(yīng)引起的渦度誤差和數(shù)值激波不穩(wěn)定性等難點問題難以解決。依靠加密計算網(wǎng)格，減少時間步長，高精度技術(shù)等傳統(tǒng)方法并不能改正這些缺陷。有時候加密網(wǎng)格還導(dǎo)致計算終止得更快。采用包含多維信息的多維黎曼解法器，是解決這類問題的可行方案。
[0005]遺憾的是，針對二維黎曼問題的解法器還沒有在任何這類問題的應(yīng)用軟件中實施。這是由于二維黎曼問題的精確求解理論上還沒有完全解決，即使使用近似解法器方法，實施過程也非常復(fù)雜?，F(xiàn)有的研究成果，基本都停留在學(xué)術(shù)層面，大部分僅僅應(yīng)用在四邊形網(wǎng)格的歐拉方法中。法國的Despres和Maire等人，分別設(shè)計了一種可以在非結(jié)構(gòu)網(wǎng)格(六邊形或者三角形)上使用的二維近似黎曼解法器，但僅僅適用于拉格朗日方法?，F(xiàn)有的任意拉格朗日歐拉方法在網(wǎng)格邊界的通量計算中，還沒有實施簡單、適用于復(fù)雜網(wǎng)格體系、健壯有效的二維黎曼解法器。

【發(fā)明內(nèi)容】

[0006]為了克服現(xiàn)有技術(shù)中的缺陷，本發(fā)明提供了一種利用二維黎曼解模擬多介質(zhì)可壓縮流動的任意拉格朗日歐拉方法，針對復(fù)雜的流體力學(xué)計算問題，研制的方法能適應(yīng)復(fù)雜的網(wǎng)格體系，具備分辨多介質(zhì)問題的能力，并取得良好的計算效果。與傳統(tǒng)方法使用局部一維的黎曼解方法不同的是，本發(fā)明設(shè)計了一個兩維的黎曼問題解法器，利用多維信息改正傳統(tǒng)方法的不足。本發(fā)明是一種健壯、簡潔的數(shù)值算法，可以明顯地減少傳統(tǒng)方法的誤差。其技術(shù)方案為，
[0007]一種基于多維黎曼解的任意拉格朗日歐拉方法，包括以下步驟:
[0008]I)網(wǎng)格生成與物理量的初始分布
[0009]所有的物理量(密度、速度、壓強、能量)均定義在網(wǎng)格的中心。
[0010]2)下一時刻網(wǎng)格的確定和生成
[0011]針對不同的網(wǎng)格策略(歐拉方法，自適應(yīng)網(wǎng)格方法和拉格朗日方法)，求解特定的網(wǎng)格方程，并形成軟件程序，主要包括橢圓型網(wǎng)格生成器，自適應(yīng)網(wǎng)格生成器等。
[0012]3) 二維黎曼解法器組件
[0013]傳統(tǒng)方法缺陷:工程應(yīng)用中的一維黎曼解算法，不具備二維黎曼解的性質(zhì)。文獻中已有的二維黎曼解算法需要利用前一個時刻的物理量，求解復(fù)雜的波系相互作用情況。由于理論不完備，這造成算法的實施非常困難。
[0014]新的發(fā)明想法:如果二維黎曼解法器僅僅提供經(jīng)過網(wǎng)格邊界的流量算法，而不計算復(fù)雜波系的相互作用，就可以大大減少問題的難度。
[0015]新的發(fā)明遇到的困難:如何從前一個時刻網(wǎng)格中心的物理量，獲得直到下一個時刻的網(wǎng)格節(jié)點和網(wǎng)格邊界的物理量。如果采用直接的插值等數(shù)值技術(shù)，由于物理量存在間斷，將會完全失效。
[0016]新的發(fā)明主要內(nèi)容包括:節(jié)點運動速度的確定方法和網(wǎng)格邊界通量的計算方法。
[0017]I)節(jié)點運動速度的確定方法
[0018]在某一網(wǎng)格V。的所有邊界上(看圖1(a))，求解沿著邊界qq+法方向的局部一維黎曼問題，獲得網(wǎng)格邊界上的運動速度。當(dāng)圍繞所要求解的節(jié)點< (看圖1(b))的所有邊界上的邊運動速度全部求出時，利用這些邊速度的法方向分量，采用最小二乘方法求出網(wǎng)格節(jié)點的運動速度 <。其具體規(guī)則是節(jié)點的運動速度在網(wǎng)格邊的投影與這些邊速度的法方
向分量的差，在最小二乘意義下最小。計算過程使用的權(quán)函數(shù)，需要特殊的設(shè)計。
[0019]2)網(wǎng)格邊界通量的計算方法
[0020]利用節(jié)點處的二維流體運動速度< ;將此速度在網(wǎng)格邊界qq+處的投影，更新為
相應(yīng)邊界的新的流動速度，再以兩邊的狀態(tài)和此流動速度，設(shè)計沿qq+法向的新的一維黎曼解。此時由于網(wǎng)格節(jié)點處的運動速度，進入了邊界通量的計算中，傳統(tǒng)一維黎曼解中的間斷跳躍條件不再成立。一種單邊的新的跳躍條件被構(gòu)造出來，一個邊界上具有了兩個流量(看圖2(b))，分別對應(yīng)不同的網(wǎng)格單元。這種方法能保證網(wǎng)格的運動方式與邊界通量完全相容。沿邊界的切向動量，也由于引進了節(jié)點速度的平均作用，導(dǎo)致計算更加穩(wěn)定。
[0021]3)邊界通量表達形式
[0022]給定了兩層網(wǎng)格節(jié)點的移動速度后，根據(jù)黎曼解的結(jié)構(gòu)，設(shè)計出邊界通量的表達形式。4)獲得下一時刻的計算結(jié)果
[0023]如果計算時間已經(jīng)到了需要的時刻，則計算終止，否則回到第2)步。
[0024]本發(fā)明的有益效果:
[0025]本發(fā)明適合復(fù)雜的網(wǎng)格，特別是無結(jié)構(gòu)網(wǎng)格及結(jié)構(gòu)和非結(jié)構(gòu)混合網(wǎng)格；清晰地分辨介質(zhì)界面；明顯地改善現(xiàn)有流體力學(xué)算法的健壯性，提高計算精度。
【專利附圖】

【附圖說明】
[0026]圖1是網(wǎng)格及其節(jié)點的幾何示意圖，其中圖1(a)是網(wǎng)格單元V。及相鄰單元，圖1(b)是網(wǎng)格節(jié)點及相鄰單元；
[0027]圖2是構(gòu)造邊界流量的算法示意圖，其中圖2(a)是傳統(tǒng)算法，圖2(b)是新提出的邊界流量算法；
[0028]圖3是新發(fā)明的局部一維黎曼解的波系分解圖，其中圖3(a)是傳統(tǒng)算法，圖3 (b)是新提出的邊界通量；
[0029]圖4是激波衍射問`題密度輪廓圖，其中圖4 (a)是采用傳統(tǒng)的HLLC —維近似黎曼解算法，圖4(b)是采用新的二維黎曼解算法；
[0030]圖5是初始網(wǎng)格一致傾斜時(Saltzmann網(wǎng)格)，其中圖5 (a)是采用傳統(tǒng)的HLLC一維近似黎曼解算法，圖5(b)是采用新的二維黎曼解算法，
[0031]圖6是初始網(wǎng)格為三角形網(wǎng)格的雙馬赫問題，其中圖6(a)是采用傳統(tǒng)的HLLC —維近似黎曼解算法，圖6(b)是采用新的二維黎曼解算法；
[0032]圖7是初始網(wǎng)格為三角形網(wǎng)格的雙馬赫問題，其中圖7(a)是為密度輪廓，圖7(b)是為網(wǎng)格圖；
[0033]圖8是初始網(wǎng)格為三角形網(wǎng)格的多介質(zhì)激波管問題，圖8(a)為初始網(wǎng)格，圖8 (b)為計算過程中的拉格朗日網(wǎng)格，圖8 (c)為計算過程中的自適應(yīng)移動網(wǎng)格；圖8 (d)是沿著橫向的密度分布，圖8(e)是沿著橫向的壓強分布。
【具體實施方式】
[0034]下面結(jié)合附圖和【具體實施方式】對本發(fā)明的技術(shù)方案作進一步詳細地說明。
[0035]使用一種基于多維黎曼解的任意拉格朗日歐拉方法，求解如下的流體力學(xué)方程組
謂 U) f)G(U)
[0036]— + -十」+ —=Hs
ο? axay
[0037]其中H表示多物理過程的源項，U，F(xiàn)，G分別是
【權(quán)利要求】
1.一種基于多維黎曼解的任意拉格朗日歐拉方法，其特征在于，包括以下步驟: 1)網(wǎng)格生成與物理量的初始分布； 2)下一時刻網(wǎng)格的確定和生成:需要求解與求解特定的網(wǎng)格方程，并形成軟件程序，主要包括橢圓型網(wǎng)格生成器，自適應(yīng)網(wǎng)格生成器； 3)黎曼解法器組件 構(gòu)造一個全新的兩維近似黎曼解法，該解法器包括節(jié)點速度的確定方法和網(wǎng)格邊界通量的計算方法，在任意網(wǎng)格下，網(wǎng)格節(jié)點的速度算法:已知網(wǎng)格單元中心的物理量，需要求解2X2的代數(shù)方程組來獲得節(jié)點處的流動速度，點處的流動速度與網(wǎng)格邊界通過一些特定方法獲得的法向速度差在最小二乘意義下最??； 構(gòu)造邊界流量的方法:利用了節(jié)點處的二維流體速度，以此速度在網(wǎng)格邊界qq+處的投影，作為相應(yīng)邊界的流動速度，再以兩邊的狀態(tài)和此流動速度；構(gòu)造沿qq+法向的一維黎曼解，由于網(wǎng)格節(jié)點處的運動速度，進入了邊界通量的計算中，這種方法能保證網(wǎng)格的運動方式與邊界通量完全相容，沿邊界的切向動量，也引進了節(jié)點速度的作用； 4)邊界通量表達形式:給定了兩層網(wǎng)格節(jié)點的移動速度后，根據(jù)黎曼解的結(jié)構(gòu)，寫出邊界通量的表達形式； 5)獲得下一時刻的計算結(jié)果 如果計算時間已經(jīng)到了需要的時刻，則計算終止，否則回到第2)步。
【文檔編號】G06F17/50GK103823916SQ201310498537
【公開日】2014年5月28日 申請日期:2013年10月23日 優(yōu)先權(quán)日:2013年10月23日
【發(fā)明者】沈智軍, 閆偉, 袁光偉 申請人:沈智軍, 閆偉, 袁光偉