2.1 쉐브론 형상 및 상관식

본 연구에서는 판형 열교환기 개발을 목적으로 쉐브론 형상 전열판에 대한 전산해석을 수행하였다. 쉐브론 형상의 개략도는 Fig. 1과 같으며, 주요 설계인자는 쉐브론 각도와 단면 형상 파라미터들이다. Fig. 1(a)는 쉐브론 형상의 평면도로 쉐브론 각도(β)를 나타내고, Fig. 1(b)는 쉐브론 단면에서 나타나는 파라미터들을 도시하였다.⁽²⁾

판형 열교환기의 열전달 및 압력손실 특성은 쉐브론 각도와 형상 파라미터들의 조합으로 결정되며 Nu수와 Fanning 마찰계수(f)로 표현된다. 관련된 무차원 수는 다음 수식과 같다. 여기서, ρ는 밀도(kg/m³), v는 속도(m/s), d는 직경(m), μ는 점성계수(N-s/m²), h는 대류열전달계수(W/m²-K), k는 열전도도(W/m-K), P는 압력(Pa), L는 전열판 길이(m)이다.

지금까지 많은 실험연구를 통하여 다양한 상관식들이 개발되어 왔으며 Ayub⁽²⁾은 이와 같은 상관식들을 총 28개로 정리하여 제시하였고, 그 중에서 대표적으로 이용되고 있는 4가지 상관식들은 Muley,⁽³⁾ Wanniarachchi,⁽⁴⁾ Kumar,⁽⁵⁾ Thonon⁽⁶⁾의 식들이다. Muley와 Wanniarachchi는 쉐브론 각도(β)와 길이비(φ)의 일반 함수식 형태로 마찰계수와 누셀트 수를 정의하였고, Kumar와 Thonon은 쉐브론 각도별로 Re 수와 Pr 수만으로 표현하였다. 하지만, Muley와 Wanniarachchi의 상관식들은 길이비(φ)에 따라 다소 차이나는 경향이 있어, Kumar 상관식이 일반적으로 통용되고 있다. 또한, 단면 파라미터에서 피치간격, 높이, 두께 등에 따라서 그 특성이 달라질 수도 있다.

2.2 전산해석 모델링 및 지배방정식

쉐브론 형상의 판형 열교환기의 채널 유로는 전열판 사이의 좁은 간격에서 복잡한 형태의 구조를 띠고 있어 이를 유동영역으로 모델링하는 데는 어려움이 따른다. 그리고 열유동해석 시에도 복잡한 유로에서 유동 및 열전달을 수렴시키기 위해서는 정밀한 해석기술이 필요하다.

전산열유동해석에 필요한 쉐브론 판형 열교환기의 유동 영역 형상을 만들기 위하여 3D 설계프로그램인 Solidworks 2017을 이용하였다. Fig. 2에 나타낸 바와 같이 우선 쉐브론 전열판 형상을 만들고, 엇갈린 두 개의 쉐브론 전열판 사이의 공간을 채널 유동 영역으로 정의하였다. Fig. 2(a)는 쉐브론 전열판의 전체적인 구조를 나타내며, (b)와 (c)는 정면도와 측면도로서 쉐브론 전열판 사이에서 형성되는 마름모 형태의 유로를 보여준다. Fig. 3은 유동방향으로 달라지는 유로의 단면 형태를 순차적으로 보여주며, 마름모꼴 단면 유로가 좌우로 반복해서 교차되는 형태로 나타난다. 이러한 형상 변화로 인하여 작동유체가 쉐브론 형상과 지속적으로 충돌하여 분리되고 합쳐지는 과정이 반복될 것을 예상할 수 있다. 이와 같은 마름모 형태 유로는 상하 전열판의 이격 정도에 따라 달라질 수 있으며 유로 최적화 시에 고려해야할 설계변수이다. 본 연구에서 사용한 쉐브론 전열판의 쉐브론 각도는 50˚이며, 전체 길이는 235 mm이고 폭은 129 mm이다. 전열판은 쉐브론 형상을 대략 16개 정도 포함하고 있어 이전 연구들에 비하여 확장된 유동 영역을 사용하였다.

판형 열교환기의 전산해석을 위해서 고온 채널 유동 영역을 tetra 요소 격자로 생성하였으며 세부적인 격자 형태는 Fig. 4와 같다. 전열판의 입출부 길이는 15 mm이며 간극은 Fig. 1의 채널 간극(b)와 같다. 그리고 적절한 격자수를 결정하기 위하여 Fig. 5와 같이 격자수(Number of cells)에 따른 입구압력의 수렴성을 확인하는 격자 테스트(Grid sensitivity test)를 수행하였다. 최적의 격자수는 대략 4백만 개 수준이다.

판형 열교환기 고온 채널에 대한 열유동 해석을 위하여 상용 CFD 해석프로그램인 FLUENT v16.2를 사용하였고, 해석을 위한 지배방정식은 다음과 같다.⁽¹⁴⁾ 여기서, ρ는 밀도, $\overrightarrow{V}$는 속도벡터(m/s), $\overrightarrow{g}$는 중력가속도(m/s²), P는 압력(Pa), i는 내부에너지(J/kg), k는 열전도도(W/m-K), T는 온도(K)를 의미한다.

2.3 경계조건 및 해석조건

쉐브론 판형 열교환기의 고온 채널 내 열유동 해석을 위한 입구 및 출구는 Fig. 6과 같으며, 입구 경계조건으로 속도는 0.1 m/s에서 1.0 m/s까지 0.1단위로 증가시켰으며 온도는 80℃를 사용하였고, 반대편 출구에서는 압력 경계조건을 사용하였다. 작동유체는 물이며 물성치 중 온도의 영향을 많이 받는 점성계수와 열전도도는 온도에 따른 값으로 정의하였지만, 좁은 채널 간극에서 밀도차에 의한 부력 효과는 미미하여 상수의 밀도값을 사용하였다. 유로의 좌우측면들은 대칭조건으로 설정하였고, 상하 전열면의 경계조건은 저온 채널과의 열전달을 모사하기 위하여 기존의 상관식으로부터 도출한 평균 열전달계수와 저온 채널의 온도를 이용하였다.

쉐브론 형상의 내부 유동에서는 유로가 복잡하여 낮은 Re 수에서 난류가 발생할 수 있기 때문에 층류 점성모델과 난류모델을 모두 사용하여 비교하였다. 난류모델은 Realizable k-ε 모델과 비평형 벽함수(Non-equilibrium wall function)을 적용하였고, 비평형 벽함수는 압력구배와 급격한 유동변화가 중요한 난류모델에서 박리, 재부착, 충돌과 같은 복잡한 유동에 적합한 것으로 알려져 있다.⁽¹⁴⁾ 그 외의 벽함수를 사용한 경우에는 수렴을 시키기 어려웠으며 이는 Jain et al.⁽¹⁰⁾이 제시한 난류모델 연구내용과 일치하는 결과이다.

쉐브론 형상 내 유로가 복잡한 구조를 가지고 있어 유동과 에너지를 동시에 수렴시키는 것은 유동조건에 따라서 어려움이 발생하였다. 그래서 유동장만을 우선 수렴시킨 후 그 뒤 에너지방정식까지 순차적으로 풀었다. 또한, 낮은 유속에서 잘 수렴시킨 해석결과를 그 다음 유속조건의 초기조건으로 설정하여 효율적으로 수렴하도록 해석을 수행하였다.

3.1 압력손실에 따른 마찰계수

쉐브론 판형 열교환기의 고온 채널부의 유동장 모사를 위하여 유속 조건별로 해석을 수행하였다. Table 1은 입구측 유속이 0.1 m/s에서 1.0 m/s까지 증가함에 따라 발생하는 압력차이가 거의 선형적으로 상승하는 결과를 나타내며, 층류 점성모델(Viscous model)보다 난류모델(Turbulence model)에서 0%에서부터 20%까지 높은 압력 차이가 발생하는 것을 보여준다. 이는 난류 운동에너지와 소산율 성분으로 높아진 난류 점성계수로 인하여 압력손실이 증가한 데 기인한 것으로 판단된다. 그리고 Fanning 마찰계수는 해석결과인 압력손실을 이용하여 계산하였고, Fig. 7(a)와 같이 기존의 상관식들과 함께 비교하여 도시하였다. 층류와 난류모델에서 얻어진 마찰계수는 기존 상관식들이 제시하는 경향과 잘 부합하는 결과를 보여주고 있어서 본 연구에서 개발한 전산해석모델이 타당함을 보여준다. 난류모델로부터 계산된 마찰계수는 낮은 Re 수에서는 Kumar 상관식과 유사하고, Re 수가 높아질수록 다른 상관식들과 거의 일치하게 됨을 확인할 수 있다. 또한, Re 수가 증가할수록 마찰계수는 점차적으로 감소하는 경향을 보인다. 그리고 층류보다 난류모델이 기존 상관식들에 더 가까운 결과를 보이는 것은 쉐브론 판형 열교환기 내부 유동이 난류에 가깝다는 것을 의미한다. 지금까지 시뮬레이션 관련 연구들도 대부분 난류모델을 사용하고 있고, 문제특성에 따라서 서로 다른 난류모델과 벽함수를 적용해왔다.

판형 열교환기의 열전달과 밀접한 Nu 수도 기존 상관식들로부터 Re 수에 따른 값들을 구할 수 있고 이를 Fig. 7(b)에 도시하였다. Nu 수는 마찰계수와 다르게 상관식마다 다소 격차가 있음을 알 수 있으며, 이는 쉐브론 형상 파라미터들을 반영하는 방식에 따라 다르게 나타난 것이다. Kumar 식과 Thonon 식은 형상 파라미터 정보가 필요하지 않지만, Muley 식과 Wanniarachchi 식은 길이비 변화에 따라 민감하게 반응하며 격차가 생기는 경향을 보인다. 본 연구에서는 일반적으로 자주 이용되고 있는 Kumar 식으로부터 계산된 Nu 수를 이용하여 표면열전달계수를 같이 도출하여 고온 채널에서 열전달 해석 시 저온 채널 표면의 경계 조건으로 사용하였다. 그리고 열유동 해석으로부터 계산한 고온 채널의 Nu 수는 Fig. 7(b)와 같으며 상관식들과 유사한 경향을 보였다.

3.2 고온 채널 내 열전달 특성

쉐브론 판형 열교환기의 고온 채널에서의 열전달 특성을 살펴보기 위하여 전산해석을 수행하고 그 결과들을 정리하였다. 복잡한 구조의 유동영역에 대한 해석에서 수렴성을 개선하기 위하여 유동장 해석을 수행하여 먼저 수렴시킨 뒤 에너지방정식을 추가하여 온도분포를 해석하였다. 그리고 상하 쉐브론 전열판에서의 열전달은 저온부 채널의 작동유체와 전열판을 통하여 이루어지는데 이를 포함한 해석은 상당한 많은 격자수가 필요하게 되어 모델링의 어려움이 발생한다. 그래서 본 연구에서는 저온 측 열전달계수와 저온부 온도 20℃를 이용하여 저온채널에서의 열전달을 경계조건으로 처리하고 고온부 채널 위주로 해석을 수행하였다.

Fig. 8은 해석결과로 얻어진 표면 열유속와 난류운동에너지 그래프이며, 두 데이터가 Re 수에 거의 선형적으로 비례하는 관계를 가진다. 저온 측 열전달계수는 Re 수에 따라 5,000 W/m²-K에서 27,000 W/m²-K 범위에 있으며 쉐브론 전열판 표면의 열유속은 Fig. 8(a)와 같이 72 kW/m²에서 330 kW/m² 까지 증가하였다. 복잡한 쉐브론 유로 내에서는 Re 수가 높아질수록 유동이 충돌 및 분리되는 과정에서 난류운동에너지(Turbulence kinetic energy)가 증가하게 되고, 이에 따라 난류강도(Turbulence intensity)도 함께 상승하며 열전달이 활발해진다. Fig. 8(b)는 Re 수에 따른 난류운동에너지와 난류강도를 그래프로 나타낸 것이며, 난류운동에너지는 지수적으로 상승하는 반면, 난류강도는 선형적으로 증가하는 경향을 보인다.

Table 2는 쉐브론 판형 열교환기 열유동 해석결과로 Re 수에 따른 출구온도, 전단응력, 평균 열유속, 전체 열유속 등을 보여준다. 전단응력은 유체가 쉐브론 형상과 충돌하며 회전하거나 박리하는 부분에서 크게 작용하게 되며 상부 전열판에 작용하는 전단응력은 Fig. 10(c)와 같다. 높은 Re 수에서 전단응력이 크게 발생하여 압력손실이 커지게 된다. 또한, 평균 표면열유속과 더불어 전열판 전체면적에서 빠져나가는 전체 열유속은 대략 2 kW에서 10 kW까지 변화한다.

판형 열교환기에서의 열전달율(Q)과 전달단위수(Number of transfer unit, NTU)는 아래 식과 같으며, 열전달율과 LMTD를 이용하여 NTU를 계산할 수 있다.⁽⁹⁾

여기에서 △T_lm는 대수평균온도차,⁽⁹⁾$\dot{m}$_h는 질량유량, C_p,h는 비열, C_min은 열용량율(Heat capacity rate), A_s는 열전달 면적, U는 총괄열전달계수(Overall heat transfer coefficient)이다.

열교환기의 유용도(Effectiveness)는 식(9)과 같이 정의되며 입출구 온도들로 표현할 수 있다.⁽⁹⁾

Fig. 9는 해석결과로부터 Re 수에 대한 출구온도와 NTU에 대한 유용도를 정리하여 나타낸 것이다. 이 유용도에 대한 결과는 대향류 열교환기(Counter flow heat exchanger)에서 알려져 있는 이론식과도 정확히 일치하는 결과이다.

여기에서 c는 열용량비(capacity rate ratio)로서C_min/C_max이다.

Fig. 10은 쉐브론 판형 열교환기의 온도 및 압력분포를 나타낸 것이며 주어진 모든 유속 조건에 대하여 거의 유사한 경향을 보인다. Fig. 10(a)는 하부 전열판의 온도와 단면 유체의 온도를 나타낸 것이다. 표면온도는 쉐브론 전열판의 굴곡에 따라 국부적으로 차이가 나지만, 유체 온도는 표면온도와는 달리 고온에서 저온으로 점차적으로 감소하는 경향을 보인다. 압력도 유동방향으로 마찰손실에 의하여 입구로부터 거의 선형적으로 감소하는 것을 알 수 있고, 벽면 전단응력은 유동이 타고 넘어가는 쉐브론 굴곡에서 집중해서 발생하는 결과를 보여준다. 전단응력이 높아지는 것은 유체가 벽면과 활발히 충돌하는 것을 의미하며, 이로 인하여 난류속도변동(Turbulence velocity fluctuation)과 난류강도가 증가하며 난류 유동을 활성화시키며 열전달 성능을 향상시키는 결과를 가져온다.

Fig. 11은 입구유속이 1.0 m/s인 경우에 쉐브론 판형 열교환기에서 발생하는 유선과 속도벡터를 나타낸 것이다. 좌측의 입구부를 통하여 작동유체가 균일하게 유입되고 쉐브론 형상 유로를 거치며 우측 출구부로 빠져나가는 결과를 보여준다. Fig. 11의 상단부분은 유동방향 단면상을 통과하는 속도벡터를 도시한 것이며 유속이 좌우로 굽이치며 전진하는 경향을 보여준다. 이는 앞서 Fig. 3에서 보인 바와 같이 상하의 쉐브론 전열판이 서로 맞물려 지그재그 형식의 유로가 발생한 것을 잘 보여준다. 쉐브론 각도와 전열판 사이 간극크기에 따라 굽이치며 돌아나가는 회전유동이 달라질 것으로 예상된다.