「シミュレーション」練習課題
1次元熱拡散方程式の初期値境界値問題：
∂u
∂t =a∂2u
∂x2, 0< x <1, t >0.
u(0, t) =u(1, t) = 0, t >0.
u(x,0) =f(x), 0< x <1.
を考える。a >0は拡散係数、f(x)は初期関数（熱の初期分布）を表す。境界条件は、ここでは ゼロ・
ディリクレ境界条件（＝境界で温度0)を課す。
講義で説明した陽的差分スキームにて、上の数値計算をせよ。
なお、Scombにサンプルファイル“heat1dex-pre.c”があるので作業フォルダにダウンロードし、利用し てよい。（差分スキーム部分を各自で書き入れること。）
サンプルでは、a= 1, ∆x= 1/N (N = 50), ∆t= 0.0001とした。(プログラム中では、メッシュサイズ はそれぞれdx, dtと表記。)
練習1.上記設定での陽解法による数値計算プログラムをつくり、数値計算を行なえ。また、計算データを gnuplotにより可視化せよ。
—————————————————————————- この陽的差分スキームでは、∆x,∆t が
a∆t
∆x2 ≤ 1 2
を満たさないと安定に計算ができないことが分かっている。(差分スキームの安定性参照。) 練習2.例えば、∆t= 0.001として計算し、数値解が振動を起こしていることを確認せよ。
追記：最初の数ステップだけ計算して、その結果を表示するようにすると、解が振動を起こして不安定化 していく様子がよく見えると思う。100ステップとか計算すると、たぶん数値的に扱えない大きい値になっ てしまい、可視化がうまくいかないと思う。
—————————————————————————-
練習3.以上の確認が終わったら、考える区間を(0,1)から(0, L) (L >0)とし、両端でディリクレ境界条 件を課した次の問題に対するプログラムを作成せよ。初期関数は各自で適当に考え、計算結果を確認せよ。
(初期関数も区間(0, L)で与えることに注意。)
∂u
∂t =a∂2u
∂x2, 0< x < L, t >0.
u(0, t) =b, u(L, t) =c, t >0.
u(x,0) =f(x), 0< x < L.
(a >0, b, c は定数。)
(プログラムはほとんど変わらないはず。)
練習4. 講義で説明した陰的差分スキームによる計算プログラムを作成せよ。毎時間ステップで連立一次方 程式Ax=bを解く必要が出てくる。行列Aと 右辺のベクトルbがどのようになるか事前に考えてからプ ログラムを作成すること。連立一次方程式の数値解法はなんでもよい。(各自連立一次方程式の数値解法等 を復習・下調べしてくること.)
1
