初めに
お味噌汁が好きな理由は様々ある．
美味しさは勿論のこと，統計学（cf. みどりの窓口の枕），熱対流を感じられることなど．
Methods
支配方程式
熱と流れをそれぞれ独立に考えたことはあるが（e.g., heat conduction in a rod, shear-driven flow in a square / cubic cavity），同時に考えたことはない．
先ず，それぞれの方程式は以下．
\begin{align}
\partial_t \rho
+ \nabla \cdot \left( \rho \boldsymbol{u} \right)
&= 0
\\
\rho
\left(
\partial_t \boldsymbol{u}
+ \left( \boldsymbol{u} \cdot \nabla \right) \boldsymbol{u}
\right)
&= - \nabla p
+ \mu \nabla^2 \boldsymbol{u}
+ \rho \boldsymbol{f}
\\
\rho c_p \left(
\partial_t T
+ \left( \boldsymbol{u} \cdot \nabla \right) T
\right)
&= \lambda \nabla^2 T
+ \rho c_p s
\end{align}
ここで，$\boldsymbol{u}$は流速$[ \mathrm{m / s} ]$，$p$は圧力$[ \mathrm{Pa} ]$，$\boldsymbol{f}$は物体力$[ \mathrm{m / s^2} ]$であり，これには重力を与える．また，運動に関するパラメータとして，$\rho$は密度$[ \mathrm{kg / m^3} ]$，$\mu$は粘度$[ \mathrm{Pa \cdot s} ]$である．
さらに，$T$は温度$[ \mathrm{K} ]$，$s$は熱源$[ \mathrm{K / s} ]$であり，ゼロとする．なお，熱に関するパラメータとして，$c_p$は比熱容量$[ \mathrm{J / (kg \cdot K)} ]$，$\lambda$は熱伝導率$[ \mathrm{W / (m \cdot K)} ]$である1．
いま，密度$\rho$を基準点$\rho_0$周りでTaylor展開する．
\begin{align}
\rho
&= \rho_0
+ \left. \frac{\partial \rho}{\partial T} \right|_{p} ( T - T_0 )
+ \left. \frac{\partial \rho}{\partial p} \right|_{T} ( p - p_0 )
+ (\text{HOT})
\\
\end{align}
但し，密度が圧力に依存して変化することはないとし，温度のみの関数として取り扱う（Boussinesq近似，$\rho \left( T, p, ... \right) \simeq \rho \left( T \right)$2）．したがって
\begin{align}
\rho
&\simeq \rho_0
+ \left. \frac{\partial \rho}{\partial T} \right|_{p} ( T - T_0 )
\\
&= \rho_0 + \rho^{\prime}
\quad \left(
\rho^{\prime}
:= \left. \frac{\partial \rho}{\partial T} \right|_{p} ( T - T_0 )
\right)
\end{align}
ここで，$\rho^{\prime}$は基準点からの変動であり，$\rho^{\prime} \ll \rho_0$．熱膨張率$\beta \ [ \mathrm{1 / K} ]$を
\begin{align}
\beta
&= - \frac{1}{\rho_0}
\left. \frac{\partial \rho}{\partial T} \right|_{p}
\end{align}
と定めると，
\begin{align}
\rho^{\prime}
&= \left. \frac{\partial \rho}{\partial T} \right|_{p} \left( T - T_0 \right)
\\
&= - \rho_0 \beta \left( T - T_0 \right)
\\
\therefore
\frac{\rho^{\prime}}{\rho_0}
&= - \beta \left( T - T_0 \right)
\end{align}
である．いま，外力$\boldsymbol{f}$に重力$\boldsymbol{g}$を与えたことを思い出せば，
\begin{align}
\frac{\rho}{\rho_0} \boldsymbol{f}
&= \left( 1 + \frac{\rho^{\prime}}{\rho_0} \right) \boldsymbol{g}
\\
&= \left( 1 - \beta \left( T - T_0 \right) \right) \boldsymbol{g}
\end{align}
であるから，
\begin{align}
\nabla \cdot \boldsymbol{u}
&= 0
\\
\partial_t \boldsymbol{u}
+ \left( \boldsymbol{u} \cdot \nabla \right) \boldsymbol{u}
&= - \frac{1}{\rho_0} \nabla p
+ \frac{\mu}{\rho_0} \nabla^2 \boldsymbol{u}
+ \left( 1 - \beta \left( T - T_0 \right) \right) \boldsymbol{g}
\\
\partial_t T
+ \left( \boldsymbol{u} \cdot \nabla \right) T
&= \kappa \nabla^2 T
\end{align}
を解くこととする．但し，$\kappa$は熱拡散率$[ \mathrm{m^2 / s} ]$である（$\kappa = \lambda / \rho c_p$）．
Results
空間微分をArakawa-B grid上で3次+2次精度で，時間微分を1次精度で近似した．離散化パラメータは$h = 2\times10^{-2} \text{ [m]}$と$\tau = 2\times10^{-3} \text{ [s]}$．
Shear-driven flow
先ずは，ただのshear-driven flowを長方形領域上で計算する．速度ノルム$| \boldsymbol{u} |$と，動圧で正規化した圧力$q = p / (1/2 \rho | \boldsymbol{u} |^2)$を示す．
| $\text{Re}$ | $\text{Results}$ |
|---|---|
| $\sim 10^2$ |  |
| $\sim 10^3$ |  |
| $\sim 10^4$ |  |
squareの自然な拡張となっており良し．
Buoyancy-driven flow
深さ$1 \ [ \mathrm{m} ]$の鍋を水で満たし，底面を沸点$100 \ [ ^{\circ}\mathrm{C} ] = 373.15 \ [ \mathrm{K} ]$で温め，表面を室温$20 \ [ ^{\circ}\mathrm{C} ] = 293.15 \ [ \mathrm{K} ]$に解放する．その他，パラメータはこちらから借用する．
- 密度 : $\rho_0 = 1.00 \times 10^{3} \ [ \mathrm{kg / m^3} ]$
- 粘性率 : $\mu = 1.00 \times 10^{-3} \ [ \mathrm{Pa \cdot s} ]$
- 動粘性率 : $\nu = \mu / \rho_0 = 1.00 \times 10^{-6} \ [ \mathrm{m^2 / s} ]$
- 熱伝導率 : $\lambda = 0.60 \ [ \mathrm{W / (m \cdot K)} ]$
- 比熱容量 : $c_p = 4.18 \times 10^{3} \ [ \mathrm{J / (kg \cdot K)} ]$
- 熱拡散率 : $\kappa = \lambda / \rho c_p = 1.44 \times 10^{-7} \ [ \mathrm{m^2 / s} ]$
- 体積膨張率: $\beta = 2.10 \times 10^{-4} \ [ \mathrm{1 / K} ]$
- 重力加速度: $g = 9.81 \ [ \mathrm{m / s^2} ]$
- Prandtl数 : $\mathrm{Pr} = \nu / \kappa = 6.97$
- Grashof数 : $\mathrm{Gr} = g \beta L^3 \Delta T / \nu^2 = 1.65 \times 10^{11}$
- Rayleigh数: $\mathrm{Ra} = \mathrm{Pr} \cdot \mathrm{Gr} = 1.15 \times 10^{12}$
先ず，以下は，
- 流速: 全ての壁面で非滑り
- 圧力: 全ての壁面で斉次Neumann
- 温度: 上下の壁面に温度差，その他全ての壁面で室温
としたときの結果．
| $\text{Dim}$ | $\text{Results}$ |
|---|---|
| $\text{2D}$ |  |
| $\text{3D}$ |  |
次いで，以下は，
- 流速: 全ての壁面で非滑り
- 圧力: 全ての壁面で斉次Neumann
- 温度: 上下の壁面に温度差，その他全ての壁面で斉次Neumann
としたときの結果．
| $\text{Dim}$ | $\text{Results}$ |
|---|---|
| $\text{2D}$ |  |
| $\text{3D}$ |  |
最後に，以下は，
- 流速: 上下の壁面で非滑り，その他全ての壁面で周期的
- 圧力: 上下の壁面で斉次Neumann，その他全ての壁面で周期的
- 温度: 上下の壁面に温度差，その他全ての壁面で周期的
としたときの結果．
| $\text{Dim}$ | $\text{Results}$ |
|---|---|
| $\text{2D}$ |  |
| $\text{3D}$ |  |
終わりに
お盆休み，僅かばかりの遊び心．
| $h$ | $1 \times 10^{-2} \text{ [m]}$ | $5 \times 10^{-3} \text{ [m]}$ | $2 \times 10^{-3} \text{ [m]}$ | $1 \times 10^{-3} \text{ [m]}$ |
|---|---|---|---|---|
| - |  |
 |
 |
 |
| - |  |
 |
 |
 |
鶏冠を立ち上げたニワトリか，湯煎されたチョコレートのよう．その昔，高校の地理で，ペルー・チリ沖で起きるという湧昇流なるものを習ったことを思い出す．また，高解像度の解を眺めていると，木星の大赤斑が思い浮かぶ．
そう言えば，密度差を陽な浮力としているから，Rayleigh-Taylor instabilityのような構造も当然生まれる（というか，普通に体積力の存在）．
| 1 | 2 | 3 | 4 |
|---|---|---|---|
 |
 |
 |
 |
Appendix
無次元量
幾つかの無次元量を導入する．
Prandtl数
動粘性と熱拡散の比．Ludwig Prandtlに由来．
\begin{align}
\mathrm{Pr}
&= \frac{\nu}{\kappa}
= \frac{\mu / \rho}{\lambda / \rho c_p}
= \frac{\mu c_p}{\lambda}
\end{align}
ここで，
- $\mu$: 粘性率$[ \mathrm{Pa \cdot s} ]$
- $\nu$: 動粘性率$[ \mathrm{m^2 / s} ]$
- $\rho$: 密度$[ \mathrm{kg / m^3} ]$
- $\kappa$: 熱拡散率$[ \mathrm{m^2 / s} ]$
- $\lambda$: 熱伝導率$[ \mathrm{W / (m \cdot K)} ]$
- $c_p$: 比熱容量$[ \mathrm{J / (kg \cdot K)} ]$
である．
Grashof数
浮力と粘性力の比．Franz Grashofに由来．
\begin{align}
\mathrm{Gr}
&= \frac{g \beta L^3 \Delta T}{\nu^2}
\end{align}
ここで，
- $g$: 重力加速度$[ \mathrm{m / s^2} ]$
- $\beta$: 体積膨張率$[ \mathrm{1 / K} ]$
- $\nu$: 動粘性率$[ \mathrm{m^2 / s} ]$
- $L$: 特徴的な長さ$[ \mathrm{m} ]$
- $\Delta T$: 特徴的な温度差$[ \mathrm{K} ]$
である．
Rayleigh数
浮力と熱拡散の比．John William Strutt, 3rd Baron Rayleighに由来．
\begin{align}
\mathrm{Ra}
&= \frac{g \beta L^3 \Delta T}{\kappa \nu}
= \mathrm{Pr} \cdot \mathrm{Gr}
\end{align}
ここで，
- $g$: 重力加速度$[ \mathrm{m / s^2} ]$
- $\beta$: 体積膨張率$[ \mathrm{1 / K} ]$
- $\kappa$: 熱拡散率$[ \mathrm{m^2 / s} ]$
- $\nu$: 動粘性率$[ \mathrm{m^2 / s} ]$
- $L$: 特徴的な長さ$[ \mathrm{m} ]$
- $\Delta T$: 特徴的な温度差$[ \mathrm{K} ]$
である．
-
「熱拡散率（thermal diffusivity）」と「熱伝導率（thermal conductivity）」と「熱伝達率（heat transfer coefficient）」は，全て異なる．熱拡散率$\kappa \ [ \mathrm{m^2 / s} ]$は，熱伝導率$\lambda \ [ \mathrm{W / (m \cdot K)} ]$，比熱容量$c_p \ [ \mathrm{J / (kg \cdot K)} ]$，密度$\rho \ [ \mathrm{kg / m^3} ]$を用いて，$\kappa = \lambda / \rho c_p$で与えられる．熱伝達率$\eta \ [ \mathrm{W / (m^2 \cdot K)} ]$は，壁面と流体などの2種類の物質間での熱エネルギーの伝わり易さである．熱移動量$Q \ [ \mathrm{W} ]$，熱流束密度$J \ [ \mathrm{W / (m^2)} ]$，電熱面積$A \ [ \mathrm{m^2} ]$，壁面の温度$T_0 \ [ \mathrm{K} ]$，媒質の温度$T \ [ \mathrm{K} ]$を用いて，$\eta = Q / A (T_0 - T) = J / (T_0 - T)$で与えられる．但し，$T_0 \gt T$． ↩
-
非圧縮性近似とは，密度$\rho$が圧力$p$に依存して変化することはない，と仮定するものである．直ちに$\rho = (\text{const.})$を言うものではない（参考）． ↩