三维有限元头颅模型参数及边界条件研究.docx
《三维有限元头颅模型参数及边界条件研究.docx》由会员分享,可在线阅读,更多相关《三维有限元头颅模型参数及边界条件研究.docx(13页珍藏版)》请在冰豆网上搜索。
三维有限元头颅模型参数及边界条件研究
三维有限元头颅模型参数及边界条件研究
一
l94一
基础研究?
生物医学工程与II缶床2005年7月第9卷第4期BME&ClinMed,July2005.Vo1.9,No.4
三维有限元头颅模型参数及边界条件研究
何黎民,卢亦成,吴建国,刘平,丁祖泉,陈学强,王保华
摘要:
目的确定三维有限元模型材料阻尼参数,实现寰枕关节模拟.方法对不同材料阻尼参数下的模型进行冲
击载荷加载,分析比较计算后不同材料阻尼下颅骨应力一时间曲线形态,确定模型采用的材料阻尼参数值.应用弹簧
单元模拟寰枕关节,比较模拟前后,模型在冲击载荷下颅骨应力,确定弹簧单元参数.结果不同颅骨材料阻尼条件下,
额部冲击区域颅骨节点的yonMises应力曲线均无第2峰值,应力峰值随阻尼增大而逐渐降低,峰值出现时间后移,枕
部颅骨节点应力曲线在冲击后期应力下降趋势更显着,在0.001—0.004之间颅骨应力曲线形态最佳.模拟寰枕关节,
弹簧单元采用颅骨材料参数,在X,l,轴方向弹簧弹性系数10N/mm,Z轴方向的为20N/mm,颅骨应力曲线后期可下
降.结论颅骨线弹性材料阻尼系数对模型颅骨应力响应有显着影响.模型脑组织应力主要受颅骨应力的影响.采用三
维弹簧单元模拟寰枕关节可降低模型颅骨应力,模拟寰枕关节对模型颅骨应力的影响较颅骨材料阻尼系数的影响小.
关键词:
有限元模型;头颅;材料参数;边界条件
中图分类号:
R68文献标识码:
A文章编号:
1009—7090(2005)04—0194—04
Studyonparametricandboundaryconditionofthefmiteelement
Yi-chen~,wuJian—guo,LIUPing,DINGZu-quan2,CHENXue—qiang2,WANGBan—hua
Neurosurgery,ChangzhengHospital,Shanghai200003,China;2.BioscienceandMedicalEngineering
Univers,Shanghai200072,China)
Correspondingauthor:
.LUYi—cheng,E-maihlychengC~
HELi—min,LU
(J.Departmentof
Academy,Tonal
Abstract:
ObjectiveTostudymaterialdampingparametersandboundaryconditioninfiniteelementmodelofthehuman
head.MethodsInthethree—dimensionalhumanfiniteelementmodel,brain,skuHanddurawereassumedtobehomogeneous
andisotropicmaterialwithlinearelasticbehavior.Aseriesofstructuraldampingparametersofskullandbrainwereadoptedto
simulateafrontalheadimpact.ThecalculatedresultsofyonMisesstressofskullandbrainwereobserved.Using2Dspring
elements,theatlantooccipitaljointwassimulatedinthree—dimensionalspace.ResultsTherewerenosecondpeakintheyon
Misesstresscurveofskull,whenskulldampingparameterWasdefined.ThehigherdampingparameterWasadopted,thelower
peakvalueofskullyonMisesstresscurvewasobserved.Itispropertodefineskullstructuraldampingfactorbetween0.001and
0.004.AftertheatlantooccipitaljointWassimulated,withdifferentelasticcoefficients(10ofspringelementindifferent
directions.10N/mminX-Ydirectionand20N/mminZdirection,theyonMisesstresscurveofskullbecamemoresmoothly.
ConclusionsTheskulldampingparameterinfluencetheyonMisesstressofskullgreatly.Thechangeofbrmnstresscurves
associatewithskullstresscurvesclosely.
Keywords:
finiteelementmodel;humanhead;materialproperties;boundaryconditions
由于头颅三维有限元模型涉及的参数较多,在模
拟计算直接冲击造成的头颅结构应力变化前,需确定
其边界条件和模型参数.在模型参数研究方面,根据
既往颅脑生物力学研究结果…,制作模型各种结构的
材料均为各向同性的线弹性材料;为提高计算效率,
作者单位:
1.上海第二军医大学长征医院脑外科,上海2o0loo3;
2_上海同济大学生命科学与医学工程研究院,上海2O0072
收稿日期:
2004—10—14;修回日期:
20o5—0l一05
作者简介:
何黎民(男,1972一,河南荥阳人),博士,主治医师,主要
从事神经外科工作和颅脑冲击伤有限元模型研究.
基金项目:
上海市科技发展基金项目(0141l9o99)
通信联系人:
卢亦成,E—mail:
lycheng~
减小计算成本,首先构建头颅模型,重点研究线弹性
材料阻尼系数对应力响应的影响,以及如何实现对寰
枕关节的模拟.
1材料与方法
1.1仪器与软件
MX8000多排螺旋CT(Philips意大利),Vi—
sion_plus1.5T核磁共振(Simenz德国).Dimension650
图形工作站(Dell美国).ANSYSY有限元软件.
1.2构建头颅有限元模型
依据正常人男性国人横轴位CT,MRI图像,半自
竺堕兰堡皇I堕7月第9卷第4期BME&ClinMed,July2005,Vo1.9,No.4
动识别重建结构,对皮肤,颅骨等实体结构采用灰度
阈值法进行图像分割,用边界提取算法进行边界提
取,再用轮廓跟踪算法和B样条曲线拟合方法,提取
边界关键点并拟合边界轮廓.对其他形态复杂的结构
进行子结构划分,人工绘制边界,确定关键点,用自编
坐标提取程序(由上海同济大学生命科学与医学工程
研究院编写)提取关键点坐标,形成有限元代码命令
流文件,导人有限元软件,对边界进行空间线性拟合.
采用实体单元,壳单元,划分网络.筛选不同结构相邻
界面上最小距离接点组合,生成粘连节点命令流文
件,在前处理中完成节点粘连.结构的头颅三维有限
元模型结构包括:
皮肤,颅骨,硬脑膜(含大脑镰,小脑
的透露三维有限元模型),大脑半球,中脑,小脑和脑
干等.模型总线数804,总面数412,总体积数24,总
单元数98744,总节点数42781.
1.3模型材料基本参数
依据既往对头颅材料的研究结果【I】和既往模型材
料参数【2-51,采用制作皮肤,颅骨,硬脑膜脑组织材料
为各向同性线弹性材料,其参数有杨氏摸量(
MPa),泊松比(和密度(p/kg/m3)取值见表1.
表1模型线弹性材料参数
1.4用二维有限单元模拟寰枕关节
在模型枕骨大孔层面的枕骨上确定代表枕骨髁
关节面的节点,图1显示选择的节点,左右两侧各19
个,记录节点序号.选用弹簧单元的参数包括:
弹性常
数,阻尼系数c,弹簧单元材料与颅骨材料相同,但
无质量.输出选定枕骨髁关节面节点的坐标,在y,
Z,yZ三个坐标轴平面生成与这些节点对应的3组
新节点,共ll4个,在对应的每两个节点间定义二维
单元,单元类型定义为上述弹簧单元.生成的二维单
元与颅骨枕骨大孔的关系如图2所示.即使用3组二
维单元模拟寰枕关节.
1.5加载冲击载荷
采用Nahum『61尸体实验的载荷条件.选择额部18
个颅骨节点,在y,Z轴方向分别加载冲击力载荷,使
合力方向等同于Nahum实验受力方向,作用力方向
平行Frankfort解剖平面,指向枕部,所选18个节点
围成面积约754mm.冲击载荷曲线如图3所示.冲
击总时程15ms,分l1个载荷步.书写有限元分析软
件载荷命令流文件.采用瞬态小变形分析模式,控制
注:
粗点所示范围为枕骨髁关节面
图l枕骨大孔层面枕骨髁关节面节点
—.
195-
注:
线条表示二维弹簧单元
图2模拟寰枕关节的单元与枕骨大孔关系
最大子载荷步数为5,选择稀疏矩阵算法进行求解.
Z
\
时间/ms
图3Nahum实验冲击载荷时间历程
2结果
O
2.1颅骨材料阻尼系数对应力响应的影响
图4显示不同颅骨材料阻尼条件下,模拟计算获
得的额部冲击区域颅骨节点的vonMises应力一时问
曲线.阻尼为O时,模拟此处颅骨应力峰值与冲击峰
值出现时间相似,应力峰值高,冲击后期出现第2峰.
设定阻尼参数后,颅骨应力峰值随阻尼系数增大而逐
渐降低,峰值出现时问后移,冲击后期无第2峰.阻
尼系数为0.004时,颅骨应力在曲线末期的值最小.
在相同载荷下,模拟计算对冲部位(枕部)颅骨节点yon
一
l96-塑垦_T程与Il缶床2005年7月第9卷第4期BblE&ClinMed,July2005,Vo1.9,No.4
0.1250.3750.6250.8751.1251.3751.625
时间/s
0.1250.3750.6250.8751.1251.3751.625
时间/s
c=0.0o8
0.1250.3750.6250.8751.1251.3751.625
时间/s
x10-2
x10-2
x10-2
8
7.2
6.4
5.6
4.8
4
3.2
2.4
1.6
0.8
0
0.1250.3750.6250.8751.1251.3751.625
时间/s
0.1250.3750.6250.8751.1251.3751.625
时间/s
c=0.01
0.1250.3750.6250.8751.1251.3751.625
时间/s
图4不同颅骨材料阻尼时额部颅骨节点yonMises应力一时间曲线
Mises应力,显示材料阻尼对颅骨应力响应有类似影
响.
模拟寰枕关节几弹簧单元参数对颅骨应力的影响
图5显示采用不同方法模拟寰枕关节时,颅骨对
冲击载荷的应力一时间曲线(A—D).A完全约束枕骨
髁平面;B仅用Z轴方向弹簧模拟寰枕关节,弹簧单
元K=10N/mm;C采用,y,Z三个轴向的弹簧单
元模拟寰枕关节,,y轴方向弹簧单元K=10N/
mm,Z轴方向弹簧单元K=20N/mm;D在C参数基
础上,增加颅骨阻尼参数c=0.00l.A的颅骨应力一时
间曲线峰值最高,末期下降不佳,存在第2峰值;B的
曲线峰值基本不变,曲线末期第2峰值较A的有改
善;C与B结果基本相同;D的颅骨应力一时间曲线有
显着改善.模拟计算结果还显示,单独设定弹簧单元
的阻尼系数或弹簧单元材料阻尼应力一时间曲线影
响情况相似,当K>500N/mm时,颅骨应力一时间曲
x10-2
x10-2
X10-2
02505OO750l0O0l2501500
时间/s
注:
A完全约束,固定枕骨髁关节面,B仅用Z轴弹簧单元模拟
寰枕关节,K=20N/mm,C用,y,Z三轴向弹簧单元模拟寰枕
关节为20N/mm,10N/mm,D在C条件下增加定义颅骨材料
阻尼c=0.001
图5同一颅骨节点在不同边界条件下的应力曲线
96308642086420
∞\∞0∞一L10∞\∞0∞一L10
∞勰加840
Bd\∞0∞差L10
8642086400
∞\∞0∞一L10
09876543210
Bd\目∞.∞羞LI.Bd\臼∞.∞差LI.
Bd\∞0∞羞L10
生物医学工程与l晦床2005年7月第9卷第4期BME&ClinMed,July2005,Vo1.9,No.4
线与边界完全约束的应力一时间曲线趋于一致.
3讨论
早期材料实验证实颅骨等生物组织具有均质线
弹性材料特征【lj,这种材料特性在模拟计算中较易实
现,但材料阻尼参数对冲击载荷下应力响应的计算结
果存在显着影响,主要表现为冲击末期结构应力存在
"震荡"表现171,应力一时间曲线出现第2,第3峰值.此
计算结果与头颅尸体冲击实验测得的颅内压曲线不
符合,Shuck的材料力学研究显示脑组织具有黏弹性
特征,因此为实现应用线弹性材料特征进行模拟,需
要设定适当的材料阻尼参数.Shuck等人测定的脑组
织阻尼系数为0.0067—0.0087,Churn与Huangl8~模型
采用的颅骨,脑阻尼系数为0.001,本实验应用线弹性
材料特征定义头颅结构时也发现上述应力响应现象,
因此分别对颅骨,脑组织线弹性材料阻尼系数进行计
算研究.结果显示颅骨阻尼系数的设定可降低颅骨
yonMises应力峰值,阻尼系数为0.004时,冲击部位
下颅骨应力曲线与冲击载荷曲线的形态最相似.不同
阻尼条件下,对冲部位颅骨应力曲线也有上述表现,
但曲线较冲击部位在时间上延迟,考虑这种情况与冲
击在颅骨上传递耗时相关.
头颅在受到冲击时,通过颈限制头颅运动,所以
要考虑这种边界条件.因为实验冲击载荷总作用时间
为10ms,在这一短暂时间内虽可完成单突触反射,
但因为头颅此时整体位移很小,颈肩部骨骼肌及其附
属结构的感受器尚未完全被激活,故颈肩肌肉对头颅
的约束尚不充分,所以在冲击时程中颈部对头颅的约
束作用主要通过寰枕关节韧带完成.因此,采用,
y,Z轴三个方向的3组二维弹簧单元与枕骨髁关节
面相连,模拟整个颈部,主要是颅颈交界的韧带和颈
椎的韧带对头颅的约束.正常情况下,寰枕关节在水
平面上可发生的一定位移,加之寰枕关节及各颈椎关
节突关节在水平面上的位移,使头颅易于在水平面上
转动.水平面即模型坐标的,y轴平面.在纵向上,
即模型坐标的Z轴方向,寰枕关节韧带和颈椎间的
纵向韧带不发生显着位移,较强地限制头颅在纵向上
的运动.所以我们对水平方向和纵向弹簧采用了不同
的值,,y轴方向弹簧单元的较小,Z轴方向的
较大,这将更好地模拟应力在颅内的分布.
计算结果显示完全约束(枕骨髁完全固定)的边
界条件下颅骨应力一时间曲线峰值最高,存在第2峰
值,末期曲线下降不佳.应用弹簧单元模拟寰枕关节
后,颅骨应力一时间曲线峰值下降,第2峰及曲线末
期改善.单用Z轴方向弹簧单元模拟寰枕关节与应
用三个轴.向弹簧单元模拟关节时颅骨应力一时间曲
—.
197-
线有差异,但无统计学意义,这可能与实验冲击载荷
作用于额部中线上,且冲击力方向平行于头颅矢状面
有关,当载荷作用部位与冲击力方向变化后,应力三
个轴向弹簧单元模拟寰枕关节将更有效.从寰枕关节
解剖角度分析,三维弹簧单元模拟更接近解剖实际.
单独设定弹簧单元阻尼或弹簧单元材料阻尼对颅骨
应力无明显影响,表明此两个参数基本不影响颅骨的
应力响应,相比较而言,颅骨材料的阻尼系数在相同
条件下对颅骨应力一时间曲线影响最显着,应力曲线
峰值显着下降,曲线无第2峰值,末期曲线下降良好,
与冲击载荷曲线形态基本一致.弹簧单元值在一
定范围内(10一l00N/mm)对颅骨应力曲线的影响较
相似,当大于500N/mm后,颅骨应力一时间曲线
与边界完全约束时曲线趋于一致,表明模拟弹簧硬
度不宜过大,否则起不到模拟关节作用.既往部分
有限元模型在冲击中假设头颈边界为自由边界【8?
91,这
对于时程较短的冲击而言,与不模拟寰枕关节所计算
的应力相似,但我们认为,为更接近实际地模拟冲击,
扩大模型模拟适用范围,应当模拟寰枕关节.
参考文献:
【1】RuanJS,KhalilT,KingAI.Humanheaddynamicresponseto
sideimpactbyfiniteelementmodelingU】.JBiomechEng,
1991,113(3):
276-283.
【2】KumaresanS,RadhakrishnanS.Importanceofpartitioningmem?
branesofthebrainandtheinfluenceoftheneckinheadin.
jurymodeling[J].MedBiolEngComput,1996,34
(1):
27-32.
【3】MesSS,ThibaultKLInfantskullandsutttreproperties:
mea?
surementsandimplicationsformechanismsofpediatricbrain
injury[J].JBiomechEng,2000,122(4):
364-371.
【4】ZhangLY,HardyW,OmoriK,eta1.Rencentadvancesinbrain
injuryresearch:
anewmodelandnewexperimentaldata[J】.
BioengineeringConference,2001,50:
833—834.
【5】KleivenS.Influenceofimpactdirectiononthehumanheadin
predictionofsubduralhematoma[J].JNeurotrauma,2003,20
(4):
365—379.
【6】NahumAM,SmithR,WardCC.Intracranialpressuredynam?
icsduringheadimpact【C】.ProceedingsofStappCarCrash
Conference,1977.339—366.
【7】ChuCS,LinMS,HuangHM,eta1.Finiteelementanalysisof
cerebralcontusion[J].JBiomech,1994,27
(2):
187-194.
【8】HuangHM,LeeMC,ChiuWT,eta1.Three—dimensionalfinite
elementanalysisofsubduralhematoma[J].JTrauma,1999,47
(3):
538-544.
【9】WillingerR,KangHS,DiawB.Three-dimensionalhumanhead
finite——elementmodelvalidationagainsttwoexperimentalim?
pacts[J].AnnBiomedEng,1999,27(3):
403--410.