Fortranを使って物理現象をコンピュータ上で再現してみよう！💻
これまでのステップでFortranの基本から応用まで学んできましたね。いよいよ集大成として、Fortranが得意とする科学技術計算の中でも代表的な「シミュレーション」に挑戦してみましょう。 ここでは、比較的シンプルな例として「熱伝導」と「粒子運動」のシミュレーションを取り上げます。
物理法則や数学モデルに基づいて、現実の現象をコンピュータ上で模擬実験することです。Fortranはその数値計算能力の高さから、古くから様々なシミュレーションに用いられてきました。
1. 熱伝導シミュレーション 🔥🧊
物体の温度が時間とともにどのように変化するかをシミュレーションします。ここでは、1次元の棒における熱の伝わり方を考えてみましょう。 簡単な方法として有限差分法を用います。これは、空間と時間を小さな区間に区切り（離散化）、各点での温度変化を計算していく手法です。
考え方（1次元・陽解法）
棒をいくつかの点（格子点）に分割し、各点の温度を配列で管理します。ある時刻 `t` における格子点 `i` の温度を `T(i, t)` とすると、次の時刻 `t+dt` における温度 `T(i, t+dt)` は、隣接する点の温度を使って以下のように近似的に計算できます（陽解法）。
T(i, t+dt) = T(i, t) + alpha * dt / dx**2 * (T(i+1, t) - 2*T(i, t) + T(i-1, t))
ここで、`alpha` は熱拡散率、`dt` は時間ステップ幅、`dx` は格子点の間隔です。これをすべての格子点について、時間を進めながら繰り返し計算します。
Fortranコード例 (1次元熱伝導 – 陽解法)
簡単な一次元熱伝導のシミュレーションコードの例です。棒の両端の温度は固定されているとします。
program heat_conduction_1d
implicit none
integer, parameter :: N = 100 ! 格子点の数
real, parameter :: L = 1.0 ! 棒の長さ
real, parameter :: D = 0.01 ! 熱拡散率 (alpha)
real, parameter :: T_end = 100.0 ! シミュレーション終了時間
real, parameter :: dx = L / real(N-1) ! 格子間隔
real, parameter :: dt = 0.1 * dx**2 / D ! 時間刻み (安定性条件から決定)
integer, parameter :: nsteps = int(T_end / dt) ! 時間ステップ数
integer :: i, n
real, dimension(0:N+1) :: T, T_new ! 温度配列 (境界条件用に+2)
real :: time
! 初期条件: 中央部のみ高温、両端は0度
T = 0.0
T(N/2) = 100.0
! --- 時間発展ループ ---
time = 0.0
do n = 1, nsteps
time = time + dt
! 内部点の温度更新 (陽解法)
do i = 1, N
T_new(i) = T(i) + D * dt / dx**2 * (T(i+1) - 2.0 * T(i) + T(i-1))
end do
! 境界条件 (両端0度固定)
T_new(0) = 0.0
T_new(N+1) = 0.0
! 温度配列の更新
T = T_new
! --- 定期的に結果を出力 (例: 100ステップごと) ---
if (mod(n, 100) == 0) then
write(*, '(F8.3, 5F10.3)') time, T(1), T(N/4), T(N/2), T(3*N/4), T(N)
end if
end do
write(*,*) 'Simulation finished.'
end program heat_conduction_1d
⚠️ 注意点
2. 粒子運動シミュレーション 🚀
質点（粒子）が力（例えば重力）を受けてどのように運動するかをシミュレーションします。ここでは、簡単な例として、重力下での投げ上げ運動を考えます。 基本的な方法としてオイラー法を用います。これは、微小時間 `dt` 後の速度と位置を、現在の速度と加速度（力から計算）を使って近似的に求める方法です。
考え方（オイラー法）
時刻 `t` における粒子の位置を `x(t)`、速度を `v(t)`、加速度を `a(t)` とします。加速度 `a(t)` は、粒子にはたらく力 `F(t)` と質量 `m` から `a(t) = F(t) / m` で求まります。 オイラー法では、微小時間 `dt` 後の速度と位置を以下のように近似します。
v(t+dt) = v(t) + a(t) * dt
x(t+dt) = x(t) + v(t) * dt
これを時間を進めながら繰り返し計算します。
Fortranコード例 (鉛直投げ上げ – オイラー法)
地面から初速度 `v0` で真上にボールを投げ上げるシミュレーションの例です。空気抵抗は無視します。
program particle_motion_1d
implicit none
real, parameter :: g = 9.8 ! 重力加速度 (m/s^2)
real, parameter :: v0 = 20.0 ! 初速度 (m/s)
real, parameter :: t_end = 5.0 ! シミュレーション終了時間 (s)
real, parameter :: dt = 0.01 ! 時間刻み (s)
integer, parameter :: nsteps = int(t_end / dt) ! 時間ステップ数
integer :: n
real :: time, y, v ! 時刻, 高さ, 速度
! 初期条件
time = 0.0
y = 0.0
v = v0
! ヘッダー出力
write(*, '(A8, A15, A15)') 'Time(s)', 'Height(m)', 'Velocity(m/s)'
write(*, '(F8.3, F15.6, F15.6)') time, y, v
! --- 時間発展ループ ---
do n = 1, nsteps
! 加速度 (一定: -g)
real :: a = -g
! 速度と位置の更新 (オイラー法)
v = v + a * dt
y = y + v * dt ! 更新前の速度 v(t) を使うことに注意 (より正確には v(t+dt)*dt や平均を使う方法もある)
! 時間更新
time = time + dt
! 地面より下に行ったら停止 (簡単な衝突判定)
if (y < 0.0) then
y = 0.0
v = 0.0
write(*, '(F8.3, F15.6, F15.6)') time, y, v
write(*,*) 'Hit the ground!'
exit ! ループを抜ける
end if
! --- 定期的に結果を出力 (例: 10ステップごと) ---
if (mod(n, 10) == 0) then
write(*, '(F8.3, F15.6, F15.6)') time, y, v
end if
end do
write(*,*) 'Simulation finished.'
end program particle_motion_1d
💡 改善のヒント
- 精度向上: オイラー法は最も単純ですが、誤差が蓄積しやすいです。より精度の高いルンゲ＝クッタ法などを使うと、より現実に近い結果が得られます。
- 多次元へ: 2次元（放物運動）や3次元に拡張することも可能です。位置、速度、力をベクトルとして扱います。
- 相互作用: 複数の粒子間の力（万有引力、クーロン力など）を考慮すると、より複雑で興味深いシミュレーション（N体シミュレーション）が可能です。
3. シミュレーションを行う上での考慮事項 🤔
- 数値計算手法の選択: 問題に応じて適切なアルゴリズム（陽解法/陰解法、オイラー法/ルンゲ＝クッタ法など）を選ぶ必要があります。精度、安定性、計算コストを考慮します。
- パラメータ設定: 時間刻み `dt` や格子間隔 `dx` は、計算の精度や安定性に大きく影響します。小さすぎると計算時間が膨大になり、大きすぎると誤差が大きくなったり不安定になったりします。
- 境界条件: シミュレーション領域の端で物理量がどう振る舞うかを設定します（例：温度固定、断熱、周期的境界など）。
- 可視化: シミュレーション結果は数値の羅列であることが多いです。Step 7で学んだGnuplotなどのツールを使ってグラフやアニメーションにすることで、結果を理解しやすくなります。
- 検証: 解析解（もしあれば）と比較したり、物理的な保存則（エネルギー保存則など）が成り立っているかを確認したりして、シミュレーションの妥当性を検証することが重要です。
4. 次のステップと参考情報 📚
今回紹介したのは非常に簡単な例ですが、Fortranによるシミュレーションの基本的な流れは掴めたでしょうか？
- 発展: 熱伝導を2次元や3次元に拡張する、粒子の相互作用（衝突、引力など）を導入するなど、モデルを複雑にしてみましょう。
- ライブラリ活用: より高度な数値計算（連立一次方程式の解法、高速フーリエ変換など）には、BLAS/LAPACKやFFTWといった最適化されたライブラリの利用も検討しましょう（Step 7参照）。
- 学習リソース:
- Fortran-lang 公式サイト: Fortranに関する最新情報やチュートリアルがあります。
https://fortran-lang.org/
https://fortran-lang.org/learn/ (学習リソース) - 数値計算の教科書: 有限差分法、有限要素法、数値積分法など、シミュレーションに使われる様々な手法について詳しく解説されています。大学の講義資料なども参考になります。
- オンラインチュートリアル: “Fortran simulation tutorial”, “Finite Difference Method Fortran”, “Particle Simulation Euler method Fortran”などで検索すると、様々な解説記事やコード例が見つかります。
- Fortran-lang 公式サイト: Fortranに関する最新情報やチュートリアルがあります。
Fortranを使ったシミュレーションは奥が深く、様々な分野で活用されています。ぜひ興味のある現象を対象に、自分だけのシミュレーションプログラムを作成してみてください！✨
コメント