１次元の非定常熱伝導方程式を表計算ソフトで解く（プチ）

１次元の非定常熱伝導方程式を表計算ソフトで解く

１次元の非定常熱伝導方程式

　　

を表計算ソフトで解けるようにしたにゃ。

ここで、ρは密度、cは比熱、λは熱伝導率、tは時刻、xは位置。

（１）は次のように変形できる。

　　

κは温度伝導率や温度拡散率などと呼ばるもので、熱ではなく、温度の伝わりやすわすさをあらわす物理量。

熱伝導率λは熱の伝わりやすさ、温度伝導率κは温度の伝わりやすさをあらわす。

熱伝導率が大きくても、密度×比熱が大きいと、温度はなかなか伝わりにくいんだにゃ。

温度Tは、時刻tと位置xの２変数関数T(x,t)になる。

粗い近似だけれど、

　　

と近似できる。

これを（２）に代入して整理すると、

　　

となる。

そして、この（３）式を元にして、（２）の偏微分方程式を解こうというわけ。

ちなみに、α=1/2のとき、

　　

という非常に簡単な計算式になる。

（４）を用いて解くことも可能。

κ=1/2とし、Δx=0.1とすると、α=1/2だから、

　　

したがって、この場合、（４）式を用いて１回計算すると、時間は0.01だけ進んだことになる。

次の問題は、この計算原理を理解するためのだケロ。

問題

　　

を、

　　

をΔx=1として、（４）式を使って解け。

【解】

（４）式を使えるのは

　　

のとき。

κ=1/2、Δx=1だから、Δt=1。

こんな計算は粗くて、この計算結果は使いものにはならないけれど、それでもこの偏微分方程式の解の特徴を比較的よくあらわしている。

最初は急激に温度が降下し、徐々に温度降下の速度が緩やかになり、両端の温度に近づいてゆく。

 

ちなみに、表計算で使った式は右の図に出ているにゃ。

C2のセルに
=(B2+D2)/2

という計算式を書き、それをD、Eにコピペする。

B列とF列は全部０だから、予め、計算をするところにすべて０をセットしておけばいい。
こうするだけで、簡単に計算できてしまう。

C３セルの計算に使っているのは、青のB2とD2セル。

一つ前の時間の隣接するセルの値をつかって、C3セルの値を計算しているわけ。

だから、連立方程式を解かずに、次々と計算してゆくことができる。

連立方程式を使って解く方法は陰解法と呼ばれる。今回紹介したのは、陽解法の中で最もシンプルなもの。
計算の精度をあげようと思えば、Δxを小さくすればいい。
いまは、Δx=1で計算しているけれど、Δx=0.1くらいにとれば、かなり精度のいい計算ができる。
ただし、κΔt/(Δx)²=1/2の制約があるから、Δt=0.01に取らなければならず、膨大な繰り返し計算をしなければいけない。
計算量は、何と、１０００倍になる！！
だから、陽解法は、計算の原理は簡単だけれど、精度よく計算しようとすると、膨大な計算量なるので、あまり使われず、多くの場合、連立方程式を解く陰解法が使われる。

そして、ネムネコは、これよりも複雑な計算ができる表計算のスプレッドシートを公開したにゃ。
Bloggerの方にそのアドレスを記しておいたので、興味のある人は見てほしいニャ。

http://nemneko.blogspot.jp/2016/11/blog-post_81.html

エクセルで計算できるように、スプレッドシートをダウンロードできるようにもしておいたにゃ。