-- Views
August 24, 26
スライド概要
流体解析プログラムあるいはスタガード格子による偏微分方程式解法ライブラリー basilisk の紹介
basilisk を使用してみた ・なぜ興味を持ったか yotakagi77さん Qiita “GerrisをUbuntuにインストール”(2021年10月23日) moritam51さん 第88回オープンCAE勉強会@関西(2022年2月5日) Basiliskによる軸対称milkCrown解析 (Openscc) #staggered grid 使用 my first CFD 任意形状への対応で挫折 VOF、IBM (Immersed Boundary Method) との親和性 詳細なウェッブサイト https://basilisk.fr 第20回OpenCAE勉強会@横浜 2026/08/22 1
Basilisk 開発者 Stéphane Popinet Centre National de la Recherche Scientifique (CNRS) Institut Jean le Rond d’Alembert、ソルボンヌ大学 有限差分法 PDEソルバー、staggered grid 前身ソフト gerris 特徴 VOFに強い 層流 ソース(C) 公開、研究者用コード 研究グループで使用したコード&開発中のコード (sandbox) UGM … 2023, 2025, 2027? 2
インストール インストールページ https://basilisk.fr/src/INSTALL インストールディレクトリー basilisk |-- src |-- ライブラリー *.h | ユーティリティー |-- example 例題 *.c |-- test サンプル *.c |-- tutorial 練習、ユーザー問題 ・ インストールされるプログラム Basilisk C; basilisk 専用プレコンパイラ 起動コマンド qcc (プレコンパイル → C) 3
解析実行のあらまし 1.個別問題向きのソース(*.c) を作成 形状設定&格子形成 初期条件 境界条件 デフォールトはノイマン&フリースリップ 実行のコントロール(出力、終了判定) 2.コンパイル(qcc、 executable; ex. bump) make を使って簡易化可能 (所要オプション付加) 3.実行 ./bump > out.ppm 2> log y 解析領域 (デフォールト) 1 0 1 x 4
チュートリアル問題
*サンブナン(浅水波方程式)の解析 2D
#include “saint-venant.h” サン-ブナン方程式の
ライブラリー
ソルバーを使用する
#define LEVEL 8
メッシュの分解能 2^8
(ほかに境界条件等)
~Events(次頁)~
解析のコントロール,出力
等など
int main() {
}
origin (-0.5, -0.5);
原点の移動 (スケール)
init_grid (1 << LEVEL);
メッシュ生成
run();
実行
ほかに、物性値等のパラメータも設定
5
チュートリアル問題 2
Event の記述
*初期値設定
タイミング
event init (t = 0) {
double a = 1., b = 200.;
foreach()
h[] = 0.1 + a*exp(- b*(x*x + y*y));
}
*名前はユーザーが付ける
*終了時の操作
event end (t = 4) {
printf ("i = %d t = %g\n", i, t);
}
6
チュートリアル問題 3 *時間ステップ毎の最小・最大値の出力 → gnuplot 等でプロット event graphs (i++) { stats s = statsf (h); h(水深)の統計量 fprintf (stderr, "%g %g %g\n", t, s.min, s.max); } *一定時間毎の画像の出力 → アニメーション作成など event images (t += 4./300.) { output_ppm (h, linear = true); scalar l[]; foreach() l[] = level; static FILE * fp = fopen ("grid.ppm", "w"); output_ppm (l, fp, min = 0, max = LEVEL); } output_ppm は ppm 形式の画像を作成する関数 1格子点⇔1ピクセル 7
チュートリアル問題 4 *mesh adaption の設定 event adapt (i++) { adapt_wavelet ({h}, (double []){4e-3}, maxlevel = LEVEL); } 水深 格子サイズ ppm → mp4 変換済 8
ライブラリーと関数の調べ方 ・ ライブラリーと関数のマニュアル reference manual 的なものはない website にソースコードと説明 サイト内検索機能 → 使用例、機能、引数、デフォールト値 srcディレクトリー grep “output_ppm” *.h ライブラリー(ソルバーとユーティリティー) saint-venant.h、 navier-stokes/centered.h two-phase.h、 two-phase-clsvof.h two-phase-levelset.h utils.h、 output.h、 view.h、lambda2.h 等 例えば、output.h で定義されている関数 output_field、output_metrix、output_ppm 等 9
例1静止液滴(無重力) test/rising.c を修正 scale SCALE 0.01、0.001、0.0001 m (直径 2., 0.2, 0.02mm) 0.2*scale y scale x 空気 密度 1.2kg/m3、 粘度 1.6E-5 Pas 水 密度1000kg/m3、粘度 0.001 Pas 表面張力係数 σ= 0.071 N/m 重力加速度 g = 0 m/s2 メッシュ分割 LEVEL 6,7,8(64,128, 256) 液滴 13, 25, 51 10
静止液滴の解析結果(VOF値) SCALE = 0.001 LEVEL 8 51/256 7 6 21/128 13/64 PPM → MP4 11
Spurious Current の変動 SCALE= 0.001, LEVEL=8 Time [0.0、0.01] 0.0001 12
Spurious Current の変動 SCALE= 0.001, LEVEL=7 Time [0.0、0.01] 0.0001 13
Spurious Current の変動 SCALE= 0.001,LEVEL=6 Time [0.0、0.01] 0.0001 14
Adaptive Mesh(VOF) LEVEL 8 Adapted(6~8) PPM → MP4 15
静止液滴の解析結果 LEVELの変動 メッシュ分割(瞬間値) 16
Spurious Current の変動 SCALE= 0.001, LEVEL=6~8 Time [0.0、0.01] 0.0001 scale;0.01 17
Spurious Current の変動 SCALE= 0.01, LEVEL=6 Time [0.0、0.1] 0.001 scale;0.1 18
Spurious Current の変動 SCALE= 0.0001, LEVEL=8 Time [0.0、0.001] 0.00001 scale;0.001 19
例2 液滴の自由落下 モデル 初期位置 scale y g 0.5x scale scale 想定落下位置 0.2*scale x 2xscale SCALE [m] 落下時間 [sec] D [mm] 0.01 0.0452 2 0.001 0.0143 0.2 g(m/s2) 0.0001 0.00452 0.02 9.81 20
任意形状物体の導入 Embeded boundary methods を実装 (embed.h) P. Schwartz, M. Barad, P. Colella, and T. Ligocki. A cartesian grid embedded boundary method for the heat equation and poisson’s equation in three dimensions. Journal of Computational Physics, 211(2):531–550, 2006 物体形状の”符号付“ 距離関数(φ、phi)を与えれば、 体積占有率、面積占有率を求める機能が用意されている。 distance.h、fractions.h φを指標にして AMR することで、 物体近傍メッシュを細分化 surface fraction (vector) volume fraction (scalar) n 21
距離関数の求め方 ・自前 ex. examples/naca2424-starting.c 格子点で関数値を与えることも可能 test/porous (porous3D) ・式を与えるとφを計算する関数 solid examples/sphere x2 + y2 –R 2 ・ stl ファイルから φ を計算 examples/Tangaroa.c、 sandbox/ysaade/f1/f1_W09.c ▲物体の境界条件はvolume fraction のみ考慮 surface fraction は計算可能 examples/distance.c embed.h との整合性 (sandbox に改良版が載っている) 22
例1 球を過ぎる流れ(Re=300) examples/rsphere.c 球 λ2 の等値面に y-z面の渦度で着色 z y x 23
例1 球を過ぎる流れ(Re=300) 2次元断面 X方向流速のコンター x 球 24
例1 球を扁平に(Re=300) 扁平 λ2の等値面に y-z面の渦度で着色 25
例2 調査船の周りの気流と波 examples/tangaroa.c Tangaroa ニュージーランドの海洋調査船 Basilisk で再現された船の形状 但し、PC上の仮想計算機 basirisk webpabe 元の stl ファイル Basilisk - src/examples/tangaroa.c 26
例2 調査船の周りの気流と波 λ2の等値面 流れ 海水面に 波高(Z座標値)で着色 27
例2 形状の再現性 Basilisk で再現された船の形状 但し、PC上の仮想計算機 元の stl ファイル Basilisk - src/examples/tangaroa.c 28
例2 調査船を球に置き換え λ2の等値面 流れ 海水面に 波高(Z座標値)で着色 29
例3 多孔質物体中の流れ test/rsphere.c ・多孔質物体の再現 解析領域中に円(球)を発生させて pore とする 中心座標と半径に乱数を利用 (一様乱数) オーバーラップしたものは和をとる 上下左右で周期的な構造 連通する流路ができる程度まで球を発生させる -percolation周期的な境界条件を設定 重力により沈降 上下、左右境界で周期的 30
例3 Porous(2D、オリジナル) 粒子数;800、最小格子幅 1/1024(210,level:10) 流速分布 気孔(流体) 固体 重力 31
例3 Porous(2D) 粒子数;800、 最小格子 Level 5→11 32
例3 Porous(2D) 粒子数;800、 最小格子 Level 5→11 33
例3 Porous(2D) 流速分布 粒子数; 800 精細度 level:10 1000 2500 海⇔島 34
纏め ・ basilisk を Ubuntu24(virtual)にインストール Example, Test の例題を実行 ・ インストールとチュートリアル以外、マニュアルはない 例題とライブラリーのソースコードが web で解説されている 関数の内容や引数の意味が簡単にわかるわけではない sandbox 検討中の拡張機能、検討中の問題 ・空気中の浮遊液滴 直径0.2mm以下で spurious current よる液滴の揺動が 顕著 自由落下の軌道に影響 35
纏め(続) ・任意形状物体への対応 距離関数が与えられれば、Embedded Boundary Method (EBM) を適用できる。 物体表面を stl ファイルで与えた場合 正式版では、EBM の適用は部分的 sandbox には、embed.h 等の拡張版がアップされている。 MPI、自由落下 36
計算機環境 全て仮想計算機(Virtual Box) ・U24.04 (win11) ・1並列 37