1.本发明涉及岩土材料领域，特别涉及一种将计算机断层扫描数据传输至有限元模型的映射方法。
背景技术：2.土材料存在很强的空间异质性，在静动力加载作用下，这种颗粒尺度的空间异质性可以引起岩土材料的不可恢复的、非弹性的非均质变形，进而促使岩土材料出现渐进失稳破坏，并进一步威胁天然岩土体以及岩土工程结构物的稳定性，甚至引起工程事故，造成人员伤亡和财产损失。从岩土工程应用角度来看，认识和调控这种材料的空间异质性，对于进一步掌握岩土材料性能变化及引起的工程问题是至关重要的。因此，有必要建立一种映射策略，将先进x射线扫描数据投射至有限元模型，从而允许在连续介质力学基础上考虑材料空间异质性初始条件，进而促进更贴近真实情况的有限元数值模拟技术的提出。
3.专利文献cn 102768699 a公开了一种基于ct图像精确重构异质材料微观有限元网格模型的方法，通过ct图像扫描得到异质材料中各组分的形状细节，直接构造得到异质材料有限元模型。只是根据异质材料中各项组分的具体形状，构建了宏观异质有限元模型，并没有在有限元模型的每个高斯积分点中引入材料点级别的材料异质性，在揭示材料特性细观结构异质性方面存在不足。这是本技术需要着重改善的地方。
技术实现要素：4.本发明所要解决的技术问题是要提供一种将计算机断层扫描数据传输至有限元模型的映射方法，在有限元模型的每个高斯积分点考虑材料初始空间异质性，在数值分析中考虑材料初始场的影响，有利于研究、解释和控制室内试验或现场观察到的各种非弹性力学行为。
5.为了解决以上的技术问题，本发明提供了一种将计算机断层扫描数据传输至有限元模型的映射方法，基于计算机断层扫描数据，确定用于传输材料点的状态变量和/或材料参数的代表性体积元尺寸；建立对应的数值计算模型，输出有限元网格单元编号以及内部高斯积分点坐标；将有限元网格与试验试样重叠，确定高斯积分点对应的试样材料点坐标，根据选定的代表性体积元尺寸和断层扫描数据，计算该材料点的状态变量和/或材料参数，并赋予对应的高斯积分点，将计算机断层扫描数据传输至有限元模型，作为有限元计算分析的初始条件，从而在有限元数值计算中考虑材料初始空间异质性，包括如下的步骤：s1：读取计算机断层扫描获取的像素数据，数据类型为1代表实体颗粒的存在，数据类型为0代表颗粒间存在的孔隙，并据此生成试样可视化剖面；具体包括：s11、利用有限元分析平台编写子程序，获得灰度图的读取函数，将计算机断层灰度图数据读取至数值矩阵中，该数值矩阵包含试样中每个像素点的灰度数据，像素点级别的数据为1时，代表颗粒实体，而像素点级别的数据为0时，代表颗粒间的孔隙；s12、利用绘图函数，将灰度矩阵中的像素点进行可视化，数据类型为1的像素点赋
予黑色，数据类型为0的像素点赋予灰色，生成试样可视化剖面。
6.s2、选取任意材料点，利用代表性体积元概念，计算该材料点附近子区域内部的状态变量或者参数；具体包括：s21、选取试样内部任意空间位置的材料点，将以该材料点为中心的子区域用圆柱形代表体积元描述；s22、通过孔隙率计算公式计算子区域内部的每个选取的材料点以及代表体积元包围的子区域内的孔隙率计算结果，孔隙率值等于子区域内以0为代表的孔隙像素点个数除以总的像素点；s23、通过粒径分布计算公式计算子区域内部每个选取的材料点以及代表体积元包围的子区域内的平均粒径计算结果；在该步骤中，首先将子区域内以1为代表的像素点团块等效为球形颗粒，球心即为颗粒重心；若球心是处在该子区域内，判断整个颗粒处于子区域内，反之判断整个颗粒位于子区域外。
7.s3、根据任一子区域状态变量和/或材料参数随代表性体积元尺寸变化的波动，获得不同代表性体积元尺寸下的子区域内的状态变量和/或材料参数，绘出状态变量和/或材料参数随代表性体积元尺寸变化的关系曲线，根据波动稳定趋势，确定用于数据传输的代表性体积元尺寸；具体包括：s31、选取试样内部分布的多个材料点，分析计算在不同代表体积元尺寸条件下所感兴趣的孔隙率和颗粒平均粒径；s32、将不同代表体积元条件下的孔隙率值和颗粒平均粒径值绘图，并得到其随代表体积元尺寸变化的波动规律，根据波动稳定趋势，确认波动曲线稳定时的尺寸，即选择大于该尺寸的代表体积元。
8.s4、建立有限元分析模型，并划分为四节点实体单元，进一步输出所有单元体编号以及内部高斯积分点坐标；具体包括：s41、借助有限元分析平台建立圆柱形有限元数值计算模型；s42、在网格模块将模型划分为四节点实体单元，每个单元具备一个高斯积分点。利用子程序将模型内所有单元体编号以及高斯积分点坐标输入文本文件内。
9.s5、确定高斯积分点对应的试样材料点：将有限元数值模型与试验试样重叠，确定与有限元模型中的高斯积分点对应的试样内部对应的材料点空间位置；具体包括：s51、根据步骤s42得到的文本文件，读取高斯积分点坐标数据，转化为试样内部材料点的坐标数据。
10.s6、利用已确定的代表性体积元尺寸以及断层扫描数据，计算已确定的所有材料点子区域内的状态变量和/或材料参数，并根据单元编号和高斯点坐标，通过子程序平台将计算结果映射回对应的有限元模型高斯积分点，作为有限元计算的初始条件。
11.s7、提交计算，通过改变需要映射的状态变量和/或材料参数，获得不同情况下的材料空间异质性。
12.本发明的优越功效在于：1）本发明的映射方法借助代表性体积元概念，利用室内单元试样获取的计算机断层扫描数据，根据材料点对应子区域内的状态变量和/或材料参数进行平均化，并赋值给有限元模型中对应空间位置的高斯积分点，从而在数值分析中考虑材料初始场的影响，有利
于研究、解释和控制室内试验或现场观察到的各种非弹性力学行为；2）有效解决了有限元模拟中无法考虑真实颗粒材料空间异质性导致的数值计算结果失真难题，从而为岩土材料塑性变形问题的连续介质数值模拟工作提供了新的研究策略，可借助本发明的映射方法，进一步针对不同的本构模型，分析不同参数，如弹性模量、孔隙率和颗粒平均粒径等在空间尺度上的不均匀性对边值问题计算结果的影响，从而为室内试验结果分析和现场工程问题提供指导。
附图说明
13.构成本技术的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：图1为本发明具体实施例的流程图；图2a-2c为本发明具体实施例根据计算机断层扫描数据生成可视化剖面；图3a-3c为本发明具体实施例代表性体积元尺寸分析；图4a-4b为本发明具体实施例计算机断层扫描数据映射至有限元模型示意图；图5为本发明具体实施例孔隙率空间分布映射结果；图6为本发明具体实施例颗粒平均粒径空间分布映射结果。
具体实施方式
14.以一个砂土一维压缩试验为例，结合附图对本发明的实施例进行详细说明。
15.图1示出了本发明具体实施例的流程图；图2a-2c示出了本发明具体实施例根据计算机断层扫描数据生成可视化剖面；图3a-3c示出了本发明具体实施例代表性体积元尺寸分析；图4a-4b示出了本发明具体实施例计算机断层扫描数据映射至有限元模型示意图；图5示出了本发明具体实施例孔隙率空间分布映射结果；图6示出了本发明具体实施例颗粒平均粒径空间分布映射结果。
16.如图1所示，本发明提供了一种将计算机断层扫描数据传输至有限元模型的映射方法，基于计算机断层扫描数据，确定用于传输材料点的状态变量和/或材料参数的代表性体积元尺寸；建立对应的数值计算模型，输出有限元网格单元编号以及内部高斯积分点坐标；将有限元网格与试验试样重叠，确定高斯积分点对应的试样材料点坐标，根据选定的代表性体积元尺寸和断层扫描数据，计算该材料点的状态变量和/或材料参数，并赋予对应的高斯积分点，将计算机断层扫描数据传输至有限元模型，作为有限元计算分析的初始条件，从而在有限元数值计算中考虑材料初始空间异质性，包括如下的步骤：s1：读取计算机断层扫描获取的像素数据，数据类型为1代表实体颗粒的存在，数据类型为0代表颗粒间存在的孔隙，并据此生成试样可视化剖面；具体包括：s11、利用有限元分析平台编写子程序，获得灰度图的读取函数，将计算机断层灰度图数据读取至数值矩阵中，该数值矩阵包含试样中每个像素点的灰度数据，像素点级别的数据为1时，代表颗粒实体，而像素点级别的数据为0时，代表颗粒间的孔隙；s12、利用绘图函数imread和imwrite，将灰度矩阵中的像素点进行可视化，数据类型为1的像素点赋予黑色，数据类型为0的像素点赋予灰色，生成试样可视化剖面；图2a-2c 给出了砂土一维压缩试验示意图以及对应的应力水平下的试样可视化剖面，在图2a中，荷
载通过对加载板施加速率恒定的位移速率0.05mm/min施加应力，试样放置于室内三轴试验装置中。在图2b中，轴向应变与轴向应力的方向与图2a中示出的荷载施加方向一致。在图2c中，
①
至
⑥
六张图为对x射线ct数据可视化后得到的试样剖面颗粒破碎动态变化图，这六张图分别对应图2b中示出的六个应力状态下得到的试样的非弹性状态。
17.s2、选取任意材料点，利用代表性体积元概念，计算该材料点附近子区域内部的状态变量或者参数；具体包括：s21、选取试样内部任意空间位置的材料点，将以该材料点为中心的子区域用圆柱形代表体积元描述，圆柱形代表体积元长径比与试样相同，即：所利用的代表体积元的高度与底面直径比值与所分析的试样的高度与底面直径比值相同；s22、通过孔隙率计算公式计算子区域内部的每个选取的材料点以及代表体积元包围的子区域内的孔隙率计算结果，孔隙率值等于子区域内以0为代表的孔隙像素点个数除以总的像素点；s23、通过粒径分布计算公式计算子区域内部每个选取的材料点以及代表体积元包围的子区域内的平均粒径计算结果；在该步骤中，首先将子区域内以1为代表的像素点团块等效为球形颗粒，球心即为颗粒重心；若球心是处在该子区域内，判断整个颗粒处于子区域内，反之判断整个颗粒位于子区域外。
18.s3、根据任一子区域状态变量和/或材料参数随代表性体积元尺寸变化的波动，确定用于数据传输的代表性体积元尺寸；通过改变代表性体积元尺寸，获得不同代表性体积元尺寸下的子区域内的状态变量和/或材料参数，绘出状态变量和/或材料参数随代表性体积元尺寸变化的关系曲线，根据波动稳定趋势，确定用于数据传输的代表性体积元尺寸；具体包括：s31、选取试样内部分布的多个材料点，分析计算在不同代表体积元尺寸条件下所感兴趣的孔隙率和颗粒平均粒径；s32、将不同代表体积元条件下的孔隙率值和颗粒平均粒径值绘图，并得到其随代表体积元尺寸变化的波动规律；图3a-图3c给出了试样内部空间分布的多个材料点以及计算得到的孔隙率和颗粒平均粒径波动曲线（图3b-图3c中的rev高度表示体积元高度），由于波动曲线在代表体积元高度大于1.2mm时稳定，因此选择大于该尺寸，即大于1.2mm时的代表体积元；图3a示出了选取的材料点的空间位置，图3b示出了孔隙率波动曲线，图3c示出了平均粒径波动曲线。
19.s4、建立有限元分析模型，并划分为四节点实体单元，进一步输出所有单元体编号以及内部高斯积分点坐标；根据所研究的试样几何条件以及所选定的代表性体积元尺寸，建立有限元分析模型，并通过网格划分模块划分为四节点实体单元，进一步利用子程序平台输出有限元网格中所有单元体编号以及内部高斯积分点坐标；具体包括：s41、借助有限元分析平台建立圆柱形有限元数值计算模型，为了避免将试样外部空白区域数据映射至有限元计算模型，模型的大小应略微小于试样的大小，如图4a所示；s42、在网格模块将模型划分为四节点实体单元，每个单元具备一个高斯积分点。利用子程序将模型内所有单元体编号以及高斯积分点坐标输入txt文件内。
20.s5、确定高斯积分点对应的试样材料点：将有限元数值模型与试验试样重叠，确定
与有限元模型中的高斯积分点对应的试样内部对应的材料点空间位置；具体包括：s51、根据上一步得到的txt文件，通过商业有限元程序的自定义单元程序中的sdv子程序进行读取，将txt读取到的高斯积分点坐标数据通过等比例变化转化为试样内部材料点的坐标数据。
21.s6、利用已确定的代表性体积元尺寸以及断层扫描数据，计算已确定的所有材料点子区域内的状态变量和/或材料参数，并根据单元编号和高斯点坐标，通过子程序平台将计算结果映射回对应的有限元模型高斯积分点，作为有限元计算的初始条件；具体包括：s61、根据s32，将代表性体积元尺寸设定为1.4mm，计算每个材料点的代表性体积元内部子区域的孔隙率和颗粒平均粒径，并再一次通过有限元软件的子程序将数据读取并赋值给有限元模型中的所有高斯积分点，图4b示出了该映射的流程示意图，在图4b中，首先，将积分点坐标投影到试验样本上，其中，每个高斯积分点都与一个rev相关联；然后，根据ct图像计算rev域内的平均材料特性和/或初始状态变量；最后，将计算得到平均值分配回相应的高斯积分点。
22.s7、提交计算，获得数据投射结果：通过改变需要映射的状态变量和/或材料参数，获得不同情况下的材料空间异质性；图5示出了基于本发明方法映射得到的孔隙率空间分布映射结果；在该图中，显示不同高斯积分点获得了不同的孔隙率值，这些空间分布的孔隙率值构成了图示的黑白云图，黑色代表孔隙率较高区域，而白色代表孔隙率较低区域。图6示出了基于本发明方法映射得到的颗粒平均粒径空间分布映射结果。在该图中，显示不同高斯积分点获得了不同的平均粒径指标参数值，这些空间分布的平均粒径指标参数值构成了图示的黑白云图，黑色代表平均粒径较大区域，而白色代表平均粒径较小区域。
23.以上所述仅为本发明的优先实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
技术特征：1.一种将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：基于计算机断层扫描数据，确定用于传输材料点的状态变量和/或材料参数的代表性体积元尺寸；建立对应的数值计算模型，输出有限元网格单元编号以及内部高斯积分点坐标；将有限元网格与试验试样重叠，确定高斯积分点对应的试样材料点坐标；根据选定的代表性体积元尺寸和断层扫描数据，计算该材料点的状态变量和/或材料参数，并赋予对应的高斯积分点，将计算机断层扫描数据传输至有限元模型，作为有限元计算分析的初始条件，从而在有限元数值计算中考虑材料初始空间异质性。2.根据权利要求1所述的将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：包括如下的步骤：s1：读取计算机断层扫描获取的像素数据，数据类型为1代表实体颗粒的存在，数据类型为0代表颗粒间存在的孔隙，并据此生成试样可视化剖面；s2、选取任意材料点，利用代表性体积元概念，计算该材料点附近子区域内部的状态变量或者参数；s3、根据任一子区域状态变量和/或材料参数随代表性体积元尺寸变化的波动，获得不同代表性体积元尺寸下的子区域内的状态变量和/或材料参数，绘出状态变量和/或材料参数随代表性体积元尺寸变化的关系曲线，根据波动稳定趋势，确定用于数据传输的代表性体积元尺寸；s4、建立有限元分析模型，并划分为四节点实体单元，进一步输出所有单元体编号以及内部高斯积分点坐标；s5、确定高斯积分点对应的试样材料点：将有限元数值模型与试验试样重叠，确定与有限元模型中的高斯积分点对应的试样内部对应的材料点空间位置；s6、利用已确定的代表性体积元尺寸以及断层扫描数据，计算已确定的所有材料点子区域内的状态变量和/或材料参数，并根据单元编号和高斯点坐标，通过子程序平台将计算结果映射回对应的有限元模型高斯积分点，作为有限元计算的初始条件；s7、提交计算，通过改变需要映射的状态变量和/或材料参数，获得不同情况下的材料空间异质性。3.根据权利要求2所述的将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：步骤s1具体包括：s11、利用有限元分析平台编写子程序，获得灰度图的读取函数，将计算机断层灰度图数据读取至数值矩阵中，该数值矩阵包含试样中每个像素点的灰度数据，像素点级别的数据为1时，代表颗粒实体，而像素点级别的数据为0时，代表颗粒间的孔隙；s12、利用绘图函数，将灰度矩阵中的像素点进行可视化，数据类型为1的像素点赋予黑色，数据类型为0的像素点赋予灰色，生成试样可视化剖面。4.根据权利要求2所述的将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：步骤s2具体包括：s21、选取试样内部任意空间位置的材料点，将以该材料点为中心的子区域用圆柱形代表体积元描述；s22、通过孔隙率计算公式计算子区域内部的每个选取的材料点以及代表体积元包围的子区域内的孔隙率计算结果，孔隙率值等于子区域内以0为代表的孔隙像素点个数除以
总的像素点；s23、通过粒径分布计算公式计算子区域内部每个选取的材料点以及代表体积元包围的子区域内的平均粒径计算结果；在该步骤中，首先将子区域内以1为代表的像素点团块等效为球形颗粒，球心即为颗粒重心；若球心是处在该子区域内，判断整个颗粒处于子区域内，反之判断整个颗粒位于子区域外。5.根据权利要求2所述的将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：步骤s3具体包括：s31、选取试样内部分布的多个材料点，分析计算在不同代表体积元尺寸条件下所感兴趣的孔隙率和颗粒平均粒径；s32、将不同代表体积元条件下的孔隙率值和颗粒平均粒径值绘图，并得到其随代表体积元尺寸变化的波动规律，根据波动稳定趋势，确认波动曲线稳定时的尺寸，即选择大于该尺寸的代表体积元。6.根据权利要求2所述的将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：步骤s4具体包括：s41、借助有限元分析平台建立圆柱形有限元数值计算模型；s42、在网格模块将模型划分为四节点实体单元，每个单元具备一个高斯积分点；利用子程序将模型内所有单元体编号以及高斯积分点坐标输入文本文件内。7.根据权利要求2所述的将计算机断层扫描数据传输至有限元模型的映射方法，其特征在于：步骤s5具体包括：s51、根据步骤s42得到的文本文件，读取高斯积分点坐标数据，转化为试样内部材料点的坐标数据。
技术总结本发明公开一种将计算机断层扫描数据传输至有限元模型的映射方法，基于计算机断层扫描数据，确定用于传输材料点的状态变量和/或材料参数的代表性体积元尺寸；建立对应的数值计算模型，输出有限元网格单元编号以及内部高斯积分点坐标；将有限元网格与试验试样重叠，确定高斯积分点对应的试样材料点坐标，根据选定的代表性体积元尺寸和断层扫描数据，计算该材料点的状态变量和/或材料参数，并赋予对应的高斯积分点，将计算机断层扫描数据传输至有限元模型，作为有限元计算分析的初始条件，从而在有限元数值计算中考虑材料初始空间异质性。本发明的优点是有效解决了有限元模拟中无法考虑真实颗粒材料空间异质性导致的数值计算结果失真难题。算结果失真难题。算结果失真难题。
技术研发人员：薛大为 吕玺琳 黄茂松
受保护的技术使用者：同济大学
技术研发日：2022.06.30
技术公布日：2022/11/1
