环形激光诱导空化泡的数值仿真分析*

章宇健，常磊鑫，曾乐淘，易健鹏，郭钟宁，邓 宇

（广东工业大学机械工程学院，广州 510006）

0 引言

激光诱导向前转移技术（Laser Induced Forward Transfer Technology,LIFT）是通过激光诱导空化形成脉冲压力来推动物体的分离。目前这个技术正应用于微电子元件打印[1－2]、生物组织工程[3－5]和聚合物材料打印[6－7]等领域，并且具有很大的发展潜力。在生物组织工程中，传统的激光诱导前向转移技术利用近似于圆形光斑诱导空化泡推动溶液形成射流[8－10]，可以实现快速、高效的溶液分离。但是，在单细胞转移技术中，由于悬浮液中细胞之间存在流体力学相互作用，单个脉冲激光能量之下，通常会有几个非目标转移细胞被一起传送[11]。这种现象大大降低了单细胞和溶液的转移精度。为了提高单细胞的转移精度，研究学者提出的改进方法可分为2 种：（1）提高悬浮液的黏性限制生物细胞的自由行程，Ozbolat[12]在细胞悬浮液中添加藻酸盐溶液，降低了单个细胞的自由程，但是由于黏度提高，溶液和细胞间的相互作用力也更加显著，同样存在多个细胞被一起传送的现象；（2）物理方法阻隔生物细胞间相互作用，王文兵等[13]在透明基底下添加了一层带有微方孔阵列的PDMS隔绝细胞，但是这种方法灵活性较低，只适合特定直径的细胞，且微结构制作困难。Yun Wang[14]利用拉曼显微镜结合激光诱导前向转移技术来对细胞进行定位，以提高转移精度，但由于这种方法在识别定位耗时上是单个脉冲传送速度的1 000倍，整体效率较低，且样本过厚易导致拉曼信号差，也会降低转移的准确性。以上的激光诱导转移技术均运用圆形激光光斑，因此在激光作用生物材料上时，部分直接作用在生物细胞上的能量也会损伤细胞转移之后的生物活性。

2012年，Faccio[15]利用飞秒脉冲激光器在液体中诱导产生了环形空化泡，并利用高速摄影拍摄了环形空泡的整个生长周期过程。其实验结果发现，环形空化泡在生长初期有向环四周膨胀的现象，溃灭之后会产生汇聚效应，并在环中心垂直方向产生冲击波，Faccio的实验结果如图1所示。为此，本文将环形空化应用在激光诱导前向转移细胞中，通过向外膨胀实现阻隔细胞或溶液间流体力的相互作用的效果，利用溃灭的汇聚冲击波推动目标溶液的转移。本文对环形空化进行有限元数值仿真，分析环形空化泡在液体中生长溃灭过程的速度场变化规律，结合响应曲面，揭示了工艺参数对阻隔效应和汇聚效应的影响。

图1 Faccio利用高速摄影拍摄的环形空化泡生长溃灭过程Fig.1 Growth and collapse process of ring－shape cavitation bubble captured by Faccio using high－speed photography

1 环形激光诱导空化泡数值仿真模型

1.1 模型建立

激光诱导前向转移技术是利于激光诱导空化泡产生的单次脉冲压力来推动材料转移，如图2所示。整个过程可分为如下3个阶段：（1）将均匀分布的激光光束通过一系列的透镜组合后聚焦到透明玻璃板上；（2）悬浮的诱导材料在激光的光、热作用下，被击穿后产生高温高压的等离子体，使得悬浮液汽化形成空化泡核，经过空化泡的生长溃灭形成垂直向下的射流并推动目标物质向下运动；（3）液滴在自有动量和重力作用下继续向下运动，最终沉积在接收板上实现转移的过程。

图2 激光诱导前向转移技术原理Fig.2 Schematic diagram of laser－induced forward transfer cells

该方法的最大特点：激光直接与目标转移物质相接触[16]；目标转移物质所获得的推动力直接来源于空泡的膨胀所产生的射流冲击波[3]；相邻细胞或颗粒物质之间存在着流体力学相互作用[11]。本文所提出的环形空化应用在激光诱导前向转移技术上，其环形空泡本身可以起到阻隔环内外物质相互作用的能力。为验证环形空化泡向外膨胀时对外部细胞悬浮液的推动效应和空泡溃灭后的能量中心汇聚效应，对其在自由流体边界中的生长模型进行有限元仿真分析。

为更好地探究环形空化泡在不受外界条件影响下的生长溃灭周期脉动特性，仿真模型忽略其壁面约束结构。其基本问题的物理描述如下：考虑在无限液相流体中有一个环形空泡，空泡初始时为圆环管形状。由于海藻酸钠溶液的密度随着压强和温度的变化仅发生微小的变化，其压缩性可忽略不计，认为液相为不可压缩的牛顿流体。空泡初始状态在溶液中近似呈中心几何对称分布，如图3所示。

图3 自由边界下环形空泡生长过程仿真模型Fig.3 Simulation model of ring－shape cavity growth process under free boundary

对模型进行网格划分，为了获得更好的计算精度，划分方法选择为“Tetrahedrons”，即四面体单元。由于空泡区域在整个流场中占比较小，需要对环形空泡区域进行加密划分，设置其划分的最大特征长度不超过1×10－5m，最后对网格划分质量进行检查使其符合计算要求，网格划分结果如图4所示。

图4 自由边界下环形空泡生长过程仿真模型网格划分Fig.4 Mesh generation for simulation model of ring－shape cavity growth process with free boundary

在计算环形空泡生长溃灭过程中引入的基本假设条件如下：（1）忽略空泡内混合物与流场之间的物质交换；（2）在计算中不考虑重力作用的影响；（3）空泡的生长溃灭是一个高速绝热过程，表面张力参与了泡壁的形变，但在自由边界生长溃灭的过程中，泡壁的运动非常剧烈，空泡溃灭收缩产生的高速动量变化是形成冲击波的主要原因[17]，在平衡状态下，维持空泡形状的表面张力对该过程的影响非常有限，为此，计算中不考虑表面张力的影响。采用VOF 模型，将该模型设置为气－液两相流，在流体相中，液相为不同溶液质量分数的海藻酸钠溶液，气相为水蒸气。选择计算模型为瞬态模型。由于空泡溃灭后推动流场中溶液产生的流速较大[18]且作用区域为微米级别，其雷诺数Re＝ρvL/μ＜20，黏性效应在此流场中的影响更加显著，模型流动状态设置为层流模型。

为了更准确地模拟自由边界状态下的空泡生长过程，理想的仿真模型中空泡的生长过程不应受到固体壁面回弹的压力冲击波的影响，因此在初始边界条件设置中将圆柱形的上下2 个面定义为压力出口，如图3 中Pressure－outlet 位置。圆柱形的曲面部分设置成无穷远边界条件，且整个液相的建模区域也要远远大于气相，液相体积约为气相的1 500倍。初始化整个模型的初始温度为300 K，且各区域初速度均为0。环形空泡的初始能量来源于激光击穿溶液之后的吸收和转化，定义初始温度和初始压力值来共同表示其初始能量，并将计算所得初始值加载到模型的气相空泡区域上。最后设置每个时间步的步长为1×10－8s开始计算。

1.2 激光诱导环形空化泡初始值计算

空泡的初始能量计算如下[19]:

式中：E0为海藻酸钠溶液在常温下的初始内能；ΔE为激光聚焦于溶液中时被吸收的内能。

由于常温下海藻酸钠溶液的初始内能E0相对于ΔE较小，因此E0可以忽略不计；ΔE为溶液吸收激光能量Ea增加后的比内能，且与单脉冲激光能量成正比。溶液在击穿区域吸收的激光能量E0之后，增加的内能ΔE可表示为：

根据质量守恒，激光击穿溶液的初始气体密度ρ0等于该处溶液的密度，即：

式中：m0为激光作用区域溶液的质量；V0为初始环形空泡的体积。

不同溶液质量分数的海藻酸钠溶液理化性质参数[20－21]如表1所示。

表1 海藻酸钠溶液参数Table.1 Parameters of sodium alginate solution

假设环形激光光斑击穿溶液生成的初始环形空泡内外径的大小相互一致，环形空化泡初始体积为：

式中：R1为内环空泡半径；R2为环形空泡外环半径。

初始压力和初始温度可以通过理想气体状态方程得到[22]：

式中：Rg为理想气体常数，Rg＝461.5 J/(kg·K－1)；i为理想气体自由数，i＝6。

1.3 模拟仿真的控制方程

为了模拟仿真空泡的生长过程，采用有限元分析法对空泡进行数值仿真，采用气－液两相流模型对空泡进行仿真。具体控制方程如下所示。

质量连续性方程：

式中：ρ为流体密度；为流体速度。

动量方程：

式中：p 为压力向量；τn为黏性应力张量；为重力；为表面张力。

能量方程：

式中：Eh为流体的焓；T为温度；keff为流场的传导系数。

2 数值模拟结果与分析

在激光能量40 μJ、环形空泡初始内外直径比值为0.8（80 μm/100 μm）、溶液质量分数为2%的条件下，由1.2节可求出环形空泡初始压强为41.51 MPa，初始温度为312.64 K。气液两相仿真结果的中心垂直截面如图5所示。

图5 环形空化泡生长周期过程的气液两相仿真结果Fig.5 Gas－liquid two phase simulation results of ring－shape cavitation bubble growth cycle

从气液两相图可以看出，图中a～e 为环形空泡经历膨胀、收缩至溃灭的第1 个生长周期，e～i 为第2 个生长周期，i～l为第3个生长周期。可知，在此激光能量大小之下，环形空化泡的周期长度约为12 μs 。在第1 个周期即将结束时，1个环形空化泡分裂成了2个形状大小一致的环形空化泡，在第2 个周期即将结束时，2 个环形空化泡各自分裂成了4 个两两几何对称的环形空化泡。另外，随着时间的推移，环形空化泡整体呈现出逐渐向中间汇聚收缩的趋势。综上两相流仿真结果与Faccio 所利用的飞秒脉冲激光诱导环形空化泡生长周期结果对比基本一致。

由于溶液在流场中被空泡所推动产生的速度更能体现其运动状态。为了进一步探究其在自由液体边界下水平方向的阻隔效应和垂直方向的汇聚效应，对其流场中的速度分布云图进行的分析。如图6所示。

图6 环形空化泡生长周期过程的速度分布Fig.6 Velocity distribution nephogram of ring－shape cavitation bubble growth cycle

在第1 个环形空泡周期中（a～e），流场中运动速度最大的流域，从膨胀、收缩至溃灭，基本均分布于环形空化泡水平方向上，环形中心垂直方向上液体被推动产生的速度大小远远小于水平方向。在第2 个环形空泡周期中（e～i），由于空化泡的汇聚收缩效应，第2个周期中的环形空泡内外径均小于第1个周期中的环形空泡。此时，在环形空泡膨胀阶段中，环形中心垂直方向上液体被推动产生的速度大小大于水平方向，且出现能量汇集效应，如图e～f。随着时间推移，环形空泡继续进行着周期振荡。在第3 个生长周期中（i～l），环形空泡的内外径变得更加小，由图5j 可得大致为初始环形空泡的一半，即直径收缩至约50 μm，且速度在垂直方向上的汇集效应更加明显。

为了利用空泡第1 周期中膨胀阶段对环外部溶液的阻隔推动作用，和第1周期溃灭后的能量中心汇集效应来推动溶液转移，对上述仿真结果中各个时间段水平和垂直2个方向上的运动速率进行了提取，如图7所示。

图7 流场中不同方向上速率随时间的变化规律Fig.7 Variation of velocity with time in different directions in flow field

结合图5和图6可知，在第1个周期时间内（1～12 μs），由于环形空泡内外径比值较大，空泡膨胀过程中的能量没有汇集在垂直方向上，而是作用在环水平位置外部区域，且在膨胀初期达到最大值Vx＝8.9 m/s。在第2、第3 个周期时间内，在13 μs 和25 μs 时，环形空泡处于膨胀初期，垂直方向上的能量得到了汇集，流场内的溶液运动速率得到指数型增长，其速度大小Vy分别为9.3 m/s和12.7 m/s。在激光诱导前向转移技术中，环中心目标转移溶液在空泡第1个周期内的膨胀至溃灭结束就已经实现了传送，因此，合理地设计环形激光光斑内外径大小、激光能量初始值和悬浮液质量分数，对溶液转移起到了关键作用。

3 各因素对水平和垂直方向最大速率的影响

3.1 激光能量

固定溶液质量分数为2%，环形光斑内外径比值为0.8（80 μm/100 μm），探讨了不同激光能量对环形空泡膨胀初期水平方向最大速率和汇聚收缩后垂直方向最大速率的影响，数值仿真模拟方法同第2节，结果如图8所示。

图8 最大速率随激光能量变化趋势Fig.8 Variation trend of maximum rate with laser energy

3.2 环形空泡内外径比值

固定激光能量为50 μJ，溶液质量分数为2%，固定初始外环空泡直径为100 μm，改变内直径大小的条件下，探讨了不同环形激光光斑内外径比值对环形空泡膨胀初期水平方向最大速率和汇聚收缩后垂直方向最大速率的影响，数值仿真模拟方法同第2节，结果如图9所示。

图9 最大速率随环内外径比值变化趋势Fig.9 Variation trend of maximum velocity with ratio of inner diameter to outer diameter of ring

3.3 溶液质量分数

固定激光能量为50 μJ，环形光斑内外径比值为0.8（80 μm/100 μm），探讨了不同溶液质量分数对环形空泡膨胀初期水平方向最大速率和汇聚收缩后垂直方向最大速率的影响，数值仿真模拟方法同第2节，结果如图10所示。

图10 最大速率随溶液质量分数变化趋势Fig.10 The trend of maximum rate with the mass fraction of solution

4 响应曲面分析

在单因素的基础上，采用响应曲面分析法中的Box－Behnken Design建立数学模型，以空泡仿真模型中水平方向最大速率vx和垂直方向最大速率vy作为评价阻隔效应和汇聚效应的指标。选择激光能量（A）、环形光斑内外径比值（B）、溶液质量分数（C）设计三因素三水平响应曲面实验，如表2所示。由1.2节可计算得出激光能量为60 μJ 时，环形空泡初始温度值为332.83 K，细胞内高温面积比可达到2.6%[23]，而要使转移之后的细胞保持良好的生物活性，环形空化泡的初始温度不得高于338 K[24－25]。实验结果如表3所示。

表2 响应曲面模型中的变量和参数值Table.2 Variable and parameter values in response surface model

对表3结果中的Vx和Vy作方差分析，计算出模型的显著性水平P 值均小于0.000 1，得到模型显著性程度较高。失拟合误差值P＝0.349 1和P＝0.134 5均大于0.05，表明失拟没有显著性，即该拟合模型均与实验相比有较高的拟合程度。计算影响水平方向最大速率Vx的方差来源激光能量A、环形光斑内外径比值B、溶液质量分数C 对应的F 值分别为142.39、8 805.19、1 445.75。计算影响垂直方向最大速率vy的方差来源激光能量A、环形光斑内外径比值B、溶液质量分数C对应的F 值分别为299.47、25 589.00、2 189.35。比较F 值可得，3 个不同因素对水平方向上的最大速率vx的影响强弱顺序依次为：环形光斑内外径比值B、溶液质量分数C、激光能量A。3个不同因素对垂直方向上的最大速率Vy的影响强弱顺序依次为：环形光斑内外径比值B、溶液质量分数C、激光能量A。

表3 Box－Behnken Design实验设计及结果Table.3 Experimental design and results of box Behnken design

图11所示为不同因素对水平和垂直两个方向最大速率综合影响的三维散点图。可以得出，在环形光斑内外径比值为0.7(70 μm/100 μm)、激光能量60 μJ、溶液质量分数为1%时达到最大，分别为vx＝14.28 m/s、vy＝11.82 m/s。其中，通过改变环形激光光斑大小调节环形空泡内外径比值所带来的影响效果最好。但是，在激光诱导前向转移技术运用中，环形激光光斑内外径大小往往是由目标转移单位的尺寸大小所决定，比如生物细胞直径或者微细颗粒直径。因此，溶液质量分数和激光能量大小的控制也尤为重要。

图11 不同因素对最大速率的影响变化三维散点图Fig.11 Three dimensional scatter plot of the influence of different factors on the maximum velocity

5 结束语

本文运用有限元分析法对环形空泡生长溃灭过程的仿真模拟计算。在激光能量为40 μJ，环形空泡初始内外直径比值为0.8（80 μm/100 μm），溶液质量分数为2%的条件下，得出空化泡水平方向上的最大速率出现在初始膨胀时期，达到了8.84 m/s；垂直方向上的最大速率出现在空泡溃灭后期，达到了9.64 m/s。水平和垂直方向上最大速率的先后作用顺序，对应了激光诱导单细胞转移中的阻隔推动作用和能量汇聚作用，进一步验证了环形激光诱导前向转移单细胞技术的可行性。通过响应面法得到不同激光能量、环形光斑内外径比值、溶液质量分数对空泡水平方向和垂直方向最大速率的影响强弱顺序，有助于环形激光诱导前向转移生物靶材实验的设计和优化。