逆动力学步态仿真髌股关节生物力学有限元分析∗

唐班鸿，阿依古丽·喀斯木，乌日开西·艾依提†，周文正

(1. 新疆大学机械工程学院，新疆乌鲁木齐 830017；2. 新疆维吾尔自治区人民医院骨科中心，新疆乌鲁木齐 830001)

0 引言

膝关节负重大且运动量大，结构复杂且是人体中伤病率最高的关节之一，针对膝关节的生物力学研究一直都是人体骨科生物力学研究的重点[1]．膝前痛是常见的膝关节损伤表现之一，又被称为髌股关节痛，患者的主要症状表现有在行走、下蹲和跳跃等状态中出现膝前疼痛，影响正常运动状态．髌骨脱位是引发膝前痛的常见诱因之一，发病率在所有膝关节损伤中占3%左右，在18岁以下人群和过量运动人群中发病率高[2-3]．髌骨脱位会加剧髌股关节面的损伤并使关节功能退变，并可能进一步导致髌股关节炎，影响患者的日常行动和生存质量[4]．髌骨脱位患者由于疼痛的影响在行走中会出现代偿性步态，表现出与正常步态的差异．对步行步态进行研究，能够精确拆分人体在步行过程中所发生的生物力学行为，从而为临床治疗人体发生病理性行走步态功能障碍提供针对性的参考依据并制定治疗方案[5]．骨骼肌肉多体动力学建模技术和逆向动力学仿真技术已经在多项步态研究中应用，Yamagata等[6]研究了各肌肉力量减少对步行步态中膝关节内、外侧接触力的影响；Asayama等[7]研究了不同深蹲动作对膝关节内侧接触力、股四头肌肌力和臀大肌肌力的影响情况；陈亚东等[8]研究了跳跃着陆动作中不同膝关节屈曲角度下人体膝关节的相关力学性能[9-11]．本研究使用逆动力学骨肌人体模型仿真输出用于驱动髌股关节运动的关节正位力和股四头肌肌力并作为有限元分析的边界条件，建立健康髌股关节、Ⅰ度髌骨脱位和Ⅱ度髌骨脱位的髌股关节有限元模型，进行步行步态支撑期内4种不同屈膝角度的仿真分析，研究髌股关节软骨间的应力分布和接触面积特征，为膝前痛患者出现代偿性步态的原因提供一定的生物力学数据支撑以及为临床指导膝前痛患者进行正确步态行走提供依据．

1 髌股关节三维数据重建与有限元模型构建

1.1 髌股关节三维模型重建

由合作单位提供了2例CT数据用于本研究，将Dicom格式数据导入医学图像数据处理软件中进行三维重建，得到股骨、髌骨、胫骨和腓骨结构，生成三维模型并以STL格式输出备用．将各部分STL格式三维模型导入曲面处理软件中转化为实体模型，以IGS格式输出备用，再根据CT图像数据上的位置、髌-股骨实体模型和解剖结构对髌骨软骨、股骨软骨、胫骨软骨和半月板进行建模．重复对骨性结构的操作优化流程，得到软骨结构的实体模型，以IGS格式输出备用．

1.2 实体模型装配及网格划分

将得到的各IGS格式实体模型导入三维建模软件中并保存为零件格式，利用模型原始数据螺旋CT的容积效应在装配体中按原点进行组装，组装完成后将装配体导入到有限元分析前处理软件中对模型进行网格划分和材料属性定义．在网格划分中，有限元计算理论上越密集的网格将减少连续体离散产生的误差．而实际中，当网格划分小到一定程度时，求解的结果精度提升并不明显，在5%内浮动，此时的网格密度是较为合理的．在本研究中对于髌股关节有限元模型，由于骨骼-软骨组织形态结构不规则，且几何尺寸较小，不能简单地以缩小网格尺寸来提高精确性，而是遵循了以下原则：首先要保证髌股关节形态的真实性，其次对于本研究关注的重点即髌股关节间软骨上的生物力学行为，针对软骨结构的网格划分相比骨骼结构网格更加精细，模型均采用四面体网格划分．

对于髌股关节内的韧带结构，利用建模功能进行绘制，韧带结构包括髌韧带、前交叉韧带、后交叉韧带、内侧副韧带、外侧副韧带、内侧髌股韧带和外侧支持带，根据韧带结构的复杂程度以及作用功效，设置一到多条LINK180单元来模拟，韧带的起始点与三维重建位置是通过解剖结构并在专业医生的指导下进行的．

本研究主要关注关节软骨间的生物力学行为，骨骼结构的形变不在研究范围之内，因此在材料属性中将骨骼结构设置为均匀各向同性材料．关节软骨的材料特性本应属于粘弹性组织，但因为本研究中一个行走步态周期（<2 s）的时间相比软骨粘弹性时间常数（1 500 s）较小，且有相关研究结果表明对于短期加载的软骨，不管是设置为弹性材料还是粘弹性材料，其力学响应没有明显变化且接近弹性，因此本文将关节软骨设置为线弹性各向同性材料，材料属性见表1．最终构建出的有限元模型见图1．

图1 髌股关节有限元分析模型

表1 髌股关节材料属性[12-14]

1.3 有限元模型验证

对模型的有效性进行验证，在股骨顶端的截面上对有限元模型施加1 150 N的力，并在股骨内外髁中点连线的中点处施加134 N的股骨后向推力．在屈曲0˚时，本研究所建模型得到的胫股关节前后、内外方向上的相对位移分别为4.36 mm、0.76 mm．在相同载荷条件下，Song等[15]得到的结果分别为4.3 mm、0.39 mm；Pe˜na等[16]得到的结果分别为4.75 mm、0.56 mm，本研究所建模型与前人研究结果相近，可以进行下一步的模拟分析．

2 骨骼-肌肉人体步态模型逆向动力学仿真

2.1 步态周期划分与屈膝角度选择

在行走步态研究中通常将一个完整的步态周期划分为支撑期与摆动期，又称着地相与摆动相，其中支撑期约占整个步态周期的60%，是下肢主要施力过程，可细分为支撑反应期、支撑早期、支撑中期和支撑末期四个阶段[17]，步态周期划分见图2．对于髌骨脱位的膝前痛患者，正常的步态周期活动范围和大角度屈曲的关节运动受到限制，同时步行步态过程中频繁的小角度屈膝变化成为髌股关节间软骨磨损加重的原因．本研究根据骨骼肌肉人体步态模型逆向动力学仿真输出的健康人体步态周期与屈膝角度变化数据，得知小角度屈膝变化主要集中在步态支撑期内，关节正位力与肌肉力量的频繁变化也主要集中在支撑期内且膝前痛患者疼痛阶段也主要存在于支撑期，因此综合仿真数据曲线和日常行走步态中屈膝角度的变化范围，选取支撑期内的关节受力峰值瞬时角度20˚、26˚和运动过程节点10˚、30˚作为研究对象，研究膝前痛患者髌股关节内部的生物力学行为变化．

图2 人体步态周期划分

2.2 行走步态模型逆向动力学仿真

本研究使用AnyBody人体骨肌建模仿真系统建立肌肉骨骼系统力学模型，包含骨关节与肌肉系统，调用AnyBody模型库中的行走步态功能模型并取两位受试者的身高、体质量平均值170 cm、66 kg作为人体模型基础参数，以右腿为研究对象，缩放基础模型并载入．为了使计算结果更符合受试者的生物力学表现特征，研究中把通过CT数据重建出的骨骼实体模型导入AnyBody仿真模型中替换原有模型，使用自体骨进行后续的逆动力学步态仿真计算．在伸膝屈膝动作中，股四头肌是主要发力的肌肉，是实现伸膝屈膝的原动肌，其中股直肌是股四头肌中唯一的双关节肌结构，是直接传力结构，因此在模型输出的有限元分析边界条件中进行筛选，选取股直肌肌力和膝关节正位力作为载荷边界条件，输出结果曲线见图3．摘取要研究的屈膝角度关键帧，输出屈膝10˚、20˚、26˚和30˚所对应的载荷条件（表2）．

图3 逆动力学仿真步态周期内关节正位力和股直肌肌力

表2 屈膝角度对应的载荷

3 有限元分析结果与讨论

3.1 有限元分析结果

根据膝关节解剖结构对有限元模型的边界条件和施加载荷进行定义，在有限元分析软件ABAQUS中的相互作用模块设置股骨与股骨软骨、髌骨与髌骨软骨为绑定接触，设置股骨软骨与髌骨软骨切向为罚函数法，摩擦系数0.002，法向硬接触．根据髌股关节正常生理运动状态，通过所建模型的韧带结构和外加边界条件约束髌股关节各个旋转自由度和冠状面、水平面的平动自由度，保留矢状面的平动自由度，对不同髌骨脱位程度、不同屈膝角度的模型设置对应的载荷，记录并研究髌股关节面上的Mises应力分布情况和接触面积变化情况．通过计算提取各组模型的应力分布云图和接触面积值，并进行趋势变化分析，髌股关节间接触应力分布见图4、图5，双侧软骨面接触应力峰值结果见表3、表4，接触面积变化结果见表5．

图4 髌股关节股骨侧应力分布

图5 髌股关节髌骨侧应力分布

表3 髌股关节股骨侧峰值应力/MPa

表4 髌股关节髌骨侧峰值应力/MPa

表5 髌股关节接触面积/mm2

3.2 讨论

髌股关节是人体许多关节中的承重关节之一，是膝关节的重要组成部分，组织结构丰富且生物力学作用环境复杂，具有独特的运动轨迹，难以直接测量，在运动过程中时刻承受着变化的压应力和张力，髌股关节中任意微小的变化都有可能对髌股关节的生物力学行为产生较大影响．本研究建立了健康髌股关节、Ⅰ度髌骨脱位和Ⅱ度髌骨脱位的髌股关节有限元模型并进行步行步态支撑期内4种不同屈膝角度的仿真分析，获得了三类髌股关节在步行步态支撑期膝关节屈曲10˚、20˚、26˚和30˚时刻下的髌股关节间接触应力分布和接触面积情况并进行对比分析．

3.2.1 健康髌股关节模型生物力学表现分析

结果表明，健康髌股关节在四种情况下的步态支撑期屈膝20˚时取得应力最大值，股骨侧为2.36 MPa，髌骨侧为3.93 MPa，这是由于屈膝20˚时刻是步态支撑期内关节力与肌肉力的峰值角度，在该情况人体处于蹬地向前的发力期；此外应力分布在软骨面上较为均匀，关节软骨间接触面积随着屈膝角度的增大而增大．可以看出健康髌股关节在步态支撑期内的应力变化并不是简单的随着屈膝角度的增大而增大，这证明了髌股关节间的应力取决于行走过程中的多因素共同作用，是综合关节正位力、肌肉力与关节间接触面积共同影响下的结果．本研究利用逆动力学骨骼肌肉模型步态仿真技术输出人体步态支撑期内不同屈膝角度时刻的对应生物力学数据来取代用同一载荷边界条件对模型进行有限元分析，得到更有针对性的分析结果，使模拟更符合真实情况，促进生物力学有限元分析的应用．

3.2.2 髌骨脱位关节模型生物力学表现分析

对比Ⅰ度髌股脱位关节、Ⅱ度髌股脱位关节与健康髌股关节的生物力学表现，可以观察到髌骨脱位模型的应力峰值明显高于健康关节，应力分布面积明显减小且主要分布在软骨外侧面，容易出现应力集中，Ⅱ度髌股脱位关节的应力集中更为严重．

结合双侧关节软骨应力峰值结果并从图6和图7可以观察到，相比健康髌股关节，屈膝10˚时Ⅰ度髌股脱位关节的关节间应力在股骨侧增加了255%，在髌骨侧增加了178%；Ⅱ度髌股脱位关节的关节间应力在股骨侧增加了282%，在髌骨侧增加了222%；屈膝20˚时Ⅰ度髌股脱位关节的关节间应力在股骨侧增加了217%，在髌骨侧增加了127%；Ⅱ度髌股脱位关节的关节间应力在股骨侧增加了325%，在髌骨侧增加了215%；屈膝26˚时Ⅰ度髌股脱位关节的关节间应力在股骨侧增加了376%，在髌骨侧增加了256%；Ⅱ度髌股脱位关节的关节间应力在股骨侧增加了828%，在髌骨侧增加了511%；屈膝30˚时Ⅰ度髌股脱位关节的关节间应力在股骨侧增加了354%，在髌骨侧增加了293%；Ⅱ度髌股脱位关节的关节间应力在股骨侧增加了647%，在髌骨侧增加了453%．

图6 髌股关节股骨侧Mises应力峰值

图7 髌股关节髌骨侧Mises应力峰值

在接触面积上，如图8所示，相比健康髌股关节，屈膝10˚时Ⅰ度髌股脱位关节间的接触面积是健康髌股关节的25.02%，Ⅱ度髌股脱位关节间的接触面积是健康髌股关节的23.22%；屈膝20˚时Ⅰ度髌股脱位关节间的接触面积是健康髌股关节的27.30%，Ⅱ度髌股脱位关节间的接触面积是健康髌股关节的20.99%；屈膝26˚时Ⅰ度髌股脱位关节间的接触面积是健康髌股关节的26.97%，Ⅱ度髌股脱位关节间的接触面积是健康髌股关节的22.52%；屈膝30˚时Ⅰ度髌股脱位关节间的接触面积是健康髌股关节的29.12%，Ⅱ度髌股脱位关节间的接触面积是健康髌股关节的25.41%．

图8 髌股关节间接触面积变化

有研究表明，若将髌骨向外侧移动5 mm就会使髌股关节的接触面积明显减小[18]．两种髌骨脱位模型的关节间接触面积同样随着屈曲角度的增大而增大，但相比健康髌股关节间接触面积明显减小；从图9和图10可以观察到两种髌骨脱位模型的关节间应力最大值时刻均出现在支撑期屈膝30˚时，这是由于正常情况下，髌骨在膝关节屈曲20˚时开始进入股骨滑车，但异常的髌骨位置导致髌骨无法正常进入滑车从而在关节面间出现异常接触且接触面积过小，导致关节间应力骤增．分析结果中过大的髌股关节应力值说明，由于髌骨脱位出现膝前痛的患者如果以正常行走状态发力，患者的髌股关节间应力将达到健康人体的一倍甚至数倍，且Ⅱ度髌股脱位关节模型的峰值应力在屈膝20˚后均大于10 MPa．在初期患者未有痛觉感受时，以正常状态行走过程中过大的髌股关节接触应力增加了关节间软骨的摩擦，有研究表明在10∼15 MPa的应力水平下软骨表层细胞会发生死亡[19]，异常的软骨间接触面也会导致软骨变性、病损．疼痛适应理论认为疼痛会改变肌肉的力学行为，从而调整人体的动作模式[20]．在患者出现痛觉后，由于疼痛导致不能以正常状态发力行走，为减轻痛感而出现代偿性步态，多表现为膝关节屈曲角度范围缩小，大角度屈曲动作受限．此时支持行走运动的肌肉无法以正常肌力运动，不利于肌肉功能的维持，将引起肌力减退甚至肌肉功能性萎缩，异常的髌股关节应力和接触面引起关节退变，最终造成髌股关节炎．

图9 髌股关节股骨侧Mises应力峰值变化

图10 髌股关节髌骨侧Mises应力峰值变化

有限元模型在人体关节生物力学的研究过程中能够提供相对准确的模拟分析，是一种有效的非侵入式生物力学研究方法[21-22]，但由于人体关节的复杂性，且缺少对有限元模型的统一评价标准，使仿真模型与实体之间仍有一定差异．本研究中将膝关节中的韧带简化为赋予材料的杆单元形式，并在加载过程中限制了髌骨及其软骨的旋转自由度，在一定程度上影响了模拟分析与真实情况的匹配程度，减少了计算量的同时可能对模型的精度带来了一定程度的影响，在后续研究中将优化计算逻辑进一步进行研究．

4 结论

通过建立健康髌股关节、Ⅰ度髌骨脱位和Ⅱ度髌骨脱位的髌股关节有限元模型，运用逆动力学骨骼肌肉模型步态仿真技术输出生物力学载荷条件并进行步行步态支撑期内4种不同屈膝角度的仿真分析，得到以下结论：

（1）步态过程中Ⅰ度髌骨脱位的关节间应力峰值在髌骨侧为10.56 MPa，股骨侧为8.67 MPa，相比健康髌股关节模型增加了293%和354%，Ⅱ度髌骨脱位的关节间应力峰值在髌骨侧为14.88 MPa，股骨侧为14.27 MPa，相比健康髌股关节模型增加了453%和647%，接触面积则小于健康髌股关节模型，在取得应力峰值的角度下Ⅰ度髌骨脱位和Ⅱ度髌骨脱位模型的关节间接触面积分别是健康关节的29.12%和25.41%，趋势与临床分析结果一致．健康髌股关节模型的应力峰值在屈膝20˚时取得，该情况为步态支撑期内关节力与肌肉力的峰值角度；由于异常的髌骨位置导致关节间接触面异常，过高的关节间应力改变了关节间生物力学环境使两种髌骨脱位模型的应力峰值角度在屈膝30˚时出现，同时增加了关节软骨病变的可能性．

（2）根据本研究的结果可知膝前痛患者出现代偿性步态的原因是由于过高的髌股关节应力导致关节软骨磨损加剧，异常的关节间接触面导致应力集中，磨损后的关节出现疼痛导致患者无法以正常行走状态发力．研究结果为膝前痛患者出现代偿性步态提供了一定的生物力学数据支撑，为临床指导膝前痛患者进行正确步态行走以及指导康复训练提供了理论依据．

（3）由研究结果可以说明恢复髌股关节间正常的应力分布和接触面积是治疗髌骨脱位膝前痛的最终目标，可通过外加康复带、肌肉锻炼康复和手术治疗等方式实现．