밑면으로부터 등가속 열속으로 가열되는 유체층에서 열적 불안정성 해석

양대준, 정웅교, 최창균

서울대학교 공과대학 응용화학부

Analysis of the Thermal Instability in a Fluid Layer Subjected to  Ramp Heating from Below

D. J. Yang, W. G. Jeong and C. K. Choi

School of Chemical Engineering, Seoul National University

서론

자연대류는 중력장 하에서 부력이나 표면장력의 차이에 의해 발생되는 물리적 현상이다. 부력에 의한 자연대류는 열 또는 물질전달이 일어나는 각종 화공장치의 설계에서 자주 볼 수 있으며, 자연대류에 의한 유동은 실제 공정에서 열 및 물질전달 특성에 영향을 미치므로 이에 대한 연구가 진행되어 왔다. 이러한 자연대류에 관한 연구는 대기의 모델링 및 제품의 균질화, 순수화를 요구하는 결정성장, 전자소재의 가공 등의 영역에서 중요한 기여를 하고 있다.

본 연구에서는 잘 알려진 선형안정성 이론 하에서 Kang 등[1]에 의해 타당성을 보이고 있는 전파이론과 Lick[2]이나 Currie[3] 등이 사용한 동결시간 모델을 이용하여 밑면으로부터 등가속 열속으로 가열되는 유체층에서의 열적 불안정성을 해석하고, 그 이론적 해석 결과와 Ahlers 등[4]의 실험결과를 비교하였다.

이론

본 연구에서는 초기에 유체층이 밑면으로부터 열속이 열속증가계수 [eqn: rm phi]로 증가하는 계를 고려하며 고찰하는 계는 Figure 1과 같다. 이러한 경우 초기 유체층의 온도분포는 비선형적으로 변하게 되며, 온도차에 의해 발생하는 부력에 의해 대류가 일어날 것이다. 이때 유체층은 비압축성 뉴튼 유체이며, 밀도가 온도에 대해 선형으로 변하고, 외력에 의한 밀도변화를 제외하고는 밀도가 일정하다는 Boussinesq가정을 도입하면 다음과 같은 운동식, 에너지식, 상태식을 적용할 수 있다.

[eqn: rm rho_i [ { partial} over {partial t } + {bold U} CDOT { TRIANGLED} ]~{bold U}=- {TRIANGLED} P+ mu TRIANGLED^2 {bold U} +rho vecg]

(1a)

[eqn: rm [ { partial} over {partial t}+ {bold U} CDOT {TRIANGLED} ]T= alpha TRIANGLED^2 T]

(1b)

[eqn: rm rho=rho_i [ 1- beta (T-T_i )]]

(1c)

여기서[eqn: rm t, `{bold U},`P, `rho_i ,`mu,`T,`T_i , `vec g , `alpha,` beta]는 각각 가열시작시간, 속도벡터, 압력, 유체층의 초기밀도, 점도, 유체층의 온도, 유체층의 초기온도, 중력가속도, 열확산계수, 부피팽창계수를 나타낸다. 위 지배방정식에 교란되지 않는 양에 미소분의 교란량을 더하고 2차이상의 교란량을 무시하는 선형안정성이론을 적용하면 교란식이 얻어진다. 이때 다음과 같은 무차원 군들이 사용된다.

[eqn: rm theta_0 = {T_0 - T_u} over DeltaT_phi ,~~DeltaT_phi = {phi L^3} over {alpha k},~theta_1 = Ra_phi T_1 over DeltaT_phi#rm Ra_phi = {g beta DeltaT_phi L^3} over {alpha nu},~ Ra = {g beta DELTA T L ^3} over { alpha nu }~~~~~~~~~~~]

(2)

위에서[eqn: rm Delta T = T_u - T_b]이다. [eqn: rm T_b]는 계의 상부와 하부의 온도가 다를 때, 하부의 온도이고 ,[eqn: rm T_u]는 상부의 온도이다. [eqn: rm DELTA T_phi]는 시간에 따라 증가하는 열속계에서 사용되는 무차원 온도차이며, 이 온도차를 일반적인 Rayleigh 수, Ra에 대입하면 무차원[eqn: rm Ra_phi]가 얻어진다. 일반적으로 Ra에서 분자는 자연대류를 유발하는 부력을 나타내고 분모는 자연대류의 발생을 억제하는 열확산과 동점도를 나타낸다. 즉 Ra수가 증가할수록 자연대류의 발생이 쉬워진다. 위에서 정의한 무차원군과 교란량을 정규섭동식으로 치환하여 □double curl□을 취한 후 z방향의 성분을 모으면 진폭함수를 얻게된다.

[eqn: rm1 over Pr {partial } over {partial tau } del^2 w_1 = del^2 del^2 w_1 + del_1 ^2 theta_1 ~~~~]

(3a)

[eqn: rm {partial theta_1 } over { partial tau } + w_1 Ra_phi { partial theta_0 } over { partial z } = Del^2 theta_1]

(3b)

여기서[eqn: rm {del ^2}_1],[eqn: rm w_1],[eqn: rm theta_1],[eqn: rm theta_0], Pr은 각각 2차원 Laplacian, z방향의 속도교란량, 온도교란량, 그리고 전도만 일어날 때의 온도, Prandtl 수를 나타낸다. 전파이론은 이러한 교란량들의 변화가 일어나는 열침투 깊이를 길이단위로 교란량을 무차원하여 안정성해석을 수행하는 것이다. 일반적인 계의 길이보다는 실제로 변화가 발생하는 열침투깊이를 이용함으로써 이론적인 배경을 갖추었다고 할 수 있다. 따라서[eqn: tau]와 z의 함수인 속도와 온도 교란량을 유사변수의 함수로 전환하기 위해 [eqn: rm [w_1 ,theta_1 ]= [ tau w^* (zeta),theta^* (zeta)]]인 관계를 도입하면 다음과 같은 초기치 문제의 지배방정식이 얻어진다.

[eqn: rm (D^2 -a^*^2 )^2 w^* = a^*^2 theta^* - 1 over Pr left[ 1 over 2 zeta D^3 w^* -{ a^*^2 } over 2 zeta Dw^* +a^*^2 w^* right]]

(4a)

[eqn: rm - zeta over 2 D theta^* + Ra_phi tau^2 w^* {partial theta_0} over {partial z} = D^2 theta^* - a^*^2 theta^*]

(4b)

계의 대류가 발생하기 전의 전도만 일어나는 상태에서의 온도분포는 이미 잘 알려져 있으며 온도분포 식은 다음과 같다[5].

[eqn: rm {partial theta_0} over {partial tau}= {partial^2 theta_0} over {partial z^2}]  [eqn: rm `theta sub 0 ~=~ gamma (z -1 )~~~~~at ~~ tau =0]  [eqn: rm ``theta sub 0 ~=~ gamma - tau ~~~~`~~~~at ~~z=0][eqn: rm ``theta sub 0 ~=~ 0 ~~~~~~~~~~``~ at ~~z=1]

(5a)

(5b)

(5c)

(5d)

여기서 [eqn: gamma]는 [eqn: rm ({T_u - T_b}) / DeltaT_phi] 로서 초기의 온도기울기를 나타내는 무차원수이다. 위의 식 (5)로부터[eqn: rm {partial theta_0 } / {partial z}]를 구하여 식 (4-b)에 대입하면 최종적으로 다음과 같은 식이 얻어진다.

[eqn: rm (D^2 -a^*^2 )^2 w^* = a^*^2 theta^* - 1 over Pr left[ 1 over 2 zeta D^3 w^* -{ a^*^2 } over 2 zeta Dw^* +a^*^2 w^* right]]

(6a)

[eqn: rm- zeta over 2 D theta^* + {Ra}^* w^* left [ {zeta e sup { ( - {zeta sup2 } over 4 )} } over root pi - LEFT ( {zeta sup2 } over 2 +1 RIGHT ) {rm erfc } LEFT ( zeta over 2 ) RIGHT ) + gamma^* right ] = D^2 theta^* - a^*^2 theta^*]

(6b)

위에서[eqn: rm Ra^*],[eqn: rm a^*],[eqn: theta^*],[eqn: rm gamma ^*],[eqn: zeta]는 각각 무차원열침투깊이와 관련이 있는 [eqn: root tau]를 이용하여 정의되었으며 각각[eqn: rm tau^3 Ra_phi],[eqn: rm root tau a],[eqn: theta_1],[eqn: rm gamma / tau],[eqn: rm z / root tau]를 나타낸다. 그리고 이때의 경계조건은 다음과 같다.

[eqn: rm w^* = Dw^* = D theta^* =0~~~~~at~~ zeta =0]

(6c)

[eqn: rm w^* = Dw^* = theta^* =0~~~~~for~~ zeta rarrow inf]

(6d)

교란량들이 시간에 따른 변화량이 거의 없다는 가정을 하여 동결시간 모델을 사용하면 식 (3)에서

[eqn: rm LEFT ( ∂sup2 over {∂z sup2 } - a sup2 RIGHT ) sup2 {w } sub1 ~=~ a sup2 {theta } sub1 ~~~~~~]

(7a)

[eqn: rm LEFT ( ∂sup2 over {∂z sup2 } - a sup2 RIGHT ) {theta } sub1 ~=~ {rm Ra sub phi } {w } sub1 { ∂theta sub0 } over { ∂ z }]

(7b)

로 되며 경계조건은 다음과 같다.

[eqn: rm w sub1 = {partial w sub1 } over { partial z }= { partial theta sub1 } over { partial z } = 0 ~``~~{ at} ~~ z =0]

(7c)

[eqn: rm `w sub1 = { partial w sub1 } over {partial z } = theta sub1 = 0 ~~~~`~~{ at } ~~ z =1]

(7d)

이상의 교란식과 경계조건을 만족하는 해를 상미분방정식의 경계치 문제로서  수치해석적으로 구하는데, 먼저 Pr,[eqn: rm Ra^*]와 [eqn: zeta]=0에서의[eqn: rm Dw^*^3],[eqn: theta^*]를 가정하고 지배방정식이 균일하므로 [eqn: rm D ^2 w^*]를 지정할 수 있다. 초기온도기울기를 지정하고 상부경계면에서 [eqn: rm w^*],[eqn: theta^*]=0이 되는 고유치[eqn: rm Ra^*],[eqn: rm a^*]를 구하여 그 최소값을 찾아야 한다. 시간이 작은 영역에서는 기본 온도 분포의 비선형성이 심하게 나타나고 온도변화가 고정 경계면 사이에 국한된다고 가정할 수 있으며 경계 조건이 모두 주어지지 않았으므로 shooting 방법이 사용될 수 있다. 이때 구하여진 [eqn: rm Ra^*],[eqn: rm a^*]와 초기 온도기울기 [eqn: gamma^*]를 이용하여 각각의 Pr에 대한[eqn: rm {Ra_phi} _c],[eqn: rm a_c]를 구할 수 있다. 그리고 대류발생시점인 [eqn: rm tau_c]를 구할 수 있다.

결과 및 토론

수치해석 결과 구해진[eqn: rm {Ra _phi}_c]를 Pr의 함수로 도시하면 Figure 2와 같고 대류에 의해 발생하는 셀의 크기와 관련된 무차원 임계파수[eqn: rm a _c]와 Pr의 관계는 Figure 3에 도시한 바와 같다. 그림들에서 알 수 있듯이 Pr가 작은 영역에서는 Pr이 작아 질수록 계가 급격히 안정해지며 Pr이 큰 영역에서는 Pr이 안정성에 큰 영향을 미치지 않음을 알 수 있다. 파수와 셀의 크기는 반비례하므로 Pr가 증가할수록 대류에 의한 셀의 크기가 커짐을 알 수 있다. Alhers 등[4]의 실험결과와 Foster[6]의 증폭이론에 따른 예측치를 위에서 계산한 값들과 비교해서 나타낸 것이 Figure 4이다. 실험적으로 대류를 관찰할 수 있는 시점이 [eqn: rm 4tau_c]가 된다는  Foster[7]의 제안을 따르면, 그림에서 보듯이 동결시간 모델보다는 전파이론으로 계산한 결과가 [eqn: rm -0.3 <~gamma~ <0]에서 거의 실험결과에 근접해 있음을 알 수 있다.

참고문헌

1.	Kang, K. H., Choi, C. K. and Hwang, I. G.: AIChE J., 46, 15(2000).

2.	Lick, W.: J. Fluid Mech., 21, 565(1965).

3.	Currie, I. G.: J. Fluid Mech., 29, 337(1967).

4.	Ahlers, G., Cross, M. C., Hohenberg, P. C. and Safran, S.: J. Fluid Mech.,  110, 297(1981).

5.	Carslaw, H. S.  and Jager, J. G.: □Conduction of Heat in Solids□, 2nd ed., Oxford Press, Clarendon(1959).

6.	Foster, T. D.: Phys. Fluids, 8, 1249(1965).

7.	Foster, T. D.: Phys. Fluids, 12, 1282(1969).

Figure 2. Critical Ra□ with respect to Pr

Figure 1. Schematic diagram of present

system

Figure 3. Critical wavenumber with respect

to Pr

Figure 4. Comparison of present results with

data of Alhers et al.'s [4]