基于Block-Controlled網(wǎng)格劃分法的兒童胸部建模及驗證*

曹立波，蔣彬輝，冒浩杰，楊金海

(1.湖南大學(xué)，汽車車身先進設(shè)計制造國家重點實驗室，長沙 410082; 2.韋恩州立大學(xué)生物工程中心，美國密西根州 底特律 48201)

前言

在美國，每年約有1 500個16歲以下的兒童死于交通事故中[1]，其中兒童胸部損傷是導(dǎo)致兒童死亡的重要原因。有限元仿真是研究人體損傷的重要手段之一。然而，受到倫理道德的限制，可用于兒童建模的數(shù)據(jù)難以獲取，因此兒童模型的開發(fā)主要采用在成人模型基礎(chǔ)上按比例縮放的方法。文獻[2]～文獻[5]中分別通過縮放50百分位的成人模型，開發(fā)了1.5、3、6、9和15歲多剛體模型，但多剛體模型只能用于動態(tài)響應(yīng)的研究，無法進行應(yīng)力分布和損傷機理的研究。文獻[6]中開發(fā)了一個6歲兒童的有限元模型，其下肢模型的幾何數(shù)據(jù)通過6歲兒童志愿者MRI掃描獲得，而上身幾何數(shù)據(jù)未說明來源。文獻[7]中通過縮放THUMS 50百分位的成人有限元模型，獲得了3歲兒童有限元模型，該模型的幾何數(shù)據(jù)、材料參數(shù)和驗證數(shù)據(jù)均是通過比例縮放獲得的，難以反映真實的兒童碰撞響應(yīng)情況。

兒童并非縮小的成人，如嬰兒胸部截面為圓形，隨著年齡的增長，其橫向半徑逐漸增大，最終成為橢圓形；此外，兒童的肋骨材料也較柔軟，易變形而不易發(fā)生骨折，但對內(nèi)臟的保護效果較差[8-9]。上述差異說明通過簡單的比例縮放無法獲得精確的兒童模型。針對目前缺乏具有詳細解剖特征的兒童胸部有限元模型，本文中采用通過兒童臨床治療過程中獲得的CT和MRI圖像數(shù)據(jù)，開發(fā)出10歲兒童胸部有限元模型；并利用臨床治療過程中獲得的兒童胸部響應(yīng)數(shù)據(jù)對模型進行驗證，使開發(fā)的兒童胸部有限元模型具有較好的生物逼真度，為兒童損傷機理研究和防護裝置的開發(fā)提供可靠的工具。

1 幾何數(shù)據(jù)和驗證數(shù)據(jù)的獲取

本課題組與美國密西根州兒童醫(yī)院合作從多名兒童臨床CT和MRI圖像數(shù)據(jù)中獲得美國10歲兒童幾何模型，該方法獲得了韋恩州立大學(xué)倫理委員會的同意和批準。獲得的幾何模型與美國10歲兒童的統(tǒng)計數(shù)據(jù)對比見表1。由表1可見，該幾何模型尺寸與美國10歲兒童平均統(tǒng)計數(shù)據(jù)較好地吻合。幾何模型包括的主要解剖結(jié)構(gòu)有:骨架、胸部器官和除大、小腸以外的其他腹部器官。

表1 10歲兒童胸部幾何模型尺寸與美國兒童統(tǒng)計數(shù)據(jù)對比

本文中采用的驗證數(shù)據(jù)，同樣是臨床兒童心肺復(fù)蘇急救(CPR)中獲得的胸部響應(yīng)數(shù)據(jù)。兒童在接受CPR時，在外力和被救者的胸部之間放入傳感器，就可以記錄被救者的胸部響應(yīng)情況，同時，該傳感器還便于在施救過程中控制施救力的大小，該方法已被廣泛用于臨床急救中。文獻[9]中就利用CPR獲得了多個年齡段兒童的胸部力-變形曲線，本文中利用其中6個平均年齡為10.5±1.75歲的兒童數(shù)據(jù)對建立的10歲兒童胸部有限元模型進行驗證。

2 有限元模型建立

2.1 網(wǎng)格劃分

人體的解剖結(jié)構(gòu)十分復(fù)雜，且構(gòu)成人體的三維實體表面往往是不規(guī)則的自由曲面。因此，在人體建模中，實體網(wǎng)格的劃分是一項十分繁瑣而又復(fù)雜的工作；而采用合理的網(wǎng)格劃分方法來控制單元網(wǎng)格質(zhì)量極具挑戰(zhàn)性。文獻[11]中通過對比研究發(fā)現(xiàn)，目前常用的網(wǎng)格劃分方法中，映射網(wǎng)格劃分法的算法簡單、計算速度快、單元質(zhì)量好、單元密度可控，既可生成結(jié)構(gòu)化網(wǎng)格又可生成非結(jié)構(gòu)化網(wǎng)格。因此，該方法在開發(fā)大型模型時，具有其他方法不可比擬的優(yōu)勢。但是，當(dāng)三維實體的表面是十分復(fù)雜的自由曲面時，須創(chuàng)建大量的輔助線和面，人工分區(qū)比較困難，不得不對外觀做較多的簡化，網(wǎng)格劃分耗時也較長。

針對映射網(wǎng)格劃分法存在的上述缺點，ANSYS ICEM CFD/HEXA在映射網(wǎng)格劃分法的基礎(chǔ)上，開發(fā)了Block-Controlled方法來減小幾何分區(qū)的難度。該軟件可直接導(dǎo)入從CT圖像中提取三角單元格式幾何數(shù)據(jù)，無須逆向建模生成線和面，導(dǎo)入后的幾何模型如圖1(a)所示。導(dǎo)入幾何模型后，通過創(chuàng)建和編輯Block完成幾何分區(qū)，如圖1(b)所示。Block是人為創(chuàng)建的規(guī)則六面體，因此軟件可自動在Block內(nèi)部生成高質(zhì)量的純六面體單元，如圖1(c)所示。將形成的規(guī)則網(wǎng)格單元最外層節(jié)點映射到原幾何面上，即完成幾何模型的網(wǎng)格劃分，如圖1(d)～圖1(f)所示。

基于映射網(wǎng)格的Block-Controlled劃分方法大大減小了分區(qū)的難度，保證了模型內(nèi)部的網(wǎng)格質(zhì)量，同時映射最外層節(jié)點到模型的幾何面上，保留了原模型的幾何特征，降低了簡化對模型精度的影響。此外，可通過控制Block的邊或整體網(wǎng)格密度，來調(diào)整整個模型的單元數(shù)量。

本文中采用ANSYS ICEM CFD/HEXA軟件對模型的骨架和內(nèi)臟等復(fù)雜的三維實體進行網(wǎng)格剖分，并利用Hypermesh生成了肌肉和皮膚層，最終模型如圖2所示。模型中還包括腹部器官、腰椎和盆腔，為對胸部加載提供支撐邊界條件。不同的結(jié)構(gòu)采用了不同的單元類型進行模擬：皮質(zhì)骨采用殼單元模擬，松質(zhì)骨和內(nèi)臟采用實體單元模擬。最終的模型包括242 266個六面體單元，188 318個殼單元，模型的單元質(zhì)量如表2所示。

2.2 材料和屬性定義

肋骨、胸軟骨和胸骨是構(gòu)成胸腔的主要結(jié)構(gòu)。由于缺少相應(yīng)的兒童材料數(shù)據(jù)，通常兒童模型中的材料參數(shù)均是通過成人材料參數(shù)比例縮放獲得，如文獻[12]中利用優(yōu)化縮放系數(shù)的辦法獲得了10歲兒童盆腔骨的材料縮放系數(shù)和骨密質(zhì)層厚度的縮放系數(shù)，如表3所示。本文中同樣采用比例縮放的辦法，獲得10歲兒童的材料參數(shù)。

表2 10歲兒童胸部模型的網(wǎng)格質(zhì)量

2.2.1 肋骨

臨床CT圖像中無法精確測量到骨密質(zhì)層的厚度，因此肋骨的骨密質(zhì)層厚度同樣采用比例縮放獲得。文獻[13]中開發(fā)的成人模型取肋骨骨密層的厚度為0.72mm。文獻[14]中測得成人肋骨骨密質(zhì)層的厚度為0.23～2.96mm,平均厚度為0.84mm。根據(jù)文獻[12]中的縮放系數(shù)，兒童胸部模型中肋骨骨密質(zhì)層的厚度取為0.57mm。

表3 10歲兒童盆腔模型的材料屬性和縮放系數(shù)

在LS-DYNA中，骨密質(zhì)和骨松質(zhì)通常都采用彈塑性材料模型來模擬。根據(jù)相關(guān)文獻記載，成人肋骨骨密質(zhì)層的彈性模量大致范圍為9～15GPa，屈服極限為64.1～88MPa，極限應(yīng)力為79～130MPa，極限應(yīng)變?yōu)?.4%～3.0%。骨松質(zhì)的彈性模量范圍為40～748.8MPa，屈服極限為2～9MPa，極限應(yīng)力為2.5～9.12MPa[13-18]。本文中對于肋骨骨密質(zhì)和骨松質(zhì)的彈性模量、屈服極限和切線模量采用了不同的縮放系數(shù)。其中，彈性模量的縮放系數(shù)采用文獻[12]中的盆骨骨密質(zhì)彈性模量的縮放系數(shù)。由于該縮放系數(shù)中并未包括屈服極限和切線模量的縮放系數(shù)，因此屈服極限則根據(jù)Currey和Butler測得的兒童與成人大腿骨極限應(yīng)力的比值作為縮放系數(shù)，該方法同樣被Mizuno等人用于3歲兒童模型的比例縮放[7,12,19]。此外，Currey和Bulter的測量結(jié)果表明，兒童與成人骨骼的切線模量并無太大的差異[19]，因此，在兒童模型中采用了與成人相同的切線模量。最終，在10歲兒童胸部模型中，肋骨骨密質(zhì)的彈性模量取為6.48GPa，屈服極限為64.6MPa，骨松質(zhì)的彈性模量為252.4MPa，屈服極限為3.54MPa。

2.2.2 胸骨

通常，胸骨的材料屬性被定義為與肋骨的完全一致[13,16]；而胸骨骨密質(zhì)層厚在本文中則通過比例縮放后取為1.6mm。

2.2.3 胸軟骨

文獻[8]中總結(jié)了0～9歲和10～19歲兒童胸軟骨極限拉伸強度和極限拉伸應(yīng)變。根據(jù)其結(jié)果計算得到10歲兒童的胸軟骨彈性模量約為1.5MPa。本文中利用線彈性材料來模擬胸軟骨，彈性模量取為1.5MPa。

2.2.4 軟組織和內(nèi)臟

文獻[17]中通過歸納分析軟組織和內(nèi)臟的彈性模量后發(fā)現(xiàn)，軟組織的材料特性隨年齡變化并不明顯，且無固定規(guī)律。因此，在老人和成人模型中并未考慮軟組織材料的差異。本模型中，軟組織和內(nèi)臟的材料屬性同樣直接采用成人模型中的材料數(shù)據(jù)[13,15]。本模型中的主要材料參數(shù)如表4所示。

表4 10歲兒童胸部模型中的主要材料參數(shù)和縮放系數(shù)

3 模型驗證

根據(jù)記錄的兒童CPR數(shù)據(jù)，用每次擠壓胸部的峰值力除以變形可求得胸部的線性剛度，對求得同一個體的所有線性剛度值求平均后，獲得每個個體的平均線性剛度值。本文中利用6個兒童的平均線性剛度值對模型進行驗證[9]。

3.1 邊界和加載條件

根據(jù)兒童CPR手冊，將傳感器模型的中心置于胸部正中間和胸骨柄與劍突的中間位置，傳感器模型的大小為：長127mm，寬62mm，形狀與圖3所示的實物保持一致。加載曲線為一正弦位移-時間曲線，加載于傳感器的質(zhì)心位置。加載過程中的最大速度約為250mm/s，與文獻[9]在實驗中測得的最大速度一致。

3.2 結(jié)果

通過輸出傳感器和胸部之間的接觸力和傳感器質(zhì)心位置的位移，可以獲得胸部的力-變形曲線，如圖4(a)所示。圖4(b)為某一兒童典型的CPR力-變形曲線[9]。比較圖4(a)和圖4(b)可以看出，仿真曲線與實驗曲線具有相似的形狀。通過仿真曲線的峰值力和位移可以求得10歲兒童胸部模型的線性剛度為7 359.2N/m，該數(shù)值位于實驗測得的平均線性剛度值范圍內(nèi)，且十分接近整體的平均值(7 363.5N/m)，如圖5所示，表明該有限元胸部模型能較好預(yù)測10歲兒童胸部在擠壓情況下的響應(yīng)。

圖6(a)和圖6(b)分別為在整個CPR加載過程中，模型預(yù)測到的肋骨最大Von-Mises應(yīng)力和第一主應(yīng)變云圖。由圖可見：第1至第7根肋骨側(cè)面具有較高的Von-Mises應(yīng)力和最大主應(yīng)變，而第8至第12根肋骨上出現(xiàn)的應(yīng)力和應(yīng)變要低于第1至第7根肋骨；其中最大的Von-Miese應(yīng)力為40.1MPa，出現(xiàn)在第2根肋骨上，最大主應(yīng)變值為0.55%，與最大Von-Miese應(yīng)力出現(xiàn)的位置基本一致。文獻[18]中測得的成人肋骨側(cè)區(qū)域的極限拉伸應(yīng)變平均值為2.62%±1.27%，該值要遠高于模型預(yù)測的最大主應(yīng)變值，結(jié)合兒童肋骨材料特性，該結(jié)果表明兒童在接受正確的CPR情況下，不易出現(xiàn)肋骨骨折。

4 結(jié)論

利用臨床治療過程中獲得的兒童CT和MRI數(shù)據(jù)開發(fā)了一個具有詳細解剖特征的10歲兒童胸部有限元模型。在兒童人體數(shù)據(jù)十分缺乏的情況下，利用臨床治療過程中獲得數(shù)據(jù)來進行兒童有限元建模和驗證，是近年來較為新穎的人體損傷生物力學(xué)研究手段之一。網(wǎng)格劃分過程中采用基于映射網(wǎng)格的Block-Controlled劃分法，該方法不但保留了映射網(wǎng)格劃分法的優(yōu)點，同時大大降低了幾何分區(qū)的難度。分區(qū)后的模型能夠再次輕易地調(diào)整網(wǎng)格密度，且模型內(nèi)部的網(wǎng)格質(zhì)量極高，是目前人體建模中比較新穎的網(wǎng)格劃分法之一[20]。

在現(xiàn)有的成人材料數(shù)據(jù)基礎(chǔ)上，通過適當(dāng)?shù)谋壤凸δ芸s放得到了10歲兒童的材料參數(shù)，并利用臨床治療過程中獲得的CPR數(shù)據(jù)對模型進行驗證。通過驗證結(jié)果可知，該兒童胸部模型與10歲兒童的CPR數(shù)據(jù)較好吻合，具有良好的生物逼真度；該模型能模擬胸部加-卸載過程中的黏性響應(yīng)。內(nèi)臟損傷通常與胸部的黏性響應(yīng)存在密切關(guān)系，因此可以利用該模型來研究兒童內(nèi)臟損傷；模型的應(yīng)力和應(yīng)變云圖預(yù)測兒童在接受正確的CPR情況下，不易出現(xiàn)肋骨骨折。

[1] Durbin D R, Weinberg S. Policy Statement-child Passnenger Safety[J/OL]. American Academy of Pediatrics,2011,127:788-793. http://pediatrics.aappublications.org/content/127/4/788.full.pdf.

[2] Liu X J, Yang J K. Development of Child Pedestrian Mathematical Models and Evaluation with Accident Reconstruction[J]. Traffic Injury Prevention,2002,3(4):321-329.

[3] Van Rooij L, Meissner M, Bhalla K, et al. A Comparative Evaluation of Pedestrian Kinematics and Injury Prediction for Adults and Children Upon Impact with a Passenger Car[C]. SAE Paper 2004-01-1606.

[4] Van Rooij L, Harkema C, de Lange R, et al. Child Poses in Child Restraint Systems Related to Injury Potential: Investigations by Virtual Testing[C]. 19th ESV Conference, Washingtong DC, USA,2005.

[5] Forbes P A, Van Rooij L, Rodarius C, et al. Child Human Model Development: a Hybrid Validation Approach[C]. Proc. 2008 International Crashworthiness Conference, Kyoto, Japan,2008.

[6] Okamoto M, Takahashi Y, Mori F, et al. Development of Finite Element Model for Child Pedestrian Protection[C]. 18th ESV Conference, Nagoya, Japan,2003.

[7] Mizuno K, Iwata K, Deguchi T, et al. Development of a Three-year-old Child FE Model[J]. Traffic Injury Prevention,2005,6(4):361-371.

[8] Yamada H. Strength of Biological Materials[M]. Baltimore: Williams & Wilkins,1970.

[9] Maltese M R, Castner T, Niles D, et al. Methods for Determining Pediatric Thoracic Force-deflection Characteristics from Cardiopulmonary Resuscitation[J]. Stapp Car Crash Journal,2008,52:83-105.

[10] Snyder R G, Schneider L W, Owings C L. Infant, Child and Teenager Anthropometry for Product Safety Design[J]. Adv Consumer Res,1978,5:499-507.

[11] Cao L, Dian B, Lu C. Development and Validation of a Finite Element Model of the Human Neck Based on CT Images[J]. Journal of Hunan University (Natural Sciences),2009,36(1):24-29.

[12] Kim J E, Li Z, Ito Y, et al. Finite Element Model Development of a Child Pelvis with Optimization-based Material Identification[J]. Journal of Biomechanics,2009,42(13):2191-2195.

[13] Ito O, Dokko Y, Ohashi K. Development of Thorax Skeletal FE Models for Adult and Elderly[J]. Transactions of Society of Automotive Engineers of Japan,2009,40(6):1457-1462.

[14] Li Z, Kindig M W, Subit D, et al. Influence of Mesh Density, Cortical Thickness and Material Properties on Human Rib Fracture Prediction[J]. Medical Engineering & Physics,2010,32(9):998-1008.

[15] Shah C S, Yang K H, Hardy W, et al. Development of a Computer Model to Predict Aortic Rupture Due to Impact Loading[J]. Stapp Car Crash Journal,2001,45:161-182.

[16] Ruan J, El-Jawahri R, Chai L, et al. Prediction and Analysis of Human Thoracic Impact Responses and Injuries in Cadaver Impacts Using a Full Human Body Finite Element Model[J]. Stapp Car Crash Journal,2003,47:299-321.

[17] Kent R, Lee S H, Darvish K, et al. Structural and Material Changes in the Aging Thorax and Their Role in Crash Protection for Older Occupants[J]. Stapp Car Crash Journal,2005,49:231-249.

[18] Kemper A R, McNally C, Pullins C A, et al. The Biomechanics of Human Ribs: Material and Structural Properties from Dynamic Tension and Bending Tests[J]. Stapp Car Crash Journal,2007,51:235-273.

[19] Currey J D, Butler G. The Mechanical Properties of Bone Tissue in Children[J]. The Journal of Bone and Joint Surgery. American volume,1975,57(6):810-814.

[20] Mao H, Gao H, Cao L, et al. Development of High-quality Hexahedral Human Brain Meshes Using Feature-based Multi-block Approach[J]. Computer Methods in Biomechanics and Biomedical Engineering,2013,16(3):271-279.