人工心臟泵流場(chǎng)的數(shù)值分析

劉漢濤，馬理強(qiáng)，李海橋

（中北大學(xué) 機(jī)械與動(dòng)力工程學(xué)院，山西 太原030051）

0 引 言

人工心臟泵采用體外灌注的方式保證衰竭心臟的血液循環(huán)，從而代替心臟的功能［1-2］.第一個(gè)人工心臟的研制距今已半個(gè)多世紀(jì)［3］，國(guó)內(nèi)外的許多學(xué)者對(duì)心臟泵的形式也進(jìn)行了不斷改進(jìn)，目前流線型葉輪葉片的離心式心臟泵在臨床實(shí)驗(yàn)中由于具有更多的優(yōu)點(diǎn)而被廣泛采用［4-6］.然而流線型人工心臟泵內(nèi)的機(jī)械磨損、溶血和血栓等現(xiàn)象仍很難消除［7-9］，這些都與血液的流體動(dòng)力學(xué)行為相關(guān)，而流動(dòng)中大尺度分離是造成升壓比下降和血細(xì)胞破壞的主要原因［10-13］.為了解設(shè)計(jì)的心臟泵內(nèi)流動(dòng)細(xì)節(jié)，本文采用計(jì)算流體動(dòng)力學(xué)（Computational Fluid Dynamics，CFD）方法［14-16］，對(duì)人工心臟泵內(nèi)的流場(chǎng)和性能進(jìn)行了模擬，得到心臟泵內(nèi)部的詳細(xì)流場(chǎng)，通過(guò)CFD數(shù)值計(jì)算結(jié)果，找出影響人工心臟泵性能的流動(dòng)因素，分析心臟泵流場(chǎng)，特別是心臟泵內(nèi)流動(dòng)分離現(xiàn)象，對(duì)造成人工心臟泵內(nèi)部流動(dòng)分離、流動(dòng)紊亂的因素進(jìn)行分析，避免或改進(jìn)設(shè)計(jì)過(guò)程中出現(xiàn)的分離現(xiàn)象.為設(shè)計(jì)更符合血液動(dòng)力學(xué)要求的人工心臟泵，對(duì)人工心臟泵改進(jìn)提出相應(yīng)措施.

1 數(shù)學(xué)模型及數(shù)值方法

文中計(jì)算采用的質(zhì)量守恒方程與Navier-Stokes方程形式為

式中：ρ為流體密度；p為流體壓強(qiáng)；u為流體速度；F為體積力；μ為粘性系數(shù).

湍流模型采用k-ε雙方程模型

式中：Gk表示由于平均速度梯度引起的湍動(dòng)能；Gb表示由于浮力影響引起的湍動(dòng)能；YM表示可壓縮湍流脈動(dòng)膨脹對(duì)總的耗散率的影響；湍流粘性系數(shù)在計(jì)算過(guò)程中，取C1ε=1.44，C2ε=1.92，C3ε=0.09，湍動(dòng)能k與耗散率ε的湍流普朗特?cái)?shù)分別為σk=1.0，σε=1.3.

在壁面采用標(biāo)準(zhǔn)壁面函數(shù)法，銜接粘性底層低雷諾數(shù)區(qū)與湍流充分發(fā)展區(qū)，即在求解近壁區(qū)的k，ε時(shí)，將k-ε方程中的對(duì)流項(xiàng)及擴(kuò)散項(xiàng)略去，得到k，ε的近似表達(dá)式

采用SIMPLE算法進(jìn)行壓強(qiáng)-速度耦合.計(jì)算區(qū)域采用非結(jié)構(gòu)化網(wǎng)格，經(jīng)過(guò)網(wǎng)格無(wú)關(guān)解的測(cè)試，最終計(jì)算采用節(jié)點(diǎn)數(shù)262 586.計(jì)算過(guò)程中采用的時(shí)間步長(zhǎng)10-4s.

2 邊界條件

設(shè)置計(jì)算邊界條件如表1所示.葉輪轉(zhuǎn)速10 000 r／min.在計(jì)算的流速范圍內(nèi)設(shè)置血液是不可壓縮的流體，密度為1.055×103kg／m3，粘度為4×10-3Pa·s.進(jìn)出口壓差為13 332.2 Pa.其它壁面定義為無(wú)滑移固壁邊界條件.

表1 計(jì)算邊界條件 Tab.1 Boundary conditions of computation

3 計(jì)算結(jié)果分析

對(duì)泵內(nèi)速度及壓力分布進(jìn)行了求解，圖1為泵內(nèi)壓力云圖，泵進(jìn)出口壓力差13 332.2 Pa，但整個(gè)泵最高壓力為18 000 Pa，出現(xiàn)在動(dòng)葉片出口處，由于局部阻力與沿程損失，在靜止導(dǎo)葉片出口處，表壓降到13 500 Pa左右，泵內(nèi)最低的壓力出現(xiàn)在動(dòng)葉片進(jìn)口位置.

圖1 心臟泵內(nèi)壓力云圖 Fig.1 Pressure cloud picture of the artificial heat pump

圖2為泵軸的中間截面流速分布情況.計(jì)算結(jié)果顯示，泵內(nèi)速度在0.2～6.9 m／s之間.

圖2 泵中間截面速度分布 Fig.2 Velocity distribution of the pump center section

從這些區(qū)域可以看到，流動(dòng)未受到葉片的干擾，沒(méi)有明顯的流動(dòng)分離、渦流、滯流等現(xiàn)象出現(xiàn).其中0.2～3 m／s的速度分布區(qū)域，6～6.9 m／s的速度分布區(qū)域分別如圖3和圖4所示.

圖3 0.2～3 m／s速度分布區(qū)域 Fig.3 Velocity domain of 0.2～3 m／s

圖4 6.0～6.9 m／s速度分布區(qū)域 Fig.4 Velocity domain of 6.0～6.9 m／s

在動(dòng)葉片區(qū)，由于葉片的高速旋轉(zhuǎn)運(yùn)動(dòng)，泵內(nèi)最高流速達(dá)到6.9 m／s，泵內(nèi)動(dòng)葉片區(qū)速度在3～6.9 m／s，流速較高的區(qū)域主要集中在動(dòng)葉片出口位置之前的一片區(qū)域.從圖5可以看出，最高流速在輪緣處.在計(jì)算過(guò)程中，對(duì)動(dòng)葉片外緣與外部套筒之間的間隙做了簡(jiǎn)化處理，但通過(guò)分析可以得出，在固定的外部套筒壁面附近，由于流體的粘性和受壁面摩擦阻力的影響，緊貼壁面處會(huì)形成速度梯度很大的邊界層，這一微小的間隙內(nèi)必定產(chǎn)生較大的剪切應(yīng)力.很容易產(chǎn)生溶血和血栓，對(duì)血液造成破壞.因此，應(yīng)正確確定間隙，既能減少軸向泄漏，又要保證適當(dāng)?shù)募羟袘?yīng)力.進(jìn)口導(dǎo)流區(qū)域及出口導(dǎo)葉區(qū)域流速較低，流動(dòng)速度在0.2～3 m／s之間.

圖6為動(dòng)葉片輪緣周向0～180°范圍內(nèi)的速度矢量，可以看出在血液由導(dǎo)流體區(qū)流入動(dòng)葉片區(qū)時(shí)，由于流動(dòng)速度和流動(dòng)方向的忽然改變，流動(dòng)出現(xiàn)了重復(fù)循環(huán)、流動(dòng)分離、渦流等不規(guī)則的流動(dòng)模式.流體流出動(dòng)葉片區(qū)進(jìn)入導(dǎo)流區(qū)時(shí)也出現(xiàn)了微小的流動(dòng)分離.為改善流體進(jìn)入動(dòng)葉片區(qū)產(chǎn)生的流動(dòng)分離，應(yīng)改變導(dǎo)流體形狀，使流出導(dǎo)流體的流動(dòng)方向與動(dòng)葉片入口處切線方向一致.在其它區(qū)域流動(dòng)速度分布比較均勻平緩，符合生理要求.

圖5 動(dòng)葉片徑向速度分布 Fig.5 Radial velocity distribution of moving blade

圖6 動(dòng)葉片區(qū)輪緣周向0～180°區(qū)流速分布Fig.6 Velocity distribution of 0～180°wheel rim

圖7為動(dòng)葉片壓力面與吸力面速度分布，可以看出流體流動(dòng)方向與流線形外形基本一致，也存在微小的流動(dòng)方向發(fā)生改變的情況，表明葉片形狀仍存在改進(jìn)的余地.

圖8為靜止導(dǎo)葉片區(qū)速度分布，可以看出，在流經(jīng)靜止導(dǎo)葉片區(qū)時(shí)，流動(dòng)方向與靜止葉片方向逐漸一致，在出口位置流動(dòng)達(dá)到比較理想的情況.

圖9中左側(cè)位置是流體入口，右側(cè)為流動(dòng)出口，從入口不同位置分別標(biāo)識(shí)了7條流體的流動(dòng)跡線，可以看出，除了進(jìn)入和流出動(dòng)葉片位置處流動(dòng)方向發(fā)了較明顯的改變外，在其它區(qū)域流動(dòng)方向與葉片形狀比較吻合.這也表明，導(dǎo)流體的外形需要做出改進(jìn).

圖7 葉片壓力面與吸力面速度分布 Fig.7 Velocity distribution of pressure surface and the suction surface

圖8 靜止導(dǎo)葉片速度分布 Fig.8 Velocity distribution of stationary guide blade

圖9 自入口到出口按速度標(biāo)識(shí)的跡線 Fig.9 Velocity trace from the entrance to the exit

4 結(jié) 論

1）從泵總體流場(chǎng)看，在導(dǎo)流區(qū)、動(dòng)葉葉片區(qū)、導(dǎo)葉區(qū)流動(dòng)比較平緩，流動(dòng)跡線與流線型外形比較一致，沒(méi)有明顯的流動(dòng)分離、渦流、滯流等現(xiàn)象出現(xiàn)；

2）在動(dòng)葉片外緣速度較高，動(dòng)葉片外緣與外部套筒之間存在邊界層區(qū)，這一區(qū)域有較大的速度梯度，產(chǎn)生較大的剪切應(yīng)力，是產(chǎn)生溶血、血栓和血液破壞的重要原因；

3）在動(dòng)葉片入口及出口位置存在流動(dòng)的突然改變，特別是動(dòng)葉片入口區(qū)，由于導(dǎo)流體外形及流向與動(dòng)葉片入口設(shè)計(jì)問(wèn)題，導(dǎo)致這一區(qū)域產(chǎn)生了重復(fù)循環(huán)、流動(dòng)分離、渦流等不規(guī)則的流動(dòng)模式；

4）動(dòng)葉片壓力面、吸力面及靜止導(dǎo)葉片部分的流動(dòng)方向與流線型葉片外型基本一致，也存在微小的流動(dòng)方向發(fā)生改變的情況，表明葉片外形的進(jìn)一步優(yōu)化可更加改善流場(chǎng)情況，更符合血液的流體動(dòng)力學(xué)行為；

5）帶動(dòng)葉片與外部套筒間隙、改進(jìn)進(jìn)口導(dǎo)流體模型等不同流動(dòng)模型的計(jì)算結(jié)果仍需進(jìn)一步研究.

［1］Tzallas A T，Katertsidis N S，Karvounis E C，et al.Modeling and simulation of speed selection on left ventricular assist devices［J］.Computers in Biology and Medicine，2014，51：128-139.

［2］Babamir S M，Dehkordi M B.Specification and simulation of behavior of the continuous infusion insulin pump system［J］.Bio-Medical Materials and Engineering，2014，24（2）：1517-1526

［3］Marsden A L，Bazilevs，Y，Long C C，et al.Recent advances in computational methodology for simulation of mechanical circulatory assist devices［J］.Wiley Interdisciplinary Reviews-Systems Biology and Medicine，2014，6（2）：169-188.

［4］阮曉東，陳松松，錢偉文，等.基于溶血性能的離心式旋轉(zhuǎn)血泵設(shè)計(jì)［J］.中國(guó)生物醫(yī)學(xué)工程學(xué)報(bào)，2011，30（3）：411-415.Ruan Xiaodong，Chen Songsong，Qian wenwei，et al.Design of centrifugal blood pump based on hemolysis estimation［J］.Chinese Journal of Biomedical Engineering，2011，30（3）：411-415.（in Chinese）

［5］Bakouri M A，Salamonsen R F，Savkin A V，et al.A sliding mode-based starling-like controller for implantable rotary blood pumps［J］.Artificial Organs，2014，38（7）：587-593.

［6］Niklas A，Hiller K A，Jaeger A，et al.In vitro optical detection of simulated blood pulse in a human tooth pulp model［J］.Clinical oral investigations，2014，18（5）：1401-1409.

［7］Cruz D O A，Pinhof T.Turbulent pipe flow prediction with a low Reynolds number type closure［J］.Journal of Non-Newtoinan Fluid Mechanics，2003，114（2／3）：149-184.

［8］劉云，李世嘉，侯擊波，等.精密鋁鑄件電磁泵負(fù)壓充型測(cè)試技術(shù)［J］.測(cè)試技術(shù)學(xué)報(bào)，2009，23（3）：213-216.Liu Yun，Li Shijia，Hou jibo，et al.Testing technique of negative pressure mold-filling of electromagnetic pump for precision aluminum casting［J］.Journal of Test and Measurement Technology，2009，23（3）：213-216.（in Chinese）

［9］Peng Y H，Wu Y Q，Tang X Y，et al.Numerical simulation and comparative analysis of flow field in axial blood pumps［J］.Computer Methods in Biomechanics and Biomedical Engineering，2014，17（7）：723-727.

［10］Wu Huachun，Wang Ziyan，LüXujun.Desigh and fluent simulation of impeller for axial maglev heart pump［J］.Applied Mechanics and Materials，2012，195：29-34.

［11］蘇波，王登偉.一種新型無(wú)閥微泵的研究［J］.測(cè)試技術(shù)學(xué)報(bào)，2011，25（2）：107-111.Su Bo，Wang Dengwei.Research of a novel valveless micropump［J］.Journal of Test and Measurement Technology，2011，25（2）：107-111.（in Chinese）

［12］Abdoli A，Dulikravich G S，Bajaj C，et al.Human heart conjugate cooling simulation：unsteady thermofluid-stress analysis［J］.International Journal for Numerical Methods in Biomedical Engineering，2014，30（11）：1372-1386.

［13］Long C C，Esmaily-Moghadam M，Marsden A L，et al.Computation of residence time in the simulation of pulsatile ventricular assist devices［J］.Computational Mechanics，2014，54（4）：911-919.

［14］Song G L，Chua L P，Lim T M.Numerical study of of bio-centrifugal blood pump with straight impeller blade profiles［J］.Artificial Organs，2009，34（3）：98-104.

［15］Lim E，Salamonsen R F，Mansouri M，et al.Hemodynamic response to exercise and head-up tilt of patients implanted with a rotary blood pump：a computational modeling study［J］.Artificial Organs，2015，39（2）：24-35.

［16］Kafagy D H，Dwyer T W，Mc Kenna K L，et al.Design of axial blood pumps for patients with dysfunctional fontan physiology：computational studies and performance testing［J］.Artificial Organs，2015，39（1）：34-42.