ABAQUS医疗健康:SMILE手术眼部受力与变形有限元仿真

1．引言

由于缺乏户外运动、用眼过度、频繁使用电子设备等原因，中国近视患者人数逐年增加。近视患者的眼睛由于无法正常聚焦，看远处物体模糊不清(巩倩文 , 2016)。随着医疗技术的进步，治疗近视的手术方法越来越丰富，近视患者可以通过手术矫正视力，不需要继续戴眼镜，小切口基质透镜取出术（SMILE）就是其中一种较为先进的屈光手术。

 

小切口基质透镜取出术也被称为“全飞秒”手术，手术原理是通过在角膜内切削出一个小透镜并取出来改变角膜表面形状，进而矫正屈光不正。手术的过程如图 1 所示，首先在角膜基质层中切削出一个小透镜，再通过微小切口将透镜取出，等待角膜自然愈合后，角膜外形改变，进而矫正屈光度，达成手术目的。

 

 

图 1. SMILE 手术过程示意图.

 

通过建立眼球有限元模型进行仿真，可以获得术后眼球变形与受力分布，对屈光手术的治疗效果与安全性进行评价。在对屈光手术的有限元建模方面，研究人员已经建立了各种不同类型的数值模型，但针对 SMILE 手术的有限元模拟较少(陈修国, 2017; 孙立科, 2021)，且缺少对去除角膜中透镜部分过程的模拟。

 

本文简化了眼球结构，针对 SMILE 手术，利用消融函数建立透镜的几何模型，选择 Ogden描述角膜力学特性的材料模型，利用生死单元法模拟去除透镜部分过程，建立有限元模型模拟SMILE 手术后角膜受力和变形情况，为评价手术效果和安全性提供依据。

 

2．模型建立

2.1 ABAQUS 二次开发

ABAQUS 是一款功能强大的有限元分析软件，随着当今用户的需求越来越多样化，要求有限元分析软件提供方便易用的二次开发界面，供用户根据实际需求进行定制。如图 2 所示。

 

 

图 2. ABAQUS 脚本接口与 ABAQUS/CAE 的关系.

 

ABAQUS/CAE 接口使用 Python 包实现。使用 Python 进行 ABAQUS 二次开发，可以代替用户进行重复的前后处理操作，根据用户输入的参数自动建立模型，提交分析方案，提取结果数据，从而节省时间，提高工作效率。在建立 SMILE 外科有限元模型的过程中，使用 Python 脚本可以方便地改变参数，快速建立一系列的模型进行分析。

 

2.2 几何模型

人眼的结构非常复杂，而屈光系统相对简单，其中角膜对屈光力的贡献高达 70%。为简化模型，角膜的内表面和外表面可以简化为一个球形表面。通过观测患者角膜前后表面中心曲率半径，在 ABAQUS/CAE 2022 中建立角膜几何模型，如图 3 所示。

 

 

图 3. 角膜几何模型.

 

SMILE 手术使用激光来切割角膜。角膜切割部分的厚度由消融函数决定。消融函数的公式表达为：

 

 

 

公式（1）中，r1为手术前角膜前表面的曲率半径，r2为手术后期望得到角膜前表面的曲率半径，p为与角膜中心轴线距离，f(p)表示在对应距离上角膜消融的厚度，Φ是视力矫正光学区直径，取6.0mm。使用公式（2），可以计算r2的值：

 

 

 

S是要矫正的近视度数，单位为。在不同矫正屈光度下，消融曲线也不相同，如图4所示，透镜厚度随着矫正屈光度的增加而增加。

 

 

图 4. 不同屈光度下的消融函数.

 

由公式（1）(2)，根据患者角膜表面曲率半径以及矫正屈光度大小得到透镜轮廓表达式，以0.03mm 为间隔生成一系列轮廓点，在草图中使用样条曲线连接轮廓点，通过扫掠方式建立透镜部分几何模型，如图 5 所示。

 

 

图 5. 透镜部分几何模型.

 

使用布尔运算功能将角膜部件和透镜部件合并，选择保留部件边界，得到已预先划分好透镜区域的角膜部件，便于之后划分网格，如图 6 所示。图中绿色区域为要切除的透镜，灰色区域为透镜剩余部分。

 

 

图 6. 角膜部件.

 

2.3 相互作用与边界条件

 

图 7. 相互作用与边界条件.

 

角膜上施加的相互作用与边界条件如图 7 所示。

 

巩膜与眼内肌相连，在 SMILE 手术中可以认为与人体相固连，因此在角膜末端与巩膜相连部分设置固定约束。角膜内部充盈着房水，在角膜内表面施加 10mmHg 的眼内压。在角膜外部切削区域，由于透镜取出后空腔内压强小于大气压，因此还会产生负压，在角膜外表面切削区域施加压强，大小与眼内压一致。

 

在透镜区域定义“模型改变”相互作用，这一方法可以使单元失效或重新生效，从而模拟对应区域的变化。在眼内压施加完成后，透镜区域消失，模拟手术去除透镜的过程。

 

2.4 材料模型

角膜的应力-应变关系呈现出非线性特征，可以使用 Ogden 二阶方程来描述，具体参数如表1 所示(孙立科, 2021)。

 

表 1. 角膜 Ogden 模型参数.

 

 

 

2.5 网格划分

划分网格时分两部分进行，角膜基质与透镜部分。由于在透镜的尖锐边缘处不适合使用六面体单元，所以在透镜部分使用四面体单元，并以自由网格的方式进行划分。角膜基质中心部分采用结构化网格，其余部分通过扫掠方式得到，均为六面体单元，如图 8 所示。

 

 

图 8. 角膜网格划分.

 

为了准确描述角膜的变形和应力分布，网格必须足够密以获得准确的数值结果；同时过密的网格也会导致计算成本的上升。为了选择计算成本小同时结果精度高的尺寸，分别使用 0.2mm至 0.8mm 共七种尺寸划分网格，分析网格密度对计算精度的影响。使用应力最大值与角膜中心点位移值作为检查网格尺寸无关性的指标，结果如图 9 所示：

 

 

图 9. 不同网格尺寸下最大应力值以及角膜中心点位移.

 

角膜中心点位移值随尺寸减小逐渐减小，在尺寸为 0.4mm 之后几乎不变。从 0.6mm 开始，最大 Mises 应力值随尺寸减小而不断增大，在尺寸为 0.3mm 之后变化很小。在随后的建模中，网格尺寸选择为 0.3mm。

 

3．结果与讨论

 

图 10. 使用 Model Change 方法去除透镜部分.

 

通过使用 Model Change 方法，在 Apply 分析步骤后，透镜部分不再参与模型计算，以模拟去除透镜的过程，如图 10 所示。

 

3.1 屈光度的影响

使用 SMILE 手术有限元模型，在眼压为 10mmHg 的情况下，分析了矫正屈光度分别为-2D、-4D、-6D、-8D 和-10D 五组手术术后角膜应力分布，结果如表 2 所示：

 

表 2. 不同屈光度下角膜上应力分布

 

 

观察云图可以发现应力值较大单元集中在角膜中心区域，这主要是由于取出透镜之后在眼内压与负压作用下角膜前后部分被挤压贴合在一起，发生较大变形，应力较大。

 

 

图 11. 矫正屈光度-最大应力曲线.

 

如图 11 所示，矫正屈光度越高，应力最大值越大，角膜中心点位移也随之升高。这主要是由于矫正屈光度越大，要切除的透镜区域越厚，变形越大。这说明高度近视患者进行 SMILE 手术之后角膜应力较大，可能有安全性上的隐患；此外由于透镜厚度较大，前后表面可能无法完整贴合，影响手术效果。

 

3.2 角膜瓣厚度的影响

角膜瓣位于是位于切削区域透镜外侧的角膜区域，其厚度对手术的结果也有影响，对屈光矫正度为从 0.09mm 到 0.13mm 共五组如图 12 所示。

 

 

图 12. 角膜瓣厚度对最大应力和角膜中心点位移的影响.

 

随着角膜瓣厚度的增加，角膜上的最大应力也随之减小，安全隐患更小，但也对角膜愈合有不利影响；角膜中心点位移则变化规律不明显，不同角膜瓣厚度下变化不大。因此，为了获得更好的手术效果，在不影响安全性的前提下，应选择尽可能小的角膜瓣厚度，更有利于愈合。

 

4．总结

本文建立了 SMILE 手术的有限元模型，对手术后眼球受力与变形进行有限元模拟，可以在术前对手术效果和安全性进行评估。使用模型模拟了不同屈光矫正度下角膜应力与位移分布，分析了不同角膜瓣厚度下角膜最大应力和角膜中心位移值，验证了有限元模型的实用性。主要结果如下：

 

1. 随着矫正屈光度的增加，角膜应力最大值逐渐增大，这表明在治疗高度近视患者时应更加注意手术的安全性。

 

2. 随着晶状体外瓣厚度的增加，角膜应力最大值逐渐减小，角膜中心位移值变化不大。

 

资料来源：达索官方