スポンサーサイト

上記の広告は1ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。

差分法による1次元の熱伝導問題 1/2

具体的な熱伝導の問題を、1次元の熱伝導方程式を使って数値的に解いてみました。

Tag : 数学 物理学 数値解析 伝熱工学
今回解くのはこのような問題
次の条件のもと、長さ10[cm]の銅の棒の温度分布を求める。

・棒の両端は常に0℃に固定しておく。(境界条件)
・熱の流出は棒の両端でのみ起こるとする。
・初期状態(t=0)での銅棒の温度分布を↓画像のように指定する。(初期条件)
初期条件
なお、銅棒の座標はこのようにとる。
銅棒の座標

温度分布は、棒の中心が一番高温(10℃)で、端にいくにつれて低温になるような形です。

銅棒の温度がどのように変化していくかを記述する熱伝導方程式
∂u/dt = α・∂2u/dx2 ……(1)
となります。
u(x,t) : 時刻t, 位置xでの温度を表す関数
α : 物体(今回は銅)の熱拡散率
なお、Wikipediaによると熱拡散率は、物体の熱伝導率k, 密度ρ, 比熱c を用いて α=k/(ρ・c) と表せ、の場合、k=401, ρ=8920, c=380 なので、α=1.18×10-4です。

以上で計算に必要な条件がそろいました。


熱伝導方程式(1)を、FTCSスキームを用いて差分化すると

{ ud(xi,tn+1)-ud(xi,tn) }/⊿t
= α・{ ud(xi+1,tn) -2ud(xi,tn) + ud(xi-1,tn) }/(⊿x2)


となります。これは前回やりました。
ここで、udはuを差分化した関数です。

この式をコンピュータで計算できるようにするために、u(xi,tn+1) = ○○ という形に変形するとこうなります。

ud(xi,tn+1) = (1-2v)・ud(xi,tn) + v・{ ud(xi-1,tn) + ud(xi+1,tn) }

( v := α・⊿t/(⊿x2) )

ここで、xiはxを離散化したときのi番目の値を、tnはtを離散化したときのn番目の値を表します。
この式と、上述の条件(を差分化した式)
初期条件 : ud(xi,0) = 4000xi(0.1-xi)
境界条件 : ud(0,tn) = ud(0.1,tn) = 0
を使えば、tnでのudの値からtn+1でのudの値を計算できます。

このように、差分方程式を作ることで微分方程式を近似的に解く方法を差分法といいます。
特に今回のような、t=0 から始めて tn の値を順次計算していく方法は陽解法と呼ばれます。


実際に計算してみました。(プログラムは自前で組みました)

計算に使ったパラメータはこんな感じ。
xのサンプリング数 : N=51
xのサンプリング間隔 : ⊿x=0.002
tのサンプリング間隔 : ⊿t=0.01
銅の熱拡散率 : α=1.18×10-4 (既出)

実は、これらの数値を決めるのにも少しばかり条件がいるのですが、それはまた今度。
キーワードは安定性です。

計算結果をすべて載せるわけにもいかないので、t=0, 1, 2, …, 30 のときの値だけをグラフ化&gifアニメ化しました。
熱伝導の計算結果
数値データからの画像生成はScilab先生に、gifアニメ化はGiamにやってもらいました。(両方ともフリーソフト)

初めのうちは速く冷めますが、時間が経つにつれて冷めるのがゆっくりになっています。
このことは、私たちの直感と一致していますね。
粗熱はすぐに取れるけれど、完全に冷めるまでは時間がかかる・・・


長くなってしまったので、記事を分割します。
後半では、この計算結果が本当に正しいかどうかを検証します。
スポンサーサイト

コメントの投稿

非公開コメント

プロフィール

fura2

Author : fura2
数学・コンピュータを中心に、考えたこと・やったことを書いていきます。

誤植等を含め、間違いはご指摘いただければ幸いです。

FC2カウンター
検索フォーム
最新記事
最新コメント
最新トラックバック
月別アーカイブ
リンク
RSSリンクの表示
上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。