理想流體力學(xué)課程設(shè)計(jì)(HessSmith方法求附加質(zhì)量)

.一、物理背景無(wú)論是船舶還是海洋平臺(tái)在海洋開(kāi)發(fā)中都起著關(guān)鍵的作用，而開(kāi)發(fā)海洋首先需要對(duì)海洋結(jié)構(gòu)物進(jìn)行深入地研究。這其中，水動(dòng)力學(xué)中的附加質(zhì)量是研究的重要方面，掌握物體附加質(zhì)量的計(jì)算無(wú)疑具有重要的意義。附加慣性力的存在使物體在理想流體中的變速運(yùn)動(dòng)相當(dāng)于物體自身質(zhì)量上增加了一個(gè)附加質(zhì)量而在真空中運(yùn)動(dòng)，換句話說(shuō)，理想流體增大了物體的慣性，使物體很難加速也難減速。計(jì)算機(jī)是求解附加質(zhì)量的重要工具，本課程設(shè)計(jì)主要依據(jù)分布源模型的面元法等知識(shí)來(lái)對(duì)圓球、橢球、圓柱、雙橢球的附加質(zhì)量進(jìn)行數(shù)值模擬計(jì)算，并進(jìn)行相關(guān)討論。二、理論依據(jù)用表示無(wú)界流中的物體表面，來(lái)流為均勻流，其未擾動(dòng)速度或無(wú)窮遠(yuǎn)處的速度為 (2.1.1)用表示定常速度勢(shì)，它在物體外部空間域中適合拉普拉斯方程，在物面上適合不可進(jìn)入條件，在無(wú)窮遠(yuǎn)處，應(yīng)該與均勻來(lái)流的速度勢(shì)吻合，即（物體外） (2.1.2)（物面上） (2.1.3)（無(wú)窮遠(yuǎn)處）其中,單位法線向量指向物體內(nèi)部。在速度勢(shì)中分出已知的均勻來(lái)流項(xiàng)，記 (2.1.4)這里的是擾動(dòng)速度勢(shì)，應(yīng)適合以下定解條件：（物體外） (2.1.5) （物面上） (2.1.6)（無(wú)窮遠(yuǎn)處） (2.1.7)易知過(guò)物面的通量為零，即所以遠(yuǎn)方條件(2.1.7)可進(jìn)一步具體化為 （） (2.1.8)用表示點(diǎn)和之間的距離，對(duì)函數(shù)和在物面外部和遠(yuǎn)方控制面的內(nèi)部之空間域內(nèi)用格林公式，當(dāng)點(diǎn)在上述空間域內(nèi)時(shí) (2.1.9)從的遠(yuǎn)方條件(2.1.8)可知，上積分趨于零，式(2.1.9)成為 (2.1.10)其中, 是物面外的任意一點(diǎn)。在物體的內(nèi)部域中構(gòu)造一個(gè)合適的內(nèi)部解，它在內(nèi)部適合拉普拉斯方程，在物面上適合某種物面條件，其具體形式將在下面給出。對(duì)于上述物體外部的點(diǎn)函數(shù)在物體內(nèi)部域中沒(méi)有奇點(diǎn)，在內(nèi)部域中對(duì)函數(shù)和用格林公式，得到 (2.1.11)式(2.1.10)和(2.1.11)中的是物體外部同一個(gè)點(diǎn)，把兩式相減，得到 (2.1.12)在物面上取適合下述兩種物面條件，得到兩種的定解條件，一種是： (2.1.13)定解問(wèn)題(2.1.13)是拉普拉斯方程的第一類邊值問(wèn)題，它的解是存在且唯一的。取式(2.1.12)中的內(nèi)部解為式(2.1.13)所決定的函數(shù)，則式(2.1.12)成為 (2.1.14)其中 (2.1.15)式(2.1.14)表示擾動(dòng)速度勢(shì)可以用物面上的分布源表示，其中分布源密度是未知函數(shù)，將由擾動(dòng)勢(shì)的物面條件（2.1.6）來(lái)決定。物體的附加質(zhì)量，表示物體沿方向運(yùn)動(dòng)引起的方向的附加質(zhì)量，公式如下: (2.2.1)利用式(2.1.14)，再結(jié)合物面條件，得到 (2.3.1)這就是分布源密度所適合的線性積分方程。把積分方程(2.3.1)轉(zhuǎn)換成線性代數(shù)方程組，即用離散量代替連續(xù)變量。把物面分成小塊，記 (2.3.2)用平面四邊形或三角形來(lái)近似代替小曲面。具體做法如下，取第小塊的四個(gè)頂點(diǎn)坐標(biāo)之算術(shù)平均值，得到中心點(diǎn)的坐標(biāo)。計(jì)算對(duì)角線向量的向量積（指向與曲面法線指向相符合），用表示該方向上的單位向量，形成以為法線且通過(guò)中心點(diǎn)的平面，再把四個(gè)頂點(diǎn)向該平面作投影，以四個(gè)投影點(diǎn)為頂點(diǎn)組成平面四邊形，用代替原來(lái)的小曲面，稱為單元。通常把小范圍內(nèi)的分布源密度作為常數(shù)，因此只要分割不太粗，可以認(rèn)為在單元上為常數(shù)，記作，從而 (2.3.3)因此物面上的積分可以用個(gè)平面四邊形（三角形）上積分之和來(lái)近似，即 (2.3.4)上式左端的未知量是連續(xù)型變量，而上式右端的未知量是個(gè)離散量。為了求解這個(gè)未知數(shù)，須要個(gè)方程。取積分方程(2.3.1)中的動(dòng)點(diǎn)為個(gè)單元的中心點(diǎn)，稱之為控制點(diǎn)，即控制物面條件使之成立的點(diǎn)。用近似式(2.3.4)代替積分方程(2.3.1)的左端，便可以寫出的階線性代數(shù)方程組： (2.3.5)其中 稱為影響系數(shù)，即第個(gè)單元上的分布源在第個(gè)控制點(diǎn)上的影響。求解線性代數(shù)方程組(2.3.5)得到的值以后，便可以得到速度勢(shì)在控制點(diǎn)處的值，即 (2.3.6) (2.3.7)另外，物面上的誘導(dǎo)速度為 (2.3.8)其中表示求和是不計(jì)這一項(xiàng)。，這里的曲面法線指向物體內(nèi)部。三、數(shù)值模型將物體表面劃分成四邊形面面元，物面為，每一個(gè)四邊形面面元為。為了簡(jiǎn)化計(jì)算，將面網(wǎng)格投影到各自對(duì)應(yīng)的平面上，使曲面網(wǎng)格變?yōu)槠矫婢W(wǎng)格。投影的方法為： 取四個(gè)頂點(diǎn)坐標(biāo)之平均值，作為中心點(diǎn)的坐標(biāo)。計(jì)算對(duì)角線連線向量的向量積并使得積的方向與流域法向相同。用表示該方向上的單位向量。設(shè)：的坐標(biāo)為：則取投影面為過(guò)并以為法向量的平面：設(shè)在該平面上的投影點(diǎn)為：而曲面四邊形某個(gè)頂點(diǎn)為：則有：因此得到由頂點(diǎn)坐標(biāo)求解投影點(diǎn)（頂點(diǎn)）坐標(biāo)的線性方程組： 由此線性方程組可解出投影點(diǎn)坐標(biāo)。 假設(shè)速度勢(shì)和分布源在上是不變的，其值為該單元中點(diǎn)（控制點(diǎn)）處的速度勢(shì)或分布源。由于所求速度勢(shì)和速度等物理量均為物面上的物理量，因此要令點(diǎn)落在物面之上。式右端分布源的法向?qū)?shù)極限由兩部分組成，一部分是P點(diǎn)附近小曲面的貢獻(xiàn)，另一部分是屋面其余部分貢獻(xiàn)。當(dāng)所趨近于的物面上的點(diǎn)作為控制點(diǎn)的單元，積分時(shí)需要考慮奇異性；其余部分為。設(shè)其中一單元為單元，其余模型為單元。對(duì)于每一個(gè)控制點(diǎn)，令循環(huán)一次求得前述方程的積分項(xiàng)（包括奇異積分）。再由可以得到組方程，進(jìn)而形成求解各個(gè)控制點(diǎn)處物理量的矩陣。點(diǎn)表示控制點(diǎn)（編號(hào)），對(duì)于每一個(gè)控制點(diǎn)的物理量，通過(guò)在和上積分得到。即：對(duì)于任意一,，有其中， 將模型控制點(diǎn)數(shù)據(jù)導(dǎo)入到程序中計(jì)和，可以得到方程組：求解上面線性代數(shù)方程組得到的值以后，便可以得到速度勢(shì)在控制點(diǎn)處的值，即 其中，每一個(gè)i點(diǎn)（控制點(diǎn)）處誘導(dǎo)速度為（是向量）：其中，在常分布單元假設(shè)條件下： 可得，于是便可求解出速度勢(shì)和物面速度。四、幾何模型4.1橢球、圓球、有限長(zhǎng)圓柱、平行橢球利用計(jì)算機(jī)編程來(lái)完成以上四類幾何模型的建立并劃分網(wǎng)格如圖4.1、圖4.2、圖4.3、圖4.4所示。其中，橢球長(zhǎng)短軸之比為5：1，有限長(zhǎng)圓柱柱體長(zhǎng)和截面直徑之比為5，平行橢球的兩個(gè)橢球相同且長(zhǎng)軸平行，間距為短半軸的3、5、7倍。圖4.1 橢球面網(wǎng)格劃分圖4.2 圓球面網(wǎng)格劃分圖4.3 有限長(zhǎng)圓柱面網(wǎng)格劃分圖4.4 平行橢球面網(wǎng)格劃分五、計(jì)算參數(shù)及結(jié)果討論5.1橢球5.1.1橢球附加質(zhì)量系數(shù)M11經(jīng)查閱資料可知，橢球附加質(zhì)量系數(shù)M11的理論值為0.059。不同面元數(shù)時(shí)，M11的計(jì)算結(jié)果見(jiàn)表5.1。表5.1橢球的附加質(zhì)量系數(shù)M11隨面元數(shù)的變化面元數(shù)橢球的附加質(zhì)量系數(shù)M113000.062044000.061385000.061176000.061017000.060908000.060749000.06064M11隨面元數(shù)變化曲線見(jiàn)圖5.1。圖5.1 M11隨面元數(shù)變化曲線5.1.2橢球附加質(zhì)量系數(shù)M33橢球附加質(zhì)量系數(shù)M33的理論值為0.894。不同面元數(shù)時(shí)，M33的計(jì)算結(jié)果見(jiàn)表5.2。表5.2橢球的附加質(zhì)量系數(shù)M33隨面元數(shù)的變化面元數(shù)橢球的附加質(zhì)量系數(shù)M333000.97065414000.96765415000.95821456000.95187417000.94532548000.94225489000.9406142M33隨面元數(shù)變化曲線見(jiàn)圖5.2。 圖5.2 M33隨面元數(shù)變化曲線5.1.3 結(jié)果討論 由計(jì)算結(jié)果可以看出，橢球在附加質(zhì)量系數(shù)M11較M33小很多，這點(diǎn)由橢球幾何形狀可以容易看出，與實(shí)際有較好的符合。隨著面元數(shù)量的增加附加質(zhì)量系數(shù)M11和M33均更接近于各自理論值，誤差都有減小。在計(jì)算過(guò)程中發(fā)現(xiàn)，面元的劃分形式存在質(zhì)量好壞的區(qū)別，對(duì)于同樣多數(shù)量的面元，面元形式的不同會(huì)造成不一樣的結(jié)果，比如同是數(shù)量為600的單元，長(zhǎng)軸方向劃分為20份，周向30份得出的M33為0.964419，而長(zhǎng)軸方向劃分為30份，周向?yàn)?0份，得出的附加質(zhì)量M33為0.9514185，更接近與理論值，誤差更小，這是由于對(duì)于本細(xì)長(zhǎng)的橢球模型而言，周方向劃分20份已經(jīng)較為緊密，而細(xì)長(zhǎng)的長(zhǎng)軸方向需要?jiǎng)澐指嗟臄?shù)量以使面元更接近于真實(shí)物面，所以在劃分面元時(shí)，需要考慮物體的實(shí)際情況，以得到更高質(zhì)量的面元。5.2圓球5.2.1圓球附加質(zhì)量系數(shù)M11經(jīng)查資料可知，圓球附加質(zhì)量系數(shù)M11的理論值為0.5。不同面元數(shù)時(shí)，M11的計(jì)算結(jié)果見(jiàn)表5.3。表5.3球的附加質(zhì)量系數(shù)M11隨面元數(shù)的變化面元數(shù)球的附加質(zhì)量系數(shù)M113240.56372983600.54678394320.53687165400.52665476480.52387287840.5142735M11隨面元數(shù)變化曲線見(jiàn)圖5.3。圖5.3 M11隨面元數(shù)變化曲線5.2.2圓球附加質(zhì)量系數(shù)M33由于對(duì)稱性，圓球附加質(zhì)量系數(shù)M33的理論值也為0.5。不同面元數(shù)時(shí)，M33的計(jì)算結(jié)果見(jiàn)表5.4。表5.4球的附加質(zhì)量系數(shù)M33隨面元數(shù)的變化面元數(shù)球的附加質(zhì)量系數(shù)M333240.53854783600.53635444320.53374525400.53145246480.52685417840.5234581M33隨面元數(shù)變化曲線見(jiàn)圖5.4。圖5.4 M33隨面元數(shù)變化曲線5.2.3 結(jié)果討論由計(jì)算結(jié)果可以看出，隨著圓球面元數(shù)的增加，計(jì)算而得的附加質(zhì)量系數(shù)逐漸接近解析解，即理論值，誤差逐漸減小，由面元法知識(shí)可知，這是由于面元數(shù)量增加，面元密度增加，單元也就更接近真實(shí)物面，從而結(jié)果更準(zhǔn)確。M11與M33計(jì)算結(jié)果基本一致，差別較小。5.3 有限長(zhǎng)圓柱5.3.1圓柱附加質(zhì)量系數(shù)M11不同面元數(shù)時(shí)，M11的計(jì)算結(jié)果見(jiàn)下表。表5.5有限長(zhǎng)圓柱的附加質(zhì)量系數(shù)M11隨面元數(shù)的變化面元數(shù)量m112000.1378352400.1391543000.1402783600.1409234000.1412154800.1416176000.1419747200.1421888000.1422879000.142383有限長(zhǎng)圓柱M11隨面元數(shù)變化曲線圖5.5圖5.5 M11隨面元數(shù)變化曲線5.3.2圓柱附加質(zhì)量系數(shù)M33不同面元數(shù)時(shí)，圓柱M33的計(jì)算結(jié)果見(jiàn)下表。表5.6有限長(zhǎng)圓柱的附加質(zhì)量系數(shù)M11隨面元數(shù)的變化面元數(shù)量m332001.0216852400.9944273000.9688393600.9526274000.9447584800.9332416000.922057200.9147698000.9111799000.907623有限長(zhǎng)圓柱M33隨面元數(shù)變化曲線圖5.6圖5.6 M33隨面元數(shù)變化曲線5.3.3 結(jié)果討論 由計(jì)算結(jié)果可以看出，對(duì)于柱體長(zhǎng)與直徑長(zhǎng)之比為5：1的有限長(zhǎng)圓柱體而言，其沿柱體方向附加質(zhì)量系數(shù)M11要小于垂直柱體方向附加質(zhì)量系數(shù)M33，并且隨著面元數(shù)量增加而增大并且有收斂趨勢(shì)，而垂直于柱體長(zhǎng)軸方向的附加質(zhì)量系數(shù)M33隨著面元數(shù)量增加而減小并且也呈現(xiàn)出收斂趨勢(shì)。5.4平行橢球5.4.1兩橢球平行時(shí)附加質(zhì)量系數(shù)M11每個(gè)橢球面元?jiǎng)澐謹(jǐn)?shù)量為400，并且劃分形式與單個(gè)橢球時(shí)一樣，以便比較二者的區(qū)別。單個(gè)橢球在以上計(jì)算中得到的附加質(zhì)量系數(shù)M11的值為0.06149，以該值作為計(jì)算兩橢球并行時(shí)相互影響的參考值。在不同軸間距時(shí)時(shí)，M11的計(jì)算結(jié)果與對(duì)比如下。表5.7兩橢球并行時(shí)M11隨長(zhǎng)軸間距的變化間距(倍數(shù))并行時(shí)附加質(zhì)量系數(shù)M11影響30.0756224.08%50.0675410.11%70.064595.10%M11隨長(zhǎng)軸間距變化曲線見(jiàn)下圖。圖5.7兩橢球并行時(shí)M11隨長(zhǎng)軸間距的變化曲線5.4.2兩橢球平行時(shí)附加質(zhì)量系數(shù)M33單個(gè)橢球在以上計(jì)算中得到的附加質(zhì)量系數(shù)M33的值為0.9675557，以該值作為計(jì)算兩橢球并行時(shí)相互影響的參考值。不同軸間距下的算結(jié)果見(jiàn)下表。表5.8兩橢球平行時(shí)M33隨長(zhǎng)軸間距變化間距并行時(shí)附加質(zhì)量系數(shù)M33影響31.11514814.90%51.0046213.87%70.9836511.59%雙橢球并行時(shí)M33隨長(zhǎng)軸間距變化時(shí)的曲線如下，圖5.8雙橢球并行時(shí)M33隨長(zhǎng)軸間距變化曲線5.4.3 結(jié)果討論由計(jì)算結(jié)果可以看出，兩個(gè)橢球并行時(shí)，二者之間會(huì)產(chǎn)生干擾，使單個(gè)橢球附加質(zhì)量系數(shù)增加，然而其影響會(huì)隨著橢球之間軸間距增大而減小，以使每個(gè)橢球附加質(zhì)量系數(shù)逐漸接近于單個(gè)橢球運(yùn)動(dòng)時(shí)的值，可以推測(cè)當(dāng)其間距足夠大時(shí)，每個(gè)橢球可以看成單獨(dú)運(yùn)動(dòng)的效果。參考文獻(xiàn)1