[发明专利]一种基于拓扑优化技术的松质骨多孔结构建模方法

摘要

本发明涉及生物医疗技术领域，一种基于拓扑优化技术的松质骨多孔结构建模方法，包括以下步骤：(1)建立松质骨模型，(2)求解松质骨模型的体积分数，(3)求解松质骨模型的力学性能，(4)提出本方法的拓扑优化模型，(5)求解多孔结构，(6)还原CAD模型。本发明考虑到真实松质骨的性能，设计得到的模型更加符合真实松质骨结构的要求，并且大大减少了模型的数据量，满足加工和后处理的要求，为进一步的骨科临床手术应用奠定了基础。

主权项

一种基于拓扑优化技术的松质骨多孔结构建模方法，其特征在于，包括以下步骤：步骤1、建立松质骨模型，为了保证松质骨模型的真实性，本发明基于医学图像建立松质骨模型，具体包括以下子步骤：(a)根据设计要求选取股骨头并对其进行物理处理，利用Micro‑CT仪器扫描得到松质骨结构的医学图像；(b)利用Mimics软件对Micro‑CT仪器扫描得到的医学图像进行处理，包括窗宽和窗位调整、阈值分割、区域增长，得到具有清晰骨小梁结构的松质骨图像，并建立松质骨结构模型；步骤2、求解松质骨模型的体积分数，将在Mimics中建立得到的松质骨模型导出为STL格式，并将其导入到Magics软件中得到松质骨模型的体积信息，并基于此体积信息对松质骨模型的体积分数Vm/V进行求解，其中Vm是松质骨实心部分的体积，V是松质骨实心部分和空心部分的总体积；步骤3、求解松质骨模型的力学性能，选用代表体元法对模型施加约束条件和外载荷，并基于ANSYS软件对松质骨模型的9个力学性能参数进行求解，具体包括以下子步骤：(a)对松质骨模型施加约束，可通过公式(1)进行描述，u(0,y,z)=0u(a,y,z)=constant=δxv(x,0,z)=0v(x,b,z)=constant=δyw(x,y,0)=0w(x,y,c)=constant=δz---(1)]]>式中，1,3,5行分别代表x＝0，y＝0和z＝0的面分别在x方向，y方向和z方向的位移为0；u，v，w分别是x，y，z三个方向的位移；a，b，c分别是代表体元三个边的长度；2,4,6行分别代表x＝a，y＝b，z＝c的面的位移分别为δx，δy，δz，它们的值由ANSYS求解得到；(b)求解松质骨模型的等效弹性模量，包括三个等效弹性模量，分别是x方向的弹性模量Ex，y方向的弹性模量Ey，z方向的弹性模量Ez，可通过公式(2)、(3)、(4)进行描述，Ex=pabcδx---(2)]]>Ey=pbacδy---(3)]]>Ez=pcabδz---(4)]]>式中，p为施加的外载荷，它的大小不影响结果；a为模型x方向的边长，b为模型y方向的边长，c为模型z方向的边长；(c)求解松质骨模型的等效泊松比，包括三个等效泊松比，分别对应于xy面的泊松比vxy，yz面的泊松比vyz和zx面的泊松比vzx，可通过公式(5)、(6)、(7)进行描述，vxy=-δxyδx---(5)]]>vyz=-δyzδy---(6)]]>vzx=-δzxδz---(7)]]>式中，δxy为在x＝a的面施加载荷p后y＝b的面的位移，δyz为在y＝b的面施加载荷p后z＝c的面的位移，δzx为在z＝c的面施加载荷p后x＝a的面的位移；(d)求解松质骨模型的等效剪切模量，包括三个等效剪切模量，分别是与xy面平行的剪切模量Gxy，与yz面平行的剪切模量Gyz，与zx面平行的剪切模量Gzx，可通过公式(8)、(9)、(10)进行描述，Gxy=Ex2(1+vxy)---(8)]]>Gyz=Ey2(1+vyz)---(9)]]>Gzx=Ez2(1+vzx)---(10)]]>步骤4、提出本方法的拓扑优化模型，其中包括目标函数、约束条件以及控制方程，具体包括以下子步骤：(a)本发明引入拓扑优化中的材料模型，选用Solid Isotropic Material WithPenalization(SIMP)方法，可通过公式(11)进行描述，E(x)＝ρ(x)pE0  (11)式中，E(x)为材料点x处的单元材料的弹性模量，E0为材料本身的弹性模量，ρ(x)为单元材料的密度，p是惩罚因子；(b)为了能够实时监测迭代过程中的材料模型的等效性能，采用均匀化方法通过公式(12)、(13)对迭代过程中的材料模型性能进行求解，∫Y(Eijkl-Eijmn∂χmkl(y)∂yn)∂vi∂yjdy=0,∀χmki(y)∈VΩ---(12)]]>EijklH=1|Y|∫Y(Eijkl-Eijmn∂χmkl(y)∂yn)dy---(13)]]>式中，Y为单胞域，是符合上式的广义位移函数；vi表示单胞中定义的连续函数；Eijkl表示单胞中材料的弹性常数；i，j,k，l表示坐标编号，二维取1，2，三维取1，2，3，Eijmn表示yn处的四阶弹性张量，m，n同样是坐标编号，二维取1，2，三维取1，2，3，yn，yj表示不同维度的坐标，VΩ表示设计域，通过公式(12)求解得到广义位移函数然后通过公式(13)求解得到等效的弹性模量；(c)基于SIMP材料插值模型和均匀化方法，本发明以松质骨结构的九个性能参数和体积分数为设计目标提出了结构建模优化列式，可通过公式(14)进行描述，findρe,e=1,2,3,...,Nf(ρe)=(Ex-Ex‾)2+(Ey-Ey‾)2+(Ez-Ez‾)2+(Gxy-Gxy‾)2+(Gyz-Gyz‾)2+(Gzx-Gzx‾)2+(vxy-vxy‾)2+(vyz-vyz‾)2+(vzx-vzx‾)2subjecttoV≤VuppE(x)=ρ(x)pE0Homogenization---(14)]]>式中，ρe表示单元密度，具体可以参考SIMP方法；e表示单元编号，N表示单元的总数目；f(ρe)表示目标函数，Ex、Ey、Ez、Gxy、Gyz、Gzx、vxy、vyz、vzx为多孔结构的九个设计参数，同时也是要设计要达到的性能，表示设计迭代过程中九个等效性能实时的值，Vupp表示设计的体积上限，V表示单元材料分布的总体积，E(x)＝ρ(x)pE0表示引入的密度插值函数，Homogenization表示均匀化方法的控制方程；步骤5、求解多孔结构，首先确定设计域，给定一个r×r×r的正方体网格设计域，建议给定敏度过滤半径rmin为1.5，建议最大迭代步数为100步，九个目标性能参数，体积约束，给定实心材料无孔隙条件下的弹性模量和泊松比，控制某一特定阈值以上的单元在ANSYS中显示，所述阈值特定值由具体情况进行选取，范围在0‑1之间；步骤6、还原CAD模型，将ANSYS中显示的多孔结构利用CAD软件建立得到符合实际的CAD模型。