简介：一份围绕 Basilisk 开源数值模拟框架整理的博士课题学习资料，面向流体力学、地球物理等方向、需要借助 C 语言与 Shell 脚本完成仿真研究的科研人员，也适合从偏微分方程数值求解入门到进阶的个人学习者。资料从有限体积法和谱方法等基础入手，串联科学计算、编译器、调试、版本控制、并行计算与可视化等环节，可帮助梳理从方程离散化到自动化模拟运行的完整链路。压缩包共 14 个文件，约 508KB，以 C 源文件、头文件、Shell 脚本和 Makefile 为主，另含 README 说明文档、LICENSE 与 .gitignore 配置。C 文件与头文件用于定义求解器并输出 VTU 格式结果，Shell 脚本负责容器环境启动与参数扫描等任务自动化，Makefile 简化构建流程。目前已有 211 人浏览学习。通过示例代码和配置脚本，可以借鉴环境搭建、编译调试、并行计算和结果分析的一套完整思路，对开展气泡轴对称/二维模拟等博士研究具有直接参考价值；内容来自网络分享，使用前请自行核对版权与可用性。
1. Basilisk模拟与我的博士研究：C_Shell方案怎么用
博士课题做到多相流阶段，我拿到一个名为“与我的博士相关的Basilisk模拟_C_Shell_下载”的压缩包，里面不是现成的可视化软件，而是一套Basilisk源码工程：两个气泡求解器文件、一个Makefile、两个容器启动脚本，以及几个后处理头文件。Basilisk是开源CFD框架，用C语言编写，配合Shell脚本完成编译、运行和结果整理。它的真正门槛在于C语言宏扩展和事件机制，第一次打开代码会感觉不像标准C，但只要理解了分层逻辑，就能很快改出自己需要的模拟场景。这篇内容面向博士研究生、科研人员以及需要用Basilisk做两相流的工程师，从框架原理讲到容器化运行和调参验证，中间会给出可以直接照抄的命令和参数表。
2. Basilisk框架与C语言实现：树状网格、VOF和事件循环
2.1 C语言宏：Basilisk的“伪面向对象”机制
Basilisk的核心不是传统CFD程序那种类封装，而是通过C预处理宏在编译期生成数据结构。以压缩包中的bubble_2D.c为例，代码里会出现scalar f[];这样的声明，这里的f不是定长数组，而是一个标量场。宏在编译展开后变成带网格单元信息的结构体，通过foreach循环遍历网格单元。这种设计让求解器代码接近纯C性能，同时允许用户在问题文件里像写脚本一样定义物理场。
这里有一段常见的Basilisk风格C代码，示意了如何在一个简单的扩散步中操作场变量：
#include "grid/octree.h" #include "navier-stokes/centered.h" scalar T[]; face vector flux[]; event diffusion_step (i++) { foreach() T[] += dt * (flux.x[] + flux.y[] - flux.x[1,0] - flux.y[0,1]) / sq(Delta); boundary ({T}); }代码中flux.x[]表示单元左侧面的通量，flux.x[1,0]表示右侧面的通量，flux.y[]和flux.y[0,1]则对应下上两个面。Delta是当前网格单元的尺寸，dt是时间步长，sq()是Basilisk提供的平方宏。foreach()遍历当前层级的所有网格单元，boundary({T})负责在自适应网格交界处同步字段值。Basilisk的索引语法看起来像标准C数组，实际上通过宏映射到了树状网格单元，理解这一点后，再去读bubble_axis_symmetric.c中的VOF输运代码就不会发怵。
2.2 事件循环：模拟节奏的控制核心
Basilisk程序通常由main()设置物理参数和初始网格，然后调用run()启动事件调度。事件用event关键字定义，括号里写触发条件。为什么不用普通的while循环？因为科研模拟往往需要在指定时间点输出、自适应加密、改变物理规则，拆成事件可以让这些逻辑互不干扰。
典型的初始化与日志事件如下：
event init (t = 0) { fraction (f, sq(x) + sq(y) - sq(0.1)); } event logfile (i++) { fprintf (stderr, "%g %g %d\n", t, dt, grid->tn); }event init在模拟开始前执行，fraction函数根据隐式方程sq(x)+sq(y)-sq(0.1)=0标记初始气泡区域。event logfile在每一个迭代步后执行，i++表示“每次迭代都触发”。grid->tn是当前网格总单元数，能实时反映自适应加密是否工作。事件机制把时间推进、网格自适应和文件输出解耦，这也是同一个求解器内核可以套用在完全不同的物理问题上的原因。
2.3 VOF方法与自适应网格的组合
气泡模拟属于典型的两相流，Basilisk的two-phase.h模块使用几何VOF方法传输体积分数f。界面处的密度和黏度由f平均，表面张力通过连续表面力模型加入动量方程。自适应加密在界面附近提高分辨率，可以大幅降低计算量。加密逻辑通常放在一个独立事件里：
event adapt (i++) { adapt_wavelet ((scalar *){f, u.x, u.y}, (double[]){1e-2, 1e-2, 1e-2}, 10, 0); }adapt_wavelet是Basilisk的小波误差控制函数，第一个参数指定要加密的字段，第二个参数指定每个字段的绝对容差，第三个参数是最大网格层级MAX_LEVEL，第四个参数是最小层级。容差越小网格越细，最大层级每增加1，理论网格量增加4倍（二维）。这里有一个常见误用：把所有字段容差设成相同值。速度场和体积分数的尺度差异很大，实际调试时建议分开设置，例如f用1e-3，速度场用1e-2，否则会出现界面附近已经很细，但速度梯度大区域仍然捕获不足的情况。
树状网格上的通量和梯度插值比结构化网格复杂，Basilisk用face vector、gradient等宏封装底层细节。对使用者来说，主要收益是可以用很少的代码完成“界面加密、远处粗化”的网格策略，这也是Basilisk在微流体和气泡模拟中流行的重要原因。
3. Shell脚本与容器化：createContainer.sh和startContainer.sh的使用逻辑
3.1 为什么科研项目要把Basilisk装进容器
Basilisk在编译时依赖特定版本的GCC、OpenMPI、GSL以及ffmpeg。不同Linux发行版自带的工具链版本不同，很容易出现编译通过但运行崩溃的情况。我在不同服务器上遇到过几次：Ubuntu 22.04上编译好的可执行文件，放到CentOS 7上直接报缺少libgomp符号。压缩包里的createContainer.sh和startContainer.sh就是为这个问题设计的，它们把工具链和运行环境固定在容器里，宿主机的系统库版本不再影响模拟结果。
常见做法是先在宿主机上准备Docker或Podman环境，然后由createContainer.sh构建一个包含固定编译工具的镜像，startContainer.sh负责创建并进入运行容器，同时把当前工作目录挂载进去。
3.2 createContainer.sh：构建可复现的编译环境
从脚本命名和博士项目习惯看，createContainer.sh负责生成镜像。一个典型实现是这样的：
#!/bin/bash set -e IMAGE=basilisk/phd:1.0 docker build -t "$IMAGE" -f Dockerfile.basilisk . echo "image $IMAGE ready"set -e是Shell脚本里的保险丝：任何一条命令返回值非零都会立即终止脚本，避免用一个损坏的镜像继续往下跑。Dockerfile.basilisk里一般会基于Debian或Ubuntu镜像，安装build-essential、git、gnuplot、libgsl-dev等依赖，然后克隆Basilisk源码并设置环境变量BASILISK。这样构建出来的镜像就是可重复的编译环境。
这里有个细节值得注意：如果docker build执行到一半因为网络原因失败，set -e会阻止脚本继续，但已经下载的镜像层会留在本地。再次执行构建时，Docker会使用缓存，速度会快很多，不是bug。
3.3 startContainer.sh：挂载、权限和并行参数
进入模拟阶段时，startContainer.sh启动交互式容器。从使用场景看，它需要把源码目录挂载进去，同时避免在容器里产生root权限的文件。一个可用模板是：
#!/bin/bash set -e IMAGE=basilisk/phd:1.0 WORKDIR=/sim docker run --rm -it \ --name basilisk_run \ -v "$(pwd)":$WORKDIR \ -w $WORKDIR \ -u "$(id -u):$(id -g)" \ -e OMP_NUM_THREADS="${OMP_NUM_THREADS:-4}" \ $IMAGE bash几个关键点：--rm保证退出后容器被删除，不占用磁盘；-v "$(pwd)":$WORKDIR把当前目录挂载为容器工作目录，模拟产生的结果文件都落在宿主机，方便后续用Paraview查看；-u把当前UID和GID写入容器，避免生成root权限的结果文件；-e OMP_NUM_THREADS控制容器内使用的OpenMP线程数。
进入容器后，就可以按通常流程执行make && ./bubble_2D。如果容器内没有编译过源码，第一次执行make会比较慢，因为Basilisk会把所有依赖头文件都扫描一遍。后续再编译时，只有改动的C文件会被重编。
3.4 Shell脚本参数速查
下面这张表总结了脚本里最常调整的参数，我平时会贴在项目目录下：
| Shell语句 | 作用 | 注意点 |
|---|---|---|
set -e | 命令失败即退出 | 管道场景下要用set -o pipefail配合 |
docker build -f Dockerfile.basilisk | 按指定文件构建镜像 | Dockerfile文件名若非默认路径必须指定 |
docker run --rm | 容器退出后自动删除 | 配合-d后台运行时不适用 |
-v "$(pwd)":/sim | 挂载当前目录 | 路径含空格时务必加引号 |
-u "$(id -u):$(id -g)" | 避免生成root文件 | Windows下要改用Docker Desktop的权限映射 |
-e OMP_NUM_THREADS | 设置OpenMP线程数 | 线程数应小于等于物理核心数 |
这张表实际使用时会反复查。一个常见问题是：在容器中运行mpirun还需要把宿主机的高速网络接口映射进去，通常需要--network host。但博士课题的模拟单机OpenMP已经够用，优先不要上MPI，因为Basilisk的MPI编译要比OpenMP更容易踩坑，而且泡模拟的单节点算力通常足够支撑十万到百万级网格。
4. 气泡模拟实战：bubble_axis_symmetric.c与bubble_2D.c编译调参
4.1 轴对称与二维模型的选型差异
压缩包里的bubble_axis_symmetric.c和bubble_2D.c代表两种不同维度的模拟方式。bubble_2D.c在二维矩形域里模拟圆形气泡的上升，计算量小，适合做参数扫描和数值实验；bubble_axis_symmetric.c通过axi.h把二维域解释为轴对称坐标，得到的流场更接近三维球形气泡。
选型时的一个基本原则：如果只关心单个气泡的终端速度和形状，优先用轴对称模型，它兼顾了二维的计算效率和三维物理特征；如果研究气泡对之间的相互作用，轴对称就无法描述非轴对称扰动，必须用二维或三维全模型。bubble_axis_symmetric.c里的重力方向通常沿轴向，半径方向由y坐标表示，所以初始气泡形状的定义要写成sq(x)+sq(y)，而不是sq(x)+sq(z)。
4.2 Makefile与编译命令
Basilisk使用自带Makefile体系。压缩包中的Makefile一般会写成：
BASILISK ?= $(HOME)/basilisk include $(BASILISK)/Makefile.defs EXECS = bubble_2D bubble_axis_symmetric bubble_2D: bubble_2D.c bubble_axis_symmetric: bubble_axis_symmetric.c include $(BASILISK)/Makefile.rules这里的BASILISK环境变量必须指向Basilisk源码目录，否则使用默认的$(HOME)/basilisk。Makefile.defs负责探测编译器和库，Makefile.rules补充通用编译规则。编译命令很直接：
make bubble_axis_symmetric OMP_NUM_THREADS=4 ./bubble_axis_symmetric > log 2>&1第一条命令生成可执行文件。第二条命令用4线程运行，并把标准和错误输出都重定向到log。Basilisk的进度信息默认写到stderr，所以如果你只重定向标准输出，tail -f log里会什么都看不到。编译时如果出现与fractions_output.h相关的错误，通常是因为当前目录下缺少Basilisk头文件，检查BASILISK路径是否设置正确。
提示：Basilisk的标准输出默认不记录进度，必须把
stderr一起重定向，否则tail -f log里永远是空的。
4.3 参数表与核心调参逻辑
下列参数在两个C文件中基本都能找到，只是初始值不同：
| 参数 | 作用 | 调节方向 |
|---|---|---|
size() | 计算域边长 | 至少为气泡直径的6倍，否则边界效应明显 |
init_grid() | 初始网格数 | 通常取1 << 6，配合后续自适应加密 |
rho1, rho2 | 液体/气相密度 | 密度比过大会导致时间步长急剧变小 |
mu1, mu2 | 液体/气相动力黏度 | 注意Basilisk默认是无量纲方程 |
f.sigma | 表面张力系数 | 减小会抑制界面变形；为零则无界面张力 |
MAX_LEVEL | 最大网格层级 | 每加1，最小网格尺寸减半 |
CFL | 时间步长安全系数 | 默认约0.5，发散时降到0.2 |
G | 重力加速度 | 控制无量纲佛洛德数 |
调参有一个常见误区：直接改MAX_LEVEL就认为会“更准确”。网格加密后，会解析出更小的涡和界面结构，但如果初始扰动或边界条件没有相应变化，结果可能更发散。我的做法是先固定其他参数，只把MAX_LEVEL从6逐步升到10，比较气泡质心速度和界面最大曲率的收敛情况。CFL参数也需要配套调整：网格变细后，同样的物理时间步需要更多迭代步，如果CFL取得过大，显式格式会在计算下一时间步时产生nan。
4.4 运行时的坑与定位手段
实际跑气泡模拟最容易遇到三种情况。
第一种是程序启动后很快输出nan，大概率是初始f.sigma偏小或CFL偏大，表面张力模型无法稳定支撑界面。可以先检查log中最后一个完整的dt值，再把CFL降到0.2重跑。第二种是网格一直加密到最大层级后程序变慢但不报错，这说明adapt_wavelet的容差设得过小，可以放宽速度场容差，例如把u.x, u.y的容差从1e-3放宽到1e-2。第三种是make编译时提示fatal error: common.h: No such file or directory，这通常是因为BASILISK环境变量没有在当前Shell会话里生效，执行：
export BASILISK=/path/to/basilisk make clean make bubble_axis_symmetricmake clean必须执行，否则已经生成的依赖文件会继续引用旧的包含路径，即使环境变量修正了，编译仍会失败。
5. 进阶：VTU输出、并行计算与网格无关性验证
5.1 把模拟结果导出为Paraview可读的VTU文件
压缩包中的output_vtu_foreach.h就是用于输出VTU的。Basilisk官方库里既有output_vtu.h也有output_vtu_foreach.h，后者更适合在OpenMP并行下写文件，因为它把输出动作放在遍历循环内，避免多个线程同时写同一个文件描述符。在C文件里加上一个事件即可：
#include "output_vtu_foreach.h" event vtu_output (t += 0.1) { char name[80]; static int fnum = 0; sprintf (name, "snap-%03d.vtu", fnum++); FILE *fp = fopen (name, "w"); output_vtu_foreach ((scalar *){f, u.x, u.y}, (char *){"f", "ux", "uy"}, fp); fclose (fp); }t += 0.1表示每增加0.1时间单位输出一个VTU文件，fnum保证文件名不重复。输出的文件可以直接用Paraview打开，通过Threshold过滤器提取气液界面。注意这里输出的是Basilisk内部的坐标，一般是无量纲的，与实验对照时必须在后处理阶段乘以特征长度。
5.2 并行计算的线程数选择
在startContainer.sh里通过环境变量设置OMP_NUM_THREADS，但Basilisk运行时读取的就是该变量。线程数不是越多越好，自适应网格的负载分布并不均匀，线程过多会导致内存带宽竞争。单机的经验阈值是8线程以内加速比接近线性，超过后收益明显下降。如果要在多节点上跑，就必须启用MPI重新编译Basilisk，这不是单机项目优先考虑的事。
5.3 网格无关性验证的快速方法
最后给一个可用的验证思路：用脚本循环改写源码中的最大网格级别，跑多个分辨率，然后比较气泡质心高度随时间的变化。可以这样组织：
for level in 8 9 10; do sed -i "s/#define MAX_LEVEL [0-9]*/#define MAX_LEVEL $level/" bubble_2D.c make bubble_2D >/dev/null 2>&1 OMP_NUM_THREADS=4 ./bubble_2D 2>level_$level.log donesed替换源码中MAX_LEVEL的值，前提是文件里已经有#define MAX_LEVEL这一行。如果没有，就在main()函数前加一行。跑完后用awk提取三个日志中气泡质心速度的稳态值，检查相邻两个级别的相对误差。如果相邻两个分辨率的稳态速度差异小于1%，就直接用较细的那个级别出结果，省下的算力足够再多跑一组参数扫描。
本文还有配套的精品资源，点击获取
