mesoScale_Mesh混凝土细观模型程序
木木最近在学习二维问题的非结构化网格划分算法，基于开源程序Mesh2D：
https://github.com/dengwirda/mesh2d
做了一系列的案例，掌握了一些网格划分时的细节问题，比如密度尺寸控制、网格优化、几何分区等实用技巧。
有了这些技能，我用 Matlab 做了一个混凝土细观模型 mesoScale_Mesh，里面包含了骨料、纤维、孔洞的投放，投放完成后，将 Mesh2D 的网格划分技术加入进来，就实现了对于内部元素自适应加密的过程，效果如下：
现介绍一下mesoScale_Mesh程序的主要功能：
随机投放各级骨料（带界面过渡区）、纤维、孔洞 模型区域支持矩形、任意多边形、圆形、半圆形 自定义网格尺寸函数，对骨料界面过渡区、纤维尖端区、孔洞周围进行自适应加密 网格剖分算法基于 meshd2d 通过修改json配置文件，进行参数修改 可视化：投放后的几何域、网格生成过程动图、网格优化过程动图、最终网格形态图，自动保存 gif 文件 导出 inp 文件在 output 文件夹下，将各级骨料、纤维、孔洞分别写入不同的单元集 合，供 Abaqus 进一步研究分析 保存高清矢量图.eps格式：程序运行完毕后，自动保存投放后的几何域和最终网格形态图在 output 文件夹下
json 文件
{
"domain":{
"type":"poly",
"points":[
[0,0],
[4,0],
[4,2.5],
[2,3],
[0,2.5]
]
},
"itz_thk":0.01,
"levels":[
{"name":"fine","r_min":0.08,"r_max":0.14,"vol_target":0.04,"safe_margin":0.006},
{"name":"medium","r_min":0.14,"r_max":0.20,"vol_target":0.06,"safe_margin":0.006},
{"name":"coarse","r_min":0.20,"r_max":0.28,"vol_target":0.08,"safe_margin":0.006}
],
"fibers":{
"N":10,
"L_min":0.20,
"L_max":0.35,
"angle_range":[0,3.141592653589793],
"safe_margin":0.01
},
"pores":{
"N":10,
"r_min":0.03,
"r_max":0.05,
"safe_margin":0.01
},
"mesh_density":{
"h_bg":0.15,
"h_fine":0.015,"w_fine":0.07,
"h_med":0.03,"w_med":0.09,
"h_coa":0.05,"w_coa":0.11,
"h_tip":0.008,
"r_tip":0.02,
"R_tip":0.06
}
}
domain
矩形区域
"domain":{
"type":"rectangle",
"param":{"width":4.0,"height":2.5}
}
多边形区域
"domain":{
"type":"polygon",
"param":{"vertices":[[0,0],[4,0],[4,2.5],[0,2.5]]}
}
圆形区域
"domain":{
"type":"circle",
"param":{"center":[0.0,0.0],"radius":4}
}
半圆区域
"domain":{
"type":"halfcircle",
"param":{"center":[2.0,1.25],"radius":1.2,"angle":[0,3.14]}
}
骨料设置
{"name":"fine","r_min":0.05,"r_max":0.07,"vol_target":0.05,"safe_margin":0.02}
name：级别名称（如 "fine", "medium", "coarse"）,可显示在 inp 文件中单元集 合的名称 r_min：本级骨料最小半径 r_max：本级骨料最大半径 vol_target：本级骨料目标体积分数（0~1），如 0.05 表示 5% safe_margin：本级骨料对任何对象的最小安全间距（中心距）
纤维设置
"fibers":{
"N":10,
"L_min":0.05,
"L_max":0.1,
"angle_range":[0,3.141592653589793],
"safe_margin":0.01
}
N：纤维根数 L_min：纤维最小长度 L_max：纤维最大长度 angle_range：纤维方向角度范围（弧度）[起始, 终止]，如 [0, 3.14] 表示仅水平方向及上半圆 safe_margin：纤维与其它所有对象（骨料/孔/纤维/边界）的最小安全距离
孔洞设置
"pores":{
"N":10,
"r_min":0.03,
"r_max":0.05,
"safe_margin":0.01
}
N：孔洞数量 r_min：孔洞最小半径 r_max：孔洞最大半径 safe_margin：孔洞与其它所有对象（骨料/孔/纤维/边界）的最小安全距离
网格密度设置
"mesh_density":{
"h_bg":0.15,
"h_fine":0.015,"w_fine":0.07,
"h_med":0.03,"w_med":0.09,
"h_coa":0.05,"w_coa":0.11,
"h_tip":0.008,
"r_tip":0.02,
"R_tip":0.06
}
h_bg：背景区域最大三角形边长 h_fine/h_med/h_coa：细/中/粗骨料附近的网格目标边长 w_fine/w_med/w_coa：细/中/粗骨料影响半径 h_tip：纤维端点极细局部网格目标边长 r_tip：纤维端点极细区域半径 R_tip：纤维端点影响区最大半径
骨料投放顺序
先投放第一级骨料
（如“fine/细”），填满后（或达到最大尝试次数/面积不足） 再投放第二级骨料
（如“medium/中”） 依次类推，直到最后一级（如“coarse/粗”）
每一级投放时，都会避开所有已经投放的前级骨料（包括ITZ安全距离），以保证“不同级别骨料之间不会重叠”。
随机生成骨料、纤维、孔洞原理
骨料投放原理
循环试探（试投法）
每一个骨料
都通过“在区域内随机生成中心坐标和半径”产生一个候选圆形 对于每一个候选骨料，都要经过一系列碰撞与距离检测，包括： 距离边界的最小安全距离（不能出界） 距离已投放骨料的最小安全距离（防止重叠、保证ITZ） 距离孔洞、纤维的安全距离（防止与其他对象碰撞）
体积分数达标
程序会持续不断地“试投”骨料，直到当前级别的总面积（体积分数）达到目标，或达到最大尝试次数（以防死循环） 每个级别都单独投放，已投放的骨料不会再被覆盖
随机性保证均匀分布
通过rand产生均匀分布的随机坐标和半径
伪代码
for 每一级骨料 level = 1 : n_level
while (当前级体积分数未达标) && (尝试次数 < 最大尝试数)
随机生成半径 r ∈ [r_min, r_max]
随机生成中心坐标 (cx, cy) ∈ 区域内部
# 检查所有碰撞条件：
if 与区域边界的距离 ≥ r + ITZ + safe_margin
且 与已投放骨料的距离 ≥ r_new + r_exist + ITZ + safe_margin
且 与孔洞、纤维的距离也满足要求
# 判定通过
录用该骨料: 存储(cx, cy), r, 级别类型
更新已投放面积
end
end
统计本级实际体积分数
end
孔洞投放原理
与骨料类似，也是循环试投 每次随机一个中心、半径 检查是否与骨料（含ITZ）、其他孔洞、区域边界有碰撞/超界，全部通过后才正式放入 投放数量为用户指定的N个孔洞，投放失败会提前终止并发出警告
伪代码
fori = 1 : N_pore
placed = false
for n = 1 : 最大尝试数
随机生成半径 r ∈ [r_min, r_max]
随机生成中心 (cx, cy) ∈ 区域内部
# 检查所有碰撞条件：
if 与区域边界的距离 ≥ r + safe_margin
且 与所有骨料的距离 ≥ (r骨料 + ITZ) + r孔 + safe_margin
且 与其他已投放孔洞的距离 ≥ r1 + r2 + safe_margin
# 判定通过
存储该孔洞: (cx, cy), r
placed = true
break
end
end
if ~placed
警告“本孔洞投放失败，已完成i-1个”
跳出
end
end
纤维投放原理
随机投放纤维线段
随机产生纤维的中心点、长度、方向角度（全部范围可控） 由中心点+长度+方向，计算出纤维两端点的坐标
重叠&越界判定
判断两端点是否都在基体区域内 判断与骨料、孔洞的最小距离（不能穿透） 判断与已投放纤维的最近距离，防止纤维重叠或靠得太近（利用线段最短距离算法）
试投限制
若连续尝试多次仍无法投放，则终止并发出警告
伪代码
fori = 1 : N_fiber
placed = false
for n = 1 : 最大尝试数
随机生成长度 L ∈ [L_min, L_max]
随机生成方向角 θ ∈ [angle_min, angle_max]
随机生成中心 (cx, cy) ∈ 区域内部
根据(cx, cy, L, θ) 计算两端点 p1, p2
# 检查所有碰撞条件：
if p1 和 p2 都在区域内部
且 与所有骨料的距离 ≥ ITZ + safe_margin
且 与所有孔洞的距离 ≥ safe_margin
且 与已投放纤维的距离 ≥ safe_margin
# 判定通过
存储: p1, p2
placed = true
break
end
end
if ~placed
警告“本纤维投放失败，已完成i-1根”
跳出
end
end
所有对象的“安全距离”处理
所有骨料、孔、纤维在“碰撞判定”时，都加上safe_margin（安全间隙），可有效防止数值重叠、网格穿透等问题 ITZ厚度会体现在骨料与其它对象的碰撞半径上（骨料实际占据空间会扩大 ITZ 厚度）
网格质量检查
可通过导出 inp 文件至 Abaqus，进行网格质量检查，从下图可以看到模型中的每一个单元质量都是合格的，没有报错或警告⚠️。
案例
目前程序支持矩形区域、任意多边形区域、圆形区域、半圆形区域，感兴趣的朋友可在此基础上开发更多几何区域，开发的过程中，重点关注如何判断内部元素不超出边界。
矩形区域
任意多边形区域
 |  |
|---|---|
 |  |
圆形区域
 |  |
|---|---|
半圆形区域
补充
由于 mesh2d 程序本身没有设置每次迭代步的网格动态变化图，我在原有 refine2 和 smooth2 文件的基础上增添一些语句，使得看到每次迭代后的网格生成过程。将原有文件保存为副本，望拿到我整套程序的朋友们知晓。
项目目录说明：
📦<YourProjectFolder>/
┣ 📂mesh2d/ % 网格生成函数
┣ 📂output/ % 所有输出结果保存目录（如图片、INP 文件等）
┣ 📜config.json % 仿真/几何配置参数（如颗粒大小、分布密度等）
┣ 📜gen_random_agg_pore_fiber.m % 随机生成骨料、孔洞、纤维的核心函数
┣ 📜hfun8.m % 网格尺寸控制函数（用于 mesh2d）
┣ 📜mesoScale_Mesh.m % 主程序入口，执行完整流程（调用生成、网格、输出等）
┣ 📜write_inp.m % 将网格/材料/边界写入 Abaqus INP 文件的脚本
后续有时间我会更新有关 mesh2D 项目的网格划分教程。
以上整套程序的源代码已经存放在知识星球中，感兴趣的小伙伴可在星球中下载。后台回复：星球即可扫码加入。
参考文献
[1]：D. Engwirda, Locally-optimal Delaunay-refinement and optimisation-based mesh generation, Ph.D. Thesis, School of Mathematics and Statistics, The University of Sydney, http://hdl.handle.net/2123/13148, 2014.
[2]：D. Engwirda, Unstructured mesh methods for the Navier-Stokes equations, Honours Thesis, School of Aerospace, Mechanical and Mechatronic Engineering, The University of Sydney, 2005.
