정지한 채로 데워지는 유체 — 자연대류와 Rayleigh 수 1708

같은 Ra에 1/3 지수를 씌우면 열유속은 약 149로 줄고, 1e9은 높이 약 0.73 m에서야 온다.

바닥을 데워도 유체는 한동안 그대로다. 열만 위로 간다. 줄무늬 소용돌이는 한꺼번에 켜진다. 그 문턱의 숫자가 Rayleigh 수 1708이다. 데워진 벽이 얼마나 식는지는 뒤의 평판 계산이 답한다.

Ra 슬라이더를 1708 아래로 내리면 입자가 멈추고 온도가 위아래로 갈린다. 그 위로 올리면 서로 반대로 도는 롤이 켜진다.

Ra = 1.0e+4

critical Ra_c = 1708 · state: steady convection rolls (Benard cells)

문턱 바로 위에서는 순환이 약하다. 슬라이더를 더 밀면 경계층이 얇아지고 셀이 세진다.

부력과 점성이 맞붙는 자리#

자연대류는 펌프 없이 밀도 차가 만드는 흐름이다. 바닥이 뜨거우면 그 위 유체가 팽창해 가벼워진다. 가벼운 덩어리는 떠오르려 한다. 이것이 부력이다.

부력만 있으면 흐름은 바로 시작된다. 점성이 떠오르는 덩어리를 마찰로 붙잡는다. 열확산은 뜨거운 덩어리의 온도를 지워, 부력이 기대는 온도차까지 없앤다.

두 제동이 부력을 이기면 유체는 정지다. 부력이 이기면 대류가 켜진다.

Boussinesq 근사 — 밀도 변화를 한 항으로#

밀도를 모든 항에서 변수로 두면 식이 무거워진다. Boussinesq 근사는 밀도를 부력 항에서만 온도의 함수로 둔다. 나머지 밀도는 상수다.

체적 열팽창계수는 온도가 1도 오를 때 밀도가 줄어드는 비율이다.

\beta \equiv -\frac{1}{\rho}\left(\frac{\partial \rho}{\partial T}\right)_P

β는 열팽창계수, ρ는 밀도, T는 온도다. 이상기체라면 β는 1/T다.

밀도 차는 온도 차로 바뀐다.

\frac{\rho_\infty - \rho}{\rho} \approx \beta\,(T - T_\infty)

ρ∞와 T∞는 멀리 떨어진 정지 유체의 밀도와 온도다. 이 줄이 수직 평판 경계층의 운동량 방정식에 부력을 심는다.

u\frac{\partial u}{\partial x} + v\frac{\partial u}{\partial y} = g\beta\,(T - T_\infty) + \nu\frac{\partial^2 u}{\partial y^2}

u와 v는 속도, g는 중력가속도, ν는 동점성계수다. 우변에서 앞 항이 부력이고 뒤 항이 점성이다. 에너지 방정식은 온도를 확산으로 끌고 간다.

u\frac{\partial T}{\partial x} + v\frac{\partial T}{\partial y} = \alpha\frac{\partial^2 T}{\partial y^2}

α는 열확산계수다. 운동량과 에너지가 부력으로 묶인 이 쌍이 자연대류의 뼈대다.

Grashof와 Rayleigh — 두 무차원수의 분업#

식을 무차원화하면 부력 항이 한 수로 뭉친다. Grashof 수다.

\mathrm{Gr} = \frac{g\beta\,(T_s - T_\infty)\,L^3}{\nu^2} \sim \frac{\text{부력}}{\text{점성력}}

Ts는 벽 온도, L은 특성 길이다. Grashof는 강제대류에서 흐름의 자리를 길이로 가르는 Reynolds 수가 하던 일을 가져간다. 부력이 점성을 얼마나 누르는지를 잰다.

적은 점성만이 아니다. 열확산도 온도차를 지운다. 그 열확산 제동은 Prandtl 수, 곧 운동량 확산과 열확산의 비 Pr = ν/α를 곱해 넣는다. 결과가 Rayleigh 수다.

\mathrm{Ra} = \mathrm{Gr}\cdot\mathrm{Pr} = \frac{g\beta\,(T_s - T_\infty)\,L^3}{\nu\,\alpha}

분모에 ν와 α가 같이 있다. 대류가 켜질지는 Gr 하나가 아니라 Ra가 정한다. 부력이 두 제동의 곱을 이겨야 흐름이 산다.

1708이라는 문턱#

아래에서 데워지는 층을 선형 안정성 분석에 넣으면 임계값이 나온다. 위아래가 고체 벽이면 Ra_c는 1708이다. 그보다 작으면 교란이 점성과 열확산에 먹혀 사라진다. 유체는 정지하고 열은 전도로만 간다.

1708을 넘으면 나란히 반대로 도는 롤, Bénard 셀이 켜진다. 셀 폭은 층 두께의 약 두 배다. 앞의 Ra 슬라이더가 그 문을 넘나든다.

Python — 수직 평판의 Nusselt를 따라가기#

벽이 얼마나 식는지는 Nusselt 수 Nu = hL/k가 답한다. 전도 대비 대류다. 수직 등온 평판에는 전 구간에서 쓰는 Churchill–Chu 식이 있다. 60도 벽을 20도 공기 속에, 높이 0.30 m로 세운다.

import numpy as np
 
def rayleigh_number(Ts, Tinf, L, beta, nu, alpha, g=9.81):
    # 부력 대 (점성 x 열확산) 의 비
    return g * beta * (Ts - Tinf) * L**3 / (nu * alpha)
 
def churchill_chu_vertical(Ra, Pr):
    # 수직 등온 평판 — 전 영역 유효
    denom = (1.0 + (0.492 / Pr)**(9 / 16))**(8 / 27)
    return (0.825 + 0.387 * Ra**(1 / 6) / denom)**2
 
# 막온도 기준 공기 물성
nu, alpha, Pr, k = 1.85e-5, 2.60e-5, 0.71, 0.027
Ts, Tinf, L = 60.0, 20.0, 0.30          # 60C 벽, 0.3 m 높이
beta = 1.0 / (0.5 * (Ts + Tinf) + 273.15)  # 이상기체 1/Tf
 
Ra = rayleigh_number(Ts, Tinf, L, beta, nu, alpha)
Nu = churchill_chu_vertical(Ra, Pr)
h = Nu * k / L
q = h * (Ts - Tinf)                     # 단위면적 열손실
 
print(f"Ra = {Ra:.2e}")
print(f"Nu = {Nu:.1f},  h = {h:.2f} W/m^2K")
print(f"q  = {q:.1f} W/m^2")
 
# Ra = 7.03e+07
# Nu = 54.9,  h = 4.94 W/m^2K
# q  = 197.6 W/m^2

Pr = 0.69Ra marker = 1e7.0

cyan: Churchill–Chu · green: laminar 0.59 Ra^1/4 · pink: turbulent 0.10 Ra^1/3

Pr 슬라이더를 액체금속 쪽으로 내리면 곡선이 내려가고, 기름 쪽으로 올리면 올라간다. Ra 마커를 1e9 너머로 밀면 기울기가 1/4에서 1/3로 바뀐다.

출력에 찍힌 수는 Ra 7.03e7, Nu 54.9, h 4.94, 열유속 197.6이다. 막온도는 40에 273.15를 더한 313.15 K이고, 팽창계수는 그 역수 약 0.003193이다. 점성과 열확산의 곱은 1.85e-5 곱하기 2.60e-5, 곧 4.81e-10이다. 높이 0.30 m의 세제곱은 0.027이다. 이 조각이 코드의 중력 9.81, 온도차 40과 만나 인쇄된 7.03e7이 된다. Nu 54.9에 0.027을 곱하고 0.30으로 나누면 4.941이라 출력의 h 4.94와 같다. 4.94 곱하기 40은 197.6이라 열유속 줄도 닫힌다. 이 Ra를 Pr 0.71로 나누면 Gr는 약 9.90e7이다. 1708과 나란히 놓을 수가 아니다.

같은 Ra에서 층류 별칭, 0.59 곱하기 Ra의 1/4승은 약 54.0이다. 출력 54.9는 그보다 약 0.9 위에 있다. 난류 별칭, 0.10 곱하기 Ra의 1/3승은 약 41.3이다. 그 Nu로 h를 다시 구하면 약 3.72, 열유속은 약 149다. 197.6의 약 4분의 3이다.

Ra가 1e9이 되려면 7.03e7의 약 14.2배가 필요하다. 높이는 세제곱으로 들어가므로 0.30 m는 약 0.73 m가 되어야 그 문에 닿는다. 그 전의 이 평판은 1/4 지수 쪽이다. 같은 공기, 같은 40도에서 Ra가 1708이 되는 길이는 약 8.7 mm다. 그 길이는 밀폐층의 두께다. 0.30 m 평판의 문으로 1708을 읽으면 기하가 다르다.

출력의 197.6은 제곱미터당 열유속이다. 폭을 높이와 같은 0.30 m로 잡으면 면적은 0.09 m^2이고, 판 하나가 잃는 열은 197.6 곱하기 0.09, 약 17.8 W다.

층류 별칭만 따라가면 Nu는 Ra의 1/4승이고 Ra는 높이의 세제곱이라, h는 높이의 -1/4승으로 움직인다. 높이를 두 배로 하면 h는 2의 -1/4승, 약 0.841배가 된다. 4.94는 약 4.15로, 열유속은 약 166으로 줄어든다. 1/3 지수에서는 그 감소가 사라지고 h가 높이에 안 남는다. 이 출력의 Ra는 아직 그 지수가 아니다.

Nusselt 상관식 한 묶음#

형상마다 식이 다르다.

수직 평판: $\mathrm{Nu} = 0.59\,\mathrm{Ra}^{1/4}$ (층류, $10^4$ – $10^9$ ), $\mathrm{Nu} = 0.10\,\mathrm{Ra}^{1/3}$ (난류, $10^9$ – $10^{13}$ ). 전 영역은 Churchill–Chu.
수평 가열판, 위로 향함: $\mathrm{Nu} = 0.54\,\mathrm{Ra}^{1/4}$ ( $10^4$ – $10^7$ ), $0.15\,\mathrm{Ra}^{1/3}$ ( $10^7$ – $10^{11}$ ).
수평 가열판, 아래로 향함: $\mathrm{Nu} = 0.27\,\mathrm{Ra}^{1/4}$ .
아래에서 데우는 밀폐 공간: $\mathrm{Ra}_c = 1708$ 을 넘으면 Bénard 셀.

수평판의 특성 길이는 면적 나누기 둘레, $L = A_s/P$ 다. 난류 지수 1/3에서는 $\mathrm{Nu}$ 가 $\mathrm{Ra}^{1/3}$ 에 비례하고 $\mathrm{Ra}$ 가 $L^3$ 에 비례하므로 $h = \mathrm{Nu}\,k/L$ 에서 $L$ 이 빠진다. 난류 자연대류의 열전달계수는 크기에 안 매달린다. 위의 0.30 m 평판은 아직 그 상쇄를 쓸 높이가 아니다.

7.03e7에 난류 지수를 씌우면 열유속이 깎인다#

1/3을 이 출력에 씌우면 Nu 54.9가 약 41.3으로 내려가고, 197.6은 약 149로 깎인다. 높이를 식에서 지우는 지수라, 0.30 m도 같이 사라진다. 층류 가지에서는 높이를 바꾸면 h가 따라 움직인다. 지수를 먼저 고르면 그 움직임이 가짜로 죽는다.

1e9은 약 0.73 m에서야 온다. 셀이 켜진 밀폐층에서 모드를 세 개만 남기면 이 열유속은 어디까지 남는지, 롤을 세 식으로 접은 Lorenz 모형이 그 질문을 들고 있다.