人类骨小梁的非线性细观有限元分析
2013-06-21 by:广州有限元分析、培训中心-1CAE.COM 来源:仿真在线
ABAQUS 的主要功能和优点
建立模拟骨骼组织机械特性的本构模型
轻松求解几何的和材料的非线性模型
支持并行求解
背景知识
骨小梁位于长骨(如股骨)的末端和立方骨(如脊骨)中,是人类骨骼中承载生物组织的主要生物组织。它的机械特性有很高的临床价值和研究价值。增进对骨小梁机械特性的了解,有助于深入研究骨骼的断裂机理,也有助于评估年龄、疾病和药物治疗的影响。骨小梁是一个充满孔洞的组织——脊骨中 85% 以上是孔洞,并有着复杂的结构,而且这两者都因不同的人和不同的解剖部位而存在很大差异(参见图 1)。因此,要用统计学的方法确定骨小梁的机械特性,就需要多个样本的机械特性数据。
图1 人类脊骨圆柱型小梁样本的透视图
细观有限元法 (μFE) 广泛应用于骨小梁机械特性的研究,包括在光谱水平和微结构水平两个方面的研究。这些模型是通过对骨小梁样本进行高分辨率成像得到的,样本被自动划分成元素为六面体的有限元网格(参见图 2)。网格中所有的单元都完全一样,一般尺寸在 50 微米。划分网格后,一个边长 5 毫米的立方体样本的 μFE 模型一般具有 50 万个自由度。与试验用样本(8 毫米直径和 15 毫米长)类似的骨骼样本的 μFE 模型则有几百万个自由度。
图 2:含有 44μm 单元的骨骼样本中一块边长2.5 毫米立方体的 μFE 网格图
在过去,这些大量的问题会使许多研究人员不得不利用自定义代码,一个单元一个单元地迭代求解。由于非线性有限元模型非常复杂,这些自定义代码只限于线弹性分析。虽然线弹性有限元模型不能模拟骨骼受损情况,但是研究人员经常利用它与试验数据校核,确定骨骼组织的弹性特性。然而,关于骨小梁非线性机械特性的许多问题还有待解决.因为 ABAQUS/Standard 能够利用并行处理能力解决大型问题,包括复杂材料模型问题,所以它非常适合这类分析。在本技术简报中,我们利用 ABAQUS/Standard 研究了几何非线性在骨小梁机械特性中的作用。我们对一个具有四百多万个自由度的模型进行线弹性分析。通过检验此分析的并行处理能力(也就是可伸缩性),我们展示了它求解大型问题的可行性。
有限元分析方法
利用显微X线断层摄影技术 (μCT 20,Scanco Medical AG, Bassersdorf, Switzerland),以 22 微米的分辨率,对容积率为 9% 的人类脊椎骨小梁样本进行成像(参见图 1)。建立了两个 μFE 模型。首先,整个圆柱形样本被划分为大小44 微米的六面体单元网格(参见图2)。然后,从圆柱体中心划出一个边长为5毫米的立体子区,建立另一个具有相同单元大小的模型。两个模型的网格数量见表 1。
表 1:μFE 模型的网格数量
模型 |
单元数量 |
节点数 |
自由度数 |
圆柱体 |
828,853 |
1,380,834 |
4,142,502 |
立方体 |
131,322 |
216,027 |
648,081 |
圆柱体模型被用来评定直接稀疏求解器的并行处理能力。在没有摩擦的情况下,在顶面和底面应用位移边界条件,模拟 1% 压缩应变。分别用 HP rx8620 计算机中 1、2和4 个 CPU 进行线弹性分析。
利用边长为 5 毫米的立方体模型进行非线性分析。这个尺寸的立方体已大得足够确定平面特性,同时又小得足够确保非线性分析的可行性。骨骼组织模型是用铸铁塑性材料制造的。铸铁塑性材料在受到拉伸和压缩时,其弹塑性状态会有不同的屈服强度和硬化,因此会产生一个非对称的单元刚度矩阵。因此,需要使用非对称存储的并行稀疏直接求解器。样本是一个弹性模量为 13.4Gpa,泊松比为 0.3 的组织。(参见参考文献 2)根据人类股骨骨小梁组织的屈服应变(参见参考文献 3),铸铁塑性模型组织拉伸时的屈服应力为 55.2MPa,压缩时为 110.6MPa。在拉伸和压缩时,使用的是相当于弹性模量 5% 的硬化斜度。在无摩擦的位移边界条件下,拉伸和压缩采用了 2% 的公称应变。在这样低的公称应变条件下,骨骼微结构的自力接触可以忽略。此外,每个模型都进行了考虑和不考虑几何非线性变形的模拟。总共进行了四个非线性分析,为了进行比较,还计算了平面屈服应变。所有对立方体的分析都是在一台 IBM Power4 计算机上进行的,使用了两个CPU。
结果和结论
利用 4 个 CPU 对圆柱体模型进行线性分析,用时不到 16 分钟,占用内存不到 11 GB(参见表 2)。表 2 还包括了平行直接求解器的计数结果;加速因数是根据求解时间得到的。对具有几何非线性的立方体 μFE 模型进行非线性分析,用时不到 7.4 小时,占用内存 4.1 GB。每个非线性分析需要大约 100 个线性方程的解,这就强调了求解器可伸缩性的重要性。骨骼结构中初始屈服点的定位使得非线性分析的收敛变得更加具有挑战性(参见图 3)。
表 2:6.4-3 版直接稀疏求解器性能
CPU 数 | 求解时间(秒) | 加速 | 总时间(秒) |
1 | 554 | 1.00 | 1348 |
2 | 295 | 1.88 | 1107 |
4 | 171 | 3.24 | 945 |
图 3:在 2% 压缩应变情况下边长为 2.5 毫米立方体的μFE 模型的骨骼结构局部应力分布图
图 4 是根据表观应变(样本长度的变化/原始样本长度)作出的表观应力(外力/横截面面积(25 平方毫米))图。初始屈服点定义为偏移量达到 0.2% 的点。与试验数据(参见参考文献 4)类似,压缩时的屈服应变比拉伸时的大。
表 3:不同模型组合的屈服应变
几何非线性 | 拉伸 | 压缩 |
考虑 | 0.61 | 0.78 |
不考虑 | 0.59 | 0.86 |
图 4:四个非线性分析的应力应变关系图几何非线性在压缩时引起软化,拉伸时引起硬化标记显示由 0.2% 偏移量方法(点划线)确定的初始屈服点
虽然组织材料在硬化,但当考虑几何非线性时,很明显地观察到有软化发生(参见图 4)。另外屈服应变与试验测量结果类似,特别是压缩的情况下(参见参考文献 4)。这些结果表明在拉伸和压缩时骨小梁组织有着不同的屈服特性,应该把 μFE 模型和几何非线性结合起来,精确地模拟骨小梁的平面屈服特性。
致谢
ABAQUS, Inc. 衷心感谢加州大学伯克利分校的 Tony M. Keaveny 教授提供骨骼样本成像数据和有限元网格。
参考文献
van Rietbergen, B.; H. Weinans; R. Huiskes; A. Odgaard, “A New Method to Determine the Trabecular Bone Elastic Properties and Loading Using Micromechanical Finite Element Models,” Journal of Biomechanics, vol. 28, pp. 69–81, 1995.
Rho, J. Y.; T. Y. Tsui; G. M. Pharr, “Elastic Properties of Human Cortical and Trabecular Lamellar Bone Measured by Nanoindentation,” Biomaterials, vol. 18, pp. 1325–1330, 1997.
Bayraktar, H. H.; E. F. Morgan; G. L. Niebur; G. E. Morris; E. K. Wong; T. M. Keaveny, “Comparison of the Elastic and Yield Properties of Human Femoral Trabecular and Cortical Bone Tissue,” Journal of Biomechanics, vol. 37, pp. 27–35, 2004.
Morgan, E. F.; and T. M. Keaveny, “Dependence of Yield Strain of Human Trabecular Bone on Anatomic Site,” Journal of Biomechanics, vol. 34, pp. 569–577, 2001.
ABAQUS 参考资料
有关本简报中提到的 ABAQUS 功能的附加信息,请参见 ABAQUS V6.4 文档中的以下内容:
Analysis User’s Manual
“Static stress analysis”,第 6.2.2 节
“Parallel execution in ABAQUS/Standard”,第 7.18.1 节
“Cast iron plasticity”,第 11.2.10 节
相关标签搜索:人类骨小梁的非线性细观有限元分析 Fluent、CFX流体分析 HFSS电磁分析 Ansys培训 Abaqus培训 Autoform培训 有限元培训 Solidworks培训 UG模具培训 PROE培训 运动仿真