実験実施:2026年4月 実験担当:チームレッド
1.はじめに
近年、Pythonの手軽さを活用して、ゲーム製作やデータ分析などの多様な分野にチャレンジする本が多数出版されています。その中で熱流体解析に関する本も何冊か出版されていました。
半導体の発展により、どのPCでも文書編集、E-mail作成などの軽い作業にはスペックをもてあますようになりましたが、普通のPCの演算能力で熱流体解析し有意味な結果を出せるのか?実際の業務に使えるレベルなのか?と疑問をもったことがきっかけで、Pythonでの熱流体解析に挑んでみました。
図1は、今回使用するPCのスペックになります。演算を行うCPUはIntel Core i5-9400、メモリは8GBです。購入当時としてはミドルレンジモデルでしたが、6世代前のCPUですので、今のエントリーモデルと大差ないと思います。

2.評価方法と目標
(1)実用性の評価
解析と並行して解析と同じサイズの炉を作成して温度測定し、熱解析データを比較します。
その結果をもとに、pythonで解析データが実用に耐えうるか評価します。
(2)解析時間(プログラムの動作時間)を評価
プログラムの実行時間として、1日の業務時間(8時間+α)以内で結果を得ることを目標とします。
3.解析モデルと炉
八光電機で製作する炉を模し、密閉空間の底面から加熱する3次元モデルを解析します。
比較実験では、ボード状の断熱材(ロックウール25mm)で6面を覆い、床面にヒーターを設置して加熱する時の自然対流を解析します。内部空間は575mm×560mm×570mmで、中心の温度をφ3.2 ガラス被覆熱電対(Kタイプ) で実測します。
ヒーターは、底面にシリコンラバーヒーター(200×550、定格200V550W)を敷き、100℃で温度調節しました。
熱流体解析においてもヒーターは同寸法とし、表面温度は一律100℃として解析を実施します。



図4が実際に制作した炉です。下部に見える板金は別の用途で使ったものですが、箱の強度を得るのに適していたため流用しています。
この上部にさらに断熱ボードでふたをし、加熱実験を行います。
4.熱流体解析の基礎方程式
・連続方程式
・運動方程式(ナビエ・ストークス方程式)
・エネルギー方程式
5.安定性条件
・流体のCFL条件
CFL<1を満たすように、今回の計算ではΔt=0.0005を設定します。
6.解析1 中央差分法による熱流体解析
6.1.流体(無次元数)パラメーター
Pr(プラントル数) : 0.71
Re(レイノルズ数) : 100
Gr(グラスホフ数) : 100
6.2.境界条件
速度 : 粘着条件反映 → 壁面のu、v、w速度を0に固定
圧力 : 境界のゴーストセル値を0に固定
温度 : 境界のゴーストセル値 = 隣のセルの値(完全断熱)
6.3.解析手法
流体は非圧縮性流体と設定し、Boussinesq近似を適用
解析手法 : MAC法
格子 : スタッガード格子
圧力方程式の収束 : SOR法(反復計算数 : 600回 、 ω=1.7)
移流項のスキーム : 中央差分法
セル:まずは動作確認のため、セル数を50×50×50個に減らして演算を軽量化し、解析を行います。
6.4.解析結果
解析結果は図5のようになりました。(図5の温度は無次元化された数値であり、実温度を表すものではありません。)
チェッカーボード現象と呼ばれる発散が発生し、正常に解析ができませんでした。原因としては、誤差を小さくできるといわれている中央差分法を採用したため、数値拡散が小さく結果的にシミュレーションの安定性が不足したと考えられます。

そのため、次のシミュレーションでは、精度と引き換えに安定性を高めるため、移流項のスキームを中央差分法から「1次精度風上差分法」に変更して解析を行うことにしました。
7.解析2 1次精度風上差分法による熱流体解析
7.1.解析方法、パラメーター
流体パラメーターと境界条件は解析1と同じです。解析手法のみ以下ののように変更しました。
解析手法 : 射影法(Projection method)
移流項のスキーム : 1次精度風上差分法
演算回数:12000Step
このほかの条件は解析1と同じです。
射影法アルゴリズム部分のメインループ部分のプログラムだけ以下に掲載します。射影法には、仮の速度場と速度場の補正を計算するという特徴があります。
for n in range(1 , time_max):
t = thermo_Solver(u, v, w, t, q , _const, heater_mask , wall_loss)
_u , _v , _w = temp_velocity_solver(u, v, w, t, _const)
p , rhs , err_list[n] , iter_list[n] = pressure_Solver(_u , _v , _w , p , _const)
u , v , w = velocity_correction(_u , _v , _w , p , _const)
7.2.解析結果
図6に解析結果を掲載します。左が初期条件(t=0),図の右が定常状態まで計算した結果です。ヒーター周辺が高温となり、上部に向かうにつれて、また中心から離れるにつれて温度が下がるという、実際の挙動に近い解析結果が得られました。チェッカーボード現象も解消され、うまく解析ができたようです。
なお、計算時間は1時間41分かかりました。

8.解析3 Multigrid法による熱流体解析
8.1.解析方法、パラメーター
ここまでは動作確認のため少なめのセル数で演算してきました。
この後、実際のチャンバーサイズ (575mm×560mm×570mm)に合わせてセル数を増やすと計算時間が急増します。そこで解析3では、演算時間を短縮するため、圧力項のソルバーをMultigrid法に変更し、計算スピードを確認します。
パラメーター、境界条件、その他条件は解析2と同じです。(セル数も50×50×50のまま計算します。)
8.2結果
計算時間が解析2の1時間41分から1時間4分に短縮され、約40分の改善が見られました。なお、解析結果自体は解析2と同様なので図は割愛します。
9.解析4 実機サイズでの最終解析
演算時間も改善が見られたため、モデルのセル数を増やして最終解析を実施します。
9.1.流体パラメーター
Pr数、Re数は解析1~3と同じです。
Gr(グラスホフ数)については、対流の影響を確認するため、弱い自然対流(Gr=100)、中程度(Gr=1000)強い自然対流(Gr=10000)まで3パターンで比較検討しました。
9.2.境界条件
速度、圧力の境界条件は試験1~3と同じです。温度の境界条件は、実験で使用する断熱材の熱伝導率を考慮し、エネルギーが外部へ逃げる放熱条件に変更します。
9.3.解析手法
実際の大きさに合わせ、セル数を58×56×58に変更します。そのほかの手法は解析3と同じです。
9.4.解析結果
セルの数の増加により1回のシミュレーションが2時間54分に増加しました。
Gr=10000のケースでは計算が発散し、結果が得られませんでした。
結果が得られたGr=100、Gr=1000の解析結果(無次元化されているデータを実温度に復元したもの)と、実機による測定データをグラフ1で比較します。

青色の線が実機からのデータ、橙・灰色が熱解析のデータになります。
加熱初期に温度が急上昇する傾向は一致し、そのあと実測では板金の温度が上昇するにつれ徐々に空気温度も上昇し、温度安定までは時間がかかっています。シミュレーションでは、板金の熱容量をモデリングしていないため、熱容量の小さい空気の温度がすぐに飽和して一定値に達します。
また、解析のほうが実測よりも低い温度で飽和しましたが、これは解析で熱損失(断熱材から外側への熱移動)を大きく見積りすぎた事が原因ではないかと考えられます。
10.結果のまとめと感想
3次元空間の熱流体解析を、Pythonと普通の事務用PCで結果を得ることができました。また、実測データとの比較で、温度上昇の傾向は推測することができました。
しかし、プログラムの作成・検証などの作業に、はじめの想定よりプログラムのスキルや熱流体解析の知識が必要ということが分かりました。
もともとは、社内に商用CFD(熱流体解析)ソフトを持っており、より手軽に使えるものを求めている中で書籍にであい、教科書を元に必要な部分だけコードを応用すれば短時間で作成できると考えてチャレンジした企画でした。しかし実際にプログラムを作成してみると、プログラムの工数が予測よりはるかに多く、CFDの理論をきちんと理解する必要がありました。また、解析モデルの妥当性も未知数で検証も容易でないため、実際の実務に活かすのは難しいと思いました。
一方、プログラムを実装したことで、CFDの歴史的変遷をはじめ、格子配置、圧力・速度の計算手法に対する理解を深める良い機会となりました。熱流体の解析業務や関連仕事をしている方であれば、解析の知識を深める貴重なチャンスであるため、是非ともお勧めしたいと思います。
11.参考文献
(1) 河村 哲也・佐々木 桃(著) Pythonによる流体解析 朝倉書店2023年
(2) 横浜国立大学教授 博 松井純(著) Pythonによるはじめての数値流体力学 森北出版 2024年
(3) 一般社団法人オープンCAE学会 編(著) OpenFOAMによる熱移動と流れの数値解析(第2版)森北出版 2021年
(4) W. L. Briggs, V. E. Henson, and S. F. McCormick, A Multigrid Tutorial, 2nd ed., SIAM, 2000年
