栗志杰，由小川，柳占立，杜智博，張 仡，楊 策，莊 茁

（1. 清華大學(xué)航天航空學(xué)院，北京 100084；2. 陸軍軍醫(yī)大學(xué)大坪醫(yī)院，重慶 400038）

現(xiàn)代局部戰(zhàn)爭(zhēng)戰(zhàn)傷分析顯示，由爆炸導(dǎo)致的單兵傷亡約占總傷亡人數(shù)的70%，其中沖擊波直接作用引起的傷害約占60%[1]。同時(shí)，國(guó)際上爆炸恐怖襲擊日益猖獗，嚴(yán)重危害民眾的生命安全。而作為人體中樞神經(jīng)系統(tǒng)的頭部，是爆炸沖擊波的重要靶器官。據(jù)相關(guān)報(bào)告顯示，美軍中受到創(chuàng)傷性腦損傷（traumatic brain injury, TBI）影響的軍人為總服役人數(shù)的20%，這其中由于爆炸沖擊波所引起的腦損傷(blast-induced traumatic brain injury, b-TBI)占比高達(dá)40%～60%[1]。爆炸沖擊波引起顱腦損傷的途徑主要包括沖擊波對(duì)頭部的直接作用、飛起破片對(duì)頭部的沖擊作用、頭部與其他物體的撞擊作用、沖擊波對(duì)頭部的電磁熱等作用[2-3]。隨著爆炸傷的日漸頻發(fā)，爆炸沖擊波作用下的顱腦致傷機(jī)理與對(duì)應(yīng)的防護(hù)策略成為研究熱點(diǎn)。關(guān)于顱腦致傷機(jī)理，當(dāng)前的研究主要集中在沖擊波對(duì)頭部直接作用所造成的顱腦損傷。

雖然，經(jīng)過(guò)大量的動(dòng)物實(shí)驗(yàn)、病理生理學(xué)分析以及數(shù)值模擬研究，建立了顱腦碰撞損傷與頭部加速度之間的聯(lián)系，并確定了相應(yīng)的損傷準(zhǔn)則[4-6]，而相關(guān)理論并不能應(yīng)用于爆炸沖擊波直接作用下的顱腦損傷分析。Cernak 等[7]在生物激波管實(shí)驗(yàn)中，將老鼠進(jìn)行固定，保證頭部在沖擊波條件下為穩(wěn)定狀態(tài)，但是實(shí)驗(yàn)結(jié)果表明，老鼠的腦組織仍然出現(xiàn)了器質(zhì)性病變。費(fèi)周等利用不同尺寸的微型球形炸藥，對(duì)大鼠進(jìn)行不同損傷程度的b-TBI 實(shí)驗(yàn)，并從神經(jīng)功能、宏觀病理學(xué)、組織病理學(xué)及各種生物標(biāo)志物表達(dá)水平方面，建立了可靠且可重復(fù)的b-TBI 模型及損傷特征指標(biāo)[8]。該實(shí)驗(yàn)將大鼠放置在泡沫墊塊上，并將球形炸藥懸空掛于大鼠頭部上面，大鼠頭部在爆炸沖擊波作用下是相對(duì)穩(wěn)定的，但是仍然出現(xiàn)了腦組織損傷。這些動(dòng)物實(shí)驗(yàn)均表明，爆炸沖擊波對(duì)顱腦的損傷可以通過(guò)直接作用產(chǎn)生。而爆炸沖擊波對(duì)頭部直接作用的時(shí)間通常在1～3 ms，是應(yīng)力波主導(dǎo)的物理現(xiàn)象。因此，基于頭部加速度建立的顱腦碰撞損傷理論無(wú)法解釋上述實(shí)驗(yàn)現(xiàn)象。

針對(duì)爆炸沖擊波直接作用下的顱腦損傷，學(xué)者們做了大量的研究工作，并提出了多種致傷機(jī)理。Cloots 等[9-11]將微觀體積單元（critical volume element, CVE）與常規(guī)有限元單元進(jìn)行結(jié)合，建立了腦組織的多尺度模型，從而得到腦組織細(xì)胞層次的力學(xué)響應(yīng)，并結(jié)合生物細(xì)胞的力學(xué)損傷實(shí)驗(yàn)，建立了相應(yīng)的損傷準(zhǔn)則。Radovitzky 等[12]、Jean 等[13]發(fā)展了腦組織的多尺度生物力學(xué)模型，研究了爆炸沖擊波作用下顱腦在組織、細(xì)胞等不同尺度的生物響應(yīng)特征，并建立了評(píng)估b-TBI 風(fēng)險(xiǎn)的動(dòng)物與人之間的尺度律。Goeller 等[14]、Salzar 等[15]提出腦脊液空穴致傷機(jī)制，通過(guò)激波管實(shí)驗(yàn)與數(shù)值模擬，對(duì)腦脊液的空穴效應(yīng)進(jìn)行了研究，研究結(jié)果表明，腦脊液在沖擊波作用下會(huì)發(fā)生空化與湮滅，形成局部點(diǎn)高壓，從而造成腦組織損傷。Franck[16]基于類腦軟材料，研究了微型空穴效應(yīng)在b-TBI 當(dāng)中的作用，并利用實(shí)時(shí)成像技術(shù)對(duì)不同應(yīng)變率下神經(jīng)元細(xì)胞的損傷模式進(jìn)行了探究。Bhattacharjee[17]提出主動(dòng)脈沖壓致傷理論，沖擊波通過(guò)壓縮主動(dòng)脈血液，從而使其形成高壓液體，并經(jīng)毛細(xì)血管對(duì)腦組織形成沖擊，造成腦損傷。同時(shí)，Courtney 等[18]設(shè)計(jì)了動(dòng)物實(shí)驗(yàn)，驗(yàn)證了沖擊波可以通過(guò)作用軀干部位間接對(duì)腦組織造成損傷。Moss 等[19]通過(guò)數(shù)值模擬研究發(fā)現(xiàn)，顱骨在爆炸沖擊波作用下會(huì)發(fā)生局部彎曲，從而使腦組織壓力出現(xiàn)較大的壓力梯度，顱骨的彎曲變形與顱內(nèi)壓力的分布密切相關(guān)。但是，顱骨局部彎曲變形與腦組織壓力之間的內(nèi)在聯(lián)系有待于進(jìn)一步確立，相應(yīng)的顱腦損傷機(jī)理仍需進(jìn)一步研究。

基于上述研究現(xiàn)狀，本文中針對(duì)顱骨局部彎曲變形這個(gè)致傷機(jī)制，開(kāi)展進(jìn)一步的深入研究。首先，建立充分反映顱腦生理結(jié)構(gòu)的三維頭部數(shù)值模型，并利用頭部碰撞試驗(yàn)數(shù)據(jù)對(duì)頭部模型的有效性進(jìn)行驗(yàn)證。然后，基于歐拉-拉格朗日耦合法（Euler-Lagrangian coupling method，CEL）建立爆炸沖擊波-頭部流固耦合模型，對(duì)頭部在正面沖擊工況下的動(dòng)態(tài)響應(yīng)過(guò)程進(jìn)行數(shù)值模擬，確立顱骨局部彎曲變形與腦組織壓力之間的內(nèi)在聯(lián)系，初步揭示腦組織的損傷機(jī)理。

1 三維頭部數(shù)值模型與驗(yàn)證

1.1 三維頭部有限元模型

在顱腦結(jié)構(gòu)中，腦組織包括大腦、小腦與腦干3 部分，其表面有3 層被膜包裹，由內(nèi)向外依次為軟腦膜、蛛網(wǎng)膜和硬腦膜。軟腦膜薄而富含血管，緊貼腦表面并深入腦的溝裂中。蛛網(wǎng)膜是半透明結(jié)締組織薄膜，它與軟腦膜之間是蛛網(wǎng)膜下腔，內(nèi)含腦脊液（cerebral spinal fluid, CSF）、蛛網(wǎng)膜小梁和較大血管。硬腦膜呈套狀包被腦，在顱腔內(nèi)分化為大腦鐮與小腦幕，兩者將顱腔分為左、右與下3 個(gè)腔體，限制腦組織相對(duì)于顱骨的過(guò)大滑移。顱骨與硬腦膜外側(cè)緊密結(jié)合，其內(nèi)、外為剛度與密度較大的骨密質(zhì)，中間為骨松質(zhì)，為典型的類三明治結(jié)構(gòu)。圖1 顯示了人體頭部的解剖結(jié)構(gòu)。

圖 1 頭部解剖結(jié)構(gòu)示意圖[20]Fig. 1 Schematic diagram for anatomical structures of head [20]

基于顱腦的核磁共振切片，建立了三維頭部有限元模型，該模型包括腦組織（大腦、小腦與腦干），膜結(jié)構(gòu)（軟腦膜、蛛網(wǎng)膜、硬腦膜、大腦鐮與小腦幕），包圍腦組織的腦脊液以及顱骨、皮膚與頸部。該頭部模型共有303 588 個(gè)單元與530 954 個(gè)節(jié)點(diǎn)，單元尺寸為2～3 mm，充分考慮了人體頭部的生理結(jié)構(gòu)與相關(guān)細(xì)節(jié)，如圖2 所示。

1.2 本構(gòu)模型與材料參數(shù)

對(duì)于大腦、小腦與腦干等腦組織，采用黏超彈性本構(gòu)模型進(jìn)行表征[21]。利用Mooney-Rivlin 模型表征其超彈性行為，對(duì)應(yīng)的應(yīng)變能函數(shù)定義為：

圖 2 三維頭部有限元模型Fig. 2 3D finite element model of human head

表 1 黏超彈性材料參數(shù)Table 1 Material properties for hyper-viscoelastic model

腦脊液與水的性質(zhì)較為接近，在蛛網(wǎng)膜下腔內(nèi)流動(dòng)；該腔是由蛛網(wǎng)膜小梁支撐的腔體，能夠傳遞剪切變形，如圖1 所示。為了便于建模，頭部模型中的腦脊液是物理腦脊液、蛛網(wǎng)膜小梁與較大血管的綜合體。對(duì)于腦脊液CSF 的本構(gòu)模型，文獻(xiàn)[22-23]中直接采用剪切模量相對(duì)很小的線彈性模型來(lái)進(jìn)行表征。為了真實(shí)反映腦脊液與蛛網(wǎng)膜小梁的力學(xué)特性，模型中采用Mie-Grüneisen 狀態(tài)方程表征腦脊液的不可壓縮特性，同時(shí)采用較小的剪切模量（22 kPa）表征蛛網(wǎng)膜小梁的剪切傳遞特性[24-26]。Mie-Grüneisen 狀態(tài)方程表示為：

表 2 線彈性材料參數(shù)Table 2 Material properties for liner elastic model

1.3 接觸關(guān)系與邊界條件

基于顱腦的生理結(jié)構(gòu)并參考相關(guān)文獻(xiàn)[23-28]，三維頭部數(shù)值模型中的各部分組織之間的接觸關(guān)系定義為綁定連結(jié)。同時(shí)，采用狀態(tài)方程的腦脊液實(shí)現(xiàn)了腦組織與顱骨之間的有限滑動(dòng)，符合顱腦的生理結(jié)構(gòu)。而對(duì)于頭部模型中頸部的邊界條件，至今沒(méi)有明確統(tǒng)一的結(jié)論。Ganpule 等[23]直接采用固定邊界條件，對(duì)頭部碰撞實(shí)驗(yàn)[29]進(jìn)行數(shù)值模擬，兩者的結(jié)果具有相對(duì)較好的吻合性。Chen等[27]對(duì)比分析了自由邊界條件與固定邊界條件兩種工況下的模擬結(jié)果，雖然均與實(shí)驗(yàn)結(jié)果存在一定的差異，但是自由邊界對(duì)應(yīng)的模擬結(jié)果更能真實(shí)反映頭部的動(dòng)態(tài)響應(yīng)過(guò)程。軀體對(duì)顱腦動(dòng)態(tài)響應(yīng)的影響是與頸部的肌肉狀態(tài)密切相關(guān)的。當(dāng)肌肉處于緊繃狀態(tài)，頸部較為“剛硬”，軀體對(duì)顱腦動(dòng)態(tài)響應(yīng)的影響較大；而當(dāng)肌肉處于放松狀態(tài)時(shí)，頸部則較為柔軟，這種影響作用相對(duì)較小。通常情況下，頸部肌肉是處于放松狀態(tài)的。本文中基于顯式有限元程序，對(duì)包含軀體的完整模型與自由邊界條件下的頭部模型在撞擊載荷作用下的動(dòng)態(tài)響應(yīng)過(guò)程進(jìn)行了數(shù)值模擬。其中，頸部采用較小的彈性模量（3.54 MPa），表征頸部肌肉的放松狀態(tài)。圖3 對(duì)比顯示了兩種數(shù)值模型對(duì)應(yīng)的前額處與枕部處腦組織的壓力曲線。由圖3 可知，采用自由邊界條件的頭部模型，其模擬結(jié)果能夠較好地吻合完整模型的結(jié)果。這表明：當(dāng)頸部處于放松狀態(tài)時(shí)，軀體對(duì)顱腦動(dòng)態(tài)響應(yīng)過(guò)程的影響在碰撞等沖擊載荷條件下（通常為幾個(gè)毫秒）基本上是可以忽略的，此時(shí)的顱腦動(dòng)態(tài)響應(yīng)過(guò)程呈現(xiàn)出明顯的波動(dòng)效應(yīng)。同時(shí)，兩種數(shù)值模型對(duì)應(yīng)的前額處腦組織壓力在t=7 ms 時(shí)開(kāi)始產(chǎn)生差異，并且呈現(xiàn)出不斷擴(kuò)大的趨勢(shì)。在這個(gè)過(guò)程中，顱腦的動(dòng)態(tài)響應(yīng)會(huì)由應(yīng)力波主導(dǎo)轉(zhuǎn)化為整體結(jié)構(gòu)響應(yīng)，受到軀體對(duì)頸部約束作用的影響。

圖 3 前額處與枕部處腦組織的壓力曲線Fig. 3 Pressure curues of brain tissue at the forehead and occiput

1.4 頭部模型驗(yàn)證與數(shù)值模擬

Nahum 等[29]進(jìn)行的頭部碰撞實(shí)驗(yàn)，已成為頭部數(shù)值模型有效性驗(yàn)證的經(jīng)典實(shí)驗(yàn)。該實(shí)驗(yàn)利用運(yùn)動(dòng)塊體撞擊前額部位，利用傳感器記錄特征部位（前額、顱頂、枕部與顱后窩處）的顱內(nèi)壓力變化過(guò)程。由于該頭部碰撞實(shí)驗(yàn)沒(méi)有提供與頭部發(fā)生碰撞塊體的材料力學(xué)參數(shù)，參考Ganpule 等[23]進(jìn)行頭部數(shù)值模型驗(yàn)證時(shí)所使用的方法：將實(shí)驗(yàn)中測(cè)量所得到的碰撞接觸力轉(zhuǎn)化為平均壓力作用于前額部位，作用面積為1 470 mm2。對(duì)應(yīng)的壓力曲線如圖4 所示?；谏鲜鲱^部碰撞模型，可獲取顱內(nèi)壓力的動(dòng)態(tài)響應(yīng)過(guò)程，與頭部碰撞實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比分析。為減小采樣位置對(duì)模擬結(jié)果與實(shí)驗(yàn)數(shù)據(jù)對(duì)比分析的影響，在前額、顱頂、枕部與顱后窩處各取3 個(gè)特征點(diǎn)（P1、P2 與P3），特征點(diǎn)位置如圖5 所示。

圖 4 碰撞壓力曲線Fig. 4 Impact pressure curve

本文中提取特征點(diǎn)處顱內(nèi)壓力的模擬結(jié)果，與頭部碰撞實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比，兩者具有較好的吻合性，如圖6 所示。由碰撞實(shí)驗(yàn)結(jié)果可知，前額與顱頂處的顱內(nèi)壓力在整個(gè)碰撞過(guò)程中基本上均為正壓，表明這兩處腦組織處于受壓縮狀態(tài)。同時(shí)，前額處的顱內(nèi)壓力要明顯高于顱頂處的相應(yīng)值，這與頭部簡(jiǎn)化模型的理論解吻合[30]。相比而言，枕部與顱后窩處的顱內(nèi)壓力在碰撞前期階段為負(fù)壓，表明腦組織處于拉伸狀態(tài)，相應(yīng)的壓力幅值較為接近。但是，這兩處的顱內(nèi)壓力在后期階段會(huì)發(fā)生“反轉(zhuǎn)”，由負(fù)壓轉(zhuǎn)變?yōu)檎龎?，正壓幅值相比于?fù)壓幅值較小?；谌S頭部數(shù)值模型所得到的模擬結(jié)果，能夠較好地預(yù)測(cè)顱內(nèi)壓力在碰撞條件下的演化過(guò)程，與實(shí)驗(yàn)數(shù)據(jù)吻合較好，從而驗(yàn)證了該頭部數(shù)值模型的有效性。

為了對(duì)比分析頭部動(dòng)態(tài)響應(yīng)在碰撞載荷與爆炸沖擊波載荷作用下的差異，提取了頭頸部在特征時(shí)刻的壓力云圖，如圖7 所示。由于顱骨與腦組織的慣性力差異，在碰撞載荷作用下，顱骨加速度明顯大于腦組織加速度。因此，前額處（撞擊處）顱骨與腦組織之間產(chǎn)生相對(duì)徑向壓縮位移，使得腦組織處于壓縮狀態(tài)，對(duì)應(yīng)的壓力為正；而枕部處（對(duì)撞處）兩者之間產(chǎn)生相對(duì)徑向拉伸位移，腦組織處于拉伸狀態(tài)，對(duì)應(yīng)的壓力為負(fù)。隨著撞擊載荷的持續(xù)作用，顱內(nèi)正、負(fù)壓的作用范圍持續(xù)擴(kuò)大，發(fā)展成為靠近前額的前半部分腦組織受壓而靠近枕部的后半部分腦組織受拉的受力狀態(tài)。隨著外部載荷的減弱和頭頸部之間相對(duì)變形的增大，頭部運(yùn)動(dòng)會(huì)在某一特定時(shí)刻發(fā)生轉(zhuǎn)向，顱骨與腦組織的慣性差異使得兩者之間的相對(duì)位移發(fā)生轉(zhuǎn)變，致使腦組織壓力的狀態(tài)發(fā)生“反轉(zhuǎn)”。但是，“反轉(zhuǎn)”后的腦組織壓力，其幅值則相對(duì)較小。因此，顱骨與腦組織間的相對(duì)徑向位移決定腦組織壓力，撞擊處前半部分腦組織主要處于壓縮狀態(tài)，而對(duì)撞處后半部分腦組織則主要處于拉伸狀態(tài)。

圖 6 顱內(nèi)壓力模擬結(jié)果與實(shí)驗(yàn)結(jié)果的對(duì)比Fig. 6 Comparison of intracranial pressure between experiments and calculations

圖 7 前額撞擊時(shí)腦組織壓力云圖Fig. 7 Nephogram of brain pressure for the forehead collision

2 爆炸沖擊波-頭部流固耦合模型

為了獲得頭部在爆炸沖擊波下的動(dòng)態(tài)響應(yīng)規(guī)律，需要建立爆炸沖擊波-頭部流固耦合模型。頭部在常規(guī)爆炸沖擊波作用下的變形相對(duì)較小，可采用拉格朗日單元進(jìn)行網(wǎng)格劃分。而空氣需要模擬爆炸沖擊波的產(chǎn)生、傳遞以及與頭部相互作用的過(guò)程，對(duì)應(yīng)區(qū)域應(yīng)當(dāng)采用歐拉網(wǎng)格，避免單元過(guò)度變形。爆炸沖擊波與頭部在接觸邊界上的相互作用，需要通過(guò)求解流固耦合方程進(jìn)行確定。通用有限元軟件ABAQUS 通過(guò)非線性瞬態(tài)程序與歐拉-拉格朗日耦合法（Euler-Lagrangian coupling method，CEL），對(duì)整體系統(tǒng)的三大守恒偏微分方程（動(dòng)量、質(zhì)量與能量）同時(shí)進(jìn)行求解。歐拉區(qū)域能夠?qū)?dǎo)致嚴(yán)重網(wǎng)格扭曲的高度動(dòng)態(tài)事件進(jìn)行模擬（如爆炸沖擊波），并為拉格朗日區(qū)域提供了壓力邊界條件；拉格朗日區(qū)域可完全或部分位于歐拉區(qū)域內(nèi)，拉格朗日表面為歐拉區(qū)域提供邊界條件，該邊界條件不允許在其表面的法線方向上有流動(dòng)產(chǎn)生。同時(shí)，增強(qiáng)浸入邊界方法提供歐拉區(qū)域和拉格朗日區(qū)域之間的耦合作用。因此，固定歐拉網(wǎng)格與基于增強(qiáng)浸入邊界法的流固界面模型，能夠模擬爆炸沖擊波與頭部之間耦合作用，獲得頭部的動(dòng)態(tài)響應(yīng)規(guī)律。

2.1 有限元離散

模擬爆炸沖擊波的歐拉區(qū)域采用正方體形狀，邊長(zhǎng)為600 mm，這個(gè)尺寸的選擇充分考慮了沖擊波與頭部的作用時(shí)間，防止沖擊波邊界反射對(duì)頭部產(chǎn)生非真實(shí)的二次沖擊作用。經(jīng)過(guò)驗(yàn)證的三維頭部數(shù)值模型，在爆炸沖擊波-頭部流固耦合模型中可以直接使用，并放置在歐拉區(qū)域的中心位置。歐拉區(qū)域使用縮減積分六面體歐拉單元（EC3D8R）進(jìn)行網(wǎng)格劃分，并采用分塊劃分單元網(wǎng)格的方法：與頭部發(fā)生接觸的區(qū)域采用較小的單元尺寸，利于模擬爆炸沖擊波與頭部之間的流固耦合作用；沒(méi)有與頭部發(fā)生作用的區(qū)域則采用相對(duì)較大的單元，利于節(jié)約計(jì)算成本。同時(shí)，針對(duì)單元尺寸對(duì)模擬結(jié)果的影響進(jìn)行參數(shù)化研究。研究表明：當(dāng)中心區(qū)域的單元尺寸為3 mm 時(shí)，相應(yīng)的模擬結(jié)果能夠收斂，與文獻(xiàn)[23]結(jié)論基本一致。因此，將單元尺寸定義為3 mm，歐拉區(qū)域共有8 524 800 個(gè)單元。歐拉區(qū)域采用理想氣體狀態(tài)方程來(lái)表征空氣，空氣密度為1.18 kg/m3，氣體常數(shù)為287.04 J/(kg·K)，比熱容為1 000 J/(kg·K)。最終，建立的爆炸沖擊波-頭部流固耦合模型如圖8 所示。

2.2 荷載與邊界條件

根據(jù)我國(guó)工程設(shè)計(jì)規(guī)范中的設(shè)計(jì)準(zhǔn)則，爆炸沖擊波的超壓值（以大氣壓力為參考零值）會(huì)隨時(shí)間呈現(xiàn)指數(shù)型衰減。其理論表達(dá)式為：

圖 8 爆炸沖擊波-頭部的流固耦合模型Fig. 8 Fluid-solid coupling model of explosive blast wave-head

基于上述理論公式，可獲得炸藥當(dāng)量為7 kg且距離爆心3.8 m 處的壓力曲線。將理論結(jié)果與實(shí)爆實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比，兩者具有較好的吻合性，如圖9 所示。

本文中，采用這個(gè)典型爆炸沖擊波壓力曲線作為歐拉區(qū)域的輸入載荷，直接作用在入口處。歐拉場(chǎng)的出口處采用無(wú)反射邊界條件，沖擊波可以直接從這里逸出，不再對(duì)歐拉場(chǎng)產(chǎn)生影響。而歐拉場(chǎng)的環(huán)向4 個(gè)邊界則采用法向速度約束邊界條件，防止沖擊波從側(cè)面逸出，模擬平面型爆炸沖擊波與頭部的相互作用，可真實(shí)反映現(xiàn)場(chǎng)環(huán)境。

圖 9 爆炸沖擊波壓力曲線的理論結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比Fig. 9 Comparison of pressure curves of explosion blast between theoretical and experimental results

2.3 接觸與初始條件

由爆炸沖擊波-頭部流固耦合模型的理論基礎(chǔ)可知，兩者之間的相互作用需要借助基于增強(qiáng)浸入邊界法的流固界面模型來(lái)實(shí)現(xiàn)。在程序中，通用接觸（general contact）采用增強(qiáng)浸入邊界法能夠較好地處理拉格朗日單元與歐拉單元之間的接觸問(wèn)題。因此，通用接觸是進(jìn)行流固耦合分析的必備條件。

對(duì)于初始條件，系統(tǒng)中有兩個(gè)方面需要進(jìn)行合理設(shè)置。（1）歐拉區(qū)域中空氣的初始溫度與參考?jí)毫?。這里設(shè)定初始時(shí)刻的空氣溫度為303 K，以標(biāo)準(zhǔn)大氣壓力作為參考值，即歐拉區(qū)域的壓力初始值為零。（2）空氣材料的初始分布。這里可以使用體積分?jǐn)?shù)工具（volume fraction tool）來(lái)設(shè)定空氣在離散歐拉場(chǎng)中的分布，保證每個(gè)位置點(diǎn)處只有一種材料。

3 結(jié)果與分析

基于上述爆炸沖擊波-頭部流固耦合模型，可以獲得頭部在爆炸沖擊波作用下的動(dòng)態(tài)響應(yīng)規(guī)律，進(jìn)而揭示顱腦的損傷機(jī)理。

3.1 流場(chǎng)分析

爆炸沖擊波從正面作用于頭部時(shí)，流場(chǎng)壓力的整個(gè)演化過(guò)程如圖10 所示。在t=0.375 ms 時(shí)，平面沖擊波與頭部即將發(fā)生接觸，此時(shí)的沖擊波壓力峰值為170 kPa。隨后，沖擊波與頭部發(fā)生相互作用，反射疊加效應(yīng)使其壓力峰值急劇上升；當(dāng)t=0.450 ms 時(shí)，前額處的壓力峰值達(dá)到600 kPa，是入射沖擊波峰值的3.5 倍。在t=0.720 ms 時(shí)，沖擊波經(jīng)過(guò)前額進(jìn)入顱頂，曲面型的顱頂結(jié)構(gòu)使其在該處發(fā)生邊界層分離現(xiàn)象。因此，沖擊波對(duì)顱頂處的壓力峰值相對(duì)較小。當(dāng)t=0.960 ms 時(shí)，由于頭部阻擋作用而分散的沖擊波會(huì)在后腦部位重新匯集，該處的壓力值迅速上升，從而對(duì)頭部造成二次沖擊。相比于面部首次沖擊的壓力峰值，二次沖擊的對(duì)應(yīng)峰值出現(xiàn)了顯著下降。

圖 10 正面沖擊時(shí)流場(chǎng)壓力分布Fig. 10 Pressure distribution of flow field in blast frontal impact

圖 11 正面沖擊時(shí)眼部的流場(chǎng)壓力分布Fig. 11 Pressure distribution of flow field around the eyes in blast frontal impact

3.2 腦組織壓力

對(duì)于正面沖擊，先提取腦組織特征點(diǎn)處的壓力時(shí)程曲線，如圖12 所示。由圖12 可知，各特征點(diǎn)處的腦組織壓力均呈現(xiàn)出明顯的波動(dòng)特性，并且隨著時(shí)間的推移，這種波動(dòng)特性會(huì)逐步弱化。同時(shí)，從前額至顱頂處，特征點(diǎn)1～4 的腦組織壓力在前期階段呈現(xiàn)出明顯的周期特性，其中前額處的頻率高達(dá)8 kHz；而特征點(diǎn)5～8 則沒(méi)有明顯的周期特性。在爆炸沖擊波下，腦組織壓力的高頻波動(dòng)特性與碰撞載荷下的腦組織響應(yīng)是完全不同的，這就決定了針對(duì)兩種不同工況需要采用不同的顱腦防護(hù)策略。結(jié)合腦組織壓力的動(dòng)態(tài)演化過(guò)程與損傷判斷準(zhǔn)則可知，多處腦組織的壓力將超過(guò)重度損傷閾值（235 kPa）[31]，造成腦組織出現(xiàn)多處重度損傷。

圖 12 正面沖擊時(shí)特征點(diǎn)處腦組織的壓力曲線Fig. 12 Pressure curves of brain tissue at feature points in blast frontal impact

基于模擬結(jié)果，可提取顱腦在特征時(shí)間點(diǎn)處的壓力分布云圖和顱骨變形云圖，如圖13 所示。在顱骨變形圖中，帶有網(wǎng)格的顱骨為初始形態(tài)，沒(méi)有網(wǎng)格的顱骨為當(dāng)前形態(tài)，采用的放大系數(shù)為500。由壓力云圖可以看出，前額處腦組織在沖擊波作用下起初處于受壓狀態(tài)（t=0.465 ms），而后轉(zhuǎn)化為受拉狀態(tài)（t=0.520 ms）；在后續(xù)過(guò)程中，腦組織壓力會(huì)在壓、拉兩種狀態(tài)之間進(jìn)行轉(zhuǎn)化。緊鄰顱骨內(nèi)表面處的腦組織，其壓力沿著顱骨內(nèi)表面出現(xiàn)正、負(fù)壓交替現(xiàn)象，形成較大的壓力梯度；而遠(yuǎn)離顱骨內(nèi)表面的腦組織，其壓力梯度則相對(duì)較小。同時(shí)，具有典型三明治結(jié)構(gòu)的頭部顱骨，利用材料阻抗不匹配特性截?cái)嗖⒎稚_擊處的能量，從而為腦組織提供防護(hù)。

從流場(chǎng)壓力的演化過(guò)程可知，爆炸沖擊波主要作用于前額等面部結(jié)構(gòu)，而顱頂?shù)炔课挥捎诖嬖谶吔鐚用撾x效應(yīng)，受到?jīng)_擊波的作用相對(duì)較小，在進(jìn)行顱骨變形分析時(shí)可以忽略這個(gè)影響。當(dāng)t=0.465 ms 時(shí)，流固耦合作用使得爆炸沖擊波達(dá)到最大值，沖擊波直接擠壓皮膚與顱骨，前額處顱骨的“平面”變形造成該處腦組織處于擠壓狀態(tài)，體積壓力為正壓。隨著前額顱骨變形的進(jìn)一步增大，當(dāng)t=0.520 ms 時(shí)，該處顱骨的“平面”變形轉(zhuǎn)變?yōu)椤皟?nèi)凹”彎曲，致使該處與顱骨內(nèi)表面緊密連結(jié)的腦脊液與腦組織處于受拉狀態(tài)，體積壓力為負(fù)壓。而緊鄰前額處的顱骨則相當(dāng)于發(fā)生了“外凸”彎曲，對(duì)該處的腦脊液與腦組織進(jìn)行擠壓，體積壓力為正壓；變形的連續(xù)性使得顱頂處顱骨產(chǎn)生較小的“內(nèi)凹”彎曲，體積壓力為負(fù)壓。當(dāng)t=0.585 ms 時(shí)，沖擊波作用的減弱與顱骨的彈性邊界使得前額處顱骨的“內(nèi)凹”彎曲在達(dá)到最大值后發(fā)生回彈，由“內(nèi)凹”彎曲轉(zhuǎn)變?yōu)椤巴馔埂睆澢⒅鲗?dǎo)著顱骨其他部位彎曲類型的改變。對(duì)比t=0.520 ms 與t=0.585 ms 時(shí)的云圖可知，各位置處的腦組織的拉壓狀態(tài)發(fā)生轉(zhuǎn)變。對(duì)比分析后續(xù)特征時(shí)刻的腦組織壓力與顱骨彎曲變形，可以得到相同的結(jié)論：理想狀態(tài)下，顱骨發(fā)生“內(nèi)凹”彎曲，該處腦組織處于受拉狀態(tài)；顱骨發(fā)生“外凸”彎曲，腦組織則為受壓狀態(tài)。而真實(shí)條件下，顱骨內(nèi)表面處腦組織在特定位置處的受力狀態(tài)同時(shí)會(huì)受到周邊顱骨變形的影響。

圖 13 正面沖擊時(shí)顱腦壓力云圖與顱骨變形云圖Fig. 13 Nephogram of brain pressure and skull displacement in blast frontal impact

本文中進(jìn)一步提取了前額處顱骨的加速度時(shí)程曲線，并將該處特征時(shí)刻的加速度標(biāo)注在該曲線上，如圖14 所示。其中實(shí)心圓點(diǎn)表示該處腦組織壓力為正時(shí)的加速度，空心圓點(diǎn)則表示壓力為負(fù)時(shí)的加速度。由圖14 可見(jiàn)，爆炸沖擊波致使前額顱骨發(fā)生局部高頻振動(dòng)，其頻率與該處腦組織壓力的波動(dòng)頻率相一致。前額顱骨最值加速度與顱骨彎曲變形模式、腦組織狀態(tài)均具有較好的相關(guān)性：顱骨受沖擊方向加速度最大時(shí)，顱骨為“外凸”彎曲模式，腦組織處于受壓狀態(tài)；而顱骨受沖擊方向加速度最?。_擊反方向最大）時(shí)，顱骨為“內(nèi)凹”彎曲模式，腦組織處于受拉狀態(tài)。

圖 14 前額處顱骨加速度曲線Fig. 14 Acceleration curves of forehead skull

對(duì)比分析腦組織壓力、顱骨局部彎曲變形與加速度的整個(gè)演化過(guò)程可知，前額處顱骨在爆炸沖擊波作用下發(fā)生彎曲變形并引起局部振動(dòng)，該處的局部變形會(huì)沿著顱骨進(jìn)行“傳播”，影響著顱骨的變化構(gòu)型，從而使腦組織在外載作用下由前額處的“單一波源”逐步發(fā)展到整個(gè)顱骨的“多處波源”，這些波源共同決定腦組織壓力的演化過(guò)程。顱骨局部彎曲致使腦組織沿著顱骨內(nèi)表面出現(xiàn)拉、壓交替狀態(tài)，形成較大的壓力梯度。同時(shí)，顱骨局部彎曲的振動(dòng)頻率較高，使得腦組織壓力具有高頻波動(dòng)特性。

4 結(jié) 論

建立了具有詳細(xì)解剖學(xué)結(jié)構(gòu)和較高生物仿真度的三維頭部數(shù)值模型，并基于Nahum 等[29]的碰撞實(shí)驗(yàn)驗(yàn)證了該頭部模型的有效性。同時(shí)，基于CEL 方法建立了爆炸沖擊波-頭部流固耦合模型，對(duì)頭部受到爆炸沖擊波正面沖擊工況進(jìn)行了數(shù)值模擬，分析了腦組織的動(dòng)態(tài)響應(yīng)規(guī)律，揭示了相應(yīng)的損傷機(jī)理。主要結(jié)論如下。

（1）在碰撞載荷下，顱骨與腦組織間的相對(duì)徑向位移決定腦組織壓力，撞擊處前半部分腦組織主要處于壓縮狀態(tài)，而對(duì)撞處后半部分腦組織則主要處于拉伸狀態(tài)。

（2）流固耦合作用使得爆炸沖擊波峰值壓力增大為入射波的3.5 倍，具有生理凹陷結(jié)構(gòu)的眼部則使該值增大到4.5 倍，是沖擊波的重要靶器官。

（3）爆炸沖擊波在顱頂處會(huì)發(fā)生邊界層分離效應(yīng)，該處的壓力峰值相對(duì)較??；沖擊波從頭部撞擊處周邊發(fā)生繞流，并在撞擊的對(duì)側(cè)處匯合，對(duì)顱腦造成相對(duì)較小的二次沖擊。

（4）正面作用時(shí)前額處顱骨的振動(dòng)頻率高達(dá)8 kHz，從而使腦組織的壓力具有高頻波動(dòng)特性。這與碰撞載荷下的腦組織動(dòng)態(tài)響應(yīng)完全不同。

（5）前額處顱骨（正面沖擊）的局部彎曲變形，會(huì)沿著顱骨進(jìn)行“傳播”，影響著顱骨的變化構(gòu)型，從而決定了腦組織壓力的演化過(guò)程。同時(shí)，顱骨最值加速度與顱骨彎曲變形模式、腦組織狀態(tài)均具有較好的相關(guān)性。