PMT水下內爆沖擊波強度與傳播特性研究

金 鍵， 侯海量， 陳鵬宇， 朱 錫， 吳林杰， 何 苗, 王天穹

(1.海軍工程大學 艦船工程系，武漢 430033；2.中國科學院 高能物理研究所，北京 100049)

中微子是組成自然界的基本粒子，它不僅是物理學家們探索宇宙形成和演化的途徑，同時，其不帶電、微質量、光速傳遞等特點使其有望作為理想信息傳遞載體在通訊、軍事等領域帶來革命性的突破。目前，世界上共有27個主要的中微子探測試驗，我國在大亞灣反應堆中微子實驗中，發現了中微子的第三種振蕩模式，首次精確地測定了中微子的混合角，被譽為基礎物理學的一大突破[1]。大亞灣實驗二期工程即江門中微子實驗(在建)旨在解決中微子質量譜排序問題[2]。

光電倍增管(Photo-Mutiplier Tube,PMT)是用于探測宇宙中微子的核心部件，2萬只PMT依附在水中的中心探測器上構成了光信號的探測構件。2001年日本超級神岡中微子探測器發生PMT內爆事故[3]，其中心探測器底部一個PMT發生內爆后，沖擊波導致7 877個PMT發生連環爆炸，爆炸直接損失超過3 000萬美元，延誤實驗時間4年。這次事故使各國大型水下中微子探測工程開始重視水下PMT的防爆必要性。

PMT是由透明玻璃吹制而成，內部接近真空，工作狀態下承受靜水壓力。PMT的意外破壞有多種原因：①玻璃表面在靜水壓下發生蠕變；②工作周期內水池反復的排水、注水致使PMT承壓狀態變化，從而導致玻殼發生疲勞；③水中微生物對玻殼的侵蝕；④不可避免的初始缺陷以及安裝過程中的磨損、擦傷。當PMT在靜水壓力下發生內爆后，產生的沖擊波在向四周傳播的過程中有可能擊破相鄰的PMT，從而導致水中的所有PMT發生連鎖內爆。Ling等[4]利用有限元仿真的方法模擬了多個PMT水下殉爆過程，說明了初始內爆引發水下PMT連環內爆的可能性。

本文針對江門中微子實驗中的PMT防爆工作，對PMT水下內爆沖擊波的產生、氣泡壓潰時間、沖擊波的強度以及沖擊波的傳播規律開展了理論分析和數值仿真研究。

1 工程背景

江門中微子中心探測器安裝于直徑43 m、深45 m的圓柱型水池中(如圖1所示)；中心探測器外層與水池之間注滿超純水，密度為1 000 kg/m3；20 000只PMT密布在圖1所示的不銹鋼網架和水池內壁上，水池內的PMT承受0.1～0.54 MPa的靜水壓力。

圖1 中心探測器結構

PMT實物圖及外形尺寸如圖2圖所示。PMT由玻璃外殼(簡稱玻殼)及內部電氣原件組成。玻殼厚度約為4～6 mm。內部電氣結構為金屬材料，工作在真空之中，內部真空度為10-5Pa以上。

(a) PMT示意圖

(b) PMT外形尺寸圖(mm)

圖2 PMT示意圖及外形尺寸圖

Fig.2 Diagram and dimension of PMT

2 球形容器水下內爆理論分析

球形容器水下內爆過程類似于水下爆炸后的氣泡脈動過程和水下空泡破裂過程。外部水壓高于氣泡內部氣體壓力，一旦平衡狀態被打破，外部高壓水迅速向中心匯涌，壓縮氣泡內部氣體，直到氣泡達到最小半徑之后，氣泡開始膨脹并向外輻射沖擊波。有所區別的是氣泡形成的方式，水下爆炸產生的氣泡是由爆轟產物的混合物形成的；典型水下空泡由相變、成核與空化產生[5]；而本例中則由容器整體破壞產生。

在對氣泡動力學特性進行分析時，通常以無限均勻溫度液體區域內的單個球形對稱氣泡為研究對象，而熱效應、浮力、液體的黏性、液體的可壓縮性和空泡表面張力等則是影響因素，在具體研究中為了簡化問題往往忽略一些因素。Rayleigh[6]首先推導并應用了不考慮表面張力和黏性的空泡動力學方程；Plesset[7]在Rayleigh工作的基礎上得到了廣義的空泡動力學R-P方程(式(1))。

(1)

式中:R為球形氣泡半徑；Pi是氣泡內初始壓力;p∞為距離氣泡無窮遠處的靜壓力;μ,ρ分別為液體的黏性和密度；σ為表面張力。

(2)

當氣體含量很小時，可近似認為氣泡半徑被壓縮為0，從而得到了氣泡由R=R0到R=0完全壓縮所需的時間為

(3)

在對氣泡回彈階段沖擊波形成的研究中，Ivany等[8]通過研究證實表面張力和黏性對沖擊波沒有重要影響。Hicking等[9]利用可壓縮流動方程的數值解探究了氣泡回彈階段沖擊波的形成過程，并指出液體的可壓縮性對氣泡動力學特性的影響是輕微的，其主要作用是在氣泡破裂后回彈階段中對沖擊波的形成，得到了壓力脈沖遠離氣泡傳播時呈現出近似幾何形狀衰減的結論。

研究球形容器水下內爆可以利用氣泡動力學的相關理論，但也有其特殊性。從能量的角度分析，由于容器薄壁的存在，容器空腔蓄積的勢能一部分需要用來打破容器壁并加速碎片向內運動，這必然會使容器水下內爆產生的沖擊波峰值和比沖量比相同尺寸、相同氣體壓力的氣泡破裂所產生的沖擊波峰值和比沖量要小。另外，Turner[10]通過將玻璃球水下內爆實驗與數值仿真結果進行對比，得出了容器結構的破壞時間對壓力峰值具有很大影響的結論。

3 數值模擬

由于水池深45 m，對整個水池進行全尺寸建模將會帶來巨大的網格數及計算量，因此只選取45 m水深處PMT的工作環境進行建模(如圖3所示)；又因為模型是軸對稱的，為進一步減少計算量，僅建立四分之一模型。有限元計算模型使用Patran進行幾何模型建模，采用kg-mm-ms單位制。

圖3 有限元計算模型示意圖

建模時，以PMT上半球最大半徑處幾何中心作為原點，在X軸方向上，0～8 000 mm劃分100個發散型網格；Y軸方向上，-8 000～0 mm、0～8 000 mm分別劃分100個發散型網格；Z軸方向上的網格劃分與X軸相同；水域外圍邊界設置為固支邊界，中間設置為對稱邊界。

用稀薄空氣來代替PMT內部真空度，假設空氣為理想氣體，其狀態方程如下

P=(γ-1)ρe

(4)

式中:P為空氣壓力，設為10-5Pa；ρ為空氣密度；e為內能，設為1.972 7×105J/kg；γ為比熱，設為1.4。

由式(4)可知，要使PMT內部真空度達到10-5Pa，對應的內部稀薄空氣密度應為1.28×10-10kg/m3。

水的狀態方程采用Mie-Gruneisen方程

(5)

(6)

式中:p為水的壓力；μ為水的壓縮度；V0為水的相對初始體積；ρ0為水的密度，取1 000 kg/m3；c為vs-vp曲線的斜率(聲速)，取1 484 m/s；γ0為格呂內森參數，取0.11；a為γ0和μ的一階體積修正量，取3；S1、S2、S3分別為曲線擬合參數，取S1=1.979、S2=S3=0。E為水的初始體積內能，由于水幾乎不可壓縮，取值設為0。

基于安全性考慮，選取水池最深45 m處為初始水壓，即0.54 MPa。通過對式(5)、(6)的計算可以得到v0=0.999 754 6。

PMT玻璃外殼采用殼單元建立(如圖4所示)，利用這些殼單元組成的Part對內部空間進行初始化，使用LS-DYNA中的*INITIAL_VOLUME_FRAC-TION_GEOMETRY關鍵字將殼單元內部包裹空間初始化成稀薄空氣。

圖4 PMT玻殼有限元模型

為了在計算開始時刪除圖4所示的殼單元以便水涌入空穴中，在LS-DYNA的材料模型中添加*MAT_ADD_EROSION關鍵字，并在FAILTM中加入失效時間。該卡片可以讓不易失效的材料在計算開始后的某個時刻強制失效，本例中將失效時間設為0.001 ms。

有限元模型中的流體全部為歐拉單元，利用LS-DYNA程序的ALE算法對模型進行計算。

4 算法驗證

為了驗證該算法的準確性，使用該算法對Turner做的玻璃球水下內爆實進行數值模擬計算。實驗工況如圖5所示，玻璃球直徑762 mm、內部壓力為0.101 3 MPa的玻璃球在6.996 MPa的靜水壓力下內爆。三個夾角120°布置的自由場壓力傳感器距球心都為101.6 mm。

圖5 玻璃球水下內爆實驗布置圖

Fig.5 Experimental arrangement of glass implosion on underwater

建立有限元模型時，水的狀態方程中的初始體積V0設為0.996 843，空氣的密度ρ設為1.29×10-9kg/m3，球形空氣半徑設為37.3 mm(玻璃球厚度0.8 mm)。數值仿真的計算結果(距離球心101.6 mm處某單元壓力、超壓比沖量時程曲線)如圖6、7所示。

圖6 壓力時程曲線

圖7 超壓比沖量時程曲線

計算開始后，從氣泡表面徑向向外輻射稀疏波，稀疏波經過的區域，其初始靜壓會突然下降，如圖6所示，大約0.04 ms時刻壓力大幅下降了3.039 MPa。隨著高壓水進一步壓縮氣泡，壓力值又緩慢上升，直到氣泡回彈階段，出現了較大的壓力峰值。數值仿真結果和實驗所測均值列于表1。

表1實驗與數值仿真結果對比

Tab.1Comparisionbetweenexperimentresultandnumericalsimulationresult

壓力峰值/MPa超壓比沖量/(Pa·s)仿真值24.47672.14實驗值26.10601.33誤差-6.24%11.75%

對比表1所示結果，就超壓比沖量而言，仿真結果要略大于實驗值，而壓力峰值的仿真結果卻小于實驗值。這是因為波陣面的厚度非常小，約為10-5～10-6cm[11]，而有限元軟件在輸出壓力結果時是以一個單元上的平均壓力輸出，單元長度又遠遠大于波陣面厚度，從而導致壓力峰值的仿真值小于實際值。另外，由于該實驗中無法準確獲得容器壓潰時間，僅對氣泡壓潰時間的仿真值與理論值進行對比。氣泡壓潰時間的仿真值為0.406 ms，根據方程4得到的半徑37.3 mm氣泡壓潰時間為0.411 ms，兩者誤差1.2%。總體而言，數值仿真計算結果與實驗值的誤差控制在一定范圍之內，兩者吻合較好。

5 PMT水下內爆數值仿真結果及分析

對PMT在0.54 MPa靜水壓力下內爆過程進行數值仿真計算時，以單元中心到原點的距離作為衡量值，選取X+,Y+,Y-三個方向上不同距離處的6個單元作為壓力輸出點，各點壓力時程曲線及比沖量時程曲線如圖8～圖13所示。

圖8 X+方向上各點壓力時程曲線

Fig.8 Pressure time history for orientation ofX+

圖9 X+方向上各點比沖量時程曲線

圖10 Y+方向上各點壓力時程曲線

Fig.10 Pressure time history for orientation ofY+

從壓力時程曲線來看，Y+方向上155 mm處單元峰值壓力發生在9.222 ms，Y-方向上136 mm處單元峰值壓力發生在9.237 ms，假設沖擊波在水中的以聲速1 484 m/s向外均勻傳播，根據兩者壓力峰值到達時間差得到PMT內爆中心應在PMT幾何中心上方20.56 mm處。這時，X+方向上141 mm處單元距內爆中心約142.5 mm，該單元壓力峰值發生在9.227 5 ms，相比Y方向上沖擊波的傳播速度，X方向上沖擊波的傳播速度比其僅大0.4%，因此，可以認為PMT內爆中心坐標為(0，20.56，0)。

根據式(4)可以計算出半徑250 mm球狀氣泡潰滅的理論時長為9.844 ms；在數值仿真計算中，距內爆中心142.5 mm處壓力峰值發生在9.227 5 ms，以聲速1 484 m/s作為沖擊波傳播速度，則內爆中心的向外傳播沖擊波的時間為9.13 ms；數值仿真得到的氣泡潰滅時長比理論時長大7.8%，兩者基本吻合。

圖11 Y+方向上各點比沖量時程曲線

Fig.11 Impulse time history for orientation ofY+

圖12 Y-方向上各點壓力時程曲線

Fig.12 Pressure time history for orientation ofY-

圖13 Y-方向上各點比沖量時程曲線

Fig.13 Impulse time history for orientation ofY-

由于PMT外形輪廓不是球形，導致其內爆產生的沖擊波并不是以球形波陣面向外傳播，沖擊波在X+、Y+、Y-三個方向上的傳播存在一定的差異。通過對三個方向上距內爆中心不同距離處的壓力峰值、超壓比沖量進行曲線擬合(如圖14、15所示)，可以較為直觀對其進行比較分析。

從圖14、15可以看出：①對于壓力峰值擬合曲線，Y方向上的兩個冪函數擬合曲線與測點吻合較好，且兩個擬合曲線的整體趨勢基本相同；而X+方向上的冪函數擬合曲線與測量點的偏差較大，這是因為X+方向所取的測點單元并不與PMT內爆中心在一條直線上，且沖擊波向外傳播的能量在空間上是不均勻的；②對于超壓比沖量擬合曲線，各數據點與擬合曲線吻合較好；③距內爆中心相同距離處，Y-方向上壓力峰值、超壓比沖量明顯比另外兩個方向上的壓力峰值要大；④PMT內爆后產生的沖擊波在遠離氣泡傳播時呈現出冪函數形式的衰減。

圖14 壓力峰值擬合曲線

圖15 超壓比沖量擬合曲線

6 結 論

通過對PMT水下內爆現象的分析，結合氣泡動力學理論、數值仿真和算例驗證，得到以下結論：

(1) 數值仿真得到的沖擊波壓力峰值和超壓比沖量與實驗算例結果吻合較好，數值仿真得到的氣泡壓潰時間與理論值吻合較好，說明本文提供的數值算法對于容器水下內爆模型具有適用性和準確性。

(2) PMT外形輪廓對沖擊波的傳播具有一定的影響，PMT水下內爆后向下方傳播的沖擊波能最大，沖擊波在遠離氣泡傳播時呈現冪函數形式的衰減。