C++程序实现元胞自动机法模拟枝晶生长：基于偏心正方算法，综合LBM对流作用的影响

原创于 2026-03-14 09:30:00 发布 · 157 阅读

4 ·

本内容遵循CC 4.0 BY-SA版权协议

GEO检测

C++程序，基于元胞自动机法模拟枝晶生长，能实现任意角度（偏心正方算法），同时采用LBM考虑了对流作用对枝晶生长的影响。

枝晶生长的模拟总带着点金属凝固的诗意。咱们今天整点硬核的，用C++同时驾驭元胞自动机和格子玻尔兹曼方法（LBM），看看带对流的枝晶怎么野蛮生长。先扔个效果图镇楼——当枝晶臂遇上湍流，那分叉轨迹比梵高的星空还带劲。

先说核心武器：偏心正方算法。传统CA模型在计算固液界面时容易受网格各向异性限制，咱这个8邻域权重矩阵能玩出花：

// 各向异性权重配置
float anisotropy_weights[8] = {
    cos(4*(theta - 0 * PI/4)),  // 0度方向
    cos(4*(theta - 1 * PI/4)),  // 45度
    //...其他6个方向
};

这个theta就是预设的生长主导角度。比如设成PI/8时，枝晶臂会像钟表指针偏转22.5度生长。实际计算固液界面推进概率时，用这个权重矩阵做卷积，相当于给不同生长方向加了"偏好滤镜"。

重点看相变条件判断的代码段：

if (cell.state == LIQUID && undercooling > threshold) {
    float growth_prob = 0;
    for (int i = 0; i < 8; ++i) {
        growth_prob += neighbor_solid_fraction[i] * anisotropy_weights[i]; 
    }
    growth_prob *= undercooling * deltaT;
    if (rand()/(float)RAND_MAX < growth_prob) {
        cell.state = SOLID;
        release_latent_heat();
    }
}

这里用邻居的固相比例和各向异性权重做概率混合，配合过冷度决定是否凝固。有意思的是这个rand()的调用——微观的随机性会在宏观上形成自然的分形结构，就像现实中的晶体缺陷导致的不可预测生长。

C++程序，基于元胞自动机法模拟枝晶生长，能实现任意角度（偏心正方算法），同时采用LBM考虑了对流作用对枝晶生长的影响。

接下来是LBM部分，处理流体对流。上经典的D2Q9模型：

struct LBMNode {
    float f[9]; // 九个方向分布函数
    float rho, ux, uy; // 宏观量
};

void stream(LBMNode& node) {
    // 传播步骤：将f[i]沿着对应方向移动到相邻节点
    for (int i=0; i<9; ++i) {
        int dx = (i%3)-1;
        int dy = (i/3)-1;
        // 这里需要处理边界条件...
    }
}

但精髓在碰撞过程的BGK近似：

void collide(LBMNode& node) {
    float feq[9];
    for (int i=0; i<9; ++i) {
        float eu = e[i][0]*node.ux + e[i][1]*node.uy;
        feq[i] = w[i] * node.rho * (1 + 3*eu + 4.5*eu*eu - 1.5*(node.ux*node.ux+node.uy*node.uy));
    }
    // BGK松弛
    for (int i=0; i<9; ++i) {
        node.f[i] += (feq[i] - node.f[i]) / tau; 
    }
}

这里藏了个骚操作：用浮力项把温度场耦合进LBM。凝固释放的潜热会影响流体密度，从而产生自然对流。看这行代码：

float buoyancy = beta * (T - T_ref); // 温度差引发的浮力
node.uy += buoyancy * deltaT; // 在碰撞步骤中加入y方向速度修正

当枝晶尖端生长时，局部温度升高产生上升流，周围的低温液体随即补充，形成复杂的涡旋结构。这种流体运动反过来影响溶质输运，让枝晶臂长得歪七扭八——就像把晶体扔进洗衣机里搅拌。

耦合的关键在物质传输方程：

// 溶质浓度更新
cell.concentration += dt * (
    diffusion * laplacian(c.conc)  // 扩散项
    - (u.dot(gradient(c.conc)))    // 对流项
    + solidification_term         // 凝固导致的溶质排出
);

这里把LBM算出的速度场u直接代入对流项，实现双向耦合。有个坑要注意：流体和固体的移动边界处理。我们采用反弹格式，遇到固体节点时让分布函数反向弹回，这样流体就会绕着枝晶结构流动。

最后来个性能优化骚操作——用多网格策略：在枝晶前沿区域用精细网格，远处用粗网格。配合CUDA并行，百万级网格也能实时跑。不过那是另一个故事了...