3岁儿童C4-C5颈椎有限元建模方法研究*

张冠军，魏　嵬，曹立波，王丹丹

(湖南大学 汽车车身先进设计制造国家重点实验室，湖南 长沙　410082)



3岁儿童C4-C5颈椎有限元建模方法研究*

张冠军，魏嵬，曹立波†，王丹丹

(湖南大学 汽车车身先进设计制造国家重点实验室，湖南 长沙410082)

摘要：使用缩放和非线性拟合的方法，获得了3岁儿童颈部材料参数；提出了一套儿童颈部韧带的力学特性缩放方法；开发了具有精确几何及解剖学结构的3岁儿童C4-C5颈椎段有限元模型，并在准静态、动态拉伸载荷下进行了验证.结果显示，准静态拉伸刚度为211.8 N/mm，动态拉伸最终失效力为759.9 N，最终失效位移为5.083 mm，均与实验值吻合良好；动态拉伸力-位移曲线与实验曲线吻合较好.验证结果表明，本模型能够较准确地反映3岁儿童C4-C5颈椎段的准静态和动态拉伸力学特性，具有较高的生物逼真度.

关键词：儿童颈椎；生物力学；拉伸；验证

在交通事故中，儿童脊柱损伤造成的死亡率高达16%～41%，并且儿童脊柱损伤约有75%发生在颈椎段[1].儿童颈部在解剖学、形态学等方面，与成人有着明显差异，如：相对躯干更大的头部质量、更松弛的颈部韧带、更纤细的颈椎骨骼等.这些差异都将增加儿童颈部损伤风险.此外，颈部作为连接头部的重要解剖学结构，其动力学响应直接影响头部响应．因此, 儿童颈部生物力学研究对儿童颈部损伤防护及头部损伤机理的研究都至关重要.

由于缺少儿童尸体样本及用于模型验证的实验数据，公开发表的文献表明，全世界范围内仅开发了四款儿童颈部有限元模型.Kumaresan等[2]，Mizuno等[3]基于线性缩放的方法获得儿童颈部模型.这些模型均无法较准确地描述儿童颈部特有的解剖学结构；采用线弹性材料模拟颈部软组织，无法较真实地描述其生物力学特性；且这些模型均未进行验证.Meyer等[4], Dong等[5]基于儿童颈部CT图片建立了具有精确解剖学结构的颈部模型.Meyer等[4]采用刚性壳单元模拟椎骨，无法模拟骨折现象；以线弹性材料模拟椎间盘，未模拟髓核、纤维环等重要解剖学结构；对模型进行了整体动力学验证，但未对颈椎段进行单独验证.Dong等[5]在拉伸载荷条件下，仅对模型的动态拉伸最终失效力和失效位移进行验证，而未对准静态和动态拉伸刚度进行验证；未对模型进行侧向弯曲、扭转等载荷条件下的验证.

本文基于某3岁儿童颈部CT图片，致力于建立具有精确几何及解剖学结构的儿童C4-C5颈椎段模型，并赋予较准确的生物材料力学参数；再对模型在准静态、动态拉伸载荷作用下的刚度、失效力、失效位移等参数进行验证.

1颈部有限元模型的建立

选取某发育正常且无颈部损伤的3岁男童颈部CT扫描图片，利用Mimics V13.0对组织结构进行识别、提取，建立3岁儿童颈椎CAD模型.在此基础上，利用Hypermesh 10.0(Altair)对几何模型进行网格划分和前处理，建立了具有较详细解剖学结构的3岁儿童C4-C5颈椎段有限元模型，并基于LS-DYNA(971, LSTC, Livermore, CA, USA)仿真环境对本模型进行了准静态和动态拉伸验证.

C4-C5颈椎段模型包括皮质骨、松质骨、纤维环、纤维环加强纤维薄膜、髓核、前纵韧带(ALL)、后纵韧带(PLL)、黄韧带(LF)、棘间韧带(ISL)、关节囊韧带(CL)、终板、终板软骨、生长板、小关节软骨、横突软骨及椎体软骨等解剖学结构，如图1所示.

图1　3岁儿童C4-C5颈椎段有限元模型

2颈部模型材料参数

为了提高模型的生物逼真度，应尽量基于生物组织样本力学实验获得模型中的材料参数.然而，由于数值模拟技术自身的局限性，很难将实验中获得的数据直接用于有限元仿真中[6].因此，选取适当的材料模型及材料参数对模型的生物逼真度有着至关重要的作用.

2.1骨骼材料参数

在正常生理状态下，骨骼(皮质骨、松质骨)近似为线弹性材料；当骨骼的应变超过屈服应变后，其弹性模量随着载荷的增加而逐渐降低.因此，本文以*MAT_POWER_LAW_PLASTICITY材料模拟骨骼.

Panzer[7]开发的成人颈部模型中，皮质骨与松质骨弹性模量分别为16 700 MPa，291 MPa.Gilsanz等[8]发现3岁儿童骨质密度仅为成人骨质密度的0.805.本文以0.805作为3岁儿童骨骼材料的缩放系数[5].结合以上成人骨骼材料参数，即可获得3岁儿童颈部骨骼材料参数，见表1.

本文将终板的强度定义为皮质骨的1/3，即可获得终板的材料参数(见表1)[7].值得注意的是，表1中皮质骨、松质骨和终板材料参数K和N分别为材料模型中的强化系数和硬化指数，即该模型对应的塑性力学参数.

表1　3岁儿童C4-C5颈椎段有限元模型材料参数

2.2椎间盘材料参数

椎间盘主要由纤维环、髓核、终板软骨及生长板构成；其中纤维环由纤维环基质及镶嵌其中的加强纤维薄膜组成，如图2 所示.

图2　3岁儿童颈部椎间盘有限元模型

2.2.1纤维环基质材料参数

采用*MAT_HILL_FOAM材料模拟纤维环基质的非线性力学特性及其在准静态条件下的可压缩性能.该材料的工程应力公式见式(1).

(1)

式中：λ1,λ2,λ3分别为该材料3个方向的伸长率；J=λ1λ2λ3为体积变化率；Cj,bj为材料参数；n为压缩性能常数；Si为工程应力.

单轴拉伸与压缩载荷条件下，S2=S3=0，λ1=λ2，且有公式(2)：

(2)

将式(2)及J代入式(1)可得工程应力公式(3)：

(3)

根据成人纤维环单轴拉伸实验数据[11]及3岁儿童纤维环材料缩放系数0.705[9]，即可获得儿童纤维环相应的应力-伸长率数据.基于最小二乘法，利用该数据，以式(3)为目标函数进行拟合(如图3所示)，并反求材料参数Cj和Bj(见表1).

伸长率

2.2.2加强纤维薄膜材料参数

加强纤维薄膜由嵌在纤维环基质中不同方向的两簇纤维组成，其所成角度沿纤维环径向由外至内从±25°逐渐变化至±45°.本模型共设置4对加强纤维薄膜，纤维所成角度由外(第一层)至内(第四层)分别为±25°,±32°,±39°,±45°, 以*MAT_FABRIC材料模拟.根据成人椎间盘的加强纤维薄膜单轴拉伸实验[12]及3岁儿童椎间盘材料缩放系数0.705[9]，即可获得3岁儿童椎间盘加强纤维薄膜工程应力-应变曲线，如图4所示.

2.2.3髓核、生长板及软骨终板的材料参数

本文以液体材料模拟髓核，并令其材料参数与成人相同[13]；以*MAT-ELASTIC材料模拟生长板及软骨终板，并分别取其弹性模量为25 MPa，21.25 MPa[14]，见表1.

在拉伸载荷下，椎间盘失效通常发生在纤维环与软骨终板的连接部位[15].因此，本文通过定义软骨终板与生长板第一主应力失效的方式模拟椎间盘破坏.根据成人椎间盘失效力571 N[16]及儿童椎间盘失效力缩放公式[5]，获得3岁儿童椎间盘失效力为163 N.将椎间盘失效力与纤维环加强膜横截面积的比值13.8 MPa定义为椎间盘失效应力[17].

应变

2.3软骨材料参数

本文以*MAT_ELASTIC材料模拟软骨组织.成人横突软骨、椎体软骨、小关节软骨弹性模量分别为25 MPa，25 MPa和10 MPa[14].根据3岁儿童软骨材料缩放系数[9]，即可获得儿童相应软骨组织弹性模量，见表1.

2.4韧带材料参数

在韧带拉伸实验中，直接获得的数据通常为韧带的力-位移曲线.由于样本差异，各样本的韧带几何尺寸存在较大差异，因此很难直接通过定义应力-应变曲线的方式模拟其力学特性.本模型采用一维离散弹簧材料模拟韧带，该材料可以直接定义其力-位移曲线，从而更准确地模拟韧带的力学特性.

韧带的力-变形特征曲线大致呈S形,该曲线由3个控制点控制其形状，如图5所示.其中，横、纵坐标分别为各韧带的失效力、失效位移正则化后的值.假设儿童韧带的力-变形特征曲线与成人完全相同,则儿童韧带失效力与失效位移更小.

ε/εmax

由以上假设，有公式(4)：

(4)

由于儿童韧带材料缩放系数等于成人韧带伸长量与儿童韧带伸长量之比[9]，由式(4)可得公式(5)：

(5)

以上公式中εi为各控制点的应变，εmax与εmaxp分别为成人与儿童韧带失效伸长应变；dmax与dmaxp分别为成人与儿童韧带失效伸长量；l与lp分别为成人与儿童韧带长度；α为韧带材料缩放系数，本文取0.744[9]；GS为儿童颈椎几何缩放系数[10].

同理可得公式(6)：

(6)

式中：Ni为模型中模拟韧带的一维弹簧单元个数；Fi为各控制点的拉力；Fmax与Fmaxp分别为成人与儿童韧带失效力，其中：

(7)

式中：A与Ap分别为成人与儿童韧带的横截面积；σmaxp为儿童韧带失效应力；E与Ep分别为成人与儿童韧带的弹性模量.将式(7)代入式(6)中可得：

(8)

综合式(5)与(8)即为儿童韧带力-变形曲线缩放公式：

(9)

结合成人韧带力-位移特征曲线控制点数据[16]及式(9)，即可获得儿童韧带各特征点的数据.利用各韧带的3个控制点进行logistic回归拟合，即可获得各韧带的力-变形曲线，如图6所示.

位移/mm

3颈部模型仿真验证

由于儿童尸体样本稀缺，目前文献中关于儿童颈部的实验数据很少.Luck[18]利用年龄范围20周的胎儿～18周岁的儿童颈部样本，进行准静态、动态拉伸实验及弯曲-伸展实验，并获得各样本的轴向拉伸刚度、失效力及失效位移、弯曲-伸展运动范围等响应.

在车辆碰撞事故离位乘员与起爆后的安全气囊相互作用导致损伤的案例中，常见的损伤类型为颈部拉伸损伤.本文基于Luck[18]获得的儿童C4-C5颈椎段准静态和动态拉伸实验数据，对本模型进行验证.由于实验样本的年龄集中于0～22个月与72～216个月两个年龄段，本文用于模型验证的数据，主要通过插值及Luck基于其实验数据进行统计分析获得的power-law函数获得.

3.1准静态拉伸验证

Luck[18]在准静态拉伸实验中，分别将C4上表面、C5下表面与实验装置的固定端、移动端进行连接.仿真中对模型的约束及加载方式与实验设置完全相同，将C4上表面完全约束，对C5下表面进行加载，加载速度为25.7 N/s，加载至120 N卸载.

由于缺少3岁儿童的C4-C5颈椎段准静态拉伸刚度实验数据，参考Luck[18]对动态拉伸实验数据进行统计分析的方法，将其获得的C4-C5段准静态拉伸刚度进行power-law函数(见(10))非线性拟合，如图7所示.利用该函数获取3岁儿童C4-C5段准静态拉伸刚度，将仿真获得的响应与该刚度值进行对比验证.

Y=A×XB+C.

(10)

式中：Y为准静态拉伸刚度；X为年龄(月)；A，B，C分别为待定参数.

拟合获得的曲线与实验数据相关系数的平方R2=0.878，式(10)中待定参数的值分别为：A=10.063，B=0.802，C=45.344.由式(10)可得3岁儿童C4-C5颈椎段准静态拉伸刚度为223.5 N/mm.

年龄/月

3.2动态拉伸验证

Luck[18]在动态拉伸实验中的约束方式与准静态实验中完全相同，但加载方式改为以位移方式输入，且加载速率为230 mm/s.实验中获得了各样本的力-位移曲线、最终失效力、最终失效位移、20%～80%及20%～50%载荷范围内的拉伸刚度等数据，并以样本年龄为自变量，获得的响应数据为因变量进行统计学分析，获得power-law回归函数.将仿真中C4-C5的响应与基于回归函数获得的3岁儿童C4-C5实验数据进行对比验证；将仿真中获得的力-位移曲线与实验中22个月的儿童C4-C5力-位移曲线、Dewit[17]对成人C4-C5颈椎段模型进行拉伸失效仿真的力-位移曲线进行对比验证.

4仿真结果及分析

4.1准静态拉伸验证

仿真获得的C4-C5颈椎段力-位移曲线如图8所示.在50%～100%载荷范围内进行线性回归，拟合所得直线的斜率即为准静态拉伸刚度[18]，如图7所示.线性回归所得的直线与仿真曲线相关系数的平方R2=0.998，其斜率为211.8，即仿真获得的准静态拉伸刚度为211.8 N/mm，与非线性插值获得的刚度值223.5 N/mm相差5.5%，吻合较好，说明本模型在准静态拉伸载荷下的力学特性具有较好的生物逼真度.

位移/mm

4.2动态拉伸验证

C4-C5颈椎段动态拉伸仿真及年龄为22个月的样本动态拉伸实验力-位移曲线如图9所示.

位移/mm

仿真中，位移为1.748 mm时，椎间盘发生损伤，此时最初失效力为397.3 N；位移为3.151 mm时，PLL断裂；位移为5.083 mm时，ALL断裂，拉伸力达到最大值，即为最终失效力795.9 N.仿真中20%～80%及20%～50%载荷范围内的动态拉伸刚度分别为128.1 N/mm和158.2 N/mm.年龄为22个月的C4-C5样本动态拉伸实验中，最初失效力约为466.8 N，最终失效力为844.8 N，最终失效位移为5.76 mm.仿真数据与实验数据对比，仿真最初失效力、最终失效力、最终失效位移与实验数据分别相差-17.5%, -6.1%和-13.3%.3岁儿童C4-C5颈椎段动态拉伸力-位移曲线与年龄为22个月的样本的实验曲线整体吻合较好.

值得注意的是，本模型的年龄为3岁，与22个月的实验样本年龄仍有一定差距.为了更好地验证该模型的生物逼真度，利用Luck[18]获得的power-law函数得到3岁儿童C4-C5颈椎段动态拉伸最终失效力、最终失效位移、20%～80%及20%～50%载荷范围内的拉伸刚度分别为698.2 N，4.94 mm，244.3 N/mm及265.3 N/mm.仿真动态拉伸最终失效力、最终失效位移与基于power-law函数获得的响应值分别相差+12.2%、+2.8%，吻合较好；而仿真动态拉伸刚度与计算所得的刚度相差较大.仿真、实验及插值计算所得的响应值见表2.

表2　C4-C5颈椎段动态拉伸验证结果对比

作者认为动态拉伸刚度相差较大的原因可能是：仿真中很难精确地模拟失效模式，当椎间盘或韧带达到失效准则时，椎间盘或韧带的所有单元几乎同时失效，从而导致仿真曲线中的拉力明显下降.然而，实验中的椎间盘与韧带失效通常需经历较长的时间历程，因此实验输出的力-位移曲线波动较小.正是由于仿真曲线中，20%～80%及20%～50%载荷范围内出现了椎间盘和PLL失效，并导致拉力明显降低，线性回归获得的直线斜率即动态拉伸刚度比计算所得动态拉伸刚度小很多.Dewit等[17]对成人C4-C5颈椎段模型进行拉伸失效仿真获得力-位移曲线，当椎间盘、PLL和ALL失效时，拉力明显下降，曲线形状与本模型动态拉伸仿真力-位移曲线形状相似，且椎间盘、PLL和ALL失效顺序基本一致.

此外作者认为，本模型动态仿真响应与Luck[18]实验响应的误差原因还可能在于：由于儿童生物力学实验数据极度匮乏，本文只能基于22个月大的儿童颈椎动态拉伸实验曲线[18]对本模型进行验证，该实验样本与本模型的目标儿童人群的年龄(3岁)仍有一定差距，并且该实验曲线仅来自于一个实验样本，该实验样本的个体特异性也可能在一定程度上导致本模型的验证误差.

综上所述，3岁儿童C4-C5颈椎段模型动态拉伸验证中，最初失效力、最终失效力、最终失效位移及力-位移曲线与实验值吻合较好，但由于仿真手段的限制等因素，仿真所得动态拉伸刚度与实验值相差较大.总体上，本模型能够较真实地反映3岁儿童C4-C5颈椎段动态拉伸载荷下的力学特性及失效模式，具有较高的生物逼真度.这也同时证明，在模型开发过程中，作者获取颈椎各组织材料参数的方法及提出的儿童韧带力学特性缩放方法是合理且可靠的.

5结论

本文基于缩放及非线性拟合的方法获得了儿童颈部的材料参数，提出了儿童韧带力学特性曲线缩放方法，并分别在准静态、动态拉伸载荷下对模型进行了验证.验证结果表明：模型与准静态实验中的拉伸刚度值吻合较好，能够较真实地反映3岁儿童C4-C5段的准静态拉伸特性；该模型的动态拉伸力学响应与实验中的最终失效力、最终失效位移及力-位移曲线吻合较好.虽然由于仿真技术的限制，动态拉伸刚度与实验值相差较大，但本模型仍能较真实地反映3岁儿童C4-C5颈椎段在动态拉伸载荷下的生物力学特性及失效模式，并具有较高的生物逼真度.

参考文献

[1]PLATZER P, JAINDL M, THALHAMMER G,etal. Cervical spine injuries in pediatric patients [J]. Journal of Trauma Injury Infection & Critical Care, 2007, 62(2): 389-396.

[2]KUMARESAN S, YOGANANDAN N, PINTAR F A,etal. Biomechanical study of pediatric human cervical spine: a finite element approach [J]. Journal of Biomechanical Engineering, 2000, 122(1): 60-71.

[3]MIZUNO K, IWATA K, DEGUCHI T,etal. Development of a three-year-old child FE model [J]. Traffic Injury Prevention, 2005, 6(4): 361-371.

[4]MEYER F, BOURDET N, ROTH S,etal. Three years old child neck FE modelling under automotive accident conditions [C]//Proceedings of IRCOBI Conference. Maastricht, Netherlands, 2007: 277-289.

[5]DONG L, LI G, MAO H,etal. Development and validation of a 10-year-old child ligamentous cervical spine finite element model [J]. Annals of Biomedical Engineering, 2013, 41(12): 2538-2552.

[6]杨济匡, 姚剑峰. 人体颈部动力学响应分析有限元模型的建立和验证 [J]. 湖南大学学报: 自然科学版, 2003, 30(4): 40-46.

YANG Ji-kuang,YAO Jian-feng. Development and validation of a human neck FE model in impact loading condition [J]. Journal of Hunan University: Natural Sciences, 2003, 30 (4): 40-46.(In Chinese)

[7]PANZER M. Numerical modelling of the human cervical spine in frontal impact [D]. Waterloo, Canada:University of Waterloo,2006.

[8]GILSANZ V, PEREZ F J, CAMPBELL P P,etal. Quantitative CT reference values for vertebral trabecular bone density in children and young adults1 [J]. Radiology, 2009, 250(1): 222-227.

[9]YOGANANDAN N, PINTAR F, KUMARESAN S,etal. Pediatric and small female neck injury scale factors and tolerance based on human spine biomechanical characteristics [C] // Proceedings of the 2000 International IRCOBI Conference on the Biomechanics of Impact.Montpellier, France:IRCOBI Board,2000:345-359.

[10]IRWIN A, MERTZ H J. Biomechanical basis for the CRABI and hybrid III child dummies [R]// SAE Technical Paper 973317. Washington, DC: SAE International, 1997.

[11]DIANE R, WAGNERLOTZ J C. Theoretical model and experimental results for the nonlinear elastic behavior of human annulus fibrosus [J]. Journal of Orthopaedic Research, 2004, 22(4): 901-909.

[12]HOLZAPFEL G A, SCHULZE-BAUER C A, FEIGL G,etal. Single lamellar mechanics of the human lumbar anulus fibrosus [J]. Biomech Model Mechanobiol, 2005, 3(3): 125-140.

[13]PANZER M B, CRONIN D S. C4-C5 segment finite element model development, validation, and load-sharing investigation [J]. J Biomech, 2009, 42(4): 480-490.

[14]YAMADA H, EVANS F G. Strength of biological materials[M]. Baltimore: Williams & Wilkins, 1970: 297.

[15]KASRA M, PARNIANPOUR M, SHIRAZI-ADL A,etal. Effect of strain rate on tensile properties of sheep disc anulus fibrosus [J]. Technol Health Care, 2004, 12(4): 333-342.

[16]YOGANANDAN N, KUMARESAN S, PINTAR F A. Biomechanics of the cervical spine Part 2 Cervical spine soft tissue responses and biomechanical modeling [J]. Clinical Biomechanics, 2001, 16(1): 1-27.

[17]DEWIT J A, CRONIN D S. Cervical spine segment finite element model for traumatic injury prediction [J]. Journal of the Mechanical Behavior of Biomedical Materials, 2012, 10: 138-150.

[18]LUCK J F, NIGHTINGALE R W, LOYD A M,etal. Tensile mechanical properties of the perinatal and pediatric PMHS osteoligamentous cervical spine [J]. Stapp Car Crash Journal, 2008, 52: 107-134.

Study on Development and Validation Method of Finite Element Models for the 3-Year Old Pediatric C4-C5 Cervical Spine

ZHANG Guan-jun, WEI Wei, CAO Li-bo†, WANG Dan-dan

(State Key Laboratory of Advanced Design and Manufacturing for Vehicle Body, Hunan Univ, Changsha, Hunan410082, China)

Abstract:Based on scaling and nonlinear-fitting methods, the material properties of the 3-year old pediatric cervical spine were acquired． A specific scaling method for the mechanical properties of the pediatric cervical ligaments was proposed, and the C4-C5 cervical segment finite element model for the 3-year old child was developed with accurate geometries and anatomical structures and validated under quasi-static and dynamic tensile loading. It was indicated that the simulation responses were consistent with those of experiments, with an quasi-static tensile stiffness of 211.8 N/mm, a dynamic tensile ultimate failure force of 759.9 N and an ultimate failure displacement of 5.083 mm in simulations. In addition, the force-displacement curve in the dynamic tensile simulation was also similar to that of the experiment. It can be concluded that this model is able to reflect the quasi-static and dynamic tensile mechanical properties of 3-year old pediatric C4-C5 cervical segment with higher biofidelity.

Key words:pediatric cervical spine; biomechanics; tensile; validation

中图分类号：U461.91

文献标识码：A

作者简介：张冠军(1981-)，男, 山东济宁人，湖南大学助理教授, 博士†通讯联系人，E-mail:hdclb@163.com

基金项目：国家自然科学基金资助项目(51205118), National Natural Science Foundation of China(51205118);中央高校基本科研业务费资助项目(531107040162)

*收稿日期：2015-01-04

文章编号：1674-2974(2016)02-0008-07