胥 义 钟彦骞 吕 娅 王丽萍 周国燕 刘宝林
(上海理工大学生物系统热科学研究所,上海 200093)
血管移植是血管与心脏修补的良好方法,但由于新鲜供体不足,应用受到很大限制。为提高有限供体的利用率,血管组织的低温保存成为人们近年来研究的重要课题。在对血管进行低温保存时,会产生较大的热应力,如果热应力过大,可能会使材料首先产生微观裂纹,微观裂纹在更大的拉应力等各种负载的作用下会扩展成可见的宏观裂纹,这在临床上是非常危险的[1]。所以,准确预测血管冻结过程中的热应力分布,对于成功低温保存血管材料是很关键的。
文献[1]最早针对血管冻结过程热应力分布进行求解,但由于当时血管的低温物理参数匮乏,计算中都是近似地取用冰或水的一些物理参数。然而,不同的生物材料由于其组分的差异,这必然会引起其低温状态下物理参数的显著不同,所以这样的计算结果与实际相差较大。虽然,文献[2-3]采用了部分实测物理参数对血管进行了热应力分析,但他们在选择计算参数时,没有综合考虑降温速率和保护剂浓度对热膨胀系数的影响,以及保护剂浓度对弹性模量的影响。而从文献[4]的研究结果表明,对血管不同的预处理方式以及在冻结过程中不同的冻结方式,都会影响血管在冻结过程中未冻水分额,从而大大影响热膨胀系数以及拉伸力学性质(如弹性模量、极限应力强度)等。因此,有必要全面了解降温速率、低温保护剂浓度对血管冻结过程中应力分布的影响。
此外,在已有的研究中,发现经过快速降温或多次冻结一复温循环后,血管就出现了肉眼可见的宏观裂纹[5],但还未有相关的热应力分析文献报道。因此,本研究还在血管实体模型中引入微裂纹,期望通过不同宽度的微裂纹来模拟多次冻结一复温循环之后的裂纹扩展情况,从而进一步讨论微裂纹的存在对血管冻结过程中应力分布的影响。
为了便于理论分析,假定血管材料在冻结过程中表现为各向同性。如果将血管看作无限长圆筒,则主动脉血管冻结过程就可以简化为沿血管径向变化的问题。其模型示意图如图1所示。
图1 血管边界条件示意(R1—血管内径;R2—血管外径;Tb—外界冷媒温度;hb—冷媒在血管内外壁表面的对流换热系数)Fig.1 Schematic diagram of boundary conditions of considered artery(R1—artery inner diameter;R2—artery outer diameter;Tb—refrigerant temperature;hb—convection heat transfer coefficient)
有关血管冻结过程的传热控制方程和热应力控制方程的描述见文献[1],这里不再赘述。本研究采用ANSYS计算软件,选择直接藕合方法,其单元类型为“SOLID5”,它的自由度包括温度、位移,其计算结果可以同时获得温度场分布、各坐标轴方向的分应力、切应力以及主应力等参数。
以兔胸主动脉血管为计算实体模型,其内径为3 mm,壁厚1 mm,对无微裂纹存在的血管实体模型进行影射划分网格后如图2(a)所示。针对血管有微裂纹的情况,由于V型切口在工程结构或材料中夹杂开裂时经常遇到,因而研究V型切口有一定现实意义[6]。再加上文献[5,7-8]的研究结果表明:血管冻存过程中出现的裂纹开口方向都是与血管轴线垂直的。因此,选择V型切口的弹性裂纹模型,并且垂直于血管轴线方向。如图2(b)是血管有穿透型微裂纹及其有限元网格划分的情况,其表面投影如图2(c)所示,其中W表示裂纹宽度,通过改变裂纹宽度,来模拟裂纹扩展情况。
文献[1]通过对血管热应力本构方程组的分析,确定影响应力分布的几个主要变量为:弹性模量E、泊松比 ν、线性膨胀系数α以及温度分布 T(r,t)等。而其中温度分布T(r,t),又与材料的热导率、表观比热容、密度等参数有关。而这些物理参数的研究目前已经比较细致和完善了,其来源见表1所示。
图2 血管实体模型及其有限元网格的划分。(a)无裂纹实体模型;(b)有微裂纹实体模型;(c)垂直于血管轴向的微裂纹(宽度为W)示意图Fig.2 Solid model and finite element mesh of artery.(a)without any crack;(b)solid model with a microcrack;(c)schematic diagram of a micro-crack with width W perpendicular to artery axis
表1 计算参数选取Tab.1 Literature source of chosen parameters
表2 不同边界条件所对应的平均降温速率Tab.2 Cooling rates corresponding to different boundary conditions
定解条件包括初始条件和边界条件。由于在实际的冻结过程中,内、外边界通常为对流边界条件,那么,通过改变对流换热系数,就可以获得内、外边界上与之相对应的平均降温速率。冷媒温度一般为低温氮气,这里取-100℃。通过大量模拟计算,针对本文研究的血管物理性质,血管外表面对流换热系数与平均降温速率存在如表2所示的对应关系。由于血管材料的冻结点温度都在0℃以下[12],因此,其初始温度为 0℃即可。
在材料力学的研究中,有多种屈服破坏准则,而最大主应力准则是其中一种最常用的准则[13]。在工程应用中,承载构件的破坏通常与构件上发生破坏一点的最大主应力有关。所以,分析结果时只是针对最大主应力来讨论。为便于讨论,选取图3所示的血管横截面,图中“A”、“B”、“C”分别表示血管壁外表面、中心以及内表面处的有限元节点。因此,这些节点处的主应力就分别代表了血管壁外表面、中心以及内表面处的主应力。
图3 血管横截面节点示意图(A—血管壁外表面节点;B—血管壁中心节点;C—血管壁内表面节点)Fig.3 Schematic illustration of finite element nodes in artery cross section(A—outer surface node;B—intermediate node;C—inner surface node)
通过大量计算,发现血管在冻结过程的温度分布中,管壁中心位置与血管壁外表面之间的温差最大。这主要是因为在冻结过程中不断释放出的相变潜热通过内外表面时,由于面积不同而引起内表面热流密度大于外表面,最终体现为外表面比内表面降温快。不同降温速率下管壁中心位置与外壁之间的差值随温度变化如图4所示。很明显,降温速率越大,这个温差越大,而且,最大差值一般都出现在-4~-5℃之间。这主要是因为血管在0~-10℃之间的热导率有一个从0.388到1.05的突变过程[11],在开始冻结时的热导率较小,随着温度的持续降低,必然在这一温区的某一个温度点温差达到最大,而后随血管中可冻水已经大部分冻结,热导率明显增大,传热能力增强,所以两者之间的温差又逐渐减小。
图4 降温过程中血管外壁与血管壁厚中心处之间的温差Fig.4 Temperature differences between the outer and intermediate nodes in artery during freezing
图5表示无低温保护剂时,降温速率为3、5、10、50℃/min时,血管壁中 A、B、C各节点处最大主应力随温度的变化规律。这几种情况下的主应力分布有几个共同点:1)血管壁内、外表面的主应力表现为压应力,中心处表现为拉应力,这与文献[2]是一致的;2)无论是拉应力还是压应力,它们都是随着温度的降低逐渐增大,达到某一峰值后快速降低,最后趋于0值;3)三者中,血管壁外表面所受的应力值最大,内表面次之,中心处的应力值略小于内表面的应力绝对值。至于第一点,主要是因为血管在冻结过程中处于膨胀状态,而血管壁内外表面首先被冻结,一旦内部冻结膨胀必然会受到表面已冻层的约束,所以,内、外表面体现为压应力。对于第二点的解释是:这主要是因为随着温度降低,热膨胀系数也是先增大后快速减小[10],而且,这个转变点温度也是应力发生转变所对应的温度,这也说明血管冻结过程的最大应力是发生在相变区间。而引起第三点的主要原因是血管冻结过程中,因为血管特殊的圆管形状使得每一层同心圆面的面积不同,从而引起每一层面的热流密度不同,最终导致不同层面的温度不同,一般是外表面温度<内表面温度<中心层面温度。所以,计算过程中同一时刻不同层面的物理参数会随温度不同而有差异,必然引起应力分布不同。从而还可以得出一个结论就是:血管冻结过程中,外表面所受应力最大,因而显得最为脆弱,可以认为血管发生裂纹断裂时是从外表面开始的。不过,这有待进一步实验观察来验证。
图5 降温速率对血管内最大主应力的影响。(a)3℃ /min;(b)5℃ /min;(c)10℃ /min;(d)50℃/minFig.5 Effects of cooling rates on maximum thermal stress distributions in artery.(a)3℃/min;(b)5℃ /min;(c)10℃ /min;(d)50℃ /min
图6给出了不同降温速率下血管壁外表面节点最大主应力对比曲线。降温速率增大时,血管冻结过程产生的热应力明显增大,特别是当降温速率达到50℃/min时(最大主应力约为2.4 MPa),其最大应力比3℃/min(最大主应力约为0.023 MPa)增大了两个数量级。从图4可知,造成这种情况是因为降温速率增大后导致血管内部各节点之间的温差增大,由此而引起各自所对应的物理参数不同,最终综合影响应力分布的结果。
图6 降温速率时外壁面最大主应力的影响Fig.6 Effects of cooling rates on maximum thermal stress distributions of outer node
图7考察了添加低温保护剂后对血管冻结过程血管壁中A、B、C各节点处热应力分布的影响。从曲线图中可以看出,其曲线特点除了2.1节所述三个共同点和原因以外,另外还有以下特点:当DMSO浓度增大到10%时,最大应力值不再出现在热膨胀系数的转折点(-32~-35℃温区),而是出现在-4℃左右。从图4可知,血管壁中最大温差一般出现在-4℃左右,这说明这种情况下温差是引起最大应力的最主要原因。究其根本原因,还是由于DMSO浓度达到10%以后,其热膨胀系数[10]和弹性模量[9]明显减小。特别是弹性模量可以说是应力分布的放大器,它本身就减小了2~3个数量级。而根据应力本构方程组中的各参数来看,热膨胀系数和弹性模量的影响显著下降,而温差就上升为主要影响因素。
图7 低温保护剂浓度对血管壁内、外表面以及中心层面最大主应力的影响(降温速率为 5℃/min)。(a)无DMSO;(b)5%DMSO;(c)10%DMSO;(d)15%DMSOFig.7 Effects of cryo-protective agent concentrations on maximum thermalstressdistributionsin arterywith cooling rate of 5℃/min.(a)Without DMSO;(b)5%DMSO;(c)10%DMSO;(d)15%DMSO
图8对比了经过不同DMSO浓度处理后,血管在冻结过程中外表面的最大主应力随温度变化的情况。很显然,经过 DMSO处理后,血管冻结过程产生的热应力显著减小。这是因为DMSO的添加使得血管冻结过程未冻水分额明显增大[13],从而大大改变了血管冻结过程的力学性质,这也说明低温保护剂的添加而避免了热应力对血管的机械损伤。这也说明,当添加低温保护剂后,为了减弱温差对热应力的影响,采用较慢的降温速率是必要的。
图8 低温保护剂浓度对血管外壁面处最大主应力的影响(降温速率为5℃/min)Fig.8 Effects of cryo-protective agent concentrations on maximum thermal stress distributions in outer node a with cooling rate of 5℃/min
在以往的研究中,经过快速降温或多次冻结-复温循环后,血管就出现了肉眼可见的宏观裂纹[5]。从强度力学的角度来看,如果血管在冻结过程中一旦有微裂纹出现,就会存在缺口敏感性问题。而材料力学的研究表明,任何材料中一旦存在一长为W缺口或裂纹,在受力情况下的力线分布如图9所示[14]。可以清楚地看出,离裂尖稍远的地方为一均匀拉伸区,裂纹两侧面划阴影的地方由于不承受拉应力,没有应变,也没有弹性应变能,所以就没有力线。此时,这一区域在无裂纹条件下承担的应力、应变和能量都集中在裂尖附近了。如果外载荷继续增大,裂尖处的应力集中一旦超过某一临界应力值,微裂纹必然扩展成为宏观裂纹。
图9 裂纹附近力线的分布Fig.9 Force line distributions closed to micro-crack
图10给出了血管以5℃/min降温时,血管微裂纹尖点附近应力场分布。显然,在冻结过程中,裂尖附近形成了应力集中现象,而且越靠近裂尖点,其应力值越大。
图10 血管裂尖附近应力场分布(局部放大)(以5℃/min降温到(45℃,W=0.5 mm)Fig.10 Stress distributions near to the artery crack tip when cooled to(45℃with cooling rate 5℃/min(W=0.5 mm)
根据2.1和2.2部分的应力分析可知,在冻结过程中,血管外壁产生的应力最大。因此,在血管外壁的裂尖节点处也是整个血管最大应力集中处。图11所示为血管以5℃/min降温过程中,不同宽度微裂纹裂尖外壁节点的最大主应力分布,其中W=0表示血管中无微裂纹。从图中曲线也可以看出,一旦血管中存在了微裂纹,其裂尖点必然出现应力集中现象;而且,裂纹宽度越小,其应力值越大。该如何解释后者呢?
图11 微裂纹宽度时对血管裂尖点最大主应力的影响(降温速率为5℃/min)Fig.11 Effects of crack width on maximum stress distributions in artery crack tip with cooling rate of 5℃/min
在对血管冻结过程建立计算模型时假设为弹性模型,那么,根据线弹性断裂力学观点来看,必然存在一个叫做应力强度因子的参数 K[15]。对于普通构件,一般形状的裂纹应力强度因子属于 KI型(即张开型)。应力强度因子KI与作用在裂纹顶端处的名义应力σ及裂纹尺寸W之间存在以下的普遍关系:
式中,Y为表征含裂纹构件几何形状的一个无因次系数,针对某一研究对象而言是常数。当外应力σ增大到一定程度时,裂纹达到失稳状态,此时,即使外力不再增加,裂纹也会迅速扩展,直到断裂。这种情况下KI达到极值KIC,这就是断裂力学研究中所称的材料断裂韧性,KIC表示的是材料的一种力学性能。
因此,对于血管材料而言,也同样存在一个 KIC值[16]。那么,在相同外界载荷下,当裂纹宽度(即 W值)越小时,其裂尖处的应力也就越大。所以,可以得出以下结论:血管中存在的微裂纹对于其冻结过程中出现的断裂现象(尤其是在断裂成核扩展初期)有着非常大的影响,特别是超微小裂纹处很可能就成为宏观裂纹出现的源头。
降温速率对裂尖点处最大主应力分布的影响如图12所示(微裂纹宽W=0.05 mm)。当降温速率达到50℃/min时,冻结过程中最大应力值达到了约12 MPa,是5℃/min降温时的约30倍。而且该值已经远远大于血管的极限拉应力强度(6 MPa左右,见文献[9]),这时必然引起裂纹的扩展,直至出现宏观裂纹现象。
图12 降温速率对裂尖点最大主应力的影响(W=0.05 mm)Fig.12 Effects of cooling rates on maximum stress distributions in artery crack tip with crack width 0.05 mm
1)血管壁内、外表面的主应力表现为压应力,中心处表现为拉应力;无论是拉应力还是压应力,它们都是随着温度的降低逐渐增大,达到某一峰值后快速降低,最后趋于0。这主要受血管冻结过程热膨胀系数和弹性模量随温度变化特点的影响。
2)血管壁外表面所受的应力值最大,内表面次之,中心处的应力值最小,因此可以认为血管冻结过程中,外表面最为脆弱,血管发生裂纹断裂时应当是从外表面开始的。
3)添加较高浓度DMSO后,由于血管的热膨胀系数和弹性模量显著减小,血管冻结过程产生的热应力也显著减小。此时,温差是引起最大应力的最主要原因,采用较慢的降温速率就可以减弱热应力的影响。
4)一旦血管中存在了微裂纹,那么其冻结过程必然出现应力集中现象。而且,裂纹宽度越小,其应力值越大,这是受血管材料断裂韧性因子制约的。这个研究结果也表明,血管中存在超微小裂纹对于血管低温保存来讲是非常危险的,而且降温速率越大,这种潜在的危险越大。因此,在对血管进行低温保存前,应尽量避免样品内缺陷的存在。
[1]Hua Zezhao,Xu Hongyan,Zhou Guoyan,et al.Analyses of thermal stress and fracture during cryopreservation of blood vessel[J].Science in China(Series E),2001,44(2):158-163.
[2]Zhang Aili,Cheng Shuxia,Gao Dayong,et al.Thermal stress study of two different artery cryopreservation methods[J].Cryoletters,2005,26(2):113-120.
[3]Zhao Gang,Liu Zhifeng,Zhang Aili,et al.Theoretical analyses of thermal stress of blood vessel during cryopreservation [J].Cryoletters,2005,26(4):239-250.
[4]VenkatasubramanianRT, Simha N, TatsutaniK, etal.Measurement of freezing induced biomechanical property changes in arteries using indentation[J].Cryobiology,2007,55(3):360-361.
[5]徐红艳.主动脉血管组织低温保存与低温断裂的初步研究[D].上海:上海理工大学,2001.
[6]匡震邦,马法尚.裂纹端部场[M].西安:西安交通大学出版社,2002.
[7]Wassenaar C,Wijsmuller EG,Van Herwerden LA,et al.Cracked in cryopreserved aortic allografts and rapid thawing[J].The Annals of Thoracic Surgery,1995,60(4):165 -167.
[8]Pegg DE,Wusteman MC,Boylan S.Fractures in cryopreserved elastic arteries[J].Cryobiology,1997,34(2):183 -192.
[9]胥义,周国燕,胡桐记,等.兔胸主动脉在冻结状态下力学性能及其影响因素的实验研究[J].上海理工大学学报,2005,27(3):198-202.
[10]Xu Yi,Hua Zechao,Sun Dawen,et al.Effects of freezing rates and dimethyl sulfoxide concentrations on thermal expansion of rabbit aorta during freezing phase change as measured by thermo mechanical analysis[J].Journal of Biomechanics,2007,40(14):3201-3206.
[11]张爱丽.血管低温损伤机理的研究[D].合肥:中国科学技术大学,2003.
[12]胥义,周国燕,胡桐记,等.用DSC测定兔主动脉血管冻结相变区间的表观比热容及其影响因素[J].制冷学报,2005,(1):38-43.
[13]胥义,华泽钊,周国燕.兔主动脉冻结过程中未冻水份额的研究[J].工程热物理学报,2006,27(3):478-480.
[14]哈宽富.断裂物理基础[M].北京:科学出版社,1999.
[15]杨卫.宏微观断裂力学 [M].北京:国防工业出版社,1995.
[16]Xu Yi,Hua Zechao,Sun Dawen et al.Experimental study and analysis of mechanical properties of frozen rabbit aorta by fracture mechanics approach [J].Journal of Biomechanics,2008,41(3):649-655.