基于有限元方法的中耳听骨链结构运动机理研究

曾新宇，殷 亮，陈向东

深圳市宝安人民医院1耳鼻喉科,2放射科，广东 深圳 518101;3深圳大学医学部，广东 深圳 518060

基于有限元方法的中耳听骨链结构运动机理研究

曾新宇1，殷 亮2，陈向东3

深圳市宝安人民医院1耳鼻喉科,2放射科，广东 深圳 518101;3深圳大学医学部，广东 深圳 518060

目的运用中耳Micro-CT扫描数据进行三维重建，采用有限元方法对听骨链运动规律及鼓膜振动特性进行生物力学研究。方法通过对比鼓膜凸、镫骨足板振幅实验数据验证模型，并进行频率响应分析及模态分析。结果在不同频率下听骨链振动与转动均不一样，三块听小骨围绕一定的转动轴进行转动和摇摆运动，接近鼓膜凸处的锤骨柄做同相位转动和平动，镫骨足板做活塞平动。在频率1000 Hz下鼓膜整体弯曲变形局部高达2.32e-006 m，出现在环韧带附近；鼓膜凸最大变形约1.0e-007 m。鼓膜环韧带在声压激励下容易发生扭曲变形。鼓膜凸附近出现应力集中，最大应力约8.33e-004 MPa。结论听骨链运动机理研究对人耳生命科学研究和临床手术均有一定的理论指导意义。

中耳；听骨链；有限元方法；三维重建；鼓膜振动

耳聋是影响人们日常生活和社会活动的常见感官功能障碍疾病之一，其临床治疗主要通过手术如：鼓膜修补术和听骨链重建术等来重建鼓室，恢复患者听力。目前手术重建的听骨链结构尚不能与患者生理状态保持一致，难以获得理想的听力恢复。随着生物力学与医学影像技术的快速发展，运用声学、医学、力学、结构学相关理论知识对人耳结构进行研究显得非常重要。利用有限元方法可以较好地仿真模拟和了解外耳、中耳、内耳的传音机理和受力形态，对手术方式选择及术后听力恢复预测均有一定意义。

自1988年起，许多国外学者开始利用实验手段、有限元方法研究中耳结构的力学行为[1-12]。随后，国内学者也通过有限元方法对人耳结构展开了相关的仿真模拟与数值分析[13-21]。本研究运用中耳听骨链结构Micro-CT扫描数据进行三维重建，通过有限元软件对三维模型进行结构动力学分析，来探索研究和分析中耳听骨链结构在低频、中频、高频下不同振型的运动机理和频率响应动力荷载下的中耳、鼓膜结构的动力学特性，以期为临床手术提供参考。

1 资料与方法

1.1 研究对象

成人男性尸头2具（由深圳大学医学部提供），切取右侧颞骨组织，包括完整的中耳、内耳、乳突及外耳道，置入50～60 ℃烤箱烘烤48 h。

1.2 研究方法

将预处理标本置入Micro-CT（SkyScan1076，比利时）检查床进行高精度扫描，层厚35 μm（图1）。扫描数据导入Mimics16.0软件（Materialise，比利时）进行阀值分割，三维重建出听骨链结构的STL模型（图2）。将模型导入至Geomagic Studio 2014（Geomagic，美国）软件对其进行体、面修复与逆向模型调整（图3）。随后将修复后的模型导入Patran/Nastran2012软件（MSC，美国）中补充鼓膜（9 mm×8 mm×0.1 mm）、韧带（锤骨前韧带、锤骨侧韧带、锤骨上悬韧带、砧骨上悬韧带、砧骨后韧带、鼓膜环韧带、镫骨环韧带）、肌肉（鼓膜张肌、镫骨肌）等中耳软组织结构，并划分有限元网格（图4）。3块听小骨、各韧带及各肌腱共计7427个节点，32 366个四节点四面体（Tet4）Solid单元；鼓膜共计361个节点，360个四节点四边形（Quad4）2D Membrane单元。

图1 颞骨Micro-CT扫描图像

图2 CT三维重建听骨链模型

参考国外相关研究的中耳听骨链、韧带、肌腱等实验及数值模拟的相关材料属性参数[3,9,18,22-23]，本研究模型所设置的材料参数如表1所示，中耳各组织结构泊松比统一取0.3，中耳模型结构阻尼系数取0.5。

在鼓膜面上直接施加90 dB（0.632 Pa）面压力均布荷载，通过有限元方法来仿真模拟100～10000 Hz之间纯音声场激励下的中耳结构运动情况。在各肌腱、各韧带的末端设置刚性固定约束（其中鼓膜环韧带、镫骨足板外侧边缘设置铰支固定约束），鼓膜与锤骨柄之间的结构接触采用多点约束MPC处理方法，确保连接可靠与保证正常的受力传递。计算过程中不考虑内耳模型及其内部运动状态，仅采取在镫骨底板处施加一定的流体质量单元约束和黏弹性约束，来考虑内耳结构对中耳的边界影响作用。

图3 听骨链3D模型修复

图4 中耳有限元网格模型

2 结果

2.1 中耳有限元模型正确性验证

在鼓膜上施加90 dB声压均布荷载，进行频率响应动力学分析，来仿真模拟100～10000 Hz之间的中耳听骨链结构的动力学特性，计算分析后提取鼓膜（TM）凸、镫骨足板（SF）中心处的振幅（Displacement）数值与2004年Gan的实验数据进行对比（图5，6）。对比结果可以看出：本研究有限元仿真数据在鼓膜凸、镫骨足板中心处在100～10000 Hz间的的振动变形情况，从变化趋势、数值量级上基本均与Gan的实验数据保持一致，且基本在实验数据的上下限制范围内，可以说明本研究在模型建立、材料设置、单元设定、边界假设，计算分析等设置是正确的。

表1 中耳听骨链、各韧带、各肌腱的材料参数

图5 鼓膜凸振幅数据对比（FEA仿真与Gan2004实验）

图6 镫骨足板中心处振幅数据对比（FEA仿真与Gan2004实验）

2.2 听骨链结构运动特性分析

对已建模的中耳听骨链结构进行模态分析，共分析了100～10000 Hz（低频、中频、高频）下10个模态情况，由于高频模态影响不太明显，因此本研究取低频、中频前6个模态振型和1个高频模态振型进行研究听骨链运动机理。取出研究用的7个不同频率下的听骨链结构振型。由图7可知，听骨链结构振型在不同频率下的振动、转动的情况都不一样，大致可总结出听骨链三块听小骨（锤骨、砧骨、镫骨）结构是围绕一定的转动轴线进行转动和摇摆运动的。主要的转动轴为锤骨侧韧带附近处与砧骨后韧带端部处的连线、砧骨后韧带端部处于镫骨肌处的连线，以及锤骨柄处的鼓膜凸与鼓膜是做同相位转动和平动，镫骨足板处主要做垂直于镫骨足板面的平动，类似机械活塞运动形式。此外，各韧带、各肌腱、锤骨柄、锤骨头、砧骨短突、砧骨长突都有一定的摇摆转动。本研究听骨链大致的运动形态与转动轴的关系整理如图8。

2.3 鼓膜结构运动特性分析

通过对中耳听骨链结构在100～10 000 Hz下的频率响应动力学分析，可以发现，在800～1000 Hz鼓膜凸振动较为激烈。因此，本研究提取1000 Hz频率下的鼓膜整体弯曲变形云图、鼓膜整体扭转变形云图、鼓膜应力云图来着重分析鼓膜的受力与运动情况（图2～9）。由图9可知，在频率1000 Hz下的鼓膜整体弯曲变形最大值为2.32 e-006 m，出现在环韧带附近，以及鼓膜凸此时的最大变形约为1.0e-007 m。鼓膜整体扭转变形最大值为7.33 e-009 m，出现在鼓膜环韧带处，即鼓膜环韧带在声压激励下容易发生扭曲变形。鼓膜凸附近有应力集中，最大应力值为8.33e-004 MPa，主要原因可能是在锤骨柄与鼓膜的连接处在动荷载力的反复激励作用下，此处容易发生应力集中。

3 讨论

图7 7个不同频率下的听骨链结构振型

目前，国内外研究中耳振动分析的研究者不乏少数，例如国内学者刘迎曦[13-15]、姚文娟[16-19]、王振龙[21]，国外学者Gan[7-11]、Wada[2-3]等人，利用有限元方法对围绕中耳手术、中耳病变、传声机理、内耳声传导等研究较多，但系统分析中耳听骨链运动机理，详细研究中耳在低频、中频、高频各阶段的声传导及中耳三块听小骨在声压激励作用下的运动形式的较少报道。本研究借助于医学影像软件Mimics16.0和工业设计专业软件Geomagic Studio2014对中耳听骨链Micro-CT扫描数据进行三维模型重建，并修复模型缺陷部分，再通过MSC.Patran2012/Nastran2012有限元软件对处理后的中耳听骨链模型进行频率响应结构动力学计算分析，将有限元计算结果与鼓膜凸、镫骨足板中心处的实验数据对比，验证模型正确性。继而，利用正确的中耳听骨链结构有限元模型，来研究分析了中耳听骨链结构100～10 000 Hz不同频率下的听骨链结构振型运动机理和中耳、鼓膜结构的动力学特性。本研究基于生物力学原理和有限元方法的研究数据将会对耳科学研究和耳手术研究提供一定的临床理论指导作用。

图8 听骨链运动机理与转动轴的关系

图9 鼓膜结构受力运动情况

本研究有限元仿真数据在鼓膜凸、镫骨足板中心处在100～10000 Hz间的的振动变形情况，从变化趋势、数值量级上基本均与Gan 2004年的实验数据保持一致[7]。基本在实验数据的上下限制范围内，可以说明本研究在模型建立、材料设置、单元设定、边界假设，计算分析等情况是设置正确的，从而可证明本模型的正确性。

听骨链结构振型在不同频率下的振动、转动的情况都不一样，大致可总结出听骨链三块听小骨（锤骨、砧骨、镫骨）结构是围绕一定的转动轴线进行转动和摇摆运动的。主要的转动轴为锤骨侧韧带附近处与砧骨后韧带端部处的连线、砧骨后韧带端部处于镫骨肌处的连线，以及锤骨柄处的鼓膜凸与鼓膜是做同相位转动和平动，镫骨足板处主要做垂直于镫骨足板面的平动，类似机械活塞运动形式。此外，各韧带、各肌腱、锤骨柄、锤骨头、砧骨短突、砧骨长突都有一定的摇摆转动。

在频率1000 Hz下的鼓膜整体弯曲变形最大值为2.32e-006 m，出现在环韧带附近，以及鼓膜凸此时的最大变形约为1.0e-007 m。鼓膜整体扭转变形最大值为7.33 e-009 m，出现在鼓膜环韧带处，即鼓膜环韧带在声压激励下容易发生扭曲变形。鼓膜凸附近有应力集中，最大应力值为8.33e-004 MPa，主要原因可能是在锤骨柄与鼓膜的连接处在动荷载力的反复激励作用下，此处容易发生应力集中。

因人耳结构非常精小、中耳实验条件不足等情况，利用有限元方法来研究与阐述人中耳结构听小骨在不同频率下的运动力学激励，以达到结合临床手术与术后患者听力恢复调查情况，后期将对中耳听骨链置换手术、鼓膜穿孔手术所涉及的假体结构进行优化研发产品提供一定的生物力学理论指导，减少对大量动物实验及人尸体实验的依赖，紧密结合临床手术的经验及有限元理论分析方法，来指导中耳听骨链假体的研发和更好地进行同类问题患者的手术，从而能更好地在术前预测中耳手术后听力恢复情况。

[1]Lesser TH, Williams KR. The tympanic membrane in cross section：a finite element analysis[J]. J Laryngol Otol, 1988, 102(3)： 209-14.

[2]Wada H, Metoki T, Kobayashi T. Analysis of dynamic behavior of human middle ear using a finite-element method[J]. J Acoust Soc Am, 1992, 92(6)： 3157-68.

[3]Wada H, Koike T, Kobayashi T. Middle ear mechanics in research and otosurgery[C]//Proceedings of the International Workshop on Middle Ear Mechanics, 1996： 76-81.

[4]Koike T, Wada H, Kobayashi T. Modeling of the human middle ear using the finite-element method[J]. J Acoust Soc Am, 2002,111(3)： 1306-17.

[5]Ferris P, Prendergast PJ. Middle-ear dynamics before and after ossicular replacement[J]. J Biomech, 2000, 33(5)： 581-90.

[6]佚名. 标题为空[C]//Proceedings of the International Workshop on Middle Ear Mechanics, Dresden, Sept.19-22, 1996, 出版年度缺失：页码范围缺失.

[7]Gan RZ, Feng B, Sun Q. Three-dimensional finite element modeling of human ear for sound transmission[J]. Ann Biomed Eng, 2004, 32(6)： 847-59.

[8]Gan RZ, Sun Q, Feng B, et al. Acoustic-structural coupled finite element analysis for sound transmission in human ear--pressure distributions[J]. Med Eng Phys, 2006, 28(5)： 395-404.

[9]Gan RZ, Reeves BP, Wang X. Modeling of sound transmission from ear canal to cochlea[J]. Ann Biomed Eng, 2007, 35(12)：2180-95.

[10]Gan RZ, Wang X. Multifield coupled finite element analysis for sound transmission in otitis media with effusion[J]. J Acoust Soc Am, 2007, 122(6)： 3527-38.

[11]Gan RZ, Cheng T, Dai C, et al. Finite element modeling of sound transmission with perforations of tympanic membrane[J]. J Acoust Soc Am, 2009, 126(1)： 243-53.

[12]Kelly DJ, Prendergast PJ, Blayney AW. The effect of prosthesis design on vibration of the reconstructed ossicular chain： a comparative finite element analysis of four prostheses[J]. Otol Neurotol, 2003, 24(1)： 11-9.

[13]刘迎曦, 李 生, 孙秀珍. 人耳鼓膜病变数值分析[J]. 医用生物力学, 2008, 23(4)： 275-8.

[14]刘迎曦, 李 生, 孙秀珍. 人耳传声数值模型[J]. 力学学报, 2008,40(1)： 107-13.

[15]Liu YX, Li S, Sun XZ. Numerical analysis of ossicular chain lesion of human ear[J]. Acta Mech Sin, 2009, 25(2)： 241-7.

[16]姚文娟, 李 武, 付黎杰, 等. 中耳结构数值模拟及传导振动分析[J]. 系统仿真学报, 2009, 21(3)： 651-4.

[17]姚文娟, 李晓青, 李 武, 等. 中耳病变及人工镫骨形体研究[J]. 医用生物力学, 2009, 24(2)： 118-22.

[18]Yao W, Li B, Huang X, et al. Restoring hearing using total ossicular replacement prostheses--analysis of 3D finite element model[J].Acta Otolaryngol, 2012, 132(2)： 152-9.

[19]Yao W, Li B, Guo C, et al. Numerical simulation and study about partial ossicular replacement prostheses with several various materials on hearing restoration[J]. Advanc Mater Res, 2011,238(3)： 2899-904.

[20]刘后广, 塔 娜, 饶柱石. 悬浮振子对中耳声传播特性影响的数值研究[J]. 力学学报, 2010, 42(1)： 109-14.

[21]王振龙, 王学林, 胡于进, 等. 基于中耳与耳蜗集成有限元模型的耳声传递模拟[J]. 中国生物医学工程学报, 2011, 30(1)： 60-6.

[22]Herrmann G, Liebowitz H. Mechanics of bone fractures[M]. new york： academic press, 1972： 772-840.

[23]Kirikae I. The structure and function of the middle ear[M]. tokyo：[s.n.], 1960： 51-7.

Structure movement mechanism for auditory ossicle chain of human middle ear based on finite element method

ZENG Xinyu1, YIN Liang2, CHEN Xiangdong31Department of Otorhinolaryngology,2Department of Radiology, The People's Hospital of Bao’an district, ShengZhen 518101, China;3Health science center, Shenzhen University, ShengZhen 518060, China

ObjectiveTo reconstruct the middle ear structure in three-dimensional by Micro-CT data, and explore the motion rule of chain ossicles and tympanic membrane vibration characteristics in biomechanical theory research by using the finite element method.MethodsThe model was verified correctly by comparing finite element simulation data with experimental data about the amplitude of tympanic membrane umbo and stapes footplate. The conclusion was drawn by frequency response analysis and modal analysis.Results(1) The vibration and rotation of ossicular chain were different under various frequencies. Three ossicles performed rotating and rocking motion around axis of rotation. Manubrium mallei close to tympanic membrane umbo performed in same phase rotation and translation. Stapes footplate performed translation like the piston. (2) At the frequency of 1000 Hz, the total bending deformation of tympanic membrane ran up to 2.32e-006 m, which appeared in annular ligament of tympanic membrane. The maximum deformation of the membrane tympani umbo was about 1.0e-007 m. Under the excitation of sound pressure, the annular ligament of tympanic membrane was prone to distortion. The stress concentration occurred in the vicinity of membrane tympani umbo.The maximum stress was 8.33e-004 MPa.ConclusionMotion of ossicular chain will provide theoretical guidance for life science research and clinical operation of human ears.

middle ear; ossicular chain; finite element method; three-dimensional reconstruction; tympanic membrane vibration

2017-01-06

深圳市科技计划项目（JCYJ20140414111512394）

曾新宇，本科，副主任医师，E-mail： szzengxinyu@163.com