비등온, 비정상상태에서의 흡착 및 탈착 공정에 관한 수치모사

최병석□이종협

서울대학교 화학공학과

Numerical Simulation on Nonisothermal and Unsteady State

Adsorption/Desorption Process

Byung-Seok Choi and Jongheop Yi

Dept. of Chemical Engineering, Seoul National University

서론

흡착과 탈착 공정은 기체 세정과 분리 등의 환경 분야에 광범위하게 사용되고 있다. 흡착 공정은 기체 흐름에서 탈수나 불순물을 제거하는 공정으로 이용되고 있으며, 경제적 가치가 높은 묽은 농도의 화학종을 재회수하는데도 주로 이용된다. 흡착 공정에 영향을 미치는 중요 인자로써 흡착탑 내부에서 유동의 분배, 흡착제의 다공도, 기공 크기, 기공 구조 및 균일도 등과 함께, 흡착 대상 화학물의 특성 등이 있다 [1]. 그동안 여러 가지 종류의 흡착 모델식이 제안되었으나 실제 계산에서 주로 사용된 것은 막확산 저항을 무시하고 단순히 선형 평형 관계식을 가정하는 것이었다 [1]. 본 연구에서는 전산 유체 역학 (Computational Fluid Dynamics: CFD) 기법을 적용하여 비등온, 비선형 흡착에 대한 수치 해석적 모델을 실제로 구현하고 흡착탑 내에서의 유동, 온도, 압력 및 농도 분포 등을 예측하고자 하였다. 또, 비정상상태 모사를 통하여 흡착 과정과 재생 과정에서의 시간에 따른 각 분포의 변화를 파악하여 조업 변수를 제시하고자 하였다. 현재, CFD 기법은 기계, 항공 분야를 비롯하여 소각로, 전기 집진기 및 사이클론 등의 환경 분야에서도 공정의 설계와 개선을 위한 중요한 도구로서 이용되고 있다 [2]. 본 연구를 통하여 CFD 기법을 이용한 흡착탑 설계의 타당성을 파악할 수 있으며, 실험 결과와의 비교를 통하여 scale-up을 위한 중요한 설계 인자도 제공할 수 있다.

이론

연속 방정식과 운동량 보존 방정식으로부터 흡착탑 내부의 유동을 표현하고 k-□ 난류 모델을 추가하여 유동 에너지의 발생과 소산을 함께 고려하였다. 유동 방정식에는 흡착층에 의하여 발생하는 운동량의 강하항이 추가되었으며 흡착 및 탈착에 의한 에너지 변화를 기상과 고체상 각각에 대하여 표현하였다. 또 기상과 흡착상의 농도 변화를 흡착등온식 함께 계산할 수 있도록 구현하였다.

각 방정식은 다음과 같이 쓸 수 있다.

연속 방정식 :

[eqn: partialrho over partialt + "div"`(rho vec v )=0]	(1)

축방향 운동량 :

[eqn: &partial over partialt (rhow) + 1 over r cdot left[partial over partialz (r rhow cdot w) + partial over partialr (r rhov cdot w) right]##=&- partialp over partialz + partial over partialz left(mu_eff partialw over partialz right)+ 1 over r partial over partialr left(r mu_eff partialv over partialz right) -S_Ergun,w]	(2)

반경방향 운동량 :

[eqn: &partial over partialt (rhov) + 1 over r cdot left[partial over partialz (r rhow cdot v) + partial over partialr (r rhov cdot v) right]##=&- partialp over partialr + partial over partialr left(mu_eff partialw over partialr right)+ 1 over r partial over partialr left(r mu_eff partialv over partialz right) - {2mu_eff v }over r^2 -S_Ergun,v]	(3)

여기서 □는 유체 밀도, w와 v는 각각 축방향과 반경 방향의 속도 성분을 의미하고 □eff는 유효 점도 계수이다. 이 때, SErgun항은 흡착층에 의한 운동량 손실을 나타내며 Kast(1994)에 의하면 다음과 같이 정의된다 [3].

[eqn: S_Ergun,w = 150 mu_eff {(1-epsilon_f )^2} over {epsilon_f^2 d_e^2} w + 1.75 rho {(1-epsilon_f ) }over {epsilon_f d_e} w |vec v |]	(4)

[eqn: S_Ergun,r = 150 mu_eff {(1-epsilon_f )}^2 over {epsilon_f ^2 d_e^2} v + 1.75 rho {(1-epsilon_f ) }over {epsilon_f d_e} v |vec v |]	(5)

이 때, □f는 흡착층의 void 분율이며, de는 유효 거리로써 흡착제 입자의 부피를 표면적으로 나눈 값을 사용한다.

기상으로부터 고체 흡착제로의 용질 흡착 속도를 나타내기 위하여 균일 고체상 확산 모델을 사용하며 막확산과 표면에서의 평형 흡착, 내부에서의 비정상상태 전달 메카니즘이 가정되었다.

기상과 흡착 고체상에 대하여 각각의 농도 보존 방정식을 사용하며 각 방정식은 다음과 같다.

기상 흡착물질 농도 :

[eqn: rho partialC_g over partialt = &delD_eff del C_g - del vec v C_g ## &+3rho(1-epsilon_f ) k_f over {epsilon_f r_p } (C_p vert _r=r_p -C_g )]	(6)

고체상 흡착물질 농도 : [eqn: partialC_s over partialt = k_i D_pore,i psi_p alpha_f (C_p -C_s )]	(7)

기상 농도식에서 Cg는 흡착 물질의 기상 농도, Deff는 기상에서의 유효 확산 계수, kf는 입자 표면에서의 물질 교환 계수, rp는 입자 반경을 의미하며 Cp는 흡착 물질의 고체상 내부 농도를 나타낸다.

고체상 농도식에서 Cs는 입자 내부의 흡착된 상의 농도를 의미하며 ki는 입자 내부에서의 물질 전달 계수, Dpore,i는 기공에서의 확산도, □f는 단위 체적당의 비표면적을 나타낸다. [eqn: psi_p `]는 입자 uptake에 대한 시간 곡선의 실제 기울기를 나타내는데 거의 1에 가까운 값이다 (Perry and Chilton, 1973). 이 식에서 고체 표면에서는 유체와 아주 빠르게 열역학적인 평형을 이루며 기공 확산도는 입자의 macropore 내부 운동에 적용되므로 유체상 내에서의 값보다 매우 작은 값을 가지는 것으로 고려되었다. 이 때, 흡착 물질의 고체상 내부 농도와 기상 농도 사이에서의 평형 관계식은 비선형인 Freundlich isotherm을 사용하며 그 식은 다음과 같다.

Freundlich isotherm :

[eqn: C_p =k(T) C_g^beta `]	(8)

여기서 흡착 상수, k(T)는 온도와 흡착열의 함수이며 그 관계는 Van't Hoff 식에 의하여 정의되고 □는 흡착시 기상 농도의 비선형 의존성을 나타내는 상수이다.

입자 내부의 온도는 균일하며 흡착층이 단열 조건이고 입자가 균일한 구형이며 고체의 물성이 온도에 의해서 변화지 않는다고 가정하면 기상과 고체상에 대한 에너지 보존 방정식은 다음과 같이 표현된다.

기상 에너지 :

[eqn: rho C_ph partialT_g over partialt = &delk_L del T_g - rho C_ph del vec v T_g ##&-left({1-epsilon_f } over epsilon_f right) 3h over r_p (T_g - T_s vert _r=r_p )]	(9)

고체상 에너지 :

[eqn: rho_p C_ph partialT_s over partialt = h {6(1-epsilon_p )} overd_p (T_g - T_s )+k_i lambda_ad (C_p -C_s )]	(10)

기상 에너지 보존 방정식에서 Cph는 열용량, Tg는 기상의 온도, kL은 열전도도이며 h는 국지 열전달 계수를 의미하고 Ts는 고체상의 온도이다. 고체상 에너지 식에서 □p는 입자의 밀도, □p는 입자의 다공도, dp는 입자의 직경을 나타내고 □ad는 흡착열을 의미한다.

난류 모델로는 표준 k-□방정식을 사용하며 난류 에너지, k와 점성 소산 속도, □ 각각에 대한 보존 방정식은 다음과 같다.

난류 에너지, k 방정식 :

[eqn: &partial over partialt (rho k) + 1 over r left[partial over partialz (r rho w k )+ partial over partialr (r rho uk )right]##& - 1 over r left[ partial over partialz left(r mu_eff over sigma_k partialk over partialz right) - partial over partialr left(r mu_eff over sigma_k partialk over partialr right)right] = G-C_D rho epsilon]	(11)

난류 소산 속도, □ 방정식 :

[eqn: &partial over partialt (rho epsilon ) + 1 over r left[partial over partialz (r rho w epsilon )+ partial over partialr (r rho v epsilon )right] ##-& 1 over r left[ partial over partialz left(r mu_eff over sigma_epsilon partialepsilon over partialz right) + partial over partialr left(r mu_eff over sigma_epsilon partialepsilon over partialr right)right] =C_1 epsilon over k G-C_2 rho epsilon^2 over k]      (12)

이 방정식에서 G는 물리적으로 난류 전단 응력에 의한 평균 흐름으로부터의 난류 에너지, k의 생성을 나타내며 다음과 같이 정의된다.

[eqn: G= mu_t left{ 2 left[ left( partialw over partialz right)^2 +left( partialv over partialr right)^2 +left( v over r right)^2 right] + left( partialw over partialr + partialv over partialz right)^2 right}]              (13)

Reynolds수가 클 때 k-□ 방정식 내의 parameter는 상수로 여겨지며 그 값들은 Table 1과 같다 [4]. 특히, 난류 에너지 및 점성 소산 속도를 계산할 때 wall function의 영향을 고려하여 국지 Reynolds 수가 11.63이하인 벽면 근처 영역에서는 viscous sublayer를 나누어 계산을 수행한다.

수치해석

반경이 잡작스럽게 확대되는 관형 흡착탑 내부에 흡착제를 충진하고 293K, 6.55bar의 흡착 조업 조건과 568K, 1.2bar의 탈착 조업 조건에 대하여 전산 유체 역학 기법을 적용하여 비정상상태 모사를 수행하였다. 독립 변수로는 축방향과 반경 방향 변수와 시간 변수가 사용되었고 종속 변수로는 두 방향의 속도 성분, 압력, 흡착상과 기상의 농도, 기상 및 고체상의 온도, 그리고 난류 에너지와 점성 소산 속도가 계산되었으며 이들 변수값에 관하여 결과를 도시하는 것이 가능하였다. 2차원, 축대칭 흐름에 관하여 각 지배방정식에 대한 유한 차분식을 유도하였으며 이들을 제어체적에 관하여 적분하고 TDMA(Tri-Diagonal Matrix Analysis) 기법을 사용하여 각 변수에 대하여 해를 구하였다. staggered grids를 이용한 계산 방법을 사용하였으며 전체 계산 영역을 축방향 20개, 반경방향 20개의 grid로 나누고 주입구 쪽과 벽면 영역으로 grid의 밀도를 조밀하게 배치하였다. 계산 시스템의 grid가 Fig. 1에 제시되어 있다. 속도-압력 중첩 알고리즘은 SIMPLE법을 사용하였으며 이산화 방정식 계수는 upwind법으로 유도하였다 [5]. 프로그램 사용 언어로는 FORTRAN이 이용되었으며 계산은 PC-Pentium Pro 200으로 수행되었다.

결과 및 토론

흡착 공정과 탈착 공정 각각에 대하여 흡착층 내부의 속도 및 압력 분포를 계산하였으며 시간 변화에 따른 온도 변화와 농도 변화 곡선을 알아냄으로써 흡착 완료 시간을 예측하였다. 이를 이용하여 breakthrough curve를 계산하고 흡, 탈착 교환 시기를 제안할 수 있었다. 또한 입자의 크기, 다공도의 변화에 따른 흡착 성능 비교 및 흡착 온도, 압력 변화에 따른 성능 비교를 통하여 최적의 조업 조건을 제시할 수 있으므로 본 연구를 통하여 CFD를 이용한 흡착탑의 수치 모사는 흡착 공정의 설계 및 개선에 중요한 방법론으로 적용될 수 있음을 파악할 수 있다.

참고문헌

1. Gouvalias, G.S. and Markatos, N.C., Mathematical modeling of heat and mass transfer in packed-bed adsorbers/regenerators, AIChE J., 39, 1799,(1993).

2. Hoffmann, Klaus A., Computational Fluid Dynamics for Engineers, Pub. of Eng. Edu. Sys., Texas, USA(1989).

3. Oliver Bey and Gerhart, Fluid flow through catalyst filled tubes, Chem. Eng. Sci., 52, 1365(1997).

4. Gosman, A.D. and Ideriah, F.J.K., A General Computer Program for Two-Dimensional, Turbulent, Recirculating Flows, Report No. FM-83-2, CA, USA(1983).

5. Patankar, S.V., Numerical Heat Transfer and Fluid Flow, McGraw-Hill, New York, USA(1980).

Table 1. Constants of k-□ model

Constant	C□	CD	C1	C2	□k	□□
Value	0.09	1.0	1.44	1.92	1.0	1.3

Fig. 1. Grid system of the calculation domain.