2.1 수직농장 수치해석 모델

본 연구에서는 CFD 기반 수직농장 증산작용 연동 모델의 개발을 위해, Fig. 1에 제시된 바와 같이 가로 × 세로 × 높이(5.84 m × 3.84 m × 3.86 m) 크기의 실내 공간 내에 3단 2열 재배 구조의 수직농장 환경을 모델링하였다.⁽¹⁰⁾ 공조 시스템은 측면 벽에 설치된 3개의 디퓨저(0.6 m × 0.2 m)를 통해 23.25℃의 공기를 공급하고, 하단에 위치한 배기구(0.4 m × 0.6 m)를 공기를 배출하는 방식으로 구성되어 있다. 디퓨저를 통해 공급된 유량은 각각 1,054.1 m³/h, 864 m³/h, 1,131.8 m³/h로 설정되었다. 작물 베드(1.35 m × 2.55 m)의 수평 면적과 0.15 m의 높이를 고려했으며, 재배 작물로는 높이 0.1 m의 상추를 적용했다. 광합성 유도를 위해 65 W LED 조명을 각 재배 단마다 8개씩, 총 6개 단에 균일하게 배치하였다. LED는 광합성 유효 파장 영역인 400－700 nm의 파장을 포함하며, 개별 LED의 방사 에너지는 Plank-Einstein의 상관관계에 의해 22.8 W로 산출되었다.⁽¹²⁾ 이를 바탕으로 산출된 작물 재배 면적을 고려한 광합성유효복사 에너지$I_{{PAR}}$는 53.0 W/m²이며, 이에 상응하는 광합성 유효광량자속 밀도 ${PPFD}$는 256 μmol/(m²․s)에 해당한다.

Fig. 1 Geometry of crop cultivation area and HVAC system in vertical farm.

2.2 수치해석 방법 및 경계 조건

Fig. 1의 수직농장 형상을 기반으로, 작물의 증산모델을 적용한 CFD 수치해석 모델의 구성 및 적용 기법을 Table 1에 정리하였다. 모델 개발은 상용 프로그램인 STAR-CCM+ 2410 version을 사용하였다. 난류 모델은 Reaynolds-Averaged Navier-Stokes을 고려했으며, 실내 열․기류 해석에 일반적으로 활용되는 Realizable $k$–$\varepsilon$ model을 적용했다.^(13,14) 연속체는 비압축성으로 고려하고, 작물의 증산 현상을 고려하기 위해 공기와 수분의 다성분 유체인 습공기를 적용했다. 작물은 다공성 매질로 설정하고 Darcy-Forchheimer 방정식 기반의 공기 유동 점성과 관성 저항을 고려했다.⁽¹⁵⁾ 대류항의 이산화에는 Second-order upwind scheme을 적용하여 해석 결과의 정밀성을 고려했으며, 압력과 속도의 연동 해석을 위해 SIMPLE 알고리즘을 적용했다. 복사열에 의한 영향을 고려하기 위해 매질의 흡수, 방출, 산란 현상을 고려할 수 있는 Discrete ordinates method(DOM)을 적용했다. 또한, 부력 효과는 수치해석의 계산 효율성을 고려하여 Boussinesq 근사 기법을 적용했다.

본 시뮬레이션에서는 Polyhedral 격자를 기본 구조로 생성했으며, 급․배기구 주변의 경계층 유동 해석 정밀도를 확보하기 위해 Prism layer mesh 기법을 적용했다. 또한, 격자 독립성(Grid independence) 검토를 위해 격자수를 단계적으로 조정하면서 해석을 수행했다. 검토한 격자 개수는 Coarse(약 1,280,000), Selected(약 1,562,000), Fine(약 1,760,000) 조건에서, 급․배기구 중심에서 0.2 m 거리에서 지면부터 천장까지의 연직 방향 속도 프로파일 3개를 기준으로 평가했다. 속도 결과 Coarse와 Fine조건의 속도 RMSE는 0.13 m/s, Selected와 Fine조건의 속도 RMSE는 0.04 m/s의 결과를 얻었다. 본 연구에서는 시뮬레이션의 정확도와 효율을 고려하여 Selected조건으로 수행했다.

Table 2에 본 연구에서 수행한 CFD 시뮬레이션에 적용된 경계 조건을 정리하였다. 급기구는 속도 입력 조건으로 설정하여 실험으로 얻어진 온도와 속도 조건을 적용했다.⁽¹⁰⁾ 이때 적용한 유입 기류의 난류 강도는 10%, 난류 길이는 0.021 m로 설정했다.⁽¹¹⁾ 시뮬레이션에는 상대습도 조건에 따른 4가지 시나리오(30%, 50%, 70%, 90%)를 구성하였으며, 증산작용의 효과를 확인하기 위해 증산모델이 적용되지 않은 Case 0을 추가했다. 각 시나리오 조건에서의 공기 및 수증기 질량 비율을 Table 3에 제시했다.

배기구는 출구 압력 조건으로 정압$p_{{s}}$가 0이 되도록 설정했다. 벽체와 LED 조명의 본체는 단열 조건으로 가정하였으며, LED 조사면은 열유속 761.7 W/m²를 갖는 Neumann 조건으로 설정했다.⁽⁹⁾ 작물 재배 구역에서의 유동 저항에 의한 운동량, 증산에 의한 현열 및 수증기 변화는 체적 소스와 싱크로 적용하였으며 이에 대응하는 관계식은 열평형에 의해 유도된 식(8)-(10)을 적용했다. 습공기에서 $Y_{{air}}$와 $Y_{{H}_{2}{O}}$는 공기와 수증기의 질량분율을 의미한다.

2.3 실내 열․기류환경 및 증산모델의 연동

작물과 주변 환경은 복사$R_{{n}{et}}$, 현열$Q_{{sen}}$, 잠열$Q_{{lat}}$의 열적 상호작용을 통해 에너지 교환이 이루어지며, 식(1)의 열평형 방정식을 만족한다.⁽¹⁵⁾ $R_{{n}{et}}$은 작물에 전달되는 총복사량, 현열은 주변 환경과 작물과의 대류 열전달, 잠열은 증산으로 발생하는 수증기량과 연관이 있다.

주변 환경으로부터 작물이 전달받는 복사량은 식(2)와 같이 표현된다.

여기서 $\alpha_{{r}}$는 작물의 반사계수, $I_{{PAR}}$는 재배 면적당 광합성유효복사 에너지로 53.0 W/m², ${CAC}$는 작물의 투영 잎 면적 대비 재배 면적 비율로 100%를 적용했다. 현열과 잠열 교환량은 식(3)과 식(4)로 표현된다.

여기서 ${LAI}$는 엽면적지수로 재배 면적에 대한 총 잎 면적의 비율을 의미한다. 현열 교환은 잎이 표면 온도$t_{{l}}$와 주변 환경의 습공기 온도$t_{{a}}$ 차이에 의해 발생하며, 잠열 교환은 잎 표면에서의 절대습도$\chi_{{l}}$와 주변 공기의 절대습도$\chi_{{a}}$ 차이에 따른 수증기 물질 전달 현상으로 나타난다. $\rho_{{a}}$와 $c_{{p}}$는 습공기의 밀도와 정압비열이다. $L_{{v}}$는 물의 기화 잠열을 나타내며, 2.45 × 106 J/kg을 적용했다. $r_{{a}}$는 공기역학적 저항으로, 잎 표면 경계층에 의해 발생하는 저항을 의미하며, $r_{{s}}$는 기공 저항으로 수증기가 기공을 통해 배출될 때 발생하는 저항을 나타낸다.

잎 표면에서의 절대습도는 포화수증기로부터 도출되며, 포화수증기 압력은 엽온에 의해 결정된다. 따라서, 식(1)-(4)에서 유일한 미지 변수는 엽온 $t_{{l}}$이며, 엽온은 식(5)을 통해 정리된다. 엽온은 반복 계산을 통해 도출되어야 하며, 연산 과정은 다음 Pseudocode에 나타냈다.

Graamans et al.⁽¹⁶⁾은 이 두 저항을 식(6)와 식(7)과 같이 제안했다.

여기서 $h$는 작물의 높이인 0.1 m, $|{u}|$는 기류 속력, ${PPFD}$는 광합성 유효광량자속 밀도로 256 μmol/(m²․s)을 적용했다.

위의 열평형 방정식에 의해 도출된 결과를 CFD에 적용하기 위해서는, 이를 체적 기반의 열전달 현상으로 변환할 필요가 있다. 이는 ${LAD}={LAI}/ h$의 상관식을 통해 식 (8)-(10)으로 정리된다. 여기서, $S_{{m}}$은 유동 저항에 의한 운동량 변화를, $S_{{e}}$는 현열에 의한 에너지 전달, $S_{{t}}$은 증산에 의한 수증기 발생량을 각각 나타낸다. $C_{{d}}$는 작물의 항력계수, ${u}$는 작물 재배 구역의 기류속도 벡터이다.

2.4 엽온 계산을 위한 알고리즘

현열에 의한 열전달 $S_{{e}}$와 증산에 의한 수증기 발생량 $S_{{t}}$를 도출하기 위해서는 엽온 $t_{{l}}$을 알아야만 한다. 이를 계산하기 위한 CFD에 적용할 알고리즘을 정리했다. 작물 재배 구역의 습공기 상태인 기온 $t_{{a}}$, 상대습도 ${RH}$, 기류속도 벡터 ${u}$, 셀의 개수$n$, 습공기의 밀도 $\rho_{{a}}$, 습공기의 정압비열 $c_{{p}}$를 입력데이터로 활용하여 엽온 $t_{{l}}$의 계산을 수행한다. 셀의 개수는 계산이 적용될 셀의 개수를 확인하고, 각 셀마다의 반복 계산을 통해 엽온 $t_{{l}}$을 도출했다. 각 셀마다 식(1)-식(4)의 관계식으로부터 열평형을 만족하는 $t_{{l}}$는 Newton-Raphson반복법에 의해 간단하게 도출되었다. 각 셀마다 계산된 엽온 $t_{{l}}$을 활용하여 식(9)와 식(10)로부터 작물의 증산에 의한 주변 공기의 냉각과 수증기 발생량을 구현했다.

3.1 실내 기류 특성

CFD에 의해 작물이 재배되고 있는 실내 기류 환경의 재현성을 검증하기 위해 기존 문헌과 비교 분석을 수행했다. Kang et al.⁽¹⁰⁾의 연구에서는 실내로 공급되는 공기의 조건으로 상대습도 78.2%에 해당하는 기온 23.25℃, 절대습도 16.4 g/m³의 습공기가 사용되었으며, 동일한 조건을 본 시뮬레이션에 적용하여 비교 검토했다. Fig. 2(a)에 나타낸 위치에서 각각 4개의 높이에 해당하는 속력의 기존 문헌의 실험 결과와 본 연구의 시뮬레이션 결과를 비교하였고, Fig. 2(b)는 그 결과를 나타낸다. $x$축은 실내 기류의 최대 속력인 2.62 m/s, $y$축은 실내 높이를 기준으로 각각 정규화하여 표현했다. 비교 결과, Line 1에서 0.050, Line 2에서 0.029의 RMSE값을 보였다. 또한, Fig. 2(c)에는 각 급기구를 통해 유입된 공기의 유선(Streamline)을 시각화했으며, 세 번째 급기구를 통한 기류가 가장 균일하게 확산됨을 확인했다.

Fig. 2 Validation of simulated velocity magnitude against experimental data reported by Kang et al.⁽¹⁰⁾.

3.2 습도 변화가 증산량과 주위 열환경에 미치는 영향

공조 기류의 상대습도는 작물의 증산작용에 결정적인 영향을 미치며, 이는 실내 열환경 변화로 직결된다. 급기 기류의 상대습도 변화가 증산작용에 미치는 영향을 분석하기 위해, 실내 평균 기온, 엽온, 작물로부터 발생한 수증기 발생량을 Table 4에 정리하였다. 시뮬레이션 결과, 상대습도가 증가할수록 증산량은 점진적으로 감소하는 경향을 보였으며, 이에 따라 실내 기온과 엽온은 모두 상승하는 양상을 나타냈다. 표준편차는 상대습도가 증가할수록 온도, 엽온, 증산량 모두 전반적으로 감소하는 결과를 보였다. 특히, 엽온의 경우, 상대습도 70% 이하 조건에서는 실내 기온보다 낮은 값으로 유지되었으나, 상대습도가 90%로 높은 조건(Case 4)에서는 오히려 실내 기온보다 엽온이 높게 형성되었다. 이러한 결과는 LED 조명으로부터 발생한 복사열이 증산에 의한 냉각 효과보다 우세하게 작용했기 때문으로 해석된다. 그러나, 최대, 최소와 같은 국부적인 조건에서는 작물의 생리반응에 의해 실내 전체 평균의 패턴과는 다른 양상을 나타내는 경향도 나타났다.

일반적으로 상대습도가 낮을수록 잎의 기공을 통한 수분 증발량이 증가하게 되며, 이는 증산에 수반되는 잠열 손실로 인해 엽온이 감소하는 결과를 초래한다. 이러한 엽온 저하는 주변 공기와의 대류 열전달을 통해 실내 기온 하강을 유도한다. 이러한 열․수분 교환 메커니즘을 바탕으로, 공조 기류의 상대습도 조건이 실내 온도 분포에 미치는 영향을 Fig. 3에 나타냈다. 시뮬레이션 결과 단면 온도 평균 및 편차는 Case 1(Section A: 23.4±0.88℃, Section B: 23.5±1.47℃), Case 2(Section A: 24.0±0.84℃, Section B: 24.1±1.43℃), Case 3(Section A: 24.5±0.88℃, Section B: 24.7±1.15℃), Case 4(Section A: 25.1±1.00℃, Section B: 25.3±1.22℃)로, LED 조명에 의한 발열은 실내 온도를 상승시키는 주요 요인으로 작용한다(Table 5). 한편, 상대습도가 30%인 공조 기류가 적용된 Case 1에서는 실내 평균 온도가 23.4℃의 결과를 보였으며, 작물 주변에서는 증산작용에 의해 상대적으로 낮은 온도 분포가 형성되는 경향이 확인되었다.

상대습도 증가에 따른 증산작용의 변화를 정량적으로 확인하기 위해, 각 시나리오에서의 실내 평균 기온을 비교하였다. 시뮬레이션 결과, Case 1에서 23.4℃, Case 2에서 23.9℃, Case 3에서 24.4℃, Case 4에서 24.8℃로 나타났다. 이는 습도가 높아질수록 증산에 따른 수분 증발량이 감소하고, 그에 따라 작물 주변에서의 냉각 효과가 저하되면서 실내 기온이 점진적으로 상승하는 경향을 나타낸 것으로 해석된다.

Fig. 4에서 공조 기류의 상대습도 변화에 따른 실내 습도 분포의 변화를 검토했다. 시뮬레이션 결과 단면 상대습도 평균 및 편차는 Case 1(Section A: 36.1±8.77%, Section B: 36.9±10.7%), Case 2(Section A: 53.3±6.49%, Section B: 53.4±7.61%), Case 3(Section A: 68.9±4.17℃, Section B: 69.1±5.27%), Case 4(Section A: 83.1±3.39%, Section B: 82.6±4.14%)로 나타났으며, 작물 재배 영역에서는 상대습도가 80–100% 수준으로 형성되었다(Table 5). 이는 개발된 CFD 알고리즘이 증산에 따른 수증기 발생 효과를 적절히 구현하고 있음을 나타낸다. 또한, 급기 기류의 초기 상대습도가 증가할수록 작물 주변의 국소 상대습도도 함께 상승하는 경향이 나타났다. 이러한 분석 결과는 제안된 CFD 기반 증산모델이 수직농장 내 미기후의 공간적 이질성을 정량적으로 해석할 수 있음을 보여주며, 향후 환경 제어 전략의 정밀 설계 및 최적화를 위한 실질적 해석 도구로서의 활용 가능성을 입증한다.

Fig. 3 Temperature distribution under varying relative humidity conditions in each case.

Fig. 4 Relative humidity distribution under varying relative humidity conditions in each case.

3.3 증산모델 적용에 따른 작물 환경 변화의 차이

증산모델 적용의 중요성을 확인하기 위해 온실 혹은 수직농장에서 일반적으로 검토되고 있는 열․기류만 적용한 결과와 비교를 수행했다. 증산작용을 제외한 Case 0과 Case 3의 온․습도 결과를 Fig. 5와 같이 비교했다. 시스템 공조 설정은 동일했지만, 시뮬레이션 결과는 다음과 같았다. 단면 온도 평균 및 편차는 Case 0(Section A: 25.6±1.21℃, Section B: 25.9±1.40℃), Case 3(Section A: 24.5±0.88℃, Section B: 24.7±1.15℃), 단면 상대습도 평균 및 편차는 Case 0(Section A: 60.4±3.98%, Section B: 59.3±1.40℃), Case 3(Section A: 68.9±4.17℃, Section B: 69.1±5.27%)의 결과로 증산모델 적용의 유무에 의해 상당한 차이를 보였다. 이러한 결과는 증산모델의 적용이 실제 작물의 생육환경 및 제어에 매우 중요한 요소임을 시사한다.

Fig. 5 Comparison of temperature and humidity with and without transpiration model.