有限元法仿真LIPUS照射骨细胞的局部声场分布

蒋文杰，孙舒心，徐 峥，他得安

(1. 复旦大学电子工程系，上海 200433；2. 同济大学声学研究所，上海 200092)

0 引 言

有文献表明，超声能够对骨骼肌组织产生各种生物学效应。特别是低强度脉冲超声(LIPUS)被认为是促进骨愈合的一种有效手段[1-4]，是用于治疗骨质疏松症的潜在手段[5]，且具有无创伤、无电离辐射、价格低廉等优点[6]。LIPUS治疗骨病的相关研究表明，LIPUS可以加速骨形成[7]，减少骨量丟失[8-9]。此外，LIPUS作为一种物理刺激，可以直接影响骨原细胞，并导致矿化结节形成[10]。

利用LIPUS治疗骨骼疾病时，其物理机制是复杂的。常见的说法有超声的声致穿孔[11]、振荡及热效应等，而声波在传播媒介中的振荡及其引起的温升是两个重要的物理效应。其中，声波在传播媒介中的振荡是影响治疗效果的主要因素，它对骨细胞(如前成骨细胞)具有直接的影响[12]。此外，温升也是影响治疗效果的因素之一，它是由超声作用下的能量沉积引起的[13]。

作为一种低强度的治疗超声，LIPUS可以有效降低对组织或细胞产生的热效应。这是因为：第一，LIPUS不是连续波，它是由脉冲重复的频率为1.0～1.5 MHz的正弦波构成；第二，LIPUS的声强较低，空间平均时间平均声强一般在100 mW/cm2以下[14]。

在LIPUS治疗骨细胞的实验中，通常将圆柱型细胞培养皿放置在非聚焦平面圆形换能器的上方，如图 1(a)所示，两者通过耦合剂连接。或直接将换能器表面浸没在培养液中，如图1(b)所示。为了研究上述方法的治疗效果及相应的物理机制，需要明确定义超声参数及换能器照射位置等实验条件。

图1 利用LIPUS治疗骨细胞的两种常见实验装置Fig.1 Two common experimental setups in LIPUS treatment of bone cells

目前，已有相关文献报道了LIPUS照射细胞培养皿的声场分布结果，研究采用仿真与实验两种方式对其进行分析及验证，得到培养皿内声压与声强分布的定量结果[15]。其中，仿真采用超声换能器及超声波传播仿真软件PZFlex，同时考虑超声的非线性传播及超声在固体材料中的纵波与横波传播。实验则采用水听器测量法，测量培养皿内部的声压分布。结果显示，实验结果与仿真数据具有一致性，验证了超声在培养皿内的声场分布。此外，由于超声的参数及照射位置会直接影响细胞治疗的效果，因此也有文献对相关问题展开了研究[16-17]。然而，当前很少有文献对超声在微观层面内的传播进行研究，尚不清楚LIPUS照射骨细胞的声场分布在局部微观区域内具有怎样的特性和规律。

为了从物理机制的角度探讨并理解超声与骨细胞间的相互作用，本文在不同位置用不同频率和激励声压的超声对培养液中的骨细胞进行照射，并对骨细胞及其周围局部区域内的声场分布进行了仿真。通过对仿真结果的分析与讨论，得到了微观角度下的声压、声强分布与超声参数、照射位置之间的关系，为实现有效的超声治疗提供重要的理论参考。

1 有限元仿真

在LIPUS治疗骨细胞的实验中，声场的分布对治疗效果起到了相当重要的作用，而不同的实验设置又会直接影响到声场分布的情况。因此，本文基于有限元原理，在图1所示实验装置的基础上，选取局部微观区域在细胞尺度上进行建模仿真，并通过设置不同的超声参数与换能器照射位置，对细胞及其周围培养液中的局部声场进行仿真计算。

1.1 仿真模型及其物理参数

本文选取培养皿内包含骨细胞在内的局部区域作为研究对象，并在该区域内建立骨细胞模型。考虑到圆柱型培养液区域的对称性，选取该圆柱形区域轴截面的一半进行仿真计算，如图2所示。求得该区域的声场分布后，沿对称轴r=0作翻转，就可以得到整个轴截面的声场分布。

图2 骨细胞-培养液模型Fig.2 Bone cells-DMEM model

该模型中，小椭圆是细胞质，大小椭圆之间的部分是细胞膜，小圆是细胞核。小圆的半径为1 µm；小椭圆的半长轴为 4.5 µm，半短轴为 4 µm；大椭圆的半长轴为 4.95 µm，半短轴为 4.4 µm[18]。本仿真中，主要研究细胞尺度下局部区域内的声场分布，旨在探究超声与骨细胞间的相互作用。因此图2中，将矩形的右侧分层设为完美匹配层，表示传播到该区域的声波将被全部吸收，而不会发生反射，这样可以有效地避免入射超声波与反射超声波的叠加效应。

1.2 软件平台及仿真参数设置

本仿真采用的COMSOL Multiphysics 5.0，是一款基于有限元方法的多物理场耦合分析软件，可用于声学、力学、光学等领域的建模仿真研究。本文使用了它的压力声学-频域模块，并对施加超声激励后的稳态声场进行仿真分析。

本仿真的求解域方程为：

在式(1)中，

在式(3)中，

在式(5)中，

上述各式中，α为声衰减系数；ρ为材料密度；c为声速。这些参数的取值均取决于材料属性。仿真中各项材料的物理参数如表1所示[15,19]。

表1 仿真模型的材料参数Table 1 Material parameters of the simulation model

1.3 网格划分

在声场仿真中，一般将网格尺寸的最大值设置在c/f/20到c/f/6之间。其中c为声速，f为超声的频率。为了得到稳定的仿真结果，本文将网格尺寸的最大值设为c/ f /20(µm)，由此得到上述模型的网格划分。

1.4 仿真思路及超声参数

本文的仿真采用非聚焦平面圆形换能器在上、下侧照射含有骨细胞的培养液区域(Dulbecco's Modified Eagle Medium，DMEM)温度为20℃。超声在培养液中传播，并对其中的骨细胞产生声学作用。本文中所有仿真均视为在室温中进行，因而把培养液的密度、声速及声吸收系数等物理参数均按照20℃的标准来设定，且这些参数值可用于不同治疗装置下的声场仿真。

仿真中，首先在矩形区域的底边施加超声，然后在上、下两边同时施加超声。超声参数如下，频率为1.0～1.5 MHz、激励声压为104～5×104Pa。最后对不同实验条件下的声场分布进行比较分析。

2 仿真结果与讨论

通过改变超声参数及辐照位置，本文得到了不同实验条件下的声场分布情况。本节将讨论各超声参数及照射方式对声压及声强分布的影响。所有仿真在频域中进行，所得声强为时间平均声强(ITA)。为了更直观地进行比较，仿真结果均采用二维平面图或一维数值曲线来展示。下面所示的各图为r-z平面内或中心轴线r=0及径向轴线上的声场分布。

2.1 仅在底边照射超声

图 3(a)、3(b)分别为仅在模型底边进行超声照射时，r-z平面上的声压与声强分布。这里选用的超声参数为：频率为1.0 MHz，激励声压为1×104Pa。由于超声是在模型的底部进行照射的，因此超声将直接通过细胞培养液进行传播。

由图3可知，当仅在底边照射超声时，随着传播距离的增加，超声会逐渐发生衰减，直至消失。此时，超声的声压与声强由近及远逐渐减小，且当超声传播至z=15 µm以后，声压与声强都将趋近于零。可以发现，上面远离声源的两个骨细胞内的声压与声强很小，而下面靠近声源的两个骨细胞内的声场数值要大一些。采用该种照射方式时，模型区域内的声场分布不均匀。图3(b)中，在下面两个骨细胞内，细胞质与细胞核区域内的声强分布发生了明显变形，这也说明骨细胞会在一定程度上影响声场的分布。

此外，在该模型区域内，超声声压与声强的传播特性也符合超声传播理论，这也说明了当超声在微观区域内传播时，其传播特性与宏观区域下超声的传播理论是一致的。

仿真结果表明，当采用单边照射的方式时，声场的分布并不均匀。在后面的仿真中尝试在模型区域的上下两边同时照射超声，并考察声场分布是否会更加均匀。

图3 底边照射超声时的声压与声强分布(1.0 MHz，104 Pa)Fig.3 Acoustic pressure and intensity distributions when ultrasound irradiates at the bottom side(1.0 MHz, 104 Pa)

2.2 上下两边同时照射超声

图 4(a)、4(b)所示分别为在模型上下两边同时照射超声时，r-z平面上的声压与声强分布。同样地，超声参数为：频率1.0 MHz，激励声压1×104Pa。超声通过细胞培养液进行传播。与 2.1节中不同的是，这里存在两个入射超声波，所以在模型的中间区域会产生波的叠加效应，这会使声场分布更为复杂。下面将对该组仿真结果进行详细讨论。

图4 上下边同时照射超声时的声压与声强分布(1.0 MHz，104 Pa)Fig.4 Acoustic pressure and intensity distributions when ultrasound irradiates at the top and bottom sides(1.0 MHz, 104 Pa)

由图4可知，当在上、下两边同时照射超声时，由于模型整体具有对称性，因此声压与声强的分布也是上下对称的。在两边入射超声波的叠加效应下，声场的分布更为均匀，进而使得骨细胞受到的声辐射力也更加均匀。由于本仿真中上下两个波源间距仅为30 µm，远小于超声波长，因此不会产生驻波。

通过比较图3(b)与图4(b)可知，相比于单边照射，双边照射使得整体的声强更大。这样，在使用相同激励声压的同时可以获得更大的输出声强。

值得注意的是，由于两边照射的超声发生了衰减，从而使得中间位置处的声压与声强要明显低于其他区域。本仿真中，细胞尺寸在微米级，超声波长为毫米级，超声波长要远大于细胞尺寸。另外在模型中，细胞质、细胞核和细胞膜的密度与声速均比培养液大。由于声阻抗Z=ρc，因此细胞内的声阻抗要大于外液体。综合细胞尺寸及细胞声阻抗来看，超声在细胞区域会发生反射、透射，但不会产生散射。仔细观察图4(b)可以发现，声强分布受到骨细胞的影响而产生变形，尤其是在细胞核区域以及细胞之间的区域。

已有文献表明，当采用声强低于30 mW/cm2(空间时间平均)的 LIPUS治疗骨细胞时，可以增加骨细胞的丝裂原活化蛋白激酶 (Mitogen-Activated Protein Kinase，MAPK)磷酸化程度、碱性磷酸酶(Alkaline Phosphatase，ALP)活性及矿化程度，因而低于30 mW/cm2的声强在实际治疗中具有重要的参考意义[20]。当仿真中采用上述超声参数及照射方式时，骨细胞区域内的声强也符合LIPUS治疗骨细胞的要求。

2.3 无骨细胞时的声强分布

目前，已有相关文献对超声治疗骨细胞时培养皿内的声场分布进行了研究，但所得到的声场分布也都是在宏观区域内的，并没有从细胞尺度来考察声场与细胞间的相互作用。

为了从微观的角度来观察骨细胞对声场造成的影响，本节将对培养液中不含骨细胞时的声强分布进行仿真计算，并与之前培养液中含有骨细胞时的情况进行比较与分析。仿真结果如图5所示，下面会结合图4(b)一起进行详细讨论。

在培养液中没有骨细胞的情况下，且同时在模型区域的上、下两边用频率 1.0 MHz、激励声压1×104Pa的超声进行照射时，r-z平面内声强分布的仿真结果如图 5(a)所示。通过对比图 4(b)与图 5(a)可以发现，在超声照射条件与材料的物理参数完全相同的情况下，有骨细胞和无骨细胞时的声场分布存在较为明显的差异。一方面，当培养液中含有骨细胞时，声强数值要比没有骨细胞时小一些；另一方面，骨细胞会在一定程度上影响到原有的声强分布，使得声强分布产生变形，尤其是在细胞核及细胞之间的区域。

图 5(b)、5(c)分别显示了在培养液中有骨细胞及没有骨细胞时，r-z平面上位于中心轴线r=0 µm及径向轴线z=8 µm 处的声强分布。由于声场分布具有对称性，因此这里只选取了中心轴线上z=0 µm到z=15 µm、径向轴线上r=0 µm 到r=12 µm 的数据进行作图。

分别对比图 5(b)、5(c)中的两条曲线可知，骨细胞的存在会影响原有声场的分布，并吸收一部分的超声能量，从而使得培养液内的声强有所减小，这也从能量的角度解释了 LIPUS治疗骨细胞的作用机制。

2.4 激励声压对中心轴线上声场分布的影响

在 LIPUS治疗骨细胞的实验中，超声的物理参数及照射位置对治疗效果具有至关重要的作用。2.3节已经对超声的照射位置进行了研究，本节着重讨论激励声压对于声场分布的影响。

图5 有/无骨细胞时的声强分布(上下边照射，1.0 MHz，104 Pa)Fig.5 Acoustic intensity distributions with or without bone cells (ultrasound irradiates at the top and bottom sides, 1.0 MHz, 104 Pa)

本节中的所有仿真均采用在模型上下两边同时照射超声的方式。由于本文研究的是LIPUS照射骨细胞时局部区域内的声场分布，所需要的声强数值较小，这里超声的激励声压取值范围为1×104～5×104Pa。

为了定量研究激励声压对r-z平面内声场分布的影响，选取中心轴线上的声压及声强分布作为研究对象，以便更直观地显示不同激励声压下的差异。

同样地，由于中心轴线上的声场分布具有上下对称性，这里也只选取了z=0 µm 到z=15 µm 区域内的数值进行作图。仿真结果如图6所示，下面将根据各项仿真结果进行分析。

图6 不同激励声压时中心轴线处的声场分布(上下边照射，1.0 MHz)Fig.6 Sound field distribution at the central axis for different acoustic pressures(ultrasound irradiates at the top and bottom sides，1.0 MHz)

当激励声压取不同数值时，中心轴线上位于z=0 µm到z=15 µm区域内的声压、声强分布如图6所示。由图 6(a)、6(b)可知，随着超声传播距离的增加，声压与声强都会不断减小。在轴线中间z=15 µm 处，声压与声强达到一个极小值。同时，图 6(a)数据显示，中心轴线上的声压随着激励声压的增大而增大，即声场中的声压与激励声压成正相关关系。另外由图6(b)可知，当激励声压增大时，中心轴线上的声强也增大，且声场中的声强与激励声压的平方成正相关关系。

因此实际治疗过程中，在保持频率不变的情况下，可以直接调节激励声压以得到理想的声压、声强分布。

2.5 频率对中心轴线上声场分布的影响

2.4节已经对超声的激励声压进行了探究，本节着重讨论超声的频率对于中心轴线上声强分布的影响。

本节中所有仿真仍然采用上下两边同时照射超声的方式。由于研究的是LIPUS作用下的声场分布，故选用的频率范围为1.0～1.5 MHz。

当频率取不同数值时，中心轴线上位于z=0 µm到z=15 µm区域内的声强分布如图7所示。由图7可以发现，随着频率的增大，中心轴线处的声强减小，两者成负相关关系。这是由于当频率增大时，媒质的声衰减系数也增大，从而加大了超声能量的衰减。由此可知，通过调节频率可以改变中心轴线处的声强分布。这也说明在LIPUS的频率范围内，超声频率的取值会对骨细胞的治疗效果起到关键作用。

图7 不同频率时中心轴线处的声强分布(上下边照射，104 Pa)Fig.7 Acoustic intensity distributions at the central axis for different frequencies(ultrasound irradiates at the top and bottom sides，104 Pa)

此外，在本文所选取的LIPUS的频率范围内，频率的改变对中心轴线上声强分布的影响较小。因此在实际治疗过程中，可以通过改变频率的方式来对声强分布实施微调。

3 结 论

本文首先对培养液内的骨细胞进行了建模，并按照实际的材料参数及仿真要求对骨细胞-培养液模型的各个参数进行设定，然后选用不同参数的超声在不同位置进行照射，最后根据各项仿真的结果进行了分析与讨论。

根据仿真结果，得到如下结论：

(1) 不同的照射位置会直接影响声场的分布。

a. 底边照射超声时，声压与声强由近及远逐渐减小，声场分布不均匀。

b. 上下两边同时照射超声时，声场的分布更加均匀。故采用该种方式进行照射时，可以让骨细胞受到的声辐射力更加均匀，从而起到更好的治疗效果。另外，相比于单边照射，双边照射使得整体的声强更大。

(2) 超声在传播过程中有一部分能量被骨细胞吸收，这也从能量的角度阐述了LIPUS治疗骨细胞的作用机制。

(3) 从细胞层次的微观角度研究超声在骨细胞及培养液中的传播时，骨细胞周围区域的声强与激励声压的平方成正相关关系，而声强与频率成负相关关系，这与超声传播理论一致。另外，在1.0～1.5 MHz的频率范围内，频率的改变对声强分布的影响相对较小。

本文对包含骨细胞的微观区域内的声场分布进行了仿真与分析，并由仿真结果得到了上述的重要结论，这些结论可以指导我们进行相关的超声治疗实验。本文仅对频域稳态下的声场进行了研究，仿真结果不能展示从施加超声到最终稳定的过程中声场分布的实时情况。因此，未来的工作可以进一步在时域内对声场分布进行探讨，并与频域中的结果进行对比。

[1]Pounder NM, Harrison AJ. Low intensity pulsed ultrasound for fracture healing: a review of the clinical evidence and the associated biological mechanism of action[J]. Ultrasonics, 2008, 48(4):330-338.

[2]Claes L, Willie B. The enhancement of bone regeneration by ultrasound[J]. Progress Biophysics & Molecular Biology, 2007, 93(1): 384-398.

[3]Dyson M, Brookes M. Stimulation of bone repair by ultrasound[J].Ultrasound Med Biol Suppl, 1983, 2(2): 61-66.

[4]Ter Haar G. Therapeutic applications of ultrasound[J]. Progr Biophys Mol Biol, 2007, 93(1-3): 111-129.

[5]Jung Y J, Kim R, Ham H, et al. Focused low-intensity pulsed ultrasound enhances bone regeneration in rat calvarial bone defect through enhancement of cell proliferation[J]. Ultrasound in Medicine & Biology, 2015, 41(4): 999-1007.

[6]李立荣, 孙舒心, 他得安, 等. 低强度脉冲超声治疗骨质疏松的可行性分析[J]. 科学通报, 2016, 61(26): 2891-2900.

LI Lirong, SUN Shuxin, TA Dean, et al. Feasibility analysis of low-intensity pulsed ultrasound treatment for osteoporosis[J]. Chinese Science Bulletin, 2016, 61(26): 2891-2900．

[7]Ganzorig K, Kuroda S, Maeda Y, et al. Low-intensity pulsed ultrasound enhances bone formation around miniscrew implants[J].Archives of Oral Biology, 2015, 60(6): 902-910.

[8]Uddin S M, Qin Y X. Dynamic acoustic radiation force retains bone structural and mechanical integrity in a functional disuse osteopenia model[J]. Bone, 2015, 75(1): 8-17.

[9]Lim D, Ko C Y, Seo D H, et al. Low-intensity ultrasound stimulation prevents osteoporotic bone loss in young adult ovariectomized mice[J]. Journal of Orthopaedic Research, 2011, 29(1): 116-125.

[10]Takayama T, Suzuki N, Ikeda K, et al. Low-intensity pulsed ultrasound stimulates osteogenic differentiation in ROS 17/2.8 cells[J].Life Sciences, 2007, 80(10): 965-971.

[11]邱媛媛, 张春兵, 屠娟, 等. 微气泡激发的声微流对细胞声孔效应的影响[J]. 声学学报, 2012, 37(1): 91-96.

QIU Yuanyuan, ZHANG Chunbing, TU Juan, et al. Effect of microbubble-induced microstreaming on the sonoporation[J]. Acta Acustica, 2012, 37(1): 91-96.

[12]Unsworth J, Kaneez S, Harris S, et al. Pulsed low intensity ultrasound enhances mineralisation in preosteoblast cells[J]. Ultrasound Med Biol, 2007, 33(9): 1468-1474.

[13]Mitragotri S. Healing sound: the use of ultrasound in drug delivery and other therapeutic applications[J]. Nat Rev Drug Discov, 2005,4(3): 255-260.

[14]Behzad Vafaeian. Computer simulation of diagnostic and therapeutic applications of ultrasound propagation in biological tissues[D]. Behzad Vafaeian, 2015.

[15]Karin Hensel, Martin P. Mienkina, Georg Schmitz. Analysis of ultrasound fields in cell culture wells for in vitro ultrasound therapy experiments[J]. Ultrasound in Med. & Biol, 2011; 37(12): 2105-2115.

[16]Jennifer, Man, Richard M. Shelton, et al. Low intensity ultrasound stimulates osteoblast migration at different frequencies[J]. Journal of Bone and Mineral Metabolism, 2012, 30(5): 602-607.

[17]Upen S Patel, Sleiman R Ghorayeb, Yuki Yamashita, et al. Ultrasound field characterization and bioeffects in multiwell culture plates[J]. Journal of Therapeutic Ultrasound, 2015, 3(1): 8.

[18]王威, 董秀珍, 付峰, 等. 基于 Comsol Multiphysics平台仿真计算分析细胞形变伴随的局部组织电阻抗变化[J]. 医疗卫生装备,2010, 31(12): 20-22.

WANG Wei, DONG Xiuzhen, FU Feng, et al. Simulated calculate and analyze the impedance changes of local tissue when the cells’shapes varied by the software of comsol multiphysics[J]. Chinese Medical Equipment Journal, 2010, 31(12): 20-22.

[19]Huang J, Holt R G, Cleveland R O, et al. Experimental validation of a tractable numerical model for focused ultrasound heating in flow-through tissue phantoms[J]. J. Acoust. Soc. Am., 2004, 116(4 Pt.1): 2451.

[20]Angle S R., Sena K, Sumner D R, et al. Osteogenic differentiation of rat bone marrow stromal cells by various intensities of lowintensity pulsed ultrasound[J]. Ultrasonics, 2011, 51(3): 281-288.