組織工程三維多孔骨支架內(nèi)部微流體流場(chǎng)研究①

趙士明 李文雷 趙靜一③ 郭 銳 高 遠(yuǎn)* 盧子帥

(* 燕山大學(xué)河北省重型機(jī)械流體動(dòng)力傳輸與控制實(shí)驗(yàn)室 秦皇島 066004)(**先進(jìn)鍛壓成形技術(shù)與科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室(燕山大學(xué)) 秦皇島 066004)(***唐山工業(yè)職業(yè)技術(shù)學(xué)院機(jī)械工程系 唐山 063299)(****河北省唐山市第二人民醫(yī)院創(chuàng)傷外科研究所 唐山 063001)

0 引 言

組織工程(tissue engineering)是利用工程學(xué)和生命科學(xué)的基本原理，開(kāi)發(fā)能恢復(fù)、維持或改善受損組織或器官功能的生物代替物。組織工程支架作為細(xì)胞和新生組織的載體要保證營(yíng)養(yǎng)液、生長(zhǎng)因子等能夠從多孔或微管進(jìn)入支架內(nèi)部，與細(xì)胞進(jìn)行物質(zhì)交換排出代謝物[1，2]。組織工程支架為細(xì)胞接種、繁殖和保持各自功能提供了保障，其可以引導(dǎo)組織再生，控制組織結(jié)構(gòu)。研究表明，流體剪應(yīng)力等各種外載力的刺激，對(duì)骨細(xì)胞的增殖分化具有重大影響。支架內(nèi)部微孔的尺寸、空間走向、結(jié)構(gòu)及布局都將影響其內(nèi)部微流體流場(chǎng)的分布。營(yíng)養(yǎng)液在管道里究竟能給多大范圍內(nèi)的細(xì)胞提供足夠的養(yǎng)分,與支架的內(nèi)部空間結(jié)構(gòu)、孔徑和微流體的流速等有密切關(guān)系[3]。

現(xiàn)有的骨組織工程支架在實(shí)驗(yàn)過(guò)程中存在局部新生組織生長(zhǎng)較快的問(wèn)題，以致于營(yíng)養(yǎng)液、生長(zhǎng)因子和代謝物的流動(dòng)被阻塞，其主要原因是未對(duì)骨支架內(nèi)部流場(chǎng)的作用進(jìn)行研究。目前，針對(duì)支架內(nèi)部的速度場(chǎng)和壓力場(chǎng)分布的研究已有報(bào)道，文獻(xiàn)[4]通過(guò)實(shí)驗(yàn)和計(jì)算流體動(dòng)力學(xué)數(shù)值模擬方法研究了骨組織工程支架內(nèi)部流場(chǎng)分布，2種方法所得結(jié)果相近，能夠相互驗(yàn)證，為量化骨組織工程支架內(nèi)部流體力學(xué)環(huán)境和細(xì)胞生長(zhǎng)、分化的關(guān)系提供了基礎(chǔ)。文獻(xiàn)[5]提出了將曲率驅(qū)動(dòng)細(xì)胞增長(zhǎng)與流體動(dòng)力學(xué)結(jié)合的模型。燕山大學(xué)趙靜一團(tuán)隊(duì)提出了組織工程中的微流體現(xiàn)象[3]，對(duì)近幾年骨組織工程支架內(nèi)微孔結(jié)構(gòu)[6]、微流體流場(chǎng)研究[7]和微流體驅(qū)動(dòng)控制研究[8]進(jìn)行了綜述，并對(duì)骨組織支架內(nèi)部微孔結(jié)構(gòu)進(jìn)行了設(shè)計(jì)和仿真研究[9]。文獻(xiàn)[10]研究了3種微管道中細(xì)胞濃度分布，并分析了微管道系統(tǒng)中細(xì)胞可能沉積的部位。生物體內(nèi)流體剪應(yīng)力、靜水壓力、機(jī)械拉伸力等機(jī)械力共同作用是骨組織生長(zhǎng)和維持的關(guān)鍵，其中流體剪應(yīng)力對(duì)細(xì)胞、組織的刺激更為突出[11]。在組織體外培養(yǎng)中各種力學(xué)刺激是必不可少的，得到越來(lái)越多的重視[12，13]。張寅星等人[14]研究了不同剪應(yīng)力水平下NO、ALP的含量，發(fā)現(xiàn)較低剪應(yīng)力條件下，成骨細(xì)胞的分化與剪應(yīng)力成一定的比例關(guān)系；而在高剪應(yīng)力水平下，兩者之間沒(méi)有明顯的比例關(guān)系。文獻(xiàn)[15]分析了不同流場(chǎng)下的剪應(yīng)力分布，在震蕩壁面邊界下，剪應(yīng)力水平明顯提高。然而，目前對(duì)骨組織工程多孔支架結(jié)構(gòu)內(nèi)部流場(chǎng)及剪應(yīng)力分布的研究較少。

本文根據(jù)工程學(xué)和生物醫(yī)學(xué)原理，設(shè)計(jì)了4種不同內(nèi)部結(jié)構(gòu)的骨組織工程支架，分析其在不同邊界條件作用下，不同支架內(nèi)部速度、壓力、流體剪應(yīng)力的大小及分布情況，通過(guò)理論分析和數(shù)值仿真發(fā)現(xiàn)，針對(duì)人體各部分組織結(jié)構(gòu)的不同，設(shè)計(jì)不同結(jié)構(gòu)形式的支架，能夠調(diào)整營(yíng)養(yǎng)液的輸送速度，通過(guò)調(diào)節(jié)流體剪應(yīng)力的范圍，可以增加骨細(xì)胞的生長(zhǎng)速度，該研究具有一定的理論和實(shí)際意義。

1 建立骨組織工程支架物理模型和數(shù)學(xué)模型

1.1 物理模型

為了確定合適的組織工程支架內(nèi)部結(jié)構(gòu)，對(duì)不同結(jié)構(gòu)內(nèi)部微流體流場(chǎng),使用CREO3.0建立了4種三維模型，并使用ICEM CFD進(jìn)行網(wǎng)格劃分如圖1～圖4所示。在人體骨組織中具有水平的浮克曼管和豎直的哈佛管，中心主管通過(guò)橫向的浮克曼管與豎直的哈佛管相連，依據(jù)這種特征設(shè)計(jì)了3種仿生結(jié)構(gòu)支架和一種螺旋支架。4種支架中圖1是最基礎(chǔ)的仿生模型，主管道直徑為2 mm，長(zhǎng)度為20 mm，支管直徑為1.2 mm。圖2和圖3是在支架1基礎(chǔ)上的改進(jìn)型，而圖4是從工程流體力學(xué)角度出發(fā)設(shè)計(jì)的內(nèi)部結(jié)構(gòu)。圖1、圖2、圖3所示3種支架內(nèi)部主管和支管尺寸相同，圖1、圖3支管與主管道垂直，圖2浮克曼管與主管道夾角呈30 °，圖3將哈佛管移至浮克曼管端部并且將每層浮克曼管間距離增大1.5倍，增加了連接哈佛管的輔助支管。圖4中螺旋管道支管直徑為1.2 mm。

網(wǎng)格劃分是計(jì)算流體力學(xué)數(shù)值模擬中重要的一環(huán)，生成一套高質(zhì)量網(wǎng)格將顯著提高計(jì)算精度和收斂速度。ICEM CFD因其友好的操作界面、完善的幾何功能、靈活的拓?fù)鋭?chuàng)建、先進(jìn)的O型網(wǎng)格技術(shù)得到快速發(fā)展。圖1、圖2、圖3三種支架因其結(jié)構(gòu)相對(duì)復(fù)雜、相貫界面較多，劃分結(jié)構(gòu)網(wǎng)格將非常復(fù)雜，消耗大量時(shí)間，所以采用非結(jié)構(gòu)網(wǎng)格。圖4結(jié)構(gòu)相對(duì)簡(jiǎn)單，使用結(jié)構(gòu)網(wǎng)格，可以大幅度提高網(wǎng)格質(zhì)量和計(jì)算時(shí)間。

圖1 仿生支架1

圖2 仿生支架2

圖3 仿生支架3

圖4 螺旋結(jié)構(gòu)支架4

1.2 數(shù)學(xué)模型

骨組織工程支架內(nèi)的微流孔道雖然達(dá)到了微米級(jí)別，在其流動(dòng)特性上呈現(xiàn)了與宏觀的不同，但是根據(jù)前人的大量研究來(lái)看，支架孔道內(nèi)的流動(dòng)仍然滿足宏觀的物理守恒定律，包括質(zhì)量守恒定律和動(dòng)量守恒定律和能量守恒定律?？變?nèi)流體流動(dòng)的質(zhì)量守恒方程即連續(xù)方程可用式(1)表示；動(dòng)量方程的本質(zhì)是滿足牛頓第二定律，可用式(2)表示；能量守恒定律本質(zhì)是熱力學(xué)第一定律，可用式(3)[7]表示。營(yíng)養(yǎng)液粘度較高，流動(dòng)速度低，由雷諾公式(4)可知Re≈12，其液在支架內(nèi)的流動(dòng)狀態(tài)為層流，其他管道與主管道相比直徑小，所以雷諾數(shù)更低。假設(shè)營(yíng)養(yǎng)液為不可壓縮流體，其運(yùn)動(dòng)微分方程可以簡(jiǎn)化為式(5)的形式[7]。骨組織工程支架內(nèi)的微流孔道達(dá)到了微米級(jí)別，在其流動(dòng)特性上呈現(xiàn)了與宏觀的不同。

(1)

式中，ρ為流體密度，單位為kg/m3。

(2)

式(2)中，τxx、τxy、τxz等是因粘性作用而產(chǎn)生的作用在微元體表面上的粘性應(yīng)力τ的分量；fx、fy、fz為x、y、z方向上的單位質(zhì)量力。

(3)

其中，E為流體微團(tuán)的總能，h為焓，keff為有效熱傳導(dǎo)系數(shù)，Sh為體積源項(xiàng)。

(4)

式(4)中，u為平均流速，單位為m/s；d為管道內(nèi)徑，單位為m；μ為動(dòng)力粘度系數(shù)，單位為kg/m·s。

(5)

式中，p為流體微元表面壓力，v為流體運(yùn)動(dòng)粘度。

流體在流動(dòng)過(guò)程中因受到管壁的摩擦阻力和管道直徑變化的影響導(dǎo)致流動(dòng)方向改變而產(chǎn)生流動(dòng)阻力，主要有沿程阻力和局部阻力。由沿程阻力引起的能量損失與流程長(zhǎng)度成正比，簡(jiǎn)稱沿程損失，用式(6)計(jì)算；過(guò)流斷面的大小、形狀、方位等發(fā)生變化時(shí)流體所受到的阻力是各種各樣的，集中在很短的流段內(nèi)，這種阻力為局部阻力，所引起的能量損失可用式(7)表示。

(6)

(7)

式中，λ為沿程阻力系數(shù)，l為流程，d為管道內(nèi)徑，v為流速，ζ為局部阻力系數(shù)。

2 通過(guò)有限元進(jìn)行數(shù)值模擬與分析

根據(jù)所建立的物理模型和數(shù)值模型，對(duì)骨組織工程支架內(nèi)部流動(dòng)的壓力與流速進(jìn)行了數(shù)值模擬分析。

2.1 壓力和速度分析

把ICEM CFD中劃分好的網(wǎng)格文件導(dǎo)入到Fluent中進(jìn)行仿真研究，營(yíng)養(yǎng)液的密度ρ為1×103kg/m3，粘度μ為0.0018 Pa·s為不可壓縮牛頓流體。假設(shè)流體與壁面的接觸邊界是無(wú)滑移壁面，選取管道入口端面為速度入口邊界條件，營(yíng)養(yǎng)液流速恒定為0.01 m/s，出口邊界壓力為0 Pa。仿生骨支架速度分布云圖，如圖5～圖8所示。

圖5 仿生支架1內(nèi)流體速度分布云圖

圖6 仿生支架2內(nèi)流體速度分布云圖

圖7(a)、7(b)為支架3在主視面和左側(cè)面的速度分布圖。由圖5、圖6、圖7速度分布云圖可知，這3種結(jié)構(gòu)支架沿著主管道從入口到出口速度變化明顯，由于沿程阻力和局部阻力以及各個(gè)支管的分流作用，速度逐漸減小。在第一層浮克曼管處的分流現(xiàn)象明顯，分流速度接近入口速度的一半。支架2和支架1、支架3相比較，其浮克曼管與主管夾角呈60 °，流體的動(dòng)量與分流支管的夾角小，流體流動(dòng)的局部阻力小，所以在支架2中支管的流速明顯高于支架1、支架3。正是因?yàn)榱黧w動(dòng)量的矢量性以及主管和哈佛管之間距離不足，無(wú)法完全消除向下的動(dòng)量分量，支架1出現(xiàn)在第一層浮克曼管處向上的速度比向下的速度低。由于支架1的浮克曼管的另一端直接出口分流作用明顯，而支架2、支架3的哈佛管與浮克曼管的另一端相連，支架2的哈佛管中流體的流速高于支架1。支架3第一層浮克曼管與哈佛管連接位置向上出口的距離小于下出口，阻力比較小，大部分流體從上出口流出。由圖8所示螺旋支架的速度分布云圖可知，其流動(dòng)速度從上到下逐漸變緩，流動(dòng)狀態(tài)均勻，脈動(dòng)小。

(a) 主視

(b) 左視

(a) 縱截面

(b) 橫截面

圖8所示支架4的最高流速在螺旋管中間達(dá)到了0.0196 m/s，支架4內(nèi)部螺旋管之間相互獨(dú)立并不存在分流現(xiàn)象。由于螺旋管不存在交叉口沒(méi)有局部阻力，但是沿程阻力會(huì)有所提高。流體在螺旋管中由于離心力的作用流速逐漸提高，阻尼隨之增大，直到阻力與離心力平衡。在4種結(jié)構(gòu)中支架4的流場(chǎng)分布更為均勻，并且可以通過(guò)減螺距、增加螺旋管的數(shù)量提高孔隙率。

仿生骨支架壓力分布云圖如圖9～圖12所示。

由于出口邊界條件為零壓，從整體趨勢(shì)上看4種支架從入口到出口壓力逐步降低到0。支架1中壓力作用范圍最小，只在入口附近存在較高的壓力。在支架2、3中哈佛管和浮克曼管交叉匯合，并且只在上下兩端存在出口，造成較高的壓力以及較大的作用范圍。對(duì)于支架3第一層交叉口處的壓力梯度高也造成了向上出口的流速高。在支架4中流體的流程遠(yuǎn)大于其他3種結(jié)構(gòu)，因此沿程阻力大，壓力損失高，會(huì)在入口處形成高壓，并且壓力梯度大，3層螺旋管中可以明顯看出從外層到內(nèi)層壓力逐步降低。

2.2 流體剪應(yīng)力分析

體外培養(yǎng)過(guò)程中流體剪應(yīng)力難以測(cè)量，實(shí)驗(yàn)復(fù)雜，數(shù)值模擬為支架內(nèi)流體剪應(yīng)力的分布及水平提供了良好的條件。牛頓流體對(duì)壁面的剪應(yīng)力可用如下公式計(jì)算：

(8)

式中，τ為流體剪應(yīng)力，μ為動(dòng)力粘度系數(shù)，vx為流體在x軸方向的速度分量，vy為流體在y軸方向的速度分量。

由式(8)可以看出剪應(yīng)力大小與流體動(dòng)力粘度、速度梯度成正比。目前，壁面流體剪應(yīng)力的調(diào)制主要有增加流體粘度、增大入口速度以及在支架內(nèi)設(shè)置可動(dòng)的芯片，但細(xì)胞間傳質(zhì)有最優(yōu)速度范圍，營(yíng)養(yǎng)液有最佳粘度等限制。通過(guò)改變支架的進(jìn)出口邊界條件，使用Fluent UDF自定義函數(shù)功能改變速度梯度，在不增加近壁面速度的條件下進(jìn)行流體剪應(yīng)力的調(diào)制。

由文獻(xiàn)[5,11]可知支架內(nèi)平均剪應(yīng)力0.232～0.304 Pa時(shí)能顯著促進(jìn)細(xì)胞分化,平均剪應(yīng)力低于0.1 Pa影響了細(xì)胞分化的速度，而剪應(yīng)力大于0.304 Pa時(shí)，細(xì)胞分化水平不再明顯變化。圖13～圖16為4種支架結(jié)構(gòu)在優(yōu)化參數(shù)后邊界條件下的剪應(yīng)力分布曲線。通過(guò)曲線對(duì)比，調(diào)制后4種結(jié)構(gòu)支架的流體剪應(yīng)力顯著提高，接近或達(dá)到最優(yōu)的范圍，支架1、2、3剪應(yīng)力提高了0.1 Pa左右，從入口到出口逐漸降低，與速度分布趨勢(shì)一致。對(duì)于支架4流體剪應(yīng)力提高了0.2 Pa，在很小的范圍內(nèi)波動(dòng)，整體上剪應(yīng)力分布均勻,支架4達(dá)到了對(duì)細(xì)胞分化最優(yōu)的剪應(yīng)力。

圖9 仿生支架1內(nèi)流體壓力分布云圖

圖10 仿生支架2內(nèi)流體壓力分布云圖

(a) 主視

(b) 左視

(a) 縱截面

(b) 橫截面

圖13 仿生支架1剪應(yīng)力調(diào)制曲線圖

圖14 仿生支架2剪應(yīng)力調(diào)制曲線圖

圖15 仿生支架3剪應(yīng)力調(diào)制曲線圖

圖16 螺旋支架4剪應(yīng)力調(diào)制曲線圖

3 結(jié) 論

使用CREO3.0通過(guò)仿生學(xué)與螺旋理論建立了4種整體結(jié)構(gòu)尺寸相近、內(nèi)部微孔道布局走向不同的骨組織工程支架三維模型；運(yùn)用N-S方程結(jié)合流體力學(xué)管路能量損失理論建立了4種支架的數(shù)學(xué)模型；使用ICEM CFD劃分網(wǎng)格，并利用Fluent進(jìn)行了數(shù)值仿真，得到了組織工程支架內(nèi)部的流動(dòng)狀態(tài)，并對(duì)骨支架內(nèi)微流體剪應(yīng)力進(jìn)行了優(yōu)化，所得結(jié)果如下。

(1) 支架2的速度場(chǎng)分布更均勻，優(yōu)于支架1和3。支架2哈佛管流速高于支架1和3，提高了營(yíng)養(yǎng)物質(zhì)的交換及代謝廢物的排出速度；螺旋支架4內(nèi)部速度場(chǎng)由于離心慣性力和沿程阻力的共同作用，流速較高，整個(gè)速度場(chǎng)分布均勻，但是流程長(zhǎng)。

(2) 4種支架的壓力流場(chǎng)中支架4更均勻，由入口到出口壓力逐步降低，支架1的壓力作用面積最小，支架4的壓力值最大，內(nèi)外層壓力梯度較高。

(3) 改變進(jìn)出口邊界條件進(jìn)行流體剪應(yīng)力的調(diào)制使其達(dá)到對(duì)細(xì)胞增殖、分化的最佳值。通過(guò)流體剪應(yīng)力曲線對(duì)比發(fā)現(xiàn)，支架1、2、3的剪應(yīng)力從入口到出口逐步減小，調(diào)制后剪應(yīng)力提高約2倍；支架4剪應(yīng)力分布均勻，調(diào)制后提高到0.3 Pa，提高了成骨細(xì)胞的生長(zhǎng)速度。

使用ANSYS Fluent對(duì)所設(shè)計(jì)的4種骨組織工程支架內(nèi)部微流體流場(chǎng)的模擬計(jì)算與分析，通過(guò)支架結(jié)構(gòu)選擇和剪應(yīng)力的優(yōu)化，可以使支架內(nèi)部的微流體流場(chǎng)分布更均勻，剪應(yīng)力更符合人體組織生長(zhǎng)需要，進(jìn)而更加適合骨細(xì)胞的增殖分化、營(yíng)養(yǎng)物質(zhì)的交換及代謝廢物的排出，可為組織工程骨支架的設(shè)計(jì)和使用提供一定的理論依據(jù)。