なぜ弱火のレシピを時短のために強火にしてはいけないか1次元お肉を熱伝導方程式で焼いてみる

はじめに

料理の失敗談として「弱火で10分のところを強火ですれば時間短縮になる!」といって「外は黒焦げ、中が生焼け」となるという話をよく(?)聞きます。これはなぜでしょうか?今回はこれについて熱伝導方程式を解いてもとめようというのが今回の記事の趣旨です。 f:id:ikaro1192:20190603134316p:plain

理論と手法

物理において何が現象の主要な支配要因になってるかを把握するためにはパラメータを最小限に抑えて、現象をもっともらしく表現できる近似をしてあげることが大切です。今回は最低限、現象を説明できる条件の1次元の熱伝導を計算します。また、厳密に解けそうな気がしますが数値計算を用いて解いていきます。

ということで1次元のお肉を焼いていきましょう。1次元の熱伝導方程式は下記のように表現できます。

  \frac{\partial T}{\partial t} = \kappa  \frac{{\partial}^2 T}{{\partial} {x}^2}

ここでTは肉内部の温度分布で時間tと位置xの関数、κは熱拡散率です。肉の熱拡散率は(杉山,2013)*1にかかれている1.30×10^-7 m2/sを使います。境界条件としてフライパンに接している側は140度と200度(IHで弱火と強火でよく使われている値)とします。とりあえずあとで数式に登場させるためこれを T_hと置きます。また、空気と接している側の境界条件ですが断熱*2とします。これは料理をしている最中に肉と比べて部屋の温度がどんどん上がっていかないという物理的条件に対応します。また熱拡散率は本来温度に依存しますが簡単のため定数とします。焼くお肉の厚さですがここは庶民ですので L=1 cmのもの*3を用意します。そして、初期条件としてお肉の温度を一様に T_0=10度とします。 以上の条件を式にすると

  T(t,x=0) = T_h

   \frac{\partial T}{\partial x}  |_{x=L}= 0

  T(t=0,x \lt L) = T_0

次に数値計算できるよう離散化をします。今回は微分方程式の解き方としては陽解法を採用します。具体的には微分を下記の差分に置き換えます。

 \frac{\partial T}{\partial t} = \frac{T(t+\Delta t) - T(t)}{\Delta t}

 \frac{{\partial}^2 T}{\partial {x}^2} = \frac{T(x+\Delta x) - 2 T(x) + T(x-\Delta x)}{( {\Delta x}^2)}

ここで \Delta t \Delta xはそれぞれを表しています。

  \frac{T(t+\Delta t) - T(t)}{\Delta t} = \kappa  \frac{T(x+\Delta x) - 2 T(x) + T(x-\Delta x)}{( {\Delta x}^2)}

整理すれば

  T(t+\Delta t, x)  = T(t, x) + \frac{\kappa \Delta t}{{\Delta x}^2}  (T(t, x+\Delta x) - 2 T(t, x) + T(t, x-\Delta x))

を得ます。この調子で断熱の境界条件も離散化しましょう。

 \frac{\partial T}{\partial x} |_{x=L} = \frac{T(L) - T(L-\Delta x)}{\Delta x} =0

  T(L)  =T(L-\Delta x)

と最後とその手前の温度が同じという条件に帰着します。

以上によって求めた計算結果を1分、2分、3分での熱の分布を140度で加熱した場合、200度で加熱した場合で比較していきます。

実装

次に実装するだけです。今回はC++で実装しました。下にソースコードを示します。 また、いくつか式には出てこない部分があるので解説をしておきます。

#include<iostream>
#include<vector>

int main(){

    int input_temprature;
    int input_time;

    std::cin >> input_temprature >> input_time ;

    const double boundary_condition = input_temprature;//境界条件 [℃]
    const double temperature_init = 10; // 初期条件 [℃]
    const double thermal_diffusivity = 0.130; // [mm^2/s]
    const double x_max = 10; // お肉の厚さ [mm]
    const double dx = 0.1;
    const int Nx = x_max/dx;

    const double t_max = input_time;//加熱時間[s]
    const double dt = 0.8 * dx*dx/(2*thermal_diffusivity); //安定条件より。0.8は念の為。解説参照
    const int Nt = int(t_max/dt);

    std::vector<double> before_temperature;
    std::vector<double> now_temperature;

    //境界条件
    before_temperature.push_back(boundary_condition);

    //初期値を入れる
    for(int i =1; i<Nx;++i){
        before_temperature.push_back(temperature_init);
        now_temperature.push_back(temperature_init);
    }

    for(int j=0; j<Nt;++j){
        //境界条件
        now_temperature[0] = boundary_condition;
        for(int i =1; i < Nx-1; ++i){
            now_temperature[i] = before_temperature[i] + ((thermal_diffusivity* dt)/(dx*dx))*(before_temperature[i+1]-2*before_temperature[i]+before_temperature[i-1]);
        }

        //境界条件
        now_temperature[Nx-2] = now_temperature[Nx-3];
        before_temperature = now_temperature; // 本文参照
    }

    double x = 0;
    for(auto temperature : now_temperature){
        x+=dx;
        std::cout << x <<" "<< temperature << std::endl;
    }

    return 0;
}

まずは安定条件について。詳しくは数値計算の教科書等を参照してほしいのですが、今回採用した方法では \Delta t \Delta xが大きいと解が安定しません。具体的には

 \Delta t \leq \frac{{(\Delta x)}^2}{2\kappa}

を満たしている必要があります。そのためdtを下記式で求めています。

const double dt = 0.8 * dx*dx/(2*thermal_diffusivity); //安定条件より。0.8は念の為。解説参照

各経過時間を配列で確保することもできるのですが、それではメモリが多量に必要になります。そこで2ステップ分の配列を用意してそれを入れ替えていくことでメモリを節約します。下記部分です。

        before_temperature = now_temperature; // 本文参照

ほんとはこんなコピー発生させてコストがーとおもいますが、今回は計算を回してみてもそんなに時間がかからなかったのでこのままにしてあります。

境界条件は怪しいですのでもしバグを見つけたらおしえてください()

計算結果

記載の便宜上、「焼きすぎ」、「適温」、「生焼け」の3つを定義します*4。100度を超えた領域を「焼きすぎ」、50度を下回る領域を「生焼け」、そして100度以下、50度以上を「適温」とします。また適温部分の長さを「可食部」とすることにします。また、140度で3分やいたとき、200度で1分焼いたとき、200度で2分20秒焼いたときの温度分布をそれぞれパターンA、パターンB、パターンCとします。

計算結果のグラフは下記となります。

f:id:ikaro1192:20190611020718p:plain
フライパンからの距離と温度の関係。140度で3分やいたとき、200度で1分焼いたとき、200度で2分20秒焼いたときのそれぞれ温度分布を赤線、青線、黒線で示してある。破線は50度(生焼けの境界温度)と100度(焼きすぎの境界温度)を示している。

パターンAでは約2.9 mmにいったところで適温となっています。これには生焼け部はありませんので、したがってパターンAの可食部は7.1 mmとなります。パターンBは焼きすぎの領域がパターンAと一致しています。パターンBは5.0 mmのところで生焼けとなります。したがって、このときの可食部は2.1 mmとなります。残りは生焼けの部分で5.0 mmです。パターンCでは4.5 mmで適温となります。したがって可食部は5.5 mmとなります。表にまとめると下記のようになります。

パターン 焼きすぎ部の長さ(mm) 生焼け部の長さ(mm) 可食部部(mm)
A 2.9 0 7.1
B 2.9 5.0 2.1
C 4.5 0 5.5

パターンAとBを比べると2.9 mmを堺にしてプライパン側ではAよりBの温度が高くなり、外側ではAよりBの温度が低くなっています。パターンAとCを比べると常にCの温度がAの温度を上回っています。

考察

どう頑張っても140度の温度分布と200度の温度分布を一致させることはできませんでした。どんな時間加熱しても200度で加熱した場合は可食部は140度3分の場合の可食部以下となってしまいます。これは加熱する温度によって温度分布は異なるものになることを示します。温度分布が異なるのですから、弱火を強火で代用することができないことがわかります。

パターンAとBの比較からは表面の焼け具合を同様にすると中が生焼けになってしまうことがわかります。また、焼きすぎの部分でも200度のほうが同じ位置なら温度が高いです。これは「外は黒焦げ、中が生焼け」という状況に対応します。パターンCのグラフの意味することは外が黒焦げとなってしまった状況に対応します。つまり200度のほうがはやく表面は焼けますが、熱が中まではやく伝わっていくかは別問題ということを意味します。これは実は数値計算をするまでもなく数式でも理解することができます。離散化した式を見てみると次の瞬間の温度というのは周囲の温度差にある定数をかけたものとよむことができます。つまり、温度差があればそれだけ早く"平均的な"温度に落ち着こうとするものの、瞬時にその温度にいけるわけではないということを示しています。そういった時間ラグがあるため肉の中を一瞬で熱は伝わらず結局、表面には早く熱が伝わり焦げ、中の方は伝わらず生焼けとなってしまうのです。

まとめ

今回は1次元の豚肉を焼くとその温度分布はどうなるかを熱伝導方程式を数値計算を用いて解くことでもとめました。結果は加熱側の温度が違えば加熱時間をどう調整しようと同じ結果にすることはできないというものです。式を考察することにより熱は一瞬で伝わることができないため加熱温度が高ければ表面には早く熱が伝わり焦げ、中の方は伝わらず生焼けとなってしまうことを明らかにしました。ということでレシピの火加減と時間はまもるようにしましょう。

*1:https://www.jstage.jst.go.jp/article/cookeryscience/46/4/46_299/_pdf

*2:結論は大きく変わりませんが本当は少しづつ空気に熱は漏れ出すの熱抵抗を使ったほうがいいというコメント頂きました

*3:実は最初は3 cmでやったのですがうまく火が通らなかったのでこうしました

*4:実際に焦げたり美味しく安全に食べれる温度とは異なります