脊柱胸腰段三维有限元模型的建立与验证

陈浩，张锦洪，贺增良

（东南大学附属第二医院骨科，江苏南京210003）

脊柱胸腰段三维有限元模型的建立与验证

陈浩，张锦洪，贺增良

（东南大学附属第二医院骨科，江苏南京210003）

目的：建立正常脊柱胸腰段三维有限元模型并验证其有效性。方法：选择一位健康成人志愿者，通过CT扫描获得胸腰段数据，Mimics 10.01软件生成T12-L2椎骨三维模型，Geomagic 11.0、UG 7.0生成三维实体模型，Hypermesh 10.0划分网格，并建立椎间盘、韧带等结构，Abaqus6.9.1加载求值，施加500 N垂直压缩力，7.5 Nm的前屈、后伸、左右侧屈及左右轴向旋转力矩，分三个载荷步（2.5、5.0和7.5 Nm）施加，观察T12-L2节段运动范围。结果：建立了正常人胸腰段三维有限元模型，整个模型共90 014个单元，158 684个节点。对模型进行前屈、后伸、左右侧屈和左右旋转作用，T12-L2节段的力矩 旋转角度曲线与文献体外实验数据基本吻合。结论：本实验建立的正常脊柱胸腰段三维有限元模型有效，可进行进一步的生物力学实验研究。

胸腰椎；有限元；生物力学

脊柱胸腰段骨折约占脊柱骨折的30%～60%，其中相当一部分患者需要手术治疗，包括前路、后路及前后路联合等手术方式［1－3］。众多学者从生物力学的角度对这些不同手术方法进行研究，然而，传统生物力学研究方法如电测法、光弹法等，很难取得全域性信息。本课题组通过建立脊柱胸腰段的有限元模型并对其有效性进行验证，以期为胸腰段骨折不同手术方式的进一步研究提供基础数据。

1 材料与方法

1.1 三维有限元模型的建立

选一名健康男性志愿者（年龄27岁，身高172mm，体质量68 kg）为研究对象，X线检查排除脊柱病变。64排螺旋CT从T12～L2进行连续扫描，扫描参数：层厚1.0mm，球管电流200mA，电压120 kV。各断层图像以Dicom 3.0格式读入Mimics 10.01软件，通过区域选取、阈值分割、蒙罩编辑等功能提取T12～L2椎骨骨骼图像，并建立STL三角网格模型；然后导入Geomagic 11.0进行修补、优化等操作，建立三维曲面模型；运用UG 7.0对曲面缝合生成实体模型，然后导入Hypermesh 10.0划分网格，并在骨性结构的基础上创建椎间盘、前纵韧带、后纵韧带、黄韧带、棘间韧带、棘上韧带、关节囊韧带等结构，髓核面积占椎间盘面积的43%［4］；皮质骨采用平均厚度为1mm的C3D8单元，松质骨采用C3D8，终板采用0.5mm厚的C3D8单元；韧带采用只有轴向平移自由度的、双节点SPRINGA弹簧单元进行划分。椎间盘（含髓核和纤维环）以及终板，采用增强沙漏控制的三维六面体减缩积分C3D8R；最后导入有限元分析软件Abaqus 6.9.1，参考文献［4－7］设置椎骨、椎间盘及韧带等各结构的材料属性（表1）；将上、下小关节面间设置为面对面接触，摩擦系数为0.2；最终生成T12～L2的完整有限元模型（图1）。

表1 胸腰段有限元模型的材料属性

1.2 边界条件和加载方式

固定L2椎体下表面，在T12椎体上表面施加载荷，压力方向垂直向下，大小500 N，均匀分布在整个T12椎体的上终板。在垂直加载500 N的同时，于T12椎体上表面施加7.5 Nm的前屈、后伸、左右侧屈及左右轴向旋转力矩，分3个载荷步（2.5，5和7.5Nm）施加。

1.3 主要观察指标

测量T12-L2节段在不同工况下的旋转角度，取中立位T12椎体正中矢状位上表面前、后缘两节点A、B，连为一直线AB，前屈或后伸时，A、B点分别移至A′、B′，A′、B′连为一直线A′B′，AB与A′B′之间的夹角即为前屈或后伸旋转角度（图2）。以同样方法测量左右侧屈或左右轴向旋转等工况下的旋转角度。

图1 T12-L2有限元模型

图2 T12-L2节段在前屈或后伸工况下的旋转角度示意图

2 结果

2.1 有限元模型的建立

建立了正常脊柱胸腰段（T12-L2）有限元模型，包括椎体皮质骨、松质骨、纤维环、髓核、小关节、椎弓根、椎板、横突、棘突、前纵韧带、后纵韧带、黄韧带、棘间韧带、棘上韧带、关节囊韧带等结构。整个模型共90 014个单元，158 684个节点。

2.2 模型的有效性验证

对T12-L2完整模型施加500 N的预载荷，并分别施加前屈、后伸、左右侧屈及左右旋转力矩7.5 Nm，测量得出T12-L2的各载荷步旋转角度，力矩－旋转角度走形，结果与Panjabi等［8］的体外实验结果基本相符合（图3）。

图3 有限元模型与体外实验数据比较

3 讨论

3.1 胸腰段三维有限元模型建立的意义

脊柱胸腰段包括下胸椎与上腰椎区域，通常为T11至L2脊柱节段。此区域为胸椎后凸与腰椎前凸的交界区，也是脊柱屈曲的顶点，所以外伤时形成的杠杆力往往集中在胸腰段，使该部位很容易受到传导暴力，造成胸腰段的损伤。该部位与硬膜囊的有效间隙相对狭窄，而且是腰膨大区域，损伤后也容易形成脊髓和神经压迫和损伤。因此对该区域的生物力学研究研究尤为重要。传统的生物力学研究常通过建立物理模型、在体模型或离体模型进行，这些方法都具有一些缺点，如物理模型无法体现生理状态；在体模型难以获取生物体尤其是内部组织的实验数据；离体模型面临标本获取、实验成本等方面的困难，也难以获取实验对象内部任意部位的生物力学信息［9］。有限元法是一种离散化的数值计算方法，其前提是建立有效的有限元模型，与传统生物力学分析方法相比，具有更多的优点：能够获得椎骨及椎间盘内部的张力、应力等信息，能够模拟在活体发生的真实现象，可模拟肌肉对脊柱力学的影响，能在没有真正使用的情况下评价新的脊柱器械等［10－11］。因此，建立脊柱胸腰段有限元模型对胸腰段损伤及手术方法的生物力学研究具有重要的意义。

3.2 胸腰段有限元模型的验证

脊柱有限元模型的验证方法通常有两种，即直接验证和间接验证。直接验证是指将有限元模型与离体模型的结果直接比较，多数情况下模型参数如几何形状和材料属性等，直接来源于离体实验样本，有限元模型与离体模型的几何形状、边界条件也基本接近，然而多数情况下离体模型仍难以复制在体的状况；间接验证是指将有限元模型与文献、临床试验或历史医学数据进行比较，如结果非常接近即可认为该模型有效，该模型也可用于进一步的生物力学研究［12］。本实验采用螺旋CT薄层扫描获得连续断面图片，以二维Dicom格式存储，然后用Mimics软件对CT图像三维重建，再运用Geomagic，UG，Hypermesh，Abaqus等一系列软件建立了胸腰段的三维有限元模型。在建模的过程中，完全采用数字化操作，并自动验证模型的准确性。而且划分的单元类型和节点多，能够精确真实地反映胸腰段复杂的实际几何形状。本模型的所有材料参数均来自于国内外其他作者的脊柱生物力学实验研究所测得的数据，数值可靠。模型的有效性验证是下一步能否用来模拟活体更复杂的情况的重要环节［12］。应力-应变曲线比较和力矩 旋转角度曲线比较是有限元模型有效性验证的常用方法［13］。本实验采用力矩 旋转角度曲线进行比较，结果基本吻合，从而验证了模型的有效性，可用于临床与实验研究。

3.3 本建模方法的优点

本建模方法具有以下优点：① 模型精度高：脊柱形态结构复杂，以往通过测量的方法所得三维数据难以表达复杂及细微的结构。Dicom文件提供了非常清晰的组织密度信息，本研究运用Mimics软件，读取Dicom格式的CT资料建立胸腰段骨性结构的几何模型，更加接近人体。② 建模速度快：Mimics软件能够自动识别、读取连续CT图像，一次可读取数百张图片，并且运行时间很短；Mimics软件可对研究对象自动建模，对边界适当编辑，建模时无须手工确定节点坐标及其他处理，建模速度大大提高。③便于手术模拟：本实验运用UG软件生成实体模型，可以根据手术实际情况，对实体进行任意切割，模拟手术过程，并且可以根据需要设计内固定装置，并直接进行装配，如椎体爆裂性骨折的后路减压椎弓根螺钉固定、前路减压椎体侧方钉棒或钉板系统固定等。

3.4 本研究存在的不足

本模型的不足：①所建立的胸腰段有限元模型的材料属性比较简单，各材料均模拟为单相，事实上各材料为双相甚至多相，所以材料多相性是需解决的问题；②材料常数均引用文献数据，未采用样本相应的材料常数，即个体化材料常数，这也是目前脊柱有限元研究的困难之一；③CT图像难以区分椎间盘髓核、纤维环、韧带等软组织的边界，本实验只能根据解剖学数据模拟纤维环、髓核及韧带等组织的几何形态，难免与特定样本的解剖存在差异而造成结构偏差。

总之，有限元法作为生物力学研究方法的一种，应用范围越来越广泛，技术也日趋成熟。在骨科生物力学研究中更是如此。本实验建立了胸腰段T12-L2有限元模型，该有限元模型前屈、后伸、左右旋转、左右侧屈等各种工况下的运动范围与文献相关资料一致，可为进一步研究胸腰段生物力学提供基础。

［1］Kifune M，Panjabi MM，Liu W，et al.Functionalmorphology of the spinal canal after endplate，wedge，and burst fractures［J］.J Spinal Disord，1997，10（6）：457－466.

［2］Dai LY，Jiang SD，Wang XY，et al.A review of the management of thoracolumbar burst fractures［J］.Surg Neurol，2007，67（3）：221－231.

［3］Radcliff K，LimthongkulW，Gruskay J，et al.Surgical planning for the treatment of thoracolumbar fractures：anterior，posterior，or combined approach？［J］.Semin Spine Surg，2012，24（4）：244－251.

［4］Polikeit A，Ferguson SJ，Nolte LP，et al.Factors influencing stresses in the lumbar spine after the insertion of intervertebral cages：finite element analysis［J］.Eur Spine J，2003，12（4）：413－420.

［5］Natarajan RN，Andersson GB.The influence of lumbar disc height and cross-sectional area on the mechanical response of the disc to physiologic loading［J］.Spine（Phila Pa 1976），1999，24（18）：1873－1881.

［6］Zhong ZC，Wei SH，Wang JP，et al.Finite element analysis of the lumbar spine with a new cage using a topology optimization method［J］.Med Eng Phys，2006，28（1）：90－98.

［7］Goel VK，Ramirez SA，KongW，etal.Cancellous bone Young′smodulus variation within the vertebral body of a ligamentous lumbar spine-application of bone adaptive remodeling concepts［J］.J Biomech Eng，1995，117（3）：266－271.

［8］Panjabi MM，Oxland TR，Lin RM，et al.Thoracolumbar burst fracture.A biomechanical investigation of its multidirectional flexibility［J］.Spine（Phila Pa 1976），1994，19（5）：578－585.

［9］王艳国，张琪，刘凯.有限元法在颈椎生物力学中应用研究分析［J］.中国中医骨伤科杂志，2013，21（7）：62－64.

［10］胡勇，谢辉，杨述华.三维有限元分析在脊柱生物力学中应用研究［J］.医用生物力学，2006，21（3）：246－250.

［11］陈浩，于晓华，华国军，等.腰椎运动节段流固耦合有限元模型的建立与验证［J］.中国组织工程研究与临床康复，2010，14（52）：9706－9709.

［12］Jones AC，Wilcox RK.Finite element analysis of the spine：towards a framework of verification，validation and sensitivity analysis［J］.Med Eng Phys，2008，30（10）：1287－1304.

［13］Guan Y，Yoganandan N，Zhang J，etal.Validation of a clinical finite elementmodel of the human lumbosacral spine［J］.Med Biol Eng Comput，2006，44（8）：633－641.

Establishment and validation of the finite elementmodel of thoracolumbar motion segment

CHEN Hao，ZHANG Jin-hong，HE Zeng-liang

（Department of Orthopedics，the Second Hospital of Nanjing，Affiliated to Southeast University，Nanjing Jiangsu 210003，China）

Objective：To establish the 3D finite elementmodel of thoracolumbarmotion segment and to verify the validity of the finite elementmodel.M ethods：Spiral CT scan was performed in the T12-L2motion segments of amAle volunteer without any spinal disease.The CT data was imported to the software of Mimics 10.01 to create the geometricmodel of the T12-L2motion segmentwhich was reconstructed to the solid model by the Geomagic 11.0 and UG 7.0 software.The intervertebral discs and ligamentswere constructed and the intactmodelwasmeshed in the Hypermesh 10.0 software.Inter-segmental ranges ofmotion were calculated after T12-L2being subjected to loads ofmoments 7.5 Nm（in three steps：2.5，5.0 7.5 Nm）and preload and L2was rigidly fixed while loadswere applied at the T12for flexion，extension，lateral bending and axial rotation.For validation of the model，their predicted inter-segmental ranges ofmotion was compared with the results in the experimental study.Results：The established final intactmodel consisted of 90，014 elements and 158，684 nodes.The curve ofmoments and angular motion were in good agreement with the published experimental results.Conclusion：The 3D finite elementmodel of thoracolumbarmotion segment can simulate the natural condition and facilitate the further biomechanical research.

thoracolumbar；finite element；biomechanics

R318.01

A

1671－7783（2014）01－0072－04

10.13312／j.issn.1671-7783.y130202

陈浩（1979—），男，江苏泰兴人，医师，硕士，主要从事脊柱外科临床研究；张锦洪（通迅作者），副主任医师，硕士，E-mail：zhangjinhong 1998＠163.com

2013－09－06 ［编辑］何承志