血管分支血液对流换热分析

解海卫 张 艳 诸 凯

（天津商业大学机械工程学院，天津市制冷技术重点实验室，天津 300134）

引言

模拟生物温度场用于辅助分析已经成为生物医学工程研究中的常用方法。血液与生物组织有复杂的传热耦合关系，准确地分析血液对流换热是模拟生物组织传热的必要前提。

许多研究者为分析血液对流换热做了大量的实验研究和理论计算。Charm进行了实验测量，并分析了血液在不锈钢毛细管（直径0.6mm）内的热交换［1］。Shah通过测量一定加热量P下对应的温升T，由已知经典的努赛尔数（Nu数）和P/T的拟合公式，确定了单管内血液的Nu数［2］。Victor对单个血管的血液流动和换热进行了模拟，计算得到了血液Nu数沿管长的变化，其中壁面条件是均匀壁温和均匀热流两种，入口为完全发展流动［3］，入口条件是均匀流速、温度两种［4］。Luisa假设血液为牛顿流体，根据动量和能量守恒，得出了单根血管局部加热时的对流换热系数公式［5］。Tungjitkusolmun假设血管壁温为常数，分析计算了热疗时大血管对生物组织温度场的影响［6］。

但是，这些研究均是以单个血管为对象，着重用不同方法、确定不同换热条件下的血液Nu数（或h），而考虑在循环系统实际分支结构的基础上，分析血液换热的文献鲜见报道。Anil研究表明，血管的分支结构对血液速度场分布有明显影响，因此分支结构对血液热交换也一定存在影响［7］。根据传热学理论，流速和管径也是影响管内对流换热的主要因素。笔者通过建立血管树简化模型，模拟计算了分支结构内血液的三维速度场、温度场和血液Nu数；从血管树的实际结构出发，分析了血管分支结构参数（半径比、分支角）、血液流速、血管半径对血液与组织之间热交换的影响。

1 血管分支结构模型

研究表明，血液和组织的热交换主要发生在动脉和小动脉处（R＞0.1mm），而这些动脉一般是以上级血管的分支形式存在，且末端为二元分支结构［8］，因此计算血液换热应在分支结构中分析。建立循环系统分支结构模型，如图1所示，模型由一个主干血管和两个分支血管组成。根据血管树生成理论，由优化原则建立血管树，其主干血管和分支血管的轴线在同一平面内［9］，即该血管树关于此平面对称，因此截取分支结构的一半作为最终建立的血管树物理模型。

选取主干血管半径为0.7mm，属于热重要血管的范围［10］，长度为20mm。分支Ⅰ与分支Ⅱ的长度为40mm；分支Ⅰ与主干血管之间的夹角（即分支角α（定义见图1），分别为90°、120°、150°、180°，分支血管Ⅰ的半径与主干血管的半径之比（即半径比）分别为0.75、0.8、0.85、0.9；分支血管Ⅱ的半径由Murray定律确定。

图1 血管分支结构Fig.1 Schematic of the furcated configuration

2 数学模型的建立

2.1 模型假设［3，11］

忽略血管的膨胀性，将流动视为稳态、层流；血液为常物性，密度、比热不随温度变化；忽略血液的代谢热；忽略近壁面处血浆层的影响；血液为不可压非牛顿流体，其非牛顿性用幂律模型描述，有

式中，γ为切变率，1/s；τ为切应力，Pa；k为常数，取0.014 67（Pa·sn）；n为非牛顿指数，取0.775 5［12］。

2.2 控制方程

计算血液与组织的对流换热，需耦合求解血液的流场和温度场，即求解以下控制方程为

质量方程（连续性方程）

(1)看皮下脂肪层的厚度。“瘦肉精”猪皮下脂肪层明显变薄，通常不足1cm;而正常猪肉在皮层和瘦肉之间会有一层脂肪，肥膘约为1～2cm厚。因此，遇到脂肪层太薄、太松软的猪肉要格外小心。

动量方程

能量方程

式中，U为流体速度，m/s；S为广义源项，Pa/m；Ψ为黏性耗散函数，1/s2；ρ为密度，kg/m3；P为压力，Pa；μ为动力黏度，Pa·s；T为温度，K；c为定压比热，J/（kg·K）；λ为导热系数，W/（m·K）。

血液各参数取值见表1。

表1 血液热物性参数取值［13］Tab.1 Properties of blood

2.3 边界条件

1）主干血管入口条件：假设主干血管入口边界条件为均匀速度和均匀温度，Tin=310K。

2）分支血管出口条件：Pout=0（相对压力），qout=0。

3）血管轴线对称面条件：对称面绝热，即q=0；流动边界条件为在垂直于对称面的方向上流速为0。

4）血管壁面条件：壁面无滑移，即Uw=0；热边界条件为均匀壁温条件，Tin=309K。

2.4 网格划分和数值模拟

图2 有限元网格。（a）直管部分网格；（b）分支处网格；（c）血管截面网格Fig.2 Finite element mesh of the vascular tree.（a）mesh on straight part of vessel；（b）mesh on crosssection；（c）mesh at the bifurcation point

利用ANSYS有限元软件中的FLOTRAN模块，分析血管分支结构内的速度场和温度场，得到血管壁面处血液与壁面的换热量，进而根据努塞尔数定义，计算得出分支血管Ⅰ管长x处（如图1所示）血管壁面任意位置ф处的血液努塞尔数Nux，Ф数，以及x处的截面平均努赛尔数Nux。

3 分支结构内的速度场和温度场

图3表示血管连接处速度和温度分布，其中（a）和（b）分别为模拟得到的主干血管血液流速为80mm/s、主干血管半径为0.7mm、分支角为120°、分支Ⅰ半径比为0.8时，血管连接区域的速度和温度分布等值线图。由图3（a）可以看出，血液从主干血管流入分支血管时，速度减小，流向改变，截面高流速区向一侧偏移，速度场由主干血管的轴对称分布变为分支血管入口处的非轴对称分布，然后逐渐恢复为轴对称分布状态；由图3（b）可以看出，分支血管入口处的温度分布也呈现与速度场类似的变化，原主干血管轴心处的低温区，随流动偏移到分支血管位置1一侧（见图2）。此外，两个分支血管的速度场和温度场的偏心程度不同，这是由于两个分支血管的分支角度和半径比不同而引起的。

从图3还可以看出，除主干血管出口区域速度场有一个较小的扰动外，两个分支血管的存在对主干血管的速度场和温度场基本没有影响，主干血管的速度和温度等值线基本呈平行状态进入血管分叉区域，分支对主干血管内流动换热的影响可以忽略，即在循环系统中，血管主要受上级来流主干血管的影响，而其下级分支血管的影响可以忽略。由此可见，在模型中，分支血管Ⅰ、Ⅱ内血液换热主要受主干血管的影响，而分支Ⅰ、Ⅱ的下级分支对其影响可忽略，所以血管树模型无需建立分支Ⅰ、Ⅱ的下级分支。

图3 血管连接处速度和温度分布。（a）速度场分布等值线（单位：m/s）；（b）温度场分布等值线（单位：℃）Fig.3 Contour plot of the numerical solutions at the bifurcation position of the vessel system.（a）contour plot of velocity distribution（Unit：m/s）；（b）contour plot of temperature distribution（Unit：℃）

4 分支血管内Nu数变化规律与分析

由上述可知，分支结构的存在会影响分支血管内的速度场和温度场，使其在入口处呈现非轴对称分布，相应地也会影响分支血管内的Nu数，使之呈现与单个直管内Nu数不同的分布规律。因此，笔者以分支血管Ⅰ为研究对象，计算了血管分支结构中分支血管Ⅰ内任意截面x（即管长x处）的平均努赛尔数Nux和径向位置1、2、3处（如图2所示）的局部努赛尔数（Nux，1、Nux，2、Nux，3），并分析了分支结构参数（半径比、分支角）、分支血液流速、分支血管半径的变化对截面平均努赛尔数Nux的影响。

4.1 分支血管内努赛尔数分布规律

图4为分支血管Ⅰ半径为0.56mm、平均流速为30mm/s，分支角为120°，半径比为0.8时，其内截面平均Nux与相同管径、流速下单个直管内截面平均Nux的对比曲线。从该图可以看出，在血管入口段，分支血管的Nux明显小于单个直管的Nux，且收敛速度大于直管。这是因为血液流经主干血管后，由入口处的速度、温度均匀的流体变为出口处速度、温度完全发展的流体，虽然流入分支后，由于方向的改变造成速度场和温度场的变化，但相对于流入单个直管的血液，更接近于热完全发展阶段，因此Nux较小，且收敛速度快。在血管出口处，分支血管和单个直管的Nux都达到收敛值，约3.75。

在上述条件下，分支Ⅰ内3个典型圆周位置处的壁面局部Nux，Ф（Ф=1，2，3）随x的变化曲线。可以看出：在分支血管入口段，不同圆周位置之间Nux，Ф差异较大，这是由于血液流入分支血管时，流向改变引起温度出现非轴对称分布，不同圆周位置处的温度梯度不同，使Nux，Ф出现差异；而在出口截面段，3个圆周位置处的Nux，Ф基本相同，约为3.75，这与分支血管内温度场由非轴对称分布恢复为轴对称分布的变化过程一致。

图4 分支血管与单个直管截面平均Nux数的对比曲线Fig.4 The comparison between the mean Nuxat daughter vessel I and at a straight vessel

图5 不同圆周位置处的局部Nux，Ф数随管长变化的曲线Fig.5 The local Nux，Ф at three different positions 1，2 and 3 in daughter vessel I versusx

4.2 分支结构参数对截面平均努赛尔数的影响

当分支角分别是90°、120°、150°、180°时，其他参数与4.1节条件相同，分支Ⅰ内Nux随x变化的曲线如图6所示；当血管半径比分别是0.75、0.8、0.85、0.9时，分支Ⅰ内Nux随x变化的曲线如图7所示。由这两图可以看出：分支角、半径比越小，入口段截面平均Nux越大；分支角、半径比变化对分支血管截面平均Nux的收敛速度（热入口段长度）没有影响，收敛值都约为3.75。这主要是由于分支角、半径比越小，分支血管入口处速度场和温度场的偏心程度越严重，相应Nux也就会越大。

4.3 分支血管管径、流速对Nux的影响

当分支血管半径为0.48、0.56、0.64、0.8mm时，或分支血管血液流速为20、30、40、50mm/s时，其他参数与4.1节条件相同，分支Ⅰ内Nux随x变化的曲线如图8和图9所示。可以看出：分支血管半径、血液流速越大，入口段截面平均Nux越大；血管半径和血液流速的变化对分支血管截面平均Nux的收敛速度也有影响，血管半径或血液流速越小，Nux的收敛速度越快，收敛值都约为3.75。这主要是由于血管半径或流速越小，血液流动中管壁对流体的黏滞力与血流的惯性力之比越大，血流从入口处的偏心流动恢复为轴对称流动的速度越快，相应热入口段的长度也越短，Nux的收敛速度越快。

图6 不同分支角度对应的Nux随x变化的曲线Fig.6 The mean Nuxwith different furcating angels

图7 不同半径比对应的Nux随x变化的曲线Fig.7 The mean Nuxwith different bifurcation ratios

图8 不同分支血管半径对应的Nux随x变化的曲线Fig.8 The mean Nuxwith different vessel radii

图9 不同分支血管流速对应的Nux随x变化的曲线Fig.9 The mean Nuxwith different blood velocities

5 结论

以反映循环系统实际结构的血管分支物理模型为平台，用有限元方法数值模拟了血液的三维速度场、温度场，得到了分支血管中典型位置处Nux，Ф和截面平均Nux，分析了血管半径、血液流速和血管树结构参数（半径比、分支角）对截面平均Nux的影响，结果表明：

1）当血液由主干血管流入分支血管时，血管截面速度场和温度场由主干血管出口的轴对称分布变为分支血管入口处的非轴对称分布，然后随着在分支血管内的流动又逐渐恢复为轴对称状态，而分支血管对主干血管内流动和换热的影响较小。

2）分支血管入口段的Nux明显小于单个直管的Nux，且收敛速度大于直管，收敛值均约为3.75。分支血管入口处血液的Nux，Ф为明显非轴对称分布，且随血液沿血管继续流动，其不对称程度逐渐变小。

3）分支角越小、半径比越小、分支血管半径越大、分支血管流速越大时，入口段截面平均Nux越大。各参数对入口段截面平均Nux的影响程度为：血液流速＞半径比，血管半径＞分支角。分支血管半径、流速越小，截面平均Nux的收敛速度越快。

［1］Charm S，Paltiel B，Kurland GS.Heat transfer coefficients in blood flow［J］.Biorheology，1968，5（1）：133-145.

［2］Shah J，Dos SI，Haemmerich D，et al.Instrument to measure the heat convection coefficient on the endothelial surface of arteries and veins［J］.Medical＆ Biological Engineering and Computing，2005，43（4）：522-527.

［3］Victor SA，Shah VL.Steady state heat transfer to blood flowing in the entrance region of a tube［J］.International Journal of Heat and Mass Transfer，1976，19（6）：777-783.

［4］Victor SA，Shah VL.Heat transfer to blood flowing in a tube［J］.Biorheology，1975，123（3）：361-368.

［5］Luisa C，Dos SI，Haemmerich D.Theoretical analysis of the heat convection coefficient in large vessels and the significance for thermal ablative therapies［J］.Physics in Medicine and Biology，2003，48（12）：4125-4134.

［6］Tungjitkusolmun S，Staelin T，Haemmerich D，et al.Threedimensional finite-element analysis for radio-frequency hepatic tumor ablation［J］.IEEE Transactions on Biomedical Engineering，2002，49（1）：3-9.

［7］Anil K，Varshney CL，Sharma GC.Computational technique for flow in blood vessels with porous effects［J］.Applied Mathematics and Mechanics，2005，26（1）：63-72.

［8］Marek K，Yan R，Johanne B，et al.Physiologically based modeling of 3-D vascular networks and CT scan angiography［J］.IEEE Transactions on Medical Imaging，2003，22（2）：248-257.

［9］张艳，张于峰，诸凯，等.三维动脉血管树的模拟与分析［J］.天津大学学报，2007，40（9）：1071-1076.

［10］Chato JC.Heat transfer to blood vessels［J］.Journal of Biomechanical Engineering，1980，102（5）：110-118.

［11］Barozzi GS，Dumas A.Convective heat transfer coefficients in the circulation［J］.Journal Biomechanical Engineering，1991，113（3）：308-313.

［12］Walburn FJ，Schneck DJ.A constitutive equation for whole human blood［J］.Biorheology，1976，13（1）：201-214.

［13］陈琦，白景峰，陈亚珠.伴行血管在高强度聚焦超声下对温度场的影响［J］.上海交通大学学报，2004，38（1）：130-134.