李曉杰，趙春風，2

（1.大連理工大學運載工程與力學學部工業(yè)裝備結構分析國家重點實驗室，遼寧 大連 116024；2.大連理工大學工程抗震研究所，遼寧 大連 116024）

由于在滑移爆轟驅動下飛板的運動狀態(tài)直接影響爆炸焊接工藝中各種參數的選擇和焊接質量，所以研究飛板的運動規(guī)律是爆炸焊接理論研究中很重要的內容。對于穩(wěn)定的滑移驅動飛板運動，將坐標系置于爆轟波頭上，爆轟產物的流場就是定常的，可按二維可壓縮定常流處理爆轟產物的流動，因此，可以采用特征線差分法對爆轟產物流場進行數值計算。以往在對這種定常爆轟流場進行計算時，一直用多方等熵方程來描述爆轟產物［1－5］，沒有使用JWL和簡化JWL等描述炸藥爆轟的專用狀態(tài)方程。這不僅影響爆轟流場的計算精確度，而且影響對爆轟驅動力學實質的分析。為此，本文中力求在使用通用狀態(tài)方程的基礎上，由二維可壓縮定常流的基本理論導出特征線方程，并針對二維滑移爆轟飛板拋擲問題重新編制計算程序，以常用的JWL狀態(tài)為例進行數值計算。沿用文獻［1］中的爆轟波后采用的弱爆轟假設來處理聲速面問題。另外，采用一種直接從物理概念出發(fā)獲得等熵無旋定常超聲速二維流動特征線方程的方法，而沒有使用通常的由控制方程出發(fā)，采用不定線法或方向導數法來推導特征線方程的方法［6］，以期更明晰地描述特征線的物理意義。

1 二維定常超聲速流的特征線方程

在均勻的超聲速流場中，處于流場中某一固定點的小擾動源，所產生的擾動如圖1所示，球對稱地以聲速c相對流體向外傳播；流體本身以速度q運動，小擾動傳播的絕對速度是聲速c與流速q的矢量迭加。對于超聲速流，有q＞c，擾動面是一個錐面，且與從擾動點出發(fā)的一串球面相切，該錐面為馬赫錐，馬赫錐的邊界面（線）稱為馬赫面（線）。半錐角α為［1，5］

通常把半錐角α稱為馬赫角，Ma稱為馬赫數。

圖1 超聲速流中的小擾動Fig.1Aperturbation in supersonic flow

對于本文中要研究的平面問題，擾動源是無限長的線源，沿擾動線一系列馬赫錐組成一個楔形。在二維平面上，如圖2所示，馬赫錐組成的楔形退化為2條馬赫線I和II，在Oxy坐標系中，這2條馬赫線的方程為［1，5］

式中：θ是流動方向角。如果馬赫波后的小擾動使流速由q變?yōu)閝＋δq，流動方向角由θ變?yōu)棣龋摩龋瑒t流動穿過馬赫波I時，沿馬赫波的平行方向的速度動量方程有

圖2 超聲速流中的特征線Fig.2Characteristic curves in supersonic flow

式中：Qm＝ρqsinα為通過馬赫波上的質量流量。展開式（3），并考慮到δq→dq，δθ→dθ，且cosδθ→1，sinδθ→dθ，并忽略二階小量，將式（1）代入式（3）可得

同理，對穿過馬赫線II的參數變化可以寫出

從數學意義上來講，馬赫線就是特征線，馬赫線I和II分別對應第1族和第2族特征線，式（2）就是這2族特征線的方程，而式（4）和（5）分別是波陣面前、后的相容關系。任何復雜的超聲速流場均可以被看成是由一系列連續(xù)的空間小擾動相互作用的結果，所以，常利用2族特征線方程和對應的相容關系來求解超聲速流場。在進行數值求解時，一般使用沿特征線的相容關系。由于在超聲速膨脹流中同族特征線不相交，一族特征線穿過另一族特征線，因此，一族特征線波陣面上的相容關系就是沿另一族特征線上的相容關系，所以可以獲得沿第1、2族特征線上的相容關系為［1，5］

從物理概念出發(fā)獲得的平面等熵定常超聲速流動的特征線方程（2）和（6），由于方程式不依賴流體物態(tài)方程相關的參數，所以式（2）和（6）是通用物態(tài)方程的特征線關系。

然后，可以利用能量方程和具體的氣體等熵方程來建立馬赫數Ma和流速q與流動方向θ之間的關系。對于定?？蓧嚎s流，用熱焓i和流速q的關系表示其伯努利方程［5，7］

式中：熱焓i＝e＋p／ρ，e為比內能，p為爆壓，ρ為炸藥密度；i0為滯止焓。微分式（7）可得

再利用Ma2＝ （q／c）2＝2 （i0－i）／c2，等熵條件di＝dp／ρ和聲速c2＝（?p／ ?ρ）S＝dp／dρ，由上面推導可以得到下面的關系

由于在等熵線上，流體的熱焓i、聲速c，包括馬赫數Ma都是密度ρ的單值函數，所以式（9）是可以積分的，故可以將式（9）的積分定義成一個函數

式中：函數ν（ρ）相當于Prandtl－Meyer（P－M）函數［6］，ρCJ為炸藥CJ狀態(tài)時的密度。

將式（10）定義為不依賴流體狀態(tài)方程的P－M函數，或稱為通用狀態(tài)方程的P－M函數。在式（6）中引入該函數寫成新的特征線相容關系，再與式（2）合并，可以寫成如下超聲速定常流特征線方程組

式中：RⅠ和RⅡ是常數，分別為第1和第2黎曼不變量。顯然，只要根據具體的流體等熵狀態(tài)方程求出熱焓i、聲速c代入式（10）求得P－M函數ν（ρ），就可以和式（11）組成封閉的求解方程組了。如對于常用的JWL炸藥狀態(tài)方程，其等熵方程如下

式中：pS為等熵壓力，eS為等熵比內能，ρ0為初始密度，A、B、C、R1、R2和ω為常數，V＝ρ0／ρ是相對比體積。根據相應的定義可以推導出JWL狀態(tài)方程的聲速和熱焓表達式

由于特征線方程（11）為常微分方程組，所包含的P－M函數ν（ρ）可以用辛普森積分等高精度積分方法求解，所以使用特征線方程可以獲得高精度的差分解。

2 滑移爆轟驅動飛板的特征線差分法

在確立了平面超聲速流的特征線方程后，可將特征線方程離散化，對圖3所示的滑移爆轟驅動飛板問題進行數值計算。

2.1 初值選取

如圖3所示，炸藥向空氣飛散在O點處可以近似用中心稀疏波處理。由于在爆轟波CJ面上的馬赫數MaCJ＝1，對應的密度ρ＝ρCJ，氣流轉角θ＝0，根據普朗特繞流理論的第2族特征線相容關系式［7］

圖3 滑移爆轟作用下的飛板運動姿態(tài)Fig.3Movement of flyer plate driven by glancing detonation

可以在ρCJ到某一較小ρn范圍內給出一系列的離散點（ρi，θi），其中i＝0，1，2，…，n；離散點的（x，y）坐標均在O點處，這些離散點就是差分的初值點。再利用下式確定滯止焓

2.2 聲速面處理

由于在爆轟產物流場中，爆轟波面后的CJ面是聲速面，從O點發(fā)出的馬赫線會在飛板壁面上垂直反射，導致無法計算差分值。具體的處理方法可以簡單地令炸藥爆轟稍向弱爆轟方向偏離CJ點［5］，只對整個流場產生微小的影響。即令＝ （1－Δ）ρCJ，使爆轟波后完全變?yōu)槌曀倭鲌?，就完全可使用特征線差分求解。這樣做引起的誤差也很小，把Δ取足夠小量，對計算結果影響不大。在實際計算中，Δ可取為10－2～10－3量級，對量綱一壓力的擾動為10－4～10－6量級。

2.3 超聲速流場域內差分

在爆轟產物的超聲速流場內部，如圖4所示，當已知點1、2的流動參數時，可以用2族特征線相交求取出未知點3上的參數。由式（11）將第1族特征線RⅠ方程由已知點1至未知點3離散成差分形式，第2族特征線RⅡ方程由已知點2至未知點3離散，經聯(lián)立整理可得

圖4 爆轟產物中差分示意圖Fig.4Scheme of differences in detonation products

式中：下標1～3分別表示點1～3。按式（16）由點1、2的參數差分解出θ3，進而用P－M反函數求得ρ3。然后，利用狀態(tài)方程，求解最后代入式（17）確定點3的坐標。

2.4 飛板邊界差分

飛板在爆轟壓力驅動下的運動微分方程為［6］

式中：R為炸藥質量比，即炸藥的厚度和密度的乘積與飛板的厚度和密度乘積的比值；＝p／pCJ，pCJ為等熵狀態(tài)的爆轟壓力；D為炸藥爆速；pk為中間變量；x、y為藥厚量綱一化的坐標。如圖5所示，當已知飛板邊界點1和流場域內點2時，可利用飛板運動方程（18），與通過點2的第1族特征線求出點3的參數。將上述邊界方程作如下離散

圖5 飛板邊界差分示意圖Fig.5Scheme of differences on the boundary of flyer plate

對通過點2、3的第1族特征線方程也離散成差分形式

將式（19）代入式（20），整理得

預估方程組（21）可以用預估－校正法獲得較高的差分精度。具體做法是，以k＊＝cot（θ2＋α2）和pk＝Rp1／（ρ0D2）作為預估初值代入式（21）解出θ3、ρ3，并用狀態(tài)方程解出對應的點3的i3、c3、Ma3和α3；再用k＊＝ ［cot（θ2＋α2）＋cot（θ3＋α3）］／2和pk＝R（p1＋p3）／（2ρ0D2）代入式（21），重新求解一遍，可獲得經過校正的點3的參數。最后，用式（19）求出點3的坐標。

2.5 實例計算

用上述通用狀態(tài)方程的差分格式和方法，編制計算程序對飛板二維拋擲問題進行數值求解。計算中TNT炸藥的JWL狀態(tài)方程參數［8］：多方指數K，2.727 6；CJ狀態(tài)的密度，1.630g／cm3；CJ狀態(tài)的爆壓，21.0GPa；CJ狀態(tài)的爆速，6.93km／s；A，373.8GPa；B，2.747GPa；C，0.734GPa；R1，4.15；R2，0.90；ω，0.35。乳化炸藥的JWL狀態(tài)方程參數［9］：多方指數K，3.008；CJ狀態(tài)的密度，1.145g／cm3；CJ狀態(tài)的爆壓，7.62GPa；CJ狀態(tài)的爆速，5.141 4km／s；A，326.42GPa；B，5.808 9GPa；C，1.032 5GPa；R1，5.8；R2，1.56；ω，0.57。圖6為所計算的 TNT和乳化炸藥滑移爆轟驅動飛板的拋擲曲線，圖7為飛板的動態(tài)彎折角與量綱一飛行距離的關系，圖8為TNT、乳化炸藥的JWL和多方方程的等熵卸載線。

由圖6～7可以看出，對于TNT炸藥，用多方狀態(tài)方程描述的對飛板的加速能力較用JWL方程描述的大；取量綱一距離飛行距離y＝1處的彎折角進行對比，分別為11.96°和10.57°，用JWL狀態(tài)方程的計算結果較用多方方程的計算結果小11.62%。對于乳化炸藥，用多方狀態(tài)方程描述的對飛板的加速能力較用JWL方程描述的小；取量綱一距離飛行距離y＝1處的彎折角對比，分別為11.20°和12.31°，用JWL狀態(tài)方程的計算結果較用多方方程的計算結果大9.91%。這與圖8等熵膨脹線所表示的炸藥膨脹做功能力完全一致，即TNT的JWL等熵膨脹線在多方方程等熵線的下方，膨脹做功能力較??；乳化炸藥的JWL等熵線居于多方方程等熵線的上方，做功驅動能力較大。由此可見，通用狀態(tài)方程特征線差分方法能夠正確地反映飛板的爆轟驅動過程。

圖6 爆轟驅動飛板的拋擲曲線Fig.6Flying curves of flyer plate driven by detonation

圖7 飛板量綱一飛行距離與彎折角的關系Fig.7Relationship of dimensionless flying distance and bending angle for flyer plate

圖8 JWL等熵線與多方等熵線的對比Fig.8Comparision of JWL and polytropic isentropic curves

3 結 論

（1）從馬赫波的物理概念出發(fā)，推導出了平面二維超聲速定常流的特征線方程；并且重新定義了以流體密度為單一自變量的Prantl－Meyer函數ν（ρ），從而構成了求解平面二維超聲速定常流的封閉方程組。由于所獲得的特征線微分方程組不依賴流體狀態(tài)方程的表達形式，所以可以用于流體各種形式狀態(tài)方程的求解，是一種通用物態(tài)方程的超聲速定常流特征線求解方法。

（2）為了檢驗通用物態(tài)方程特征線解法，針對滑移爆轟驅動飛板運動問題構建了爆轟產物流場內部和飛板邊界特征線差分法格式，對爆轟波頭聲速面的處理采用了向超聲速的弱爆轟攝動的處理方法。

（3）對TNT炸藥和乳化炸藥采用JWL狀態(tài)方程與多方方程進行對比計算的結果表明，這種通用狀態(tài)方程特征線差分方法可以準確地反映炸藥狀態(tài)方程對飛板拋擲的作用，飛板的爆轟驅動過程與其狀態(tài)方程所表述的做功驅動能力一致。

最后應該說明的是，特征線差分方法所求解的方程組本質上是常微分方程組，所以完全可以參照常微分方程組的數值求解方法（如歐拉預估－校正法）獲得高精度特征線差分格式，從而提高定常爆轟驅動問題的求解精度。另外，由于特征線差分方法計算簡單，并且可以很好地反映滑移爆轟驅動問題的實質，如果將通用狀態(tài)方程特征線差分程序作為參數優(yōu)化程序的內嵌，可以構成精確簡單的炸藥狀態(tài)方程參數擬合程序，簡化炸藥狀態(tài)方程的參數擬合計算過程。

［1］李曉杰.滑移爆轟下帶覆蓋物的飛板拋擲及對切割的影響［D］.大連：大連理工大學，1987.

［2］Besshaposhnikov Y P.Motion of a plate thrown by aglancing detonation wave［J］.Combustion，Explosion，and Shock Waves，1999，35（1）：109－112.

［3］Ivanov A G，Kochkin L I，Ogorodnikov V A，et al.Characteristics of the acceleration of plates by aglancing detonation wave in the presence of an additional or concentrated mass［J］.Combustion，Explosion，and Shock Waves，1990，26（5）：612－614.

［4］邵丙璜，張凱.爆炸焊接原理及其工程應用［M］.大連：大連工學院出版社，1987.

［5］呂洪生，曾新吾.連續(xù)介質力學（中冊）：流體力學與爆炸力學［M］.長沙：國防科技大學出版社，1999.

［6］Courant R，F(xiàn)riedrichs K O.Supersonic flow and shock waves［M］.5th Edition.New York：Springer－Verlag，1976：37－75.

［7］趙春風.基于通用狀態(tài)方程的爆炸焊接特征線法研究［D］.大連：大連理工大學，2010.

［8］徐錫申，張萬箱.實用物態(tài)方程理論導引［M］.北京：科學出版社，1986：402－403.

［9］宋錦泉.乳化炸藥爆轟特性研究［D］.北京：北京科技大學，2000.