血管中血液的粘彈性流動分析

田長留,花少震

(河南工學院 機械工程學院,河南 新鄉 453003)

0 引言

血液是一種復雜的懸浮液體,主要由紅細胞、白細胞、血小板、電解質、水等組成[1]。人體中的血液在營養輸送、維持體內環境平衡和防御功能等方面起著重要作用對血液的研究已經成為醫學科學的一個重要組成部分,其中血液在心血管中的流動問題是當前生物力學研究中較為活躍的分支之一。隨著計算機技術的不斷提高,越來越多的科研工作者采用數值手段模擬血管中血液的流動[2-4],目的是根據血液動力學因素(壓力、剪切速率等)和血液流變學性質(粘度、粘彈性)為一些現象(觸變性、血小板激活、凝血塊)以及血管類疾病(動脈硬化、動脈瘤等)的研究提供參考和借鑒。現有的血液數值計算中常將血液視做理想流體或者僅考慮粘性的非牛頓流體[5],然而血漿中含有7%—10%的大分子蛋白質,具有明顯的粘彈性(viscoelastic)特征。粘彈性流體同時表現出粘性和彈性,其流變特征在很多方面都表現出和粘性流體很不同的規律,有時候甚至出現截然相反的流動現象[6]。雖然血液具有粘彈性性質,但是由于粘彈性數值求解的限制,目前對血液流動的數值模仿真中一般忽略血液的彈性,僅將血液視作理想流體或粘性流體。當血液在一個直徑比較大的血管(如主動脈等)中流動時,此時血液流動的雷諾數Re比較大,可將血液視作理想流體[7];研究血液和血管流固耦合仿真時也將血液視作理想流體[8];當血液在直徑較小的血管(如毛細血管、小動脈、小靜脈等)中流動時,此時剪切速率較小,剪切變稀效應明顯,可將血液視作非牛頓流體[9]。另一種理論認為,當血液流動產生的剪切速率小于100s-1時,血液的非牛頓流體性質明顯;當剪切速率大于100s-1時,可將血液視作牛頓流體。但是,相同位置處的血液在一個循環周期內產生的剪切速率始終是變化的。以大動脈為例,在一個循環周期內,血液流動的剪切速率可從0增長到1000s-1,因此究竟將血液視作牛頓流體還是非牛頓流體,要視情況而定[10]。

盡管對血液在血管中流動仿真的關注日益增強,但是仍然缺少關于血液粘彈性流動的數值模擬研究。本文擬在考慮血液流動過程中的粘彈性流變行為基礎上,構建粘彈性算法,模擬血液在血管中的粘彈性流動,探尋血液在血管中流動的粘彈性質。

1 血液在血管中流動粘彈性數學模型

假設血液在血管中流動時為不可壓縮流體,忽略其溫度的變化,根據血液在血管中流動的質量守恒和動量守恒,則可得到控制方程:

·u=0

(1)

(2)

式中,ρ為血液密度,u、p分別為速度和壓力,ηs(u+uT)代表粘性所貢獻偏應力,τ通過求解粘彈本構方程獲得。采用Giesekus本構模型表征血液的粘彈流變性質,則有:

(3)

(4)

2 數值算法

采用有限體積法離散控制方程(1)、(2)和粘彈性本構方程(3)、(4)。對于任一時間步長Δt,由于速度、壓力及速度、應力之間的高耦合性,采用PISO(Pressure-Implicit with Splitting of Operators)算求解法速度-壓力,然后將t+Δt時刻的速度作為已知量,進而求解應力。首先,可將式(3)、(4)整理為含瞬態項、對流項和源項的偏微分方程:

(5)

(6)

式中,下標f為對應控制體的f面上的物理量,Vp為控制體體積。求解式(6)即可得到t+Δt時刻的偏應力張量τ。黏彈性流動求解的難點在于算法的不穩定性,為了得到高穩定性算法,本文借鑒DEVSS(Discrete Elastic Viscous Split-Stress)算法,將該算法從有限元算法擴展到有限體積法,則(2)式調整為:

(7)

式中,等號右邊為源項,φ(ut)=-p+·τ+ρg-·(ηsut)。將(7)式采用有限體積法離散,則有:

(8)

式中,n代表單位外法向量,A為控制體的面的面積大小。聯立(1)式和(8)式,即可用PISO算法求得t+Δt時刻血液在血管中流動過程中速度u和壓力p的大小。

3 血液在血管中流動數值模擬

建立人體血管的三維簡化模型,如圖1所示,其中血管的直徑為2mm,入口到A截面、C截面到出口為直管,ABC為等曲率彎管,其曲率半徑為2mm。血液所對應的粘彈性參數為:η=0.0025Pa·s,ηs=0.0025Pa·s,λ=0.015s,α=0.3。血液的流動速度為0.13m/s。采用上文所發展粘彈性算法模擬血液在圖1所示血管中的粘彈性流動,并取A、B兩個截面為研究對象,分析其粘彈性流動。由于C截面為彎管的末端,直管的始端,其流動特征和A截面差別較小,在此不單獨分析C截面的粘彈性流動。

圖1 血管三維簡化模型

圖2(a)和圖2(b)分別為A、B截面上的二次流,圖中圓截面上部和下部分別代表A、B截面在彎管的內側和外側。A截面為彎管始端,速度矢量為繞截面流動的環流,速度大小由截面外層到圓心逐漸減小。數值結果表明,A截面上環流強度約為主流的百萬分之一。B截面位于ABC彎管的中部,如圖2(a)所示,其截面上速度矢量表現出兩個關于垂直于彎曲方向且兩邊對稱的渦,渦的強度約為主流的千分之一。渦的方向總是由內側中心點沿圓心附近區域向外側中心點流動。隨著流動距離的增加,由于內摩擦等因素,剪切速率不斷降低,慣性力減少而粘性力逐漸增大,在驅動力不變情況下,速度先增大后減小,且離圓心越近速度越大。然后,渦流沿血管外側中心點沿截面外層回流至內側中心點。由于靠近血管壁處的剪切速率比較大,粘彈性流體的第二法向應力差占主導作用,速度不斷增大。

(a)B截面上的渦流

圖2(a)顯示渦由內側中心點沿圓心附近區域向外側中心點流動,由于流體總是從壓力高的區域流向壓力低的區域,所以彎曲血管內側的血液流動壓力高于血管外側的血液流動壓力。圖3為B截面上壓力的分布云圖,其單位為kPa。由于A截面為等直血管的末端,其受直管中流動科式力的影響較大,內外側壓力、速度分布差別較小,本文不作展示和分析。考慮血液的粘彈性質后,從B截面壓力分布可以發現,由于渦關于內外側直徑對稱,所以對于壓力的分布也關于內外側直徑對稱。彎曲血管內側到外側的壓力逐漸減小,這和傳統黏性流體在彎管中流動時內外側壓力分布規律相反。黏性流體在彎管中的流動時外管側的壓力大于內管側的壓力[11]。本文中,考慮血液的彈性后,由于第二法向應力差作用,使得彎曲血管內側壓力大,彎管內側最大壓力為15.25kPa,外側最小壓力為14.25kPa,內外側最大壓力差達到0.8kPa。對于直線型血管A截面,考慮血液流動后,橫截面上壓力差別較小。

圖3 B截面上壓力分布(圓截面上部為彎管內側,單位kPa)

圖4 A、B截面上由內側到外側直徑上的壓力P同平均壓力的比值分布

上文分析了血液在彎曲血管和等截面直管中的粘彈性流動,數值結果研究表明血液的彈性對血液在等截面直血管中的流動影響較小,但是血液的彈性對于血液在彎曲血管中的流動非常重要。可以發現由于血液的粘彈性,彎曲血管內側的血壓要大于外側的血壓,即血液在流動過程中對彎曲血管內側的沖擊更強,也即血管內側的血管壁所需要承受的壓力更大。而且由于內側血液流速緩慢,當血液流動的慣性作用小于黏滯力時,彎曲血液內側中的蛋白質、脂肪等大分子即有可能停止流動,附于血管內壁,導致一些疾病。

圖5 B截面上由內側到外側直徑上的速度大小

4 結語

本文基于Giesekus本構關系建立了血液的等溫不可壓縮的粘彈性流動模型,采用有限體積法離散控制方程,發展了血液的粘彈性流動算法。在對等直血管和彎曲血管的血液流動仿真中發現血液的粘彈性對彎曲血管的作用更為明顯。和粘性流體在彎管中流動時內側壓力低、外側壓力高的特點相反,血液的粘彈性使其在流動過程中對彎曲血管內側的沖擊要強于外側。本文算法能有效模擬血液在血管內的粘彈性流動,可以為一些血液類疾病的病因分析與治療以及相關的生物醫學工程提供借鑒。