劉德平, 鄭 凱, 李冬梅
(1.鄭州大學(xué) 機械與動力工程學(xué)院,河南 鄭州 450001;2.廣東順德工業(yè)設(shè)計研究院,廣東 佛山 528300)
氧合器是ECMO(體外膜肺氧合)系統(tǒng)的重要組成部分,可在一定時間內(nèi)完全代替人體肺的功能。氧合器內(nèi)部血流的運動特性對其性能有著重要影響,但其結(jié)構(gòu)復(fù)雜,難以通過實驗直接觀察內(nèi)部血液流動特性[1]。計算流體力學(xué)(CFD)能夠分析各類流體機械內(nèi)部運動特性,利用CFD對氧合器進行分析,可以使其內(nèi)部流體運動狀態(tài)可視化,減少設(shè)計中的性能實驗,節(jié)省研發(fā)成本[2]。
目前,已有學(xué)者對氧合器開展了相關(guān)研究,如Consolo等[3]比較了5種不同換熱室壁面結(jié)構(gòu)對其傳熱性能的影響;Kirsch等[4]開發(fā)了一種氣血傳輸模型來預(yù)測氧合器內(nèi)的傳質(zhì)性能;Wickramasinghe等[5]通過實驗驗證了氧合器傳質(zhì)性能和摩擦因數(shù)間的相關(guān)性;Taskin等[6]建立了單纖維和多纖維模型研究其內(nèi)部氧分壓分布;Zhang等[7]研究了3種纖維膜排列方式下的流體切應(yīng)力分布情況。可見,學(xué)者們關(guān)注更多的是氧合器的傳熱、傳質(zhì)性能和血液在纖維膜內(nèi)的局部流動,而針對其整體血液流動特性的研究相對較少。
本研究通過實驗與數(shù)值模擬相結(jié)合的方式,對一款典型氧合器內(nèi)部整體血液流動特性進行分析,同時建立微尺度模型輔助驗證,獲得了氧合器流場內(nèi)速度、壓力和壁面剪切力等的分布規(guī)律;由溶血數(shù)值預(yù)估模型得到了其在工作流量范圍內(nèi)的標準溶血指數(shù)值,找到了血液在其內(nèi)部最易受到損傷的區(qū)域,為氧合器后續(xù)的設(shè)計和改進提供了參考。
圖1所示為氧合器的三維模型圖。其主要由內(nèi)腔、中腔、外腔以及各流體的進出口等組成。內(nèi)腔中空,起支撐作用;中腔和外腔被管狀纖維膜填充,分別形成變溫室和氧合室。氧合器在工作時,如圖1(b)所示,血液由入口管道流入,首先到達位于中腔的變溫室,隨后血液在纖維膜管外從上向下流動,熱交換水在纖維膜管內(nèi)從下向上流動,完成對血液的溫度調(diào)節(jié)。血液再由底部中腔與外腔的間隙流入外腔,隨后血液在纖維膜管外從下向上流動,氧氣從纖維膜管內(nèi)自上而下通過,血液和氧氣利用纖維膜管上的微孔通過彌散作用完成二氧化碳和氧氣的交換,最終血液由出口流出。
圖1 氧合器三維模型圖Figure 1 Three dimensional model figures of oxygenator
多孔介質(zhì)模型能夠?qū)ρ鹾掀鞯牧鲌鲂阅苓M行有效分析[8-9]。它將氧合器內(nèi)部的纖維膜假設(shè)為多孔介質(zhì),血液視為不可壓縮牛頓流體,設(shè)定該區(qū)域的孔隙度和流體的黏性阻力與慣性阻力系數(shù)進行計算。通過動量方程和連續(xù)性方程來描述笛卡爾坐標系中不可壓縮黏性流體的穩(wěn)態(tài)運動:
(1)
(2)
式中:ρ為密度,kg/m3;ui、uj為流體速度,m/s;p為壓力,Pa;τij為剪切應(yīng)力張量;gi為重力加速度,m/s2;Si為源項。
多孔介質(zhì)模型在動量方程的右側(cè)增加了一個Si源項,其由黏性損失和慣性損失組成,用來模擬計算域中的動量損失[10],即
(3)
在低流速下(Q<6.00 L/min),慣性損失可以忽略不計[11]。故多孔介質(zhì)的動量損失僅由達西定律計算獲得的黏性阻力決定,即
(4)
式中:Q為流量,L/min;A為流量通過纖維束的橫截面面積,m2;Δp為纖維束上的壓降,Pa;μ為流體黏度,Pa·s;L為中空纖維束長度,m。
所用血液損傷預(yù)估模型考慮了流場整體平均效應(yīng),通過對整個穩(wěn)定流場的速度域和剪切應(yīng)力域進行線性場積分計算[12]。血液中流體剪切應(yīng)力的張量表達式為
(5)
(6)
式中:μ為湍流黏度;k為湍流動能;δij為Kronecker數(shù)。通過米勒斯屈服準則,剪切應(yīng)力的張量形式可簡化為等效標量形式:
(7)
Garon等[13]以雙曲型輸運方程為基礎(chǔ)研究得到一種快速溶血預(yù)估模型,即雙曲型輸運方程式:
(8)
Dt=D1/0.785。
(9)
式中:v為速度,m/s;Dt為線性溶血指數(shù),%;D為溶血值,%。
Giersiepen等[14]的研究發(fā)現(xiàn),溶血百分數(shù)、剪切應(yīng)力τ以及暴露時間texp之間有如下關(guān)系:
D=ΔHb/Hb×100%=
(10)
式中:ΔHb為受損傷血紅蛋白濃度,g/L;Hb為總的血紅蛋白濃度,g/L。
從而可得到σ的表達式:
σ=(3.62×10-7)1/0.785τ2.416/0.785。
(11)
在計算獲得一個穩(wěn)定流場后,考慮流場整體平均溶血特性,可計算得
(12)
由式(9)指數(shù)換算可得溶血值D:
(13)
進而可以求得標準溶血指數(shù)值NIH:
NIH=Hb·D×100。
(14)
壓降實驗裝置如圖2所示。采用某型號膜式氧合器進行實驗,質(zhì)量分數(shù)為33%的甘油水溶液在室溫25 ℃下與37 ℃的人體血液中具有相同黏度,可以代替血液進行壓降實驗。通過離心泵將甘油水溶液抽入管道,導(dǎo)入氧合器,利用流量傳感器和壓力傳感器監(jiān)控氧合器的進、出口流量和壓強的變化。
圖2 壓降實驗裝置圖Figure 2 Diagram of pressure drop experimental device
應(yīng)用Pelosi等[11]的研究方法,在氧合室和變溫室之間的界面創(chuàng)建一個檢修孔,分別在流體進、出口和檢修孔處連接壓力傳感器,測量流體在3個位置的壓力值。經(jīng)過簡單的計算即可得到所需壓降,依據(jù)壓降-流量擬合計算纖維束的黏性阻力。值得注意的是,該方法測量出的氧合室和變溫室壓降值均包含了進、出口部分,盡管進出口部分的壓降占比很小,但測量結(jié)果與實際結(jié)果有一定偏差。
如圖3所示,通過線性插值法擬合得到變溫室和氧合室的Δp/Q值。將其分別代入式(4)求得兩室的黏性阻力:變溫室黏性阻力1/α=6.781×107m-2;氧合室黏性阻力1/α=7.181×108m-2。多孔介質(zhì)被認為是各向同性的,因此每個方向的黏性阻力相同。
圖3 變溫室和氧合室壓降-流量擬合曲線Figure 3 Pressure drop-flow fitting curve of the variable greenhouse and oxygenation chamber
如圖4所示,建立氧合器血液流道模型進行仿真計算。該模型包括多孔介質(zhì)和外部流體域。將由纖維束填充的變溫室和氧合室視為多孔介質(zhì),剩余區(qū)域(血液出口管道及其他間隙)為外部流體域。使用速度進口邊界,由入口流量和截面積可得進口速度為0.147~0.737 m/s。使用壓力出口邊界,設(shè)出口壓力為0 Pa,環(huán)境壓力默認為標準大氣壓。設(shè)置血液密度為1 055 kg/m3,黏度為2.36 Pa·s。所有壁面設(shè)為無滑移固壁界面,采用標準k-ε湍流模型,設(shè)置為SIMPLE算法、二階迎風離散格式。
圖4 血液流道簡化模型與徑向中間截面剖面圖Figure 4 Simplified model of blood flow channel and sectional view of middle plane
利用ICEM劃分流場網(wǎng)格,在入口流量為1.00~5.00 L/min時,計算流體介質(zhì)進出口的壓降,分別以3種不同的網(wǎng)格數(shù)(1.60×106、2.70×106、4.49×106)進行計算以驗證網(wǎng)格無關(guān)性與模擬計算準確性。將模擬結(jié)果與實驗結(jié)果進行比較,如圖5所示,模擬值與實驗值相近,證明了氧合器血液流道仿真模型的有效性,并且計算結(jié)果受網(wǎng)格數(shù)變化影響不大。
圖5 網(wǎng)格無關(guān)性驗證Figure 5 Grid independence verification
多孔介質(zhì)模型能夠從宏觀角度很好地模擬氧合器內(nèi)部流場,但卻不能觀察纖維膜上的實際微觀流動現(xiàn)象。因此,本文在以上工作基礎(chǔ)上,根據(jù)纖維膜的實際排列情況,建立了一個微尺度模型,來觀察流體在纖維束中的運動狀態(tài)。微觀三維模型如圖6所示,為六面體形狀,流體從六面體上方流入,前后兩排纖維束分別與六面體上下兩底面呈21°和159°交叉排列,具體尺寸如表1所示。利用ICEM構(gòu)建流場網(wǎng)格劃分,考慮到湍流的影響,模型采用RNGk-ε模型,所有外壁使用周期性邊界條件,采用SIMPLE算法和二階迎風格式,進行壓力-速度耦合求解。
圖6 微尺度三維模型Figure 6 Microscale 3D model
表1 微觀模型結(jié)構(gòu)參數(shù)Table 1 Structural parameters of microscopic model
通過仿真計算得到了氧合器內(nèi)部的速度、壓力及壁面剪切應(yīng)力等云圖。圖7所示為氧合器進口流量在3.0 L/min時的橫向截面速度云圖。其中,截面A為流體進口的橫向截面;截面B為進口流體域和變溫室相接處的橫向截面;截面C為流體出口的橫向截面;截面D為變溫室和氧合室相接處的橫向截面。各截面上流體域的整體速度都比較大,多孔介質(zhì)域內(nèi)速度則明顯偏小。
圖7 橫向截面速度云圖Figure 7 Velocity cloud diagram of cross-section
如圖8所示為徑向中間截面速度云圖,速度分布情況與橫向截面速度云圖相對應(yīng),流體速度在進出口時較高,在多孔介質(zhì)區(qū)域時有明顯下降。較大的流速變化往往可能產(chǎn)生漩渦,會對血液造成較大的損傷。如圖9所示給出了圖8中的流體域A、B、C處的局部放大速度矢量圖。
圖8 徑向中間截面速度云圖Figure 8 Velocity cloud diagram of the radial middle section
圖9 局部放大速度矢量圖Figure 9 Partial magnified velocity vector diagram
圖10為流量在3.0 L/min時的徑向中間截面壓力云圖。在出口流體域轉(zhuǎn)角處出現(xiàn)了一處負壓區(qū),負壓區(qū)可能會導(dǎo)致該處出現(xiàn)漩渦或回流等現(xiàn)象,造成嚴重的血液損傷,對膜式氧合器的性能產(chǎn)生不利影響,這與速度矢量圖9中的回流區(qū)相對應(yīng)。表2為不同流量下進、出口壓降模擬均值與實驗均值,壓降的模擬均值與實驗均值相對誤差起始較小。隨著流量的增加,兩者間的誤差也越大,這是因為實驗測量值包含了進、出口部分的壓力,表明進、出口部分的壓力受流量變化的影響較大。
圖10 氧合器壓力云圖Figure 10 Pressure cloud diagram of oxygenator
表2 模擬值與實驗值壓降比較Table 2 Comparison of pressure drops between simulated and experimental values
從微尺度數(shù)值模擬中獲得的速度和標量應(yīng)力分布云圖如圖11、12所示,最小速度位于纖維壁面處,最大速度位于兩纖維壁面的中間位置。最大應(yīng)力位于纖維壁面,最小應(yīng)力在遠離壁面的位置。在微尺度模型內(nèi),得到單位長度的平均壓降為45 000 Pa/m,根據(jù)纖維束的總長度計算總壓降,結(jié)果與宏觀分析所得壓降相近,誤差在15%以內(nèi)。誤差的產(chǎn)生可能是由于實驗得到的黏性阻力與實際值有偏差和使用了各向同性多孔介質(zhì)模型。
圖11 微尺度橫向和縱向截面速度分布云圖Figure 11 Microscale velocity distribution cloud diagram in transverse and longitudinal sections
圖12 標量應(yīng)力分布云圖Figure 12 Scalar stress distribution cloud diagram
通過用戶自定義函數(shù)將血液損傷預(yù)估模型編譯到FLUENT中,在獲得穩(wěn)定的流場后計算氧合器的標準溶血指數(shù)NIH,同時將壁面剪切應(yīng)力作為輔助參數(shù)分析血液損傷程度[15]。壁面剪切應(yīng)力分布云圖如圖13所示,出口管道與氧合室相連接處的壁面剪切應(yīng)力有最大值,為30.9 Pa。氧合器內(nèi)部流場的等效標量剪切應(yīng)力極值和平均值如表3所示。由表3可知,大剪切應(yīng)力分布于血液的進、出口管道區(qū)域,對應(yīng)了圖7、8所示的高流速區(qū)。在氧合器纖維束區(qū)域,血液流速和壁面剪切應(yīng)力都較小,對血細胞的破壞能力有限。血液進、出口管道由于有高流速和大壁面剪切應(yīng)力,最容易引起血液損傷。
圖13 壁面剪切應(yīng)力分布云圖Figure 13 Contour of wall shear stress distribution cloud diagram
表3 等效標量剪切應(yīng)力分布Table 3 Equivalent scalar shear stress distribution
在進口流量為1.0~5.0 L/min時,利用血液損傷預(yù)估模型對氧合器內(nèi)部流場進行血液損傷估算,結(jié)果如圖14所示,標準溶血指數(shù)在0.004 8~0.049 2 g/100 L,與文獻[15]所求范圍相近,小于人體生理最大的許用值0.1 g/100 L[16],該氧合器的溶血指數(shù)滿足使用要求。
圖14 標準溶血指標隨流量的變化Figure 14 Change of standard hemolysis index with flow
(1)在低流量范圍內(nèi),各向同性多孔介質(zhì)模型較為準確地模擬了氧合器內(nèi)部流體流動,模擬值與實驗值偏差較小,但隨著流量增加,兩者的偏差增大,表明與多孔介質(zhì)域相比,進、出口流體域的壓力受流量變化的影響較大。
(2)氧合器內(nèi)部纖維束區(qū)域的壓力分布呈同心均勻梯度分布狀態(tài)并且壓力變化與流量大小呈正相關(guān),纖維束區(qū)域是壓力損失的主要位置。
(3)在血液的進、出口管道處出現(xiàn)了漩渦和回流,同時對應(yīng)著高壁面剪切應(yīng)力值,表明這些位置是造成血液損傷的主要區(qū)域,可以針對性地對其進行結(jié)構(gòu)優(yōu)化以降低血液損傷。