ニュートンの冷却の法則を使って温度上昇を予測してみよう

tumblrに2018年に上げた記事の引っ越しです。

ニュートンの冷却の法則を使って温度上昇を予測してみよう

https://ja.wikipedia.org/wiki/ニュートンの冷却の法則

タイトルの意味が分からない人はここを参照していただければ。

簡単に言うと、熱いモノをその辺に放置しておくと、おおよそ以下のような温度変化での冷却が起こることが経験的に知られています。

$$\theta{} = \theta{}_0 e^{ -kt }$$

ただし、記号の意味は以下となります。

• $\theta{}$ : 現在のモノと周囲の温度差(K)
• $\theta{}_0$ : 初期のモノと周囲の温度差(K)
• $k$ : 冷却係数
• $t$ : 時間(s)

そもそも何故そんな事を考えたの？

以前働いていたお仕事の関係で、コイルを300℃に温める必要があったのですが、

• コイルを温める手段が2つあり、それらを同時に使う
  1. オーブンみたいな加温機(温度は可変、初期温度は適当で、後から調節)
  2. コイルに電流(固定)を流した際に発生するジュール熱
• 温度を計測する方法がコイルに付けた熱電対温度計しかない

という状況から温度上昇を予測してみようと思った次第です。
この時に、ニュートンの冷却の法則を「マイナス方向の冷却(＝加熱)」でも出来るのではないか、と考えたわけです。
なので、以下の状況を想定しなければなりません。

• 初期の温度差、つまり何℃に向かって温度が変化するかは不明
• 冷却係数なんか分かるわけがない
• 計測可能なデータは現在の温度と実験開始からの経過時間のみ

結構エグい条件ですが、参ります。

計測する上での前提条件

計算を楽にするために、以下の通りに設定します。

$$\begin{aligned} t_1 &\coloneqq{} 任意のコイル温度計測時間(秒)(既知) \\ \theta{}_0 &\coloneqq{} 初期のコイル自体の温度差(K)(不明) \\ \theta{}_1 &\coloneqq{} t_1 時のコイル自体の温度差(K)(不明) \\ \theta{}_2 &\coloneqq{} 2 \times t_1 時のコイル自体の温度差(K)(不明) \\ \end{aligned}$$

この時、以下の事が容易に導けると思います。

$$\begin{aligned} | ( 測定開始時の温度 ) - ( t_1 秒後の温度 ) | &= \theta{}_0 - \theta{}_1 \\ | ( t_1 秒後の温度 ) - ( 2 × t_1 秒後の温度 ) | &= \theta{}_1 - \theta{}_2 \\ \theta{}_1 &= \theta{}_0 e^{ -k t_1 } \cdots{} \textcircled{1} \\ \theta{}_2 &= \theta{}_0 e^{ -2 k t_1 } \cdots{} \textcircled{2} \\ \end{aligned}$$

いざ導出

まず、$\textcircled{ 1 }$と$\textcircled{ 2 }$から、以下が導出できます。

$$\begin{aligned} \theta{}_1 - \theta{}_2 &= \theta{}_0 e^{ -k t_1 } - \theta{}_0 e^{ -2 k t_1 } \\ &= \theta{}_0 e^{ -k t_1 } ( 1 - e^{ -k t_1 } ) \\ &= \theta{}_1 ( 1 - e^{ -k t_1 } ) \\ 1 - ( \frac{ \theta{}_2 }{ \theta{}_1 } ) &= 1 - e^{ -k t_1 } \\ \frac{ \theta{}_2 }{ \theta{}_1 } &= e^{ -k t_1 } \cdots{} \textcircled{3} \\ \end{aligned}$$

この時、$\textcircled{ 3 }$を変形すると$\theta{}_2 = \theta{}_1 e^{ -k t_1 }$と奇しくも$\textcircled{ 1 }$や$\textcircled{ 2 }$に近い形になります。
さて、$\textcircled{ 1 }$から$\textcircled{ 2 }$を引くと

$$\frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 } = e^{ -k t_1 } - e^{ -2 k t_1 }$$

という式が導出できますが、そこに$\textcircled{ 2 }$と$\textcircled{ 3 }$の式を代入すると

$$\begin{aligned} \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 } &= \frac{ \theta{}_2 }{ \theta{}_1 } - \frac{ \theta{}_2 }{ \theta{}_0 } \\ &= \frac{ \theta{}_2( \theta{}_0 - \theta{}_1 ) }{ \theta{}_0 \theta{}_1 } \\ \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } &= \frac{ \theta{}_2 }{ \theta{}_1 } \\ &= e^{ -k t_1 } \\ -k t_1 &= \ln{ ( \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } ) } \\ -k &= \frac{ \ln{ ( \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } ) } }{ t_1 } \cdots{} \textcircled{4} \\ \end{aligned}$$

これでさっき「分かるかー！」と言っていた冷却係数が計算出来ることが判明しました。やったぜ！
あとは、$\textcircled{ 1 }$と$\textcircled{ 2 }$から導いた式に冷却係数を代入してやれば超いい感じになるのでは？

$$\begin{aligned} \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 } &= e^{ -k t_1 } - e^{ -2 k t_1 } \\ &= \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } - ( \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } )^2 \\ &= \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 }( 1 - \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } ) \\ &= \frac{ \theta{}_1 - \theta{}_2 }{ \theta{}_0 - \theta{}_1 } \frac{ ( \theta{}_0 - \theta{}_1 ) - ( \theta{}_1 - \theta{}_2 ) }{ \theta{}_0 - \theta{}_1 } \\ \frac{ 1 }{ \theta{}_0 } &= \frac{ ( \theta{}_0 - \theta{}_1 ) - ( \theta{}_1 - \theta{}_2 ) }{ ( \theta{}_0 - \theta{}_1 )^2 } \\ \theta{}_0 &= \frac{ ( \theta{}_0 - \theta{}_1 )^2 }{ ( \theta{}_0 - \theta{}_1 ) - ( \theta{}_1 - \theta{}_2 ) } \cdots{} \textcircled{5} \\ \end{aligned}$$

これでめでたく、計測した温度の差さえ分かれば初期温度差を求めることが出来るようになりましたとさ。パンパカパーン！
あとは求めた$\theta{}_0$を、場合によって測定開始時の温度に足したり引いたりしてやれば到達温度がだいたい予想できる、という具合です。

到達温度に達するまでどれくらいかかる？

到達温度が予想できるようになって次に気になってくるのは、「じゃあ到達温度に達するまでどれくらいかかるの？」という疑問です。
早速今求めた値($\textcircled{ 4 }と$$\textcircled{ 5 }$)を利用して推測してみましょう。
まず、使用する温度計の検出限界となる最小の単位を$\alpha{}$℃とおきます。
例えば、0.1℃単位で測れる温度計ならば、$\alpha{} = 0.1$です。
この時、以下の不等式が立式できると思います。

$$\alpha{} > \theta{}_0 e^{ -k t_C }$$

この式で出てきた$t_C$は到達温度に達するのに必要な時間(秒)です。
この式を変形すると、

$$\begin{aligned} \frac{ \alpha{} }{ \theta{}_0 } &> e^{ -k t_C } \\ \ln{ ( \frac{ \alpha{} }{ \theta{}_0 } ) } &> -k t_C \\ t_C &> \frac{ \ln{ ( \frac{ \alpha{} }{ \theta{}_0 } ) } }{ -k } \cdots{} \textcircled{6} \\ \end{aligned}$$

おしまい。案外と単純です。イェイ。

触れるべき欠点

この手法の最も重大な問題は、ある程度余裕を持った$t_1$の設定をしないと(場合によりますが数分は欲しいかなぁ)、予測$\theta{}_0$がありえない大きさにすっ飛んで行ってしまうことです。
また、測定しているモノの相が変化しない(蒸発や凝固を起こさない)ことが話の大前提なので、そこが崩れると一切使えない方法論になってしまいます(特に到達温度に達する時間の予測)。
最後に、そもそも論として、今回利用したニュートンの冷却の法則は経験則による近似であるため、必ず正しい値が出る保証は一切出来ません。なので厳密な予測が必要な場合には全く使えません。
(一応、必要に迫られてこの方法論を作った際には、多少の温度のズレはご愛嬌レベルの話だったのですが、意外と実測値とのズレは少なく、十分実用に耐えうるレベルでしたよ、とは言えます。)