nanoVoxel 2000

高分辨显微CT

本研究中4D CT技术采用天津三英精密nanoVoxel 2000型，使用CT技术测试高温前后样本的孔隙，量化孔隙数量、位置、孔隙体积百分比，建立孔隙网络模型，将CT测试后的孔隙数据纳入修正后的本构模型，对构件承载力进行模拟分析，实现从 “材料损伤” 到 “构件性能” 的量化预测。通过无损检测（超声波）– 微观结构量化（CT）– 数值模拟” 技术流程，可直接应用于火灾后梁、柱等受损结构构件的损伤评估。

本实验所用水泥为P·O 42.5 普通硅酸盐水泥，其化学成分见表 1。细骨料采用天然河砂，黏土含量（质量分数）＜1.5%，粒径范围为 0.15 mm～4.75 mm，颗粒大小分布见表 2。减水剂采用第三代聚羧酸系高效减水剂。

表1. 水泥化学成分

表2. 河沙粒径分布

按照表 3 所示配合比，制备了两种砂浆试件（M30 和 M60）。样品尺寸分为两种：20 mm × 20 mm × 40 mm和40 mm × 40 mm × 160 mm。砂浆试件在模具中浇筑成型，自然养护1天后脱模，随后置于标准养护环境（温度 20±2 ℃，相对湿度＞95%）中养护 28 天。

表3. 砂浆配比

本实验采用天津三英精密仪器有限公司生产的高分辨率 X 射线三维显微 CT （nanoVoxel-2000型），对尺寸为 20 mm×20 mm×40 mm 的砂浆样品进行 CT 扫描。扫描参数：加速电压 150 kV，加速电流 60 μA，试件绕竖轴旋转 360°，曝光时间 0.55 s。

本研究采用基于原位 4D CT 技术扫描，研究高温对样本内部孔隙结构的影响。为保证结果有效，制备了同一批次的 5 个样本。每个样本均进行两次 CT扫描（常温 25 ℃下及高温作用后），以实现原位对比分析，避免材料不均匀对实验数据产生影响。首先对同一样本高温作用前后分别进行 CT 扫描，保证扫描位置与初始状态完全一致，实现相同区域的原位对比。具体步骤如下：对样品 SA 进行 X 射线 CT 扫描，获得其原始状态的二维投影图像；随后对样品 SA 进行高温处理，得到试件 SAH，再次进行 X 射线 CT 扫描，获得高温作用后状态的二维投影图像。完成数据采集，采用滤波反投影算法将二维投影数据重建为三维体数据。最后利用 Avizo 软件对三维体数据进行处理，包括边缘裁剪、平滑去噪、阈值分割，得到样本高温作用前后的三维重建结果，用于原位对比分析。CT 数据处理流程如图 1 所示。

图1. CT数据处理工作流程

图像对齐步骤：采用带有侧向限位块与中心十字标记的专用亚克力基板。初次扫描后，记录基板与扫描部件间的精确距离；后续扫描通过限位块保持该间距，将物理定位偏差控制在 ±0.1 mm 以内。在样本底部与基板中心均设置精密十字标记，确保每次扫描前精准对齐。在样本侧面对称粘贴泡棉胶，保证样本与基板紧密贴合，防止扫描过程中发生松动或位移。所有扫描参数（管电压、管电流、层厚、重建间隔等）均严格保持一致，避免因参数变化引起坐标系偏差。采用4D CT技术研究时，不仅可以运用本文的研究方法，还可以使用天津三英精密仪器所配备的具有原位加载装置的4D CT。该装置不仅可以模拟高低温场，还可以提供压力/拉力、充/放电等环境下样本内部结构观测。使用该原位加载装置4D CT检测时，样品无需从设备取出，可实现全程原位检测，有效保障检测过程的真实性与准确性，最大程度还原样品在实际加载状态下的特性。

边缘裁剪：在三维模型重建过程中，为避免边缘硬化效应对数据分析产生影响，需对 CT 图像边缘部分进行裁剪。边缘裁剪后，沿 Z 轴获得 1600 张二维 CT 切片，沿 X 轴和 Y 轴各获得 800 张二维 CT 切片。分析区域实际尺寸为 19.944 mm×19.944 mm×39.888 mm。

平滑去噪：采用高斯滤波进行降噪与图像平滑处理，以降低图像对比度、柔化材料边缘。

阈值分割：以 600 ℃作用后的砂浆为例，利用灰度直方图对二维图像进行二值化分割，如图 2(a) 所示。通过两条直线拟合灰度直方图分布曲线确定分界点，该点的横坐标（灰度值）即为二值化阈值 。计算得到阈值为 93.914，灰度值低于该阈值的区域为孔隙区，高于该阈值的区域为砂浆基体区，所得二值化图像如图 2(b) 所示。其他温度下的灰度直方图见附录 A。

图2. 砂浆 CT 二值化分割图像

三维重建与可视化：采用体绘制法对二值化后的二维图像进行三维重建与可视化。该方法假定数据体内的每个点均可发光与吸光，通过包含透明度的颜色传递函数，由标量数据确定发光量、颜色及吸光程度，从而实现可视化效果。

超声波检测时，采用的超声波频率大于 20000 Hz 的弹性波。该波段在砂浆中的传播特性受孔隙与微裂缝影响显著。超声波波速是表征材料损伤程度的高灵敏度指标。超声波波速（v）计算公式如下：

v=L/t

公式中：v 为超声波波速（km/s）；L 为砂浆试件长度（mm）。

为提高超声波对微观结构变化的敏感性，采用对测法测量超声波波速,频率范围0.6 MHz～1.2 MHz  。实验过程中，为改善样本与换能器之间的耦合效果，在样本表面与换能器表面均匀涂抹凡士林。样本尺寸为 40 mm × 40 mm × 160 mm。每种情况选取3 个砂浆试件，如图 3 所示，以棱柱体的正方形端面为检测面，超声波测点位于样本端面中心位置。

图3. 砂浆表面的超声波检测点

高温加热实验采用 一体化智能马弗炉。炉膛尺寸为 400 mm × 250 mm × 160 mm，最高设定温度可达 1000 ℃。

为减缓因升温速率过快导致的砂浆试件爆裂剥落，将升温速率控制在 5 ℃/min。炉温达到目标温度后，恒温保温 3 h。马弗炉内部升温曲线如图 4所示。

图4. 加热曲线

保温结束后，打开炉门将样本取出，置于室内自然冷却（温度：25 ± 2 ℃，相对湿度：50 ± 5 %）。

本研究采用微机控制电子万能实验机，进行单轴压缩试验，测定样本的轴心抗压强度。实验采用位移控制方式，加载速率为 0.5 mm/min。为获取样本的应力-应变曲线，分别在距样本上、下表面 10 mm 处布置两个位移传感器（LVDT1 和 LVDT2），两传感器之间的距离为 20 mm。实验装置示意图如图 5 所示。

图5. 单轴压缩实验

采用原位 4D CT 技术对样本热处理前后进行观测，以探究温度变化对其内部孔隙结构的影响。利用 Avizo 软件对经历不同温度作用的 M30 和 M60 试样进行三维展示，结果分别如图 6和图 7 所示。矩形框为分析区域，绿色区域代表孔隙。绿色区域密度越高，砂浆中的孔隙含量越大，表明内部结构损伤程度越高。 如图 6 和图 7所示，在 25 ℃时，砂浆中孔隙数量极少，且 M30 中绿色孔隙区域的密度高于 M60，说明 M60 的结构更为密实。当温度升至 200 ℃时，砂浆内部绿色区域与 25 ℃时相比无明显变化。在 400 ℃时，砂浆中产生大量新孔隙，且孔隙密度显著增大。 体孔隙率（P）定义为孔隙体积与砂浆分析区域总体积之比，以百分数表示，如公式（1）所示

V_pore 为孔隙体积（μm³），V_total  为分析区域的总体积（μm³）。

砂浆是一种非均质多相复合材料，其孔隙同时具有空间非均质性与随机分布特征。在相同试验条件和步骤下制备的样品，其孔隙率也存在离散性。为减小由材料非均质性和孔隙随机分布带来的误差，并更系统地研究温度对砂浆孔隙结构的影响，对砂浆试件在高温作用前后进行了原位扫描。 体孔隙增长率 (L_P) 定义为高温处理后的体积孔隙率与 25 ℃时的体孔隙率之比，如公式（2）所示

P_T 为高温处理后砂浆试件的孔隙率（%），P₂₅ 为 25 ℃时砂浆试件的孔隙率（%）。采用P和L_P方法可以评估温度对砂浆内部损伤的影响，不同温度下砂浆的P和L_P结果如表4所示。

表4. 不同温度条件下的孔隙率

图6. M30在不同温度下的孔隙三维结构

图7. M60在不同温度下的孔隙三维结构

体孔隙增长率随温度的变化关系如图 8a所示，该曲线图充分体现了高温作用对砂浆体孔隙率的影响。图中实线为实验结果，虚线为基于实验数据拟合的曲线。由图 8a可知，在 400 ℃以内，M60 的体孔隙增长率高于 M30；当温度超过 600 ℃时，M30 的体孔隙增长率高于 M60。

图8(a). 孔隙率随温度的变化

如图 8a所示，M30 和 M60 的体孔隙增长率随温度呈非线性增长。分别采用玻尔兹曼函数对 M30 和 M60 的体孔隙增长率与温度之间的关系进行拟合，结果如公式 (3) 和 (4) 所示，其相对系数（R²）分别为 0.952 和 0.976。

图 8b 给出了砂浆相对抗压强度与体孔隙率之间的关系，其中数据点为实验测试结果，虚线为基于 Schiller提出的函数拟合的实验结果曲线。拟合函数分别如公式 (5) 和公式 (6) 所示，相关系数（R²）分别为 0.963 和 0.981。研究表明，随着体孔隙率的增大，相对抗压强度逐渐降低，说明抗压强度与体孔隙率之间呈负相关。

图8(b). 砂浆相对抗压强度与体积孔隙率之间的关系

根据方程 (3) 和式 (4) 推导可得，砂浆的体孔隙增长率随温度升高呈上升趋势。方程 (5) 和式 (6) 表明，砂浆的抗压强度随体孔隙率的增加而降低。通过方程 (1)、(2) 与方程 (5)、(6) 的关系式相结合，建立了表征抗压强度、孔隙率与温度三者相互关系的热损伤模型。M30 和 M60 对应的热损伤模型分别由式 (7) 和式 (8) 给出

基于“温度(T)→孔隙增长速率（LP）→相对抗压强度（Δf）→高温抗压强度（f_{c,T）”构建本损伤模型，该方法严格遵守损伤力学的基本原理，即“细观缺陷的演化决定宏观性能的劣化”。为量化前述各独立方程误差的累积效应，采用高斯误差分布定律 (GEP)对模型开展系统性不确定性分析。该定律依托 “小误差假设” 与 “随机误差正态分布”，通过函数的线性近似（一阶泰勒展开），计算自变量误差对因变量的影响。} 

基于 600 ℃温度下各参数（温度、孔隙、强度）的误差敏感性最高、叠加效应最显著，选择该温度进行误差分布分析。详细分析过程见附录 B，误差分析结果如下： 在热损伤模型的级联传递过程中，总误差并非通过 “各步骤误差直接累加” 得到，而是采用均方根叠加法计算，以避免对误差的高估或低估。总相对误差 (δtotal) 主要由孔隙增长率相对误差 (δLp)、相对抗压强度相对误差 (δΔf) 及常温下相对误差 (δf25) 构成。经计算，总相对误差 δtotal 为4.85%，与直接由高温抗压强度 f_c,T 绝对误差得到的 4.49% 相对误差高度吻合。这验证了误差传递路径的正确性，同时表明总误差处于可接受范围，进而证明了热损伤模型的有效性。为进一步验证热损伤模型的有效性，将 M30 与 M60 砂浆的平均残余抗压强度实测值 (f_c) 与热损伤模型计算值 (f_c,T ) 进行对比分析，结果如表 5所示。M30 与 M60 的平均相对误差 (Δ) 分别为 7.33% 和 8.14%。结果表明，所建立的砂浆热损伤模型能够有效表征砂浆强度、孔隙率与温度三者之间的关系。

表5. 实验值与计算值对比

使用 Avizo 软件，可将 CT 三维体素数据直接转换为相单元。这种“体素直接转单元”的方法是计算材料学中处理复杂多孔介质的标准方法。其核心优势在于能够真实、无损地重构材料内部包括孔隙与微裂纹在内的复杂几何形态。 由于模型内部结构复杂、孔隙分布不规则，网格划分时选用四面体单元；相较于线性单元，其拟合能力与计算精度更优，可更精准地捕捉应力集中现象。对相单元执行四面体自动网格划分。网格划分过程中未生成点、线、面、体等几何特征，从而避免了网格划分失效。同时自动网格程序在孔隙、裂纹等复杂特征区域生成更细密的网格，最小网格尺寸为 0.874 μm。 图 9为常温 25 ℃下 M30 砂浆的三维网格模型生成结果。网格划分完成后，该模型单元数量约为104.8 万个；其他温度下模型的单元数量在96.8 万～138.5 万个之间（图10）。

图10. 不同网格尺寸下的应力 - 应变曲线

将三维网格模型导入 ABAQUS 软件，并指定单元类型为四面体单元。在三维网格模型上表面中心位置设置参考点 RP1，在参考点 RP1 与砂浆上表面之间施加耦合约束。同理，在模型下表面中心设置另一参考点 RP2，并在 RP2 与砂浆下表面之间施加耦合约束。采用位移控制方式设定边界条件：约束 RP1 在 X、Y 方向的平动自由度（U1=U2=0），并沿 Z 轴施加 0.5 mm 的位移荷载（U3=0.5 mm）；约束 RP2 的全部平动自由度（U1=U2=U3=0）。 由于网格模型无明显几何特征，采用网格角度（20°）判据选取上下表面所有节点。为使型分析过程中的收敛，在网格模块的“单元指派类型” 命令中，将四面体单元的几何阶次从线性单元修改为二次单元。

本研究采用软件对砂浆样品进行自由网格划分。该方法的优势在于能够在孔隙、裂缝等复杂特征区域自动生成加密网格。以 25 ℃条件下 M30 砂浆的模拟为例，设计了三种不同的最小网格尺寸（0.351 μm、0.874 μm 和 1.563 μm）。以峰值应力和峰值应变为评价指标，验证结果如下：

精度分析表明，随着最小网格单元尺寸增大，25 ℃下 M30 砂浆的模拟误差显著上升。最小网格尺寸为 0.351 μm 时，误差仅为 1.87%；最小网格尺寸为 0.874 μm 时，误差升至 2.97%，仍处于工程精度可接受范围内；当最小网格尺寸增大至 1.563 μm 时，误差骤增至 13.42%，应力–应变曲线与实验值偏离明显，无法准确反映材料的力学行为。

计算效率分析：网格单元尺寸越小，所需计算单元数量呈几何级数增长。最小网格尺寸为 0.351 μm 时，模型需要 3 159 383 个四面体单元，单次计算耗时约 13 h，计算时间与资源消耗大幅增加。相比之下，最小网格尺寸为 0.874 μm 时，计算耗时可降至约 5 h。该方案在保证计算精度完全满足工程实际需求的同时，有效降低了计算成本，实现了 “精度达标、效率可控” 的目标。

基于砂浆经 200 ℃高温实验后孔隙率与抗压强度变化极小，而 1000 ℃时砂浆内部结构受损严重、强度基本丧失，因此不对这两个温度开展模拟分析。本次模拟选取的温度为 25 ℃、400 ℃、600 ℃和 800 ℃。

不同温度下砂浆的峰值应力差异显著，以 M30 为例：25 ℃时峰值应力约为 47 MPa，而 1000 ℃时仅约 2.5 MPa，相差近 19 倍。为最大程度凸显各温度下的应力梯度与损伤集中区域，本研究根据各温度下样本的实际应力范围优化后，采用独立的动态色标。图 11 和图 12 分别为 M30 和 M60 在不同峰值应变下的应力云图演变规律。为便于观察模型变形，对砂浆模型的变形进行了放大处理。

图11. 不同温度下、不同峰值应变时 M30 的应力图

图12. 不同温度下、不同峰值应变时 M60 的应力图

观察结果表明：（1）模型所承受的最大应力高于实验测得的峰值应力，原因是模型中存在大量孔隙结构，导致压缩模拟过程中出现应力集中。（2）以 εₚ表示峰值应变。样本在 0.5εₚ之前应力快速增长；在 0.5εₚ～1.0εₚ之间，模型进入塑性变形阶段，应力增长速率较前一阶段放缓；当应变达到 1.5εₚ时，样本所能承受的最大应力下降，此时样本中部应力明显低于上下部位，表明样本中部结构已发生损伤，这与实验破坏模式高度吻合。（3）试件的极限抗压强度随温度升高呈单调下降趋势，与实验结果一致。

模拟得到的峰值应力与峰值应变同实验实测值存在一定偏差，主要原因如下：（1）本构关系的理想化导致模型简化，忽略了材料不均匀性与内部含水率的影响。（2）网格划分存在单元过渡不连续问题，可能无法精准捕捉应力梯度。（3）混凝土塑性损伤（CDP）模型中材料本构参数的理论假设，以及对内部孔隙缺陷的表征不够精确。 因此，模拟结果与试验结果之间存在一定误差，且该误差不可避免。 

图 13 为不同温度下单轴压缩实验与有限元数值模拟得到的砂浆应力–应变曲线对比。可以看出，有限元模拟得到的砂浆单轴压缩应力–应变曲线变化趋势与试验结果基本一致。在 25 ℃、400 ℃、600 ℃、800 ℃条件下，模拟得到的 M30 峰值应力分别为 45.612 MPa、33.011 MPa、15.959 MPa、5.333 MPa；M60 则分别为 79.501 MPa、50.101 MPa、41.599 MPa、10.467 MPa。与实验结果相比，M30 峰值应力的偏差依次为−2.974%、−0.246%、−2.564%、7.972%，平均偏差为 2.033%；M60 则分别为 6.251%、7.676%、4.714%、6.545%，平均偏差为 6.296%。模拟峰值应力与实验值十分接近，误差较小。

模拟得到的峰值应变：M30 分别为 0.323%、0.554%、0.868%、0.923%；M60 分别为 0.372%、0.484%、0.742%、0.877%。 与实验结果相比，M30 在 25 ℃、400 ℃、600 ℃、800 ℃ 下的峰值应变偏差分别为 -9.524%、-0.727%、-5.372%、-1.988%，平均偏差为 0.109%；M60 对应偏差分别为 -1.639%、0.820%、4.381%、4.466%，平均偏差为 2.007%。 峰值应变模拟值与实验值同样较为接近。由此可见，不同温度下砂浆的有限元模拟应力–应变曲线与实验应力–应变曲线具有良好的一致性。 

图13. 不同温度下砂浆单轴压缩试验与数值模拟得到的应力 - 应变曲线

图 14 和图 15 分别展示了 M30 与 M60 砂浆在不同温度、不同峰值应变下的损伤劣化规律。图 14 和图 15 中的模型同样进行了变形放大处理。孔隙是砂浆内部的薄弱部位。 

图14. 不同温度、不同峰值应变下 M30 的损伤劣化行为

图15. 不同温度、不同峰值应变下 M60 的损伤劣化行为

在达到 0.5εₚ之前，仅 800 ℃试件表面出现轻微损伤。在 0.5εₚ～1.0εₚ阶段，较大孔隙首先出现大范围损伤，随着位移荷载增大，细小孔隙也开始产生损伤。在 25 ℃条件下，砂浆仅在 1.0εₚ时在孔隙及其附近小范围区域出现损伤。随着温度升高，砂浆的损伤程度逐渐加剧，在 1.0εₚ时呈现出损伤扩展、贯通与劣化的过程。当位移荷载达到 1.5εₚ时，样本中部区域几乎完全发生损伤。 研究表明，砂浆样本在荷载作用下的损伤分布随温度升高逐渐加剧。高温作用前后，M30 与 M60 砂浆在 1.5εₚ时的损伤均主要集中在样本中段，这与高温处理后砂浆样本单轴压缩试验中观察到的裂缝发展规律一致。 与 M30 结果相比，相同温度下，当应变小于 1.0εₚ时，M60 砂浆的损伤范围更小。这是由于强度等级更高的砂浆内部结构更密实，具备更强的抗变形能力。当应变超过 1.0ε后，M60 砂浆的损伤分布更为显著，这是因为高强度砂浆的脆性特征表现得更为突出。