基于多体动力学和有限元方法的三级跳运动人体膝关节冲击损伤分析研究

第一作者刘同众男，副教授，1968年2月生

通信作者朱林男，博士后，副教授，硕士生导师，1971年10月生

基于多体动力学和有限元方法的三级跳运动人体膝关节冲击损伤分析研究

刘同众1,朱林2,3,程曦2,彭双双2,戚银银2,秦义先3,尹成龙2,章文峰2,朱德泉2

(1. 安徽农业大学体育系,合肥230036； 2. 安徽农业大学工学院机械系,合肥230036；3. 美国纽约州立大学生物工程系,纽约11794-5281)

摘要：为了研究三级跳训练、竞赛对运动员膝关节冲击损伤影响程度，以1名健康男性运动员左膝关节为研究对象，利用膝部CT扫描和磁共振(MRI)相结合的方法，建立了包括股骨、胫骨、腓骨、髌骨以及膝部各主要软骨、韧带在内完整的三维膝关节立体模型；采用多体动力学(MDA)与有限元分析(FEA)相结合的方法，计算出三级跳起跑、跳跃、落地全过程运动员膝关节处冲击状态下的应力、应变及位移变化状况。计算显示：在三级跳运动过程中，膝关节胫骨外侧接触区域受到的载荷，发生的位移最大，并随着跳跃和落地的冲击程度增强而变大。这加剧了胫骨部位的磨损程度，因而易诱使运动员膝关节发生骨损伤，计算结果与国外相关文献相比较为一致。结果表明：借助有限元与多体动力学相结合的方法，能够较好地研究三级跳训练和竞赛过程中人体膝关节的冲击振动行为，为膝关节运动损伤开辟了一条新的研究方法，也为膝关节运动损伤的防治提供了一定的参考依据。

关键词：膝关节；三级跳；冲击行为；多体动力学；有限元

基金项目：教育部第46批留学回国人员科研启动基金项目；国家留学基金资助；教育部科学研究项目(12yjazh075)

收稿日期：2014-07-10修改稿收到日期：2014-09-01

中图分类号:TJ818

文献标志码：A

DOI:10.13465/j.cnki.jvs.2015.17.008

Abstract:In order to better understand the effects of triple-jump movement on jumpers, impact behaviors of left knee joint of a healthy man athlete were studied. Based on CT and MRI scanning pictures of the jumper’s knee joint, a three-dimensional model was built, including complete femur, tibia, fibular, patellar, cartilages and main ligaments of the knee joint. Then, a method combining finite element analysis (FEA) and multi-body dynamic analysis (MDA) was used to simulate a real human triple-jump. Finally, the loading data, such as, muscle force and joint reaction obtained from an MDA were impulted into a finite element model to compute stress, strain and displacement variations of the knee joint under impact conditions. The numerical simulation results showed that the stress, stain and the displacement of medical cartilage are larger than those of lateral cartilage; the contact force and contact areas increase with increase in load, consequently, this can incur severer wear contours of the knee joint and badly injure the jumper. The computation results matched reasonably the research results in the published literatures. The results demonstrated the importance of performing MDA as a preliminary step to FEA in biomechanical fields, and the proposed method had a potential to generate a better understanding of injury mechanisms of jumpers.

Impact behaviors of human knee joint during triple-jump based on multi-body dynamics and finite element analysis

LIUTong-zhong1,ZHULin2,3,CHENGXi2,PENGShuang-shuang,QIYi-yin2,QINYi-xian3,YINCheng-long2,ZHANGWen-feng2,ZHUDe-quan2(1. Department of Physical Education, An Hui Agricultural University, Hefei 230036，China；2. Mechanical Engineering Department, An Hui Agricultural University, Hefei 230036，China；3. Bioengineeing Department, The State University of New York at Stony Brook, NY 11794-5281, USA)

Key words:knee joint; triple-jump; impact behavior; MDA; FEA

膝关节是人体最大、最复杂的关节，位居髋关节和踝关节之间，是下肢活动的枢纽(见图1 )。它的任何一个主要组成部分的损坏都会引起膝关节的反常运动，进而造成软骨和半月板磨损、变性，甚至出现运动损伤[1]。

图1　膝关节结构示意图 Fig.1 Anatomy of the knee

三级跳是一项技能性要求较高的体育项目，它是人体经过快速助跑后，由单腿跳、跨步跳和跳跃所组成的连续三次腾跃水平距离的运动(见图2)，因而对膝关节的功能要求极高，对膝关节的冲击损伤也较为严重，目前是所有运动项目损伤之首，约占总伤病的50%[2]。

图2　运动员三级跳过程示意图 Fig.2 Human triple jump for track and field sports

膝关节不仅结构复杂，而且其承受冲击受载的状况也较为复杂。相关研究结果表明，无论是运动损伤还是关节累积负重使关节软骨退化所造成的骨关节病，其生物力学根源主要是与关节运动冲击和关节内部应力分布不合理有关[3]。因此，研究三级跳过程膝关节部位的运动特征、冲击状况下应力应变变化，对揭示其发病机理，更好地制定预防和康复计划，减少运动员的运动损伤非常重要。然而，由于直接测量膝关节内部的应力应变、关节的接触面积以及半月板、韧带变形与受力较为困难，因此目前很多学者采用有限元方法研究膝关节部位的冲击和生物力学行为，并取得了一些较为重要的研究成果[4]。但值得注意的是，在利用有限元方法进行求解时，由于常常直接选取三级跳整个运动过程中最大的载荷进行加载，因而无法考虑膝关节累积疲劳损伤效应，这样，所得到计算结果与膝关节实际损伤机理相脱节，致使人们无法更好地理解膝关节部位生物力学行为。

基于此，本文在考虑膝关节运动和累积疲劳损伤事实的基础上，引入机械动力学设计思想，采用有限元分析技术与多体动力学相结合的方法，利用牛顿拉夫森迭代算法模拟三级跳起跑、跳跃和落地全过程，然后提取上述过程膝关节受力最大载荷，通过有限元分析手段对其刚度和强度进行更深一步的数值模拟，并将计算结果与国外相关文献相比较，具有较好的一致性，为膝关节运动损伤开辟了一条新的研究方法，也给膝关节运动损伤的预防和治疗提供了一定的参考依据。

1膝关节三维模型的建立

研究对象：选取一名20岁健康男性志愿者，身高175cm，体重70kg，无膝关节痛病史、外伤史，且X线检查膝关节也未见畸形及退变。

数据采集：采用螺旋CT和核磁共振成像(MRI)相结合的方法对受试者左膝关节处进行扫描，采集该部位软组织和骨组织几何数据；扫描主要集中在股骨中下段及胫腓骨中上段部位；扫描时，实验者取仰卧位，双膝关节位于扫描视野中心，保持纵轴方向不动，尽量保持扫描断面与身体长轴垂直；扫描图像层厚0.5mm，512×512像素，共获得图像序列194张(见图3)。

图3　膝关节部分CT断层图像序列示意图 Fig.3 Picture of the knee joint using CT

建模步骤：首先，将上述获取的膝关节部位几何数据以DICOM格式存储在计算机中；然后，将这些数据导入医学图像处理平台Mimics (Materialise’s Interactive Medical Image Control Sytem, Belgium)，利用断层扫描图像三维建模技术构建膝关节的三维实体模型；最后，将上述得到的模型导入机械三维建模软件SolidWorks中，借助反求工程法和基于特征的零件实体建模技术进一步修改、完善膝关节的三维实体模型，见图4。

图4　在(A)Mimics, (B)SolidWorks中构建的膝关节模型 Fig.4 3D models of the knee joint in (A)Mimics, (B) SolidWorks

2膝关节动力学建模

根据文献[3]，可以将人体视为一个完整的多刚体系统，能借助系统动力学方法研究其运动规律，即将人体运动等效成多刚体运动，刚体之间有相对运动，且通过一定形式的约束连接，如铰链、轴承和弹簧等。

由此，本文建立了包括股骨、胫骨、腓骨、髌骨以及膝部各主要软骨、韧带在内的完整的三级跳膝关节多体动力学计算模型(见图5)；并依据机械多体动力学思想分别定义了关节软骨、胫骨、肌肉和韧带等组织的材料特性[5]。

图5　三级跳膝关节多体动力学计算模型框图 Fig.5 Block diagram showing the calculation of the knee joint dynamics during triple jump

图5中:q为三级跳运动中人体的位移变量，lM为膝关节部位肌肉变化长度，vM为肌肉运动速度，fM为肌肉产生的载荷，ta为膝关节处主动关节受到的力矩，tp为被动关节受到的力矩，d为关节处转矩力臂。

3膝关节运动的动态仿真

3.1基于DDM的仿真计算流程

本文基于牛顿拉夫森迭代算法，利用内嵌于SolidWorks之上的DDM(Dynamic Designer Motion ，美国MDI公司研制的全功能动力学仿真软件)模拟三级跳助跑、起跳、腾空和落地运动全过程，对膝关节进行冲击动力学特征分析，图6为求解计算流程图。

图6　求解运算流程图 Fig.6 The flowing chart using DDM

3.2初始运动状态的建立

整个动力学模拟计算依据三级跳运动特点，共涉及三个阶段四个步骤，即助跑、起跳、腾空、落地，见图7，重点分析三级跳运动过程中人体膝关节受到地面冲击载荷的变化状况。

图7　本文计算的三级跳运动过程三个阶段示意图 Fig.7 Three contact phases of the analyzed triple jump

设置各仿真计算初始参数如下：

(1)计算距离：14.0 m。

(2)膝关节骨组织参数：依据骨组织材料非均匀性和各向异性的特点，采用均分方法，利用Mimics软件，并根据前述膝关节CT扫描图像不同的灰度值，建立图像灰度值与其材料属性函数关系，计算出骨组织材质属性及力学特性参数如表1所示。

(3)膝关节韧带特性参数：考虑到韧带材料非线性特点，设置各韧带弹性模量为430 MPa，泊松比0.45[6]，其它相应参数如表2所示。

表1　膝关节骨组织材质属性

表2　膝关节各条韧带的材料参数 [6]

表2中:MCL为内侧副韧带，LCL为外侧副韧带，ACL为前交叉韧带，PCL为后交叉韧带；C1为Neo-Hookean常数，D为体积弹性模量k的倒数，C2为指数应力系数，C3为胶原纤维拉直率，C4为拉直纤维的弹性模量，λ*为纤维被拉直后与原长的长度比。

(4)膝关节韧带初始应变：生物软组织材料在体内的连续生长、重建或破坏时， 往往具有一定的残余应力，因此为了获得精确膝关节的运动状态，设置膝关节各条韧带的初始应变见表3[6]。

表3　各条韧带的初始应变

3.3结果分析

图8　三级跳三个阶段运动过程中左膝关节承受地面垂直和水平载荷随时间变化示意图 (R/G表示地面作用载荷与人体体重的比值) Fig.8 MDA simulation results of the first, second and third phrases: the vertical and horizontal reactions from the ground (normalized with respect to the bodyweight, R/G)

图9　三级跳三个阶段运动过程中左膝关节受到的总的肌肉扭矩随时间变化示意图 Fig.9 MDA simulation results of the first, second and third phrases: the resultant muscle torques in the supporting left knee joints

由图9可以看出，在上述三级跳助跑、起跳、腾空和落地三个阶段运动过程中，人体左膝关节部位受到总的最大肌肉扭矩也出现在第三个阶段。这是由于人体在腾空、落地运动过程中，关节部位肌肉、韧带伸展、收缩程度加剧和加快造成的。此时如果运动员膝关节处韧带、肌肉承受不了这样大的载荷和伸展幅度，易引发该处韧带和肌肉出现疲劳运动损伤。事实证明，三级跳过程中韧带、肌肉拉伤约占整体损伤的40%和28%[7]。

由此从图8～9可以看出：在三级跳第三个阶段中的t=2.35s时，膝关节部位承受的冲击载荷最大，进一步对其加载进行有限元分析计算，以获得三级跳运动过程中膝关节处冲击、生物力学行为，为膝关节运动损伤的预防和治疗提供参考依据。

4膝关节有限元分析

4.1膝关节有限元模型的建立

大量研究表明，对于膝关节，无论是运动损伤还是关节累计负重使关节软骨退化所造成的骨关节病，其生物力学根源在于关节冲击运动造成的内部不均匀的应力分布[3]；临床统计分析进一步表明，胫骨应力骨折是膝关节骨损伤发病最高区域，约占膝关节整体损伤的57%[7-10]。因此，利用有限元分析方法，研究三级跳运动过程中胫骨部位的应力应变分布变化非常必要。

对膝关节进行有限元分析之前，首先将前述计算三级跳运动全过程的膝关节三维实体模型导入图像分割软件AMIRA中进行网格划分；然后利用有限元分析软件ANSYS，依据膝关节结构特点，建立其三维有限元模型。模型包括：股骨中下段，胫骨和腓骨中上段(采用SOLID185单元)，ACL、PCL、MCL和LCL(采用三维仅受拉杆单元BEAM188单元)；共生成31889个单元，17668个节点(见图10(a))。

图10　左膝关节 (a)有限元分析模型及(b)约束加载示意图 Fig.10 Finite element (a) and constraint & loaded (b) models of the total knee joint

根据膝关节运动实际状况，设置计算必须的约束和加载条件，其中关节间的连接考虑为有滑动的接触算法，见图10(b)；为了获得精确的计算结果，依据膝关节结构的差异性，对其网格的智能划分进行了一定的处理；在进行静力分析时，选用了FFE(Fast Finite Element) 方法，该方法的特点是：采用了迭代算法，并利用了数学中的稀疏矩阵及先进的数据管理技术，因而方程组收敛的速度较快，所得到的结果精度也较高。

4.2结果分析

由图11～14可以看出：在三级跳助跑、起跳、腾空和落地整个运动过程中，胫骨内外两侧均产生较为明显的变形(见图11)，并出现了向下的塌陷位移，最大位移为3.856mm(见图12)；胫骨外侧平台处等效应力最大，最大莫尔应力为385.12MPa(见图13)；XY方向上的最大剪切应力也分布在胫骨的外侧平台上，最大剪切应力为44.51MPa(见图14)。较大的压力、应力易促使胫骨部位出现磨损，诱使胫骨发生疲劳损伤。

为了验证本文计算方法的正确性，我们将预测分析的结果与国外权威期刊发表的研究结果进行了对比。

图11 胫骨关节变形状况示意图Fig.11Deformationpatternofarticularcartilage图12 胫骨关节整体位移示意图Fig.12Displacementpatternofarticularcartilage图13 胫骨关节莫尔应力示意图Fig.13VonMisesstressdistributiononarticularcartilage

图14 胫骨关节剪切应力示意图Fig.14Shearstressdistributiononarticularcartilage图15 左膝关节胫骨处磨损状况示意图Fig.15Theweardepthcontoursonthearticularcartilageoftheleftkneejoint图16 左膝关节胫骨处接触压力分布状况示意图Fig.15Contactpressurecontourcomparisononthearticularcartilageoftheleftkneejoint

图15～16是美国Denver大学Rullkoetter课题组对受试者经过500万步态训练后，采用数值模拟、实验测量方法获得的其左膝关节胫骨部位接触压力分布及磨损状况示意图[10-11]。

从图上可以看出，胫骨上受损最严重的、接触压力最大的区域恰恰正是本文计算预测的胫骨受力变形、发生位移及应力应变最大的部位。但需要说明的是，虽然图13～14的研究结果不是取自三级跳运动过程，但从田径训练统计结果来看，三级跳运动过程中胫骨处磨损是目前所有运动骨损伤之首，因而可以断定，三级跳运动膝关节胫骨部位出现的磨损损伤肯定比文献[10-11]更为严重。大量临床统计结果也表明，膝关节胫骨部位运动损伤和累计疲劳损伤是致使膝关节发生骨损伤的主要原因之一[1]。

5结论

本文针对三级跳运动冲击载荷对膝关节损伤极为严重，为所有运动损伤之首的事实，以一名身高175 cm、体重70 kg，20周岁男性健康志愿者为研究对象，首先利用CT、磁共振扫描技术和基于特征零件实体建模技术构建膝关节三维实体模型，然后引入机械动力学设计理念，运用多体动力学(MDA)和有限元分析技术(FEA)相结合的方法模拟三级跳运动全过程，并进一步完成膝关节骨损伤最严重部位-胫骨的有限元分析，结果表明：

(1)在三级跳运动助跑、起跳、腾空和落地三个阶段过程中，胫骨内外两侧均产生较为明显的变形，出现向下的塌陷位移，产生较大压力、应力应变，说明胫骨内外部位在膝关节运动过程中易发生磨损且磨损程度最为严重，是运动员出现骨损伤的关键。

(2)虽然本文所得到的计算结果与国外相关文献有些差别，但预测的膝关节骨损伤严重区域与国外相关研究结果较为吻合，与临床实际也较为一致，这说明利用多体动力学和有限元分析相结合的方法研究三级跳膝关节生物力学行为可行。

(3)多体动力学和有限元分析相结合的方法，与直接通过施加最大载荷的有限元方法相比，能更好地考虑膝关节的累积疲劳损伤效应，为人们理解膝关节实际损伤机理、膝关节生物运动力学行为提供了一种全新的研究方法，也为膝关节损伤治疗和康复训练、临床医学观察与手术方案的设计提供一定的参考依据。

参考文献

[1]宋德海.对我国优秀跳远、三级跳运动员的运动损伤流行病学调查与分析[J].四川体育科学,2011(6):47-50.

SONG De-hai. Epidemiological investigation and analysis on elite of China athletes of long jump and triple jump of sport injuries [J]. Sichuan Sport Sciences, 2011(6):47-50.

[2]吴迪.基于CT、MRI的人体膝关节三维有限元模型的建立及应用[D].大连:大连医科大学,2007.

[3]Allen S J. Optimization and performance in the triple jump using computer simulation[D]. Ph.D. thesis, Loughborough University, Loughborough, UK, 2009.

[4]郝智秀,冷慧杰,曲传咏,等.步态变化对人体膝关节接触生物力学特性的影响[J].固体力学学报,2010, 31(12):603-612.

HAO Zhi-xiu, LENG Hui-jie, QU Chuan-yong, et al. Biomechanics of the bone and the knee joint [J].Chinese Journal of Solid Mechanics, 2010, 31(12):603-612.

[5]Curtis N, Kupczik K, O’Higgins P. Predicting skull loading: applying multibody dynamics analysis to a Macaque skull[J]. Anat Rec,2008(291): 491-501.

[7]赵文凯，李庆华，巩磊.对三级跳远运动员运动损伤的调查及原因分析[J].贵州体育科技，2008(18):59-61.

ZHAO Wen-kai, LI Qin-hua, GONG Lei. Investigation and analysis on triple jump of sport injuries [J].Guozhou Sports Science and Technology，2008(18):59-61.

[8]Sarakkala S, Julkunen P, iviranta P, et al. Depth-wise progression of osteoarthritis in human articular cartilage: investigation of composition, structure and biomechanics[J]. Osteoarthritis and Cartilage, 2010(18):73-81.

[9]倪虹，甄永乾，汪延成，等，蝗虫起跳运动分析及仿生机构研究[J].振动与冲击，2014,33(2):7-13.

NI Hong, ZHEN Yong-qian, WANG Yan-cheng, et al. Bio-jumping mechanism inspired from locust take-off jumping [J]. Journal of Vibration and Shock, 2014,33(2):7-13.

[10]Halloran J P, Easley S K, Petrella A J, et al. Comparison of deformable and elastic foundation finite element simulations for predicting knee replacement mechanics[J]. Journal of Biomechanical Engineering, 2005(127):813-818.

[11]Halloran J P, Petrella A J,Rullkoetter P J. Explicit finite element modeling of total knee replacement mechanics [J]. Journal of Biomechanical Engineering, 2005(38):323-331.