吳 偉,楊曉京,袁銳波,楊向紅,王 瑞,李宗睿

(1.昆明理工大學 機電工程學院,昆明 650504;2.昆明醫(yī)科大學 附屬延安醫(yī)院,昆明 650000)

0 引言

阻塞性睡眠呼吸暫停低通氣綜合征(obstructive sleep apnea hypopnea syndrome,OSAHS)[1]是一種常見的呼吸睡眠疾病,在夜間睡眠過程中,上氣道出現(xiàn)反復坍塌導致患者呼吸不通暢,患者常出現(xiàn)睡眠時鼾聲大、白天過度嗜睡等癥狀。上氣道局部坍塌是造成睡眠時呼吸暫停的主要因素,研究表明[2]上氣道是否坍塌取決于維持上氣道擴張的力是否大于呼吸氣體時所產(chǎn)生的壓力。目前,OSAHS的治療手段主要采用手術治療與非手術治療,手術治療是通過切除上氣道中阻塞的部位,如扁桃體切除手術、腭咽成型手術以及鼻中隔偏曲矯正手術等。口戴矯治器[3]已經(jīng)成為OSAHS非手術治療的有效手段之一,其中下頜前伸式矯治器在臨床治療中效果最佳。隨著醫(yī)工結(jié)合的快速發(fā)展,計算流體力學(CFD)[4]的方法被廣泛應用于口腔醫(yī)學領域中,數(shù)值模擬計算已經(jīng)成為研究人體上氣道流體力學的重要方法。陳柳潔等[5]在不同的條件下,通過流體力學的方法研究OSAHS患者正常睡眠時的呼吸特性,探索影響病人呼吸的因素,郭宇峰等[6]利用流體力學仿真對比OSAHS兒童與正常兒童上氣道流場的差異。本文利用計算流體動力學分析下頜前伸式止鼾器對OSAHS患者上氣道內(nèi)流體的影響,可用于評估下頜前伸式止鼾器的治療療效[7],為OSAHS疾病的研究和止鼾器的設計提供理論依據(jù)。

1 下頜前伸式止鼾器的結(jié)構(gòu)及原理

下頜前伸式止鼾器屬于口腔矯治器的一種,主要是由上、下牙以及中間的連接裝置組成,如圖1所示為延安醫(yī)院[8]制作的下頜前伸式止鼾器,上下牙套部分是由壓膜塑料制作,中間的不銹鋼金屬連桿起固定支撐與調(diào)節(jié)前伸距離的作用。首先通過口內(nèi)測量尺測量OSAHS患者上下頜之間的最大前伸距離,調(diào)節(jié)中間連桿使下頜前伸距離為測量距離的75%[9],同時,保持上下壓尖之間垂直開口高度為4～5 mm。下頜前伸不僅能帶動舌體向前移動,使舌咽部氣道暢通,還可以增大舌肌張力,以減小上氣道壁坍塌的風險。

圖1 下頜前伸式止鼾器實物圖

2 上氣道三維模型的建立

選取1例在昆明醫(yī)科大學附屬延安醫(yī)院口腔科就診的輕中度OASHS患者為志愿者(男,25歲),口腔檢查全口牙齒完好,無其他口腔疾病,無上氣道手術史且顳下頜關節(jié)正常。

采用卡瓦公司的I-CAT型號的CBCT對患者進行頭部掃描,掃描參數(shù):電壓:120 KV,電流:5 mA,掃描時間:25 s?；颊弑３终藙?雙目正視前方并保持眶耳平面與地面平行,分別對患者進行2次CT掃描,患者未佩戴止鼾器以及習慣佩戴止鼾器后進行掃描,單次掃描獲得CT影像數(shù)據(jù)528張,以標準DICOM格式導出并保存。

將上述CT數(shù)據(jù)導入Mimics 19.0軟件中,將蒙版閾值范圍[10]設置為-1 024～-434,以此提取上氣道三維模型,利用蒙版分割、區(qū)域增長等指令對上氣道蒙版[11]進行處理,以確保重建上氣道模型的準確性,圖2為重建后的上氣道蒙版。

圖2 上氣道蒙版

運用Calculate3D指令將上氣道蒙版計算生成上氣道三維模型并以STL文件格式導出。在Geomagic studio 2014逆向工程軟件中對上氣道三維模型進行表面優(yōu)化處理[12],去除模型中額竇、篩竇等多余部分。利用去除特征命令優(yōu)化模型表面的釘狀物,如圖3所示。使用精確曲面功能重構(gòu)出上氣道三維實體模型,最終以X-T格式導入UG11.0軟件中,圖4分別為佩戴止鼾器前后的上氣道三維模型,在UG中檢查上氣道三維模型表面是否存在曲面未擬合等缺陷。

圖3 優(yōu)化后的上氣道三維模型

圖4 UG中的上氣道三維模型

3 上氣道流體力學數(shù)值計算

3.1 邊界條件與控制方程的設置

人體的上氣道[13]主要由鼻、咽、喉3部分組成,鼻和喉位于上氣道結(jié)構(gòu)中的上下兩端,由軟骨等組織構(gòu)成,因而具有一定的張力。而中間的口咽部分是由肌肉構(gòu)成的軟性管道,由于其缺少骨性結(jié)構(gòu)的支撐,所以容易發(fā)生上氣道壁坍塌的現(xiàn)象。

相關研究[14]表明,人在呼吸過程中,呼氣階段的軟腭向后運動的趨勢大于吸氣階段,與健康人相比,OSAHS患者的上氣道更狹小,軟腭向后移動會增大腭咽部位發(fā)生氣道阻塞的風險。因此,采用正常呼氣方式來設置上氣道的邊界條件,如圖5所示,設置喉咽下表面A為速度入口inlet,2個前鼻孔表面B和C分別為壓力出口outlet1和outlet2。考慮到人體前鼻孔與正常大氣相通,所以出口壓力值設為0,進氣口的流速[15]計算如下所示:

圖5 上氣道模型邊界條件的設置示意圖

v=Q/A

(1)

式中:Q為入口處流體體積流量,L/s;A為入口處的面積,m2;單位體積流量[16]選取600 mL/s。

選擇上氣道中的流體介質(zhì)為空氣,空氣的質(zhì)量密度ρ=1.225 kg/m3,動力黏性系數(shù)μ=1.789 4×10-5kg/(m·s)-1。

M=v/c

(2)

式中:v為進氣口流速,m/s;c為音速,默認為340 m/s。

上氣道中流體的流況依據(jù)馬赫數(shù)來判斷,由于式(1)中計算出上氣道入口空氣流速相對音速較小,式(2)中馬赫數(shù)M的計算結(jié)果小于0.3[15],所以上氣道內(nèi)的流體視為不可壓縮流體。上氣道中流體的流動方式可依據(jù)雷諾數(shù)[17]的大小判斷,通常管道內(nèi)雷諾數(shù)小于2 100的流體為層流,雷諾數(shù)大于4 000為湍流,雷諾數(shù)計算如下所示:

(3)

式中:ρ為空氣密度,kg/m3;V為流體流速,m/s;μ為空氣動力黏性系數(shù),kg/(m·s)-1;D為水力直徑,m。其中D按如下公式計算:

(4)

式中:S為橫截面積,m2;L為截面圓周長度,m。

在Mimics軟件測得喉咽部下邊緣平面面積S為409.748 mm2,根據(jù)式(3)中雷諾數(shù)的計算結(jié)果大于4 000,選擇上氣道內(nèi)流體的流動形式為湍流。仿真采用標準k-ε湍流[18]計算模型,假設上氣道壁光滑無彈性,視整個上氣道壁面為剛性,使用SIMPLE算法[19]求解壓力與速度。當計算湍流流動時,除了考慮連續(xù)性方程和動量方程外,標準的k-ε模型需要對湍動能k及其耗散率ε方程進行計算。

湍流動能方程:

Gk+Gb-ρε-YM

(5)

湍流耗散方程:

(6)

其中湍流黏性系數(shù)μt定義如下:

(7)

式中:Gk為平均速度梯度引起的湍動能;Gb為浮力造成的湍動能;YM為可壓縮對總耗散率的影響;C1ε=1.44,C2ε=1.92,C3ε=0.09為經(jīng)驗常數(shù);湍動能k的普朗克常數(shù)為σk=1.0耗散率ε的湍流普朗克常數(shù)項為σε=1.3。

3.2 網(wǎng)格劃分與無關性檢驗

由于上氣道三維模型的結(jié)構(gòu)復雜且不規(guī)則,在Fluent軟件中通過Tetrahedrons方法[20]對三維模型進行網(wǎng)格的生成。為了確定最優(yōu)的網(wǎng)格數(shù)量,因此對網(wǎng)格的數(shù)量進行無關性校驗,其結(jié)果如表1所示。

表1 網(wǎng)格無關性檢驗結(jié)果

考慮壓力與速度之間的偏差,選擇110萬～170萬之間的網(wǎng)格數(shù)量。上氣道三維模型網(wǎng)格劃分如下,佩戴止鼾器前的模型網(wǎng)格數(shù)量為1 615 507,節(jié)點數(shù)量為37 166,佩戴止鼾器后的模型網(wǎng)格數(shù)量為166 527,節(jié)點數(shù)量為36 047,圖6為佩戴止鼾器前后模型的網(wǎng)格劃分示意圖。

圖6 上氣道網(wǎng)格劃分示意圖

3.3 上氣道三維模型截面的選取

運用CFD-post對模擬計算結(jié)果后處理,首先對上氣道三維模型進行截面選取。OSAHS患者與正常人相比,上氣道結(jié)構(gòu)中口咽部分較狹窄,其中腭咽[21]是口咽部分最容易發(fā)生阻塞的位置,所以除了鼻咽和鼻腔外,將整個上氣道劃分為11個截面,如圖7所示。其中截面a至截面g為腭咽部分,截面g至截面k為舌咽和喉咽部分。

圖7 上氣道中選取的截面示意圖

4 仿真結(jié)果分析

4.1 上氣道結(jié)構(gòu)形態(tài)改變

如圖4所示,在OSAHS患者在佩戴下頜前伸止鼾器后,整個上氣道體積較之前得到增加,在腭咽段體積增加較為明顯,這與趙明莉等[22]的研究一致,表明下頜前伸式止鼾器能夠擴大OSAHS患者的上氣道體積。

4.2 上氣道流速與壓力變化

在OSAHS患者呼吸的過程中,壓力的變化與上氣道坍塌的發(fā)生有著密切的關聯(lián),流速的變化則影響患者的呼吸。利用CFD-post模塊對上氣道中每個截面進行數(shù)據(jù)統(tǒng)計并截取上氣道正中位的速度與壓力云圖,在Origin軟件中對仿真所得的數(shù)據(jù)進行描繪處理,截面最大流速圖能反映上氣道中不同截面流速的變化,而截面的最小壓力值與上氣道壁坍塌形變有關,負壓[23]能使上氣道壁向內(nèi)形變。圖8和圖9為止鼾器治療前后上氣道中截面的最大流速和截面最小壓力曲線。

圖8 止鼾器治療前后截面大流速曲線

圖9 止鼾器治療前后截面最小壓力曲線

由圖8可知,OSAHS患者在佩戴止鼾器后,腭咽各平面的流速相對之前減小,喉咽平面的流速略有增加且變化很小。由圖9可知,OSAHS患者未佩戴止鼾器前,腭咽處的壓力為負壓,而喉咽處的壓力值為正壓,負壓值可以使上氣道壁向內(nèi)變形,能夠加大腭咽處發(fā)生氣道壁坍塌的風險,這與臨床中腭咽是上氣道最容易塌陷的結(jié)論相符。經(jīng)止鼾器的治療后,OSAHS患者腭咽處的壓力值開始由負壓變?yōu)檎龎?正壓可以使上氣道壁向外變形,能夠增加維持上氣道壁的張力,同時上氣道中流速的降低,兩者都有利于患者的呼吸暢通。

分析OSAHS患者上氣道的流速與壓力云圖,可以了解上氣道中流速與壓力的分布情況。圖10和圖11為止鼾器治療前后正中平面速度與壓力云圖。

圖10 止鼾器治療前后正中平面速度云圖

圖11 止鼾器治療前后正中平面壓力云圖

如圖10和圖11所示,上氣道中的最大流速區(qū)域集中于腭咽段,最大壓力區(qū)域位于喉咽下端。在OSAHS患者佩戴矯治器后,鼻腔中的流速與之前相比略有增大,腭咽上部處出現(xiàn)較小壓力值。這是因為腭咽上端出現(xiàn)渦流[24],渦流區(qū)會出現(xiàn)局部壓力最小值且渦流區(qū)還會使過流面積變小。結(jié)合圖8和圖9所示,上氣道中流速最大和壓力最小的位置為平面f,即腭咽最小平面。經(jīng)下頜前伸式止鼾器治療后,上氣道中流速最大值由 16.163 m/s降低到12.57 8m/s,最小負壓值由-86.21 Pa 變?yōu)?50.675 Pa,腭咽最小處壓力與流速的改善可以表明下頜前伸止鼾器對OSAHS患者治療的有效性。

5 結(jié)論

1) 通過Mimics軟件建立的醫(yī)學模型能夠準確反映人體生理結(jié)構(gòu),本文所建立的模型與上氣道的解剖形態(tài)一致,能夠幫助醫(yī)生更全面地了解OSAHS患者病情,做出更加合理的診斷,同時,該模型也可以轉(zhuǎn)化為通用三維模型格式,亦可作為有限元分析的導入模型。

2) OSAHS患者在佩戴下頜前伸式止鼾器后,上氣道內(nèi)空氣的流動速度下降,腭咽處壓力值也從負壓變?yōu)檎龎?。流速的降低可以減小患者在睡眠過程中空氣對其上氣道的沖擊,在一定程度上降低鼾聲。上氣道中壓力的變化可以降低上氣道壁發(fā)生坍塌的風險,有利于OSAHS患者病情改善。

3) 運用計算流體動力學方法評估止鼾器的療效,在口腔醫(yī)學中運用有限元分析理論將仿真結(jié)果與臨床治療相結(jié)合,可為OSAHS患者治療方案的正確制訂以及止鼾器的合理設計提供新思路。