Original Article

Journal of Computational Fluids Engineering. 31 March 2026. 65-78
https://doi.org/10.6112/kscfe.2026.31.1.065

ABSTRACT


MAIN

  • 1. 서 론

  • 2. 연구 방법

  •   2.1 지배 방정식 및 유동 해석자

  •   2.2 벽면 거리 계산 방법

  • 3. 해석 결과

  •   3.1 HB-2 형상

  •   3.2 ONERA M6 날개

  •   3.3 NASA CRM(Common Research Model)

  • 4. 결 론

1. 서 론

최근 컴퓨팅 성능이 비약적으로 향상됨에 따라 항공우주공학 분야에서도 전산유체역학(Computational Fluid Dynamics, CFD)을 활용한 비행체의 공력 해석 및 설계가 활발히 수행되고 있다. 일반적인 비행체 또는 발사체는 매우 빠른 속도로 비행하므로, 그 주위에는 압축성 효과를 동반한 난류 유동이 발달한다. 이러한 난류 유동을 효과적으로 해석하기 위해 다양한 모델링 및 해석 기법이 개발되어 왔으며, 대표적인 접근 방법으로는 Reynolds-Averaged Navier-Stokes(RANS) 기법, Large Eddy Simulation (LES), 그리고 Direct Numerical Simulation(DNS)이 있다. 이 가운데 RANS 기법은 유동 박리가 심한 영역에서는 LES나 DNS에 비해 정확도가 제한될 수 있으나, 상대적으로 낮은 계산 비용으로 순항 영역 주위의 정상 유동 특성을 예측할 수 있다는 장점이 있어[1,2] 비행체 공력 해석 및 설계에 널리 활용되고 있다.

RANS 방정식은 시간 평균 과정에서 발생하는 레이놀즈 응력 항으로 인해 방정식 계가 닫히지 않는 closure 문제가 존재하며, 이를 보완하기 위해 다양한 난류 모델이 제안되었다. 대부분의 RANS 기반 난류 모델은 벽면 근처에서의 난류 감쇠 현상과 유동 특성을 모사하기 위해 벽면 거리를 주요 입력 변수로 사용한다. 벽면 거리는 난류 모델에서 소산항의 모델링이나 특성 길이 스케일의 정의에 직접적으로 사용되므로, RANS 해석에서 난류 모델링의 정확도를 좌우하는 중요한 매개변수이며, 벽면 거리의 정확도는 CFD 해석에 직접적인 영향을 미친다[3]. 일반적으로 벽면 거리는 유동장 격자점에서 가장 가까운 벽면까지의 최소 수직 거리로 정의된다. 단순한 격자점 간의 점–점 거리 계산은 실제 벽면과의 수직 거리를 정확히 반영하지 못하므로, 신뢰성 있는 벽면 거리 계산을 위해서는 점–면 또는 점–요소 기반의 기하학적 고려가 필요하다.

한편, 유동 해석에 사용되는 격자수가 증가함에 따라 벽면 거리 계산에 요구되는 연산량 또한 급격히 증가하며, 이는 대규모 격자계를 이용한 RANS 해석에서 전체 계산 효율을 저해하는 주요 요인으로 작용한다. 이러한 문제를 완화하기 위해 벽면 거리 계산의 정확성을 유지하면서도 계산 비용을 절감하고자 하는 다양한 기법들이 제안되었으며, 최근까지도 효율적인 알고리즘에 대한 연구가 지속되고 있다. Zhang[4]은 벽면 거리 계산을 두 단계로 분리하여, 먼저 벽면 격자점과 유동장 격자점 간의 꼭지점–점 거리를 계산한 후, 미리 정의된 거리 범위 내에 위치한 점들에 대해서만 벽면까지의 최소 거리를 재계산함으로써 계산 효율과 정확도를 동시에 향상시켰다. Roget 등[5]은 구형 복셀화(sphere voxelization) 기법을 기반으로 계산 효율을 개선하고 대규모 병렬 환경에서 기존 방법 대비 우수한 성능을 보고하였으나, 복셀 크기 결정이 경험식에 의존하며 알고리즘 구조가 비교적 복잡하다는 한계를 가진다. 한편, 벽면 거리 계산을 최근접 탐색 문제가 아닌 물리적 전파 문제로 재해석하여 Eikonal 및 Poisson 방정식과 같은 미분 방정식을 이용하는 접근법도 제안되었다. Tucker 등[6]은 Eikonal, Hamilton–Jacobi, Poisson 방정식 기반 기법들의 계산 효율성과 구현 특성을 비교 분석하였으며, Eikonal 및 Hamilton–Jacobi 방법은 구현 난이도가 상대적으로 높지만 적절한 수치 이산화가 적용될 경우 비교적 정확한 거리 계산이 가능한 반면, Poisson 방법은 구현이 용이하고 계산 효율이 우수하나 벽에서 멀어진 영역이나 복잡한 형상 주변에서 거리 정확도가 저하되는 경향을 보인다고 보고하였다.

본 연구에서는 KD-트리 기반의 최근접 탐색을 두 단계로 분리한 새로운 벽면 거리 계산 기법을 제안한다. 제안된 방법은 KD-트리 자료구조를 활용하여 벽면 요소에 대한 탐색 범위를 효과적으로 축소함으로써, 대규모 격자계에서도 효율적이고 안정적인 벽면 거리 계산이 가능하도록 설계되었다. 기존의 전통적인 접근법들이 모든 유동 격자점에 대해 전체 벽면 요소를 대상으로 거리 계산을 수행함에 따라 계산 비용이 급격히 증가하는 한계를 가지는 반면, 본 기법은 탐색 대상 후보를 사전에 제한함으로써 계산 비용을 크게 감소시킬 수 있다. 제안된 알고리즘의 성능을 검증하기 위해 HB-2 형상, ONERA M6 날개, 그리고 NASA에서 제안한 CRM(Common Research Model) 형상에 대해 벽면 거리 계산을 수행하였으며, 벽면 거리 정확도 분석, 계산 시간 비교, 그리고 유동장 해석을 통해 본 기법의 유효성을 평가하였다.

2. 연구 방법

2.1 지배 방정식 및 유동 해석자

2.1.1 지배방정식 및 난류 모델

본 연구에서는 압축성 난류 유동 해석을 위해 Reynolds-Averaged Navier-Stokes(RANS) 방정식을 지배방정식으로 사용하였다. RANS 난류 모델로서 항공우주공학에서 널리 사용되는 Spalart-Allmaras(SA) 모델[7]과 Menter의 Shear Stress Transport(SST) 모델[8]을 적용하였다. 두 난류 모델은 모두 벽면 인근 유동 특성을 모사하는 과정에서 격자점으로부터 벽면까지의 거리 정보를 주요 입력 변수로 사용하며, 난류 점성의 감쇠 및 벽면 효과의 공간적 범위를 결정하는 데 중요한 역할을 한다. 본 논문에서는 난류 모델의 세부 물리적 특성보다는, SA 및 SST 모델에서 공통적으로 요구되는 벽면 거리 계산의 정확성과 효율성에 초점을 두고, 대규모 격자계에서도 적용 가능한 효율적인 벽면 거리 계산 기법을 제안한다.

2.1.2 유동 해석자

본 연구에서는 Python 기반의 압축성 난류 유동 해석 오픈소스 프로그램인 pyBaram[9]을 사용하여 전산 해석을 수행하였다. pyBaram은 비정렬 격자계에 대해 유한체적법(Finite Volume Method, FVM)을 기반으로 공간 이산화를 수행하며, 압축성 유동의 지배방정식인 Euler 방정식, Navier-Stokes 방정식, 그리고 Reynolds-Averaged Navier-Stokes (RANS) 방정식의 해석을 지원한다. 또한 충격파 주위 유동을 안정적으로 포착하면서도 2차 공간 정확도를 확보하기 위한 MLP-u 기울기 제한자[10]와, Rotated-RoeM[11]을 포함한 다양한 강건한 비점성 수치 플럭스 기법을 제공한다. 정상 유동 해석을 위한 내재적 시간 전진 기법으로는 LU-SGS 방법을 포함한 여러 시간 적분 기법을 사용할 수 있다. LU-SGS 등을 포함한 희소행렬 연산의 계산 효율 향상을 위해 행렬의 대역폭을 감소시키는 Reverse Cuthill-McKee(RCM) 재배열을 적용하며, 이는 격자 읽기 또는 격자 분할 과정에서 수행된다.

pyBaram의 소스코드는 모두 Python으로 작성되었으나, 대규모 3차원 유동 해석을 위한 계산 가속 및 병렬화를 위해 고성능 Python JIT(Just-In-Time) 컴파일러인 Numba와 MPI 기반 병렬 통신을 지원하는 mpi4py 패키지를 사용하였다. 그 결과 도메인 분할 기법을 이용한 MPI(Message Passing Interface) 기반 병렬 계산과 OpenMP 기반의 멀티 스레딩 계산을 모두 지원한다.

2.2 벽면 거리 계산 방법

RANS 기반 난류 모델을 이용한 유동 해석에서는 계산 초기 단계에서 모든 유동장 격자점에 대해 벽면까지의 거리를 산출하는 과정이 요구된다. 특히 3차원 형상의 경우, 곡면으로 구성된 비행체 표면의 기하학적 특성을 고려하여 벽면 거리를 정밀하게 계산하는 것이 중요하다. 일반적인 접근법에서는 각 유동장 격자점에 대해 모든 표면 격자 요소와의 거리를 계산한 후 최솟값을 선택하게 되며, 이로 인해 계산 격자수가 증가할수록 벽면 거리 계산에 소요되는 연산량과 계산 시간이 급격히 증가하게 된다.

현재 pyBaram에서는 모든 표면 격자 노드점과 유동장 격자 노드점 사이의 거리를 직접 탐색하는 가장 기본적인 전수 탐색(Brute-force) 방식을 사용하고 있다. 이러한 방법은 구현이 단순하다는 장점이 있으나, 계산 비용이 매우 크며, 곡률이 크거나 종횡비가 큰 격자에서는 실제 벽면과의 최소 수직 거리를 정확히 반영하기 어렵다는 한계를 가진다. 본 연구에서는 이러한 문제를 개선하기 위해 KD-트리 기반의 두 단계 벽면 거리 계산 기법을 제안하며, 다음 절에서는 제안된 방법의 절차와 특징을 상세히 설명한다.

2.2.1 노드점 기반 전수 탐색 방식

노드점 기반 전수 탐색 방식은 각 유동 격자 요소를 구성하는 노드점에서 모든 벽면 격자 노드점까지의 유클리드 거리(Euclidean distance)를 계산한 후, 최소 거리를 기반으로 격자 중심점에서의 벽면 거리를 근사적으로 산출하는 방법이다. 일반적으로 격자 요소를 구성하는 노드점에서 계산된 최소 거리의 평균값을 해당 격자 요소 중심에서의 벽면 거리로 정의한다. 이 방법은 알고리즘 구조가 단순하여 구현이 용이하다는 장점이 있으나, 계산 비용과 정확도 측면에서 한계를 가진다.

계산 비용 측면에서 볼 때, 유동장 격자점의 개수를 Nf, 벽면 격자 노드점의 개수를 Nw라 하면, 전체 연산량은 O(Nf×Nw)에 비례하여 증가한다. 따라서 격자수가 수백만 개 이상으로 증가하는 대규모 3차원 공력 해석 문제에서는 벽면 거리 계산에 소요되는 계산 비용이 비약적으로 증가하여, 대규모 격자계에서는 벽면 거리 계산만으로도 상당한 계산 시간이 요구될 수 있다. 이러한 이유로 노드점 기반 전수 탐색 방식은 소규모 문제나 알고리즘 검증 용도로는 활용 가능하나, 대규모 격자계 해석에 직접 적용하기에는 실질적인 제약이 따른다.

정확도 측면에서도 한계가 존재한다. 노드점 기반 전수 탐색 방식은 실제 벽면과의 최소 수직 거리를 직접적으로 계산하지 않고, 노드점 간의 거리 정보를 평균화하여 격자 중심점에서의 벽면 거리를 근사한다. 이로 인해 곡면 형상이나 격자 직교성이 낮은 영역, 또는 종횡비가 큰 격자에서는 실제 벽면 거리와의 오차가 증가할 수 있다.

2.2.2 KD-트리 기반 2단계 벽면 거리 계산 방법

KD-트리는 k차원 공간에 분포한 점들을 효율적으로 관리하기 위한 공간 분할 자료구조로, 최근접 탐색 과정에서 불필요한 거리 계산을 줄임으로써 탐색 비용을 효과적으로 감소시킬 수 있다. 본 연구에서는 KD-트리를 활용하여 벽면 거리 계산을 두 단계로 분리함으로써, 대규모 격자계에서도 계산 효율과 안정성을 동시에 확보할 수 있는 벽면 거리 계산 기법을 제안한다. 제안된 방법은 첫 번째 단계에서 거리 계산 대상이 되는 벽면 요소의 후보군을 효과적으로 축소한 후, 두 번째 단계에서 제한된 후보에 대해서만 기하학적으로 정확한 벽면 거리 계산을 수행하는 구조로 설계되었다. 제안된 KD-트리 기반 2단계 벽면 거리 계산 기법의 전체 절차를 Fig. 1에 개략적인 순서도로 나타내었다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F1.jpg
Fig. 1.

Flowchart of the KD-tree based two-stage wall-distance computation method

첫 번째 단계에서는 벽면 격자를 삼각형 요소로 분할하여 처리한다. 사각형 격자의 경우 네 개의 삼각형으로 분할한 뒤 각 삼각형의 중점을 계산하며, 이 벽면 삼각형 중점들을 기준으로 KD-트리를 구성한다. 삼각형 요소는 점–면 거리 계산에 있어 기하학적으로 명확한 정의가 가능하므로, 이후 단계에서의 정밀한 벽면 거리 계산을 고려하여 채택하였다. 이후 유동 격자 요소의 중심점을 질의점(query point)으로 사용하여 KD-트리 기반 최근접 탐색을 수행함으로써, 전체 벽면 요소 중 유동 격자와 기하학적으로 인접한 후보들을 선별한다. 이때 탐색 범위는 비행체 대각 길이의 약 5% 이내로 제한하였으며, 이는 RANS 난류 모델에서 벽면 인근 영역의 유동 특성이 지배적으로 나타나고 벽면 거리의 정확도가 특히 중요하다는 점을 고려한 설정이다.

두 번째 단계에서는 첫 번째 단계에서 선별된 벽면 요소에 대해 보다 정밀한 벽면 거리 계산을 수행한다. 이 단계에서는 1단계와 달리 KD-트리를 이용하여 최근접 이웃을 탐색한 후, 유동 격자 중심점과 이들 벽면 삼각형 요소 사이의 수직 거리를 계산한다. 그러나 곡률이 큰 영역이나 날개의 뒷전(trailing edge)과 같이 기하학적으로 복잡한 구간에서는 실제 최소 수직 거리를 가지는 벽면 요소가 단순 최근접 탐색에서 제외될 가능성이 존재한다. 이를 보완하기 위해 단일 최근접 요소 대신 n개의 최근접 이웃을 후보로 선정하였다. 최근접 이웃의 개수 n은 전체 표면 삼각형 수에 비례하도록 설정하였으며, 기본값은 전체 삼각형 수의 약 0.1%로 하되 최소 50개의 후보를 보장하도록 하였다. 해당 파라미터는 다양한 격자 규모에 대한 사전 테스트 및 후보 수에 대한 민감도 분석을 통해 결정되었으며, 그 결과 0.1% 수준에서도 최소 수직 거리를 충분한 정확도로 산정할 수 있음을 확인하였다. 이후 선별된 n개의 벽면 삼각형 후보에 대해서만 David Eberly의 point-to-triangle 기하학적 방법[12]을 적용하여 벽면 거리를 정밀하게 계산하였다.

이와 같은 두 단계 접근법을 통해 전수 탐색 방식에 비해 벽면 거리 계산에 요구되는 연산량을 크게 감소시킬 수 있으며, 격자 규모가 증가할수록 그 효과는 더욱 두드러진다. 또한 본 연구에서는 SciPy 에서 제공하는 KD-트리 구현을 활용하였으며, 해당 패키지들은 OpenMP 기반 멀티 스레딩을 지원하여 KD-트리 구성 및 최근접 탐색 연산을 병렬화할 수 있다[13]. 이를 통해 KD-트리의 생성 및 탐색 과정이 효율적으로 수행되며, 반복 계산이 집중되는 기하학적 거리 계산 구간에는 Numba를 적용하여 효과적으로 연산을 가속하였다.

2.2.3. 삼각형 요소에 대한 최소 벽면 거리 계산

선별된 n개의 벽면 삼각형 후보에 대해서는 삼각형 요소로 표현된 벽면과 유동 격자점 간의 최소 거리를 계산하기 위해 David Eberly의 point-to-triangle 기하학적 거리 계산 알고리즘[12]을 적용하였다. 본 알고리즘은 삼각형 평면 상의 최근접 점을 barycentric 좌표를 이용하여 표현하고, 점의 투영 위치가 삼각형 내부 또는 외부에 존재하는지를 판별함으로써 최소 거리를 기하학적으로 일관되게 계산하는 방법이다.

삼각형의 꼭짓점을 (v0,v1,v2) 유동 격자점을 P라 할 때, 삼각형 평면 위의 최근접 점 q는 barycentric 좌표 (s,t)를 이용하여 다음과 같이 표현된다.

(1)
Q=v0+s(v1-v0)+t(v2-v0)

이때 Fig. 2와 같이 s0,t0,s+t1 조건이 만족하는지 여부에 따라 점의 투영 위치가 삼각형 내부 또는 외부에 존재하는지를 판별한다. 투영 위치가 삼각형 내부에 존재하는 경우, 벽면 거리는 유동 격자점 P에서 삼각형 평면으로의 수선 거리인 d=P-Q로 계산하며, 계산되며 더 이상의 후보 삼각형에 대한 탐색을 수행하지 않고 종료한다.

반면, 투영 위치가 삼각형 외부에 존재하는 경우에는 삼각형의 변 또는 꼭짓점에 대한 점–선분 또는 점–점 거리 계산을 수행하여 최소 거리를 산출한다. 이러한 경우 추가적인 기하학적 연산이 요구되며, 특히 skewness가 크거나 직교성이 저하된 격자에서는 더 많은 후보 삼각형에 대한 평가가 필요할 수 있다. 이러한 기하학적 처리 방식을 통해 삼각형 요소에 대한 최소 벽면 거리를 안정적이고 일관되게 계산할 수 있으며, 본 연구에서 제안한 KD-트리 기반 후보 축소 기법과 결합하여 대규모 격자계에서도 효율적이면서 정밀한 벽면 거리 계산이 가능하다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F2.jpg
Fig. 2.

Schematic of the point-to-triangle distance calculation: (a) nearest-point projection outside the triangle, and (b) nearest-point projection inside the triangle

3. 해석 결과

본 연구에서는 제안한 KD-트리 기반 벽면 거리 계산 기법의 계산 효율을 기존 노드점 기반 전수 탐색(Brute-force) 방식과 비교·분석한다. 특히 격자 분할 방식에 따른 병렬화 특성을 고려하여, 두 방법에 대해 MPI 및 OpenMP 기반 병렬 계산 성능을 중심으로 평가하였다. 성능 비교를 위해 서로 다른 특성을 갖는 세 가지 비정렬 격자에 대해 전산 해석을 수행하였으며, 벽면 거리 계산에 소요되는 계산 시간과 병렬 성능을 비교하였다.

또한 계산된 벽면 거리를 이용하여 실제 유동 해석을 수행함으로써, 벽면 거리 계산 방식의 차이가 유동 해석 결과에 미치는 영향을 함께 검토하였다. 유동 해석에는 pyBaram 해석자를 사용하였으며, 공간 이산화 기법으로 MLP-u2 기울기 제한자와 Roe의 비점성 수치 플럭스 기법을 적용하였다. 난류 모델은 SST 또는 SA 모델을 사용하였다. 계산은 AMD EPYC 9554(64코어) 프로세서가 장착된 로컬 클러스터 환경에서 수행되었으며, 64코어 프로세서 1개에서 최대 8개까지를 사용하여 총 64–512 코어 규모의 계산을 수행하였다. 병렬 계산 방식으로는 코어 단위의 MPI 병렬 계산과, CPU 내부에서는 OpenMP 멀티스레딩을 적용하고 CPU 간에는 MPI 통신을 사용하는 MPI/OpenMP 하이브리드 병렬 계산을 각각 적용하여 성능을 비교하였다.

3.1 HB-2 형상

미국 AEDC의 Von Karman Gas Dynamics Facility에서 초음속 및 극초음속 시험에 사용된 HB-2(Hypervelocity Ballistic Models-2) 표준 모델[14]을 대상으로 수치 해석을 수행하였다. HB-2 형상은 노즈 콘 반각 25°를 갖는 원뿔형 노즈, 원통형 동체 그리고 Flare로 구성된 축대칭 형상이다. 해석 조건으로 마하수는 2.0, 동체 직경 (D) 기준 레이놀즈수는 1.7×106인 초음속 유동 이다. 격자는 Table 1에 정리한 바와 같이 약 110만 개의 혼합 요소로 구성된 비정렬 격자를 사용하였다. 벽면 근처 난류 유동을 정확히 해석하기 위해 벽면에서 첫 번째 격자 높이는 y+1를 만족하도록 설정하였다.

Table 1.

Grid configuration for the HB-2 model

Parameter Value
Number of nodes 547,647
Number of cells 1,096,280
Hex Elements 413,070
Prism Elements 30,113
Pyramid Elements 50,495
Tet Elements 602,602
First layer height 2.0×10-5 D
Boundary layer growth rate 1.2

Fig. 3은 벽면 거리 계산 시간을 노드점 기반 전수 탐색 방식과 제안한 KD-트리 기반 2단계 탐색 방법에 대해 비교한 결과를 나타낸다. 본 계산에서는 격자 규모가 비교적 작기 때문에, 64코어 프로세서 1개를 사용한 조건에서 MPI 기반 도메인 분할(domain decomposition) 병렬화와 OpenMP 기반 멀티스레딩(multi-threading) 병렬화를 각각 적용하여 계산 시간을 비교하였다. 제안한 KD-트리 기반 방법은 기존 전수 탐색 방식에 비해 MPI 기반 도메인 분할 병렬 환경에서 벽면 거리 계산 시간을 약 81.3% 감소시켰으며, OpenMP 기반 멀티스레딩 병렬 환경에서도 약 54.7%의 계산 시간 감소를 보였다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F3.jpg
Fig. 3.

Performance comparison of wall-distance computation methods for the HB-2 grid

Fig. 4는 제안한 방법을 이용하여 계산된 벽면 거리의 정확도를 검증하기 위해, 벽면으로부터 첫 번째 격자 요소의 중심점까지의 거리를 비교한 결과를 나타낸다. 모든 표면에 대해 계산된 벽면 거리는 격자 생성 시 부여한 첫 번째 격자 요소 높이 2×10-6 의 절반에 해당하는 1×10-6 수준으로 정확하게 계산됨을 확인할 수 있었다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F4.jpg
Fig. 4.

Distribution of the first-layer height in the HB-2 grid

이와 같이 계산된 벽면 거리를 사용하여 난류 유동 해석을 수행하였다. Fig. 5(a)는 받음각 0° 조건에서의 마하수 분포를 나타내며, 선두부 형상과 후방 플레어 각도로 인해 형성되는 충격파 구조가 적절히 포착됨을 확인할 수 있다. Fig. 5(b)는 받음각을 −1°에서 13°까지 변화시키며 계산한, 노즈로부터 1.95 D 떨어진 지점에서의 피칭 모멘트 계수를 풍동시험 결과와 비교한 것이다. 전산 해석 결과는 풍동시험 결과와 전반적으로 부합한다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F5.jpg
Fig. 5.

HB-2 flow simulation results

3.2 ONERA M6 날개

ONERA M6 날개 주위의 천음속 유동은 날개 상면에서 형성되는 3차원 λ-shock 구조와 압력 분포 특성으로 인해, 천음속 유동 해석 기법의 성능을 검증하기 위한 대표적인 벤치마크 문제로 널리 활용되어 왔다[15]. 본 연구에서는 이러한 ONERA M6 날개 형상을 대상으로 천음속 조건에서의 유동 해석을 수행하였다. 해석 조건으로는 마하수 0.84, 평균 공력 시위(mean aerodynamic chord, cmac) 기준 레이놀즈수 11.72×106, 그리고 받음각 3.06°이다.

격자는 Table 2에 정리한 바와 같이 약 370만 개의 혼합 요소로 구성된 비정렬 격자를 사용하였다. 벽면 근처 유동과 천음속 충격파 구조를 적절히 해석하기 위해 벽면에서 첫 번째 격자 높이는 y+≈1을 만족하도록 설정하였다.

Table 2.

Grid configuration for ONERA M6 wing

Parameter Value
Number of nodes 2,420,785
Number of cells 3,674,867
Hex Elements 2,057,850
Prism Elements 45,328
Pyramid Elements 252,576
Tet Elements 1,319,113
First layer height 1.5×10-6 cmac
Boundary layer growth rate 1.2

Fig. 6은 ONERA M6 날개 형상에 대해 노드점 기반 전수 탐색 방식과 제안한 KD-트리 기반 벽면 거리 계산 방법의 계산 성능을 비교한 결과를 나타낸다. 본 계산에서는 격자 규모가 상대적으로 크기 때문에, 64, 128, 256 코어 조건에서 MPI 병렬 계산과 MPI/OpenMP 하이브리드 병렬 계산 환경을 적용하여 벽면 거리 계산 시간을 비교하였다. 이때 하이브리드 병렬 계산은 64코어 조건에서는 단일 프로세서 내 OpenMP 기반 멀티스레딩을 적용하였으며, 128코어 이상 조건에서는 프로세서 간 MPI 기반 도메인 분할과 프로세서 내부 OpenMP 멀티스레딩을 결합한 방식으로 수행하였다. 벽면 거리 계산 방법으로는 KD-트리 기반 방법과 전수 탐색 방법을 각각 사용하였다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F6.jpg
Fig. 6.

Performance comparison of wall-distance computation methods for ONERA M6 wing

계산 성능 비교 결과, ONERA M6 형상에서는 모든 코어 조건에서 KD-트리 기반 방법이 전수 탐색 방법에 비해 더 짧은 계산 시간을 보였다. 특히 HB-2 형상에 비해 격자수가 약 4배 증가함에 따라, 두 방법 간 계산 시간 차이는 더욱 뚜렷하게 나타났다. 64코어 조건에서 전수 탐색 방식의 경우 벽면 거리 계산 시간이 HB-2 형상 대비 10배 이상 증가한 반면, KD-트리 기반 방법은 HB-2 형상과 ONERA M6 형상 모두에서 벽면 거리 계산이 수 초 이내에 수행되어 격자수 증가에 따른 절대 계산 시간 증가는 상대적으로 제한적인 것으로 나타났다. 또한 코어 수를 증가시킬수록 전수 탐색 방식에서는 계산 시간이 점진적으로 감소하는 경향이 관찰되었으며, KD-트리 기반 방법에서도 유사한 감소 경향이 나타났으나 그 감소 폭은 비교적 크지 않았다.

Fig. 7은 ONERA M6 날개 격자에서 계산된 벽면 거리의 정확도를 검증하기 위해, 벽면으로부터 첫 번째 격자 요소 중심점까지의 거리를 비교한 결과를 나타낸다. 격자 생성 시 부여한 첫 번째 격자 요소 높이는 1.5×10-6 이며, 계산된 벽면 거리는 이 값의 절반에 해당하는 수준으로 도출되었다. 이는 벽면 거리가 벽면과 격자 중심점 간의 거리로 정의된 본 연구의 계산 방식과 격자 생성 조건에 부합하는 결과이다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F7.jpg
Fig. 7.

Distribution of the first-layer height in the ONERA M6 wing

Fig. 8은 ONERA M6 날개에 대해 계산된 압력 계수 분포를 풍동 실험 결과와 비교한 것이다. 본 연구에서는 AGARD AR-138에 수록된 ONERA M6 날개의 풍동 실험 데이터를 기준으로 전산 해석 결과를 검증하였다. 압력 계수 분포는 세미스팬 길이 b로 정규화한 위치 y/b=0.650.95에서의 단면을 기준으로 비교하였으며, 여기서 y는 날개 루트로부터의 거리이다. 해당 단면 위치는 Fig. 7에 실선으로 표시하였다. 전산 해석 결과는 두 단면 위치 모두에서 실험 데이터와 전반적으로 유사한 분포 경향을 보였다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F8.jpg
Fig. 8.

Pressure coefficient distribution on the ONERA M6 wing

3.3 NASA CRM(Common Research Model)

NASA CRM(Common Research Model) 형상은 AIAA Drag Prediction Workshop(DPW)에서 천음속 조건의 공력 성능 예측을 검증하기 위해 제안된 대표적인 벤치마크 형상으로, 다양한 수치 해석 기법의 성능 비교에 널리 활용되고 있다. 본 연구에서는 제7회 Drag Prediction Workshop(DPW-7)[16]에서 사용된 JAXA의 비정렬 3차원 격자를 대상으로 NASA CRM 형상 주위의 유동 해석을 수행하였다. 해당 격자는 MEGG3D를 이용하여 생성되었으며, 실제 항공기 형상의 기하학적 특징을 충실히 반영한 비정렬 격자 구조를 갖는다. 본 연구에서는 JAXA에서 제작한 두 종류의 격자를 사용하였으며, 이들의 상세한 구성은 Table 3에 정리하였다. 사용된 격자는 Tiny 및 Coarse 격자로 구분되며, 각각 약 2천5백만 개와 7천6백만 개의 혼합 격자 요소로 구성되어 있다.

Table 3.

Grid configurations for NASA CRM

Parameter Tiny Coarse
Number of nodes 8,698,930 26,891,512
Number of cells 25,294,690 76,058,884
Hex Elements 2,739,352 9,084,908
Prism Elements 6,105,626 18,760,276
Pyramid Elements 119,327 275,934
Tet Elements 16,330,385 47,937,766
First layer height 8.455×10-7 cmac 5.638×10-7 cmac
Boundary layer growth rate 1.323 1.205

해석 조건으로는 마하수 0.85, 평균 공력 시위 기준 레이놀즈수 2.0×107이다. DPW-7에서는 받음각 2.5~3도 사이에서 양력 계수가 CL=0.58 조건을 탐색하였으며, 본 연구에서는 해당 범위의 대표 조건으로서 받음각 2.75°를 선택하여 해석을 수행하였다. NASA CRM 형상은 앞선 두 형상에 비해 날개–동체 결합부가 존재하는 등 기하학적 복잡성이 높은 경우에 해당하며, 이러한 조건에서 벽면 거리 계산의 정확성과 계산 효율을 함께 확인하였다.

Fig. 9는 NASA CRM 형상에 대해 Tiny 격자와 Coarse 격자를 대상으로, 벽면 거리 계산 방법에 따른 계산 성능을 비교한 결과를 나타낸다. 두 격자 모두 64–256 코어(Tiny) 및 512 코어(Coarse) 조건에서 MPI 병렬 계산과 MPI/OpenMP 하이브리드 병렬 계산 환경을 적용하여 성능을 평가하였다.

Tiny 격자에 대한 계산 결과, MPI 병렬 환경에서는 64 및 128 코어 조건에서 KD-트리 기반 방법이 전수 탐색 방법 대비 각각 약 45.7% 및 39.5%의 계산 시간 감소를 보였으나, 256 코어 조건에서는 전수 탐색 방법이 약 15.4% 더 짧은 계산 시간을 나타내는 경향이 관찰되었다. 이는 MPI 도메인 분할 과정에서 일부 rank가 벽면에서 멀리 떨어진 격자들로만 구성된 영역을 담당하게 되면서, 해당 영역에 대한 1단계 KD-tree 구성 시간이 증가하고 계산 부하가 불균형하게 분포되기 때문이다. 이로 인해 병렬 효율이 저하되며, 결과적으로 전체 계산 성능 저하로 이어졌다.. 반면, MPI/OpenMP 하이브리드 병렬 환경에서는 모든 코어 조건에서 KD-트리 기반 방법이 전수 탐색 방법 대비 일관되게 우수한 성능을 보였으며, 계산 시간 감소율은 약 90% 내외로 나타났다. 특히 하이브리드 환경에서 KD-트리 기반 방법의 계산 시간은 64 코어 조건에서도 50초 이하로 유지되어, MPI 환경에서 수백 초가 소요된 경우와 비교할 때 현저한 차이를 보였다.

CRM Coarse 격자에 대한 512 코어 계산 결과에서도 유사한 경향이 확인되었다. MPI 병렬 환경에서 KD-트리 기반 방법은 전수 탐색 방법 대비 약 18.9%의 계산 시간 감소를 보였으며, MPI/OpenMP 하이브리드 병렬 환경에서는 계산 시간이 약 30.3% 감소하였다. 격자 규모가 크게 증가한 조건에서도 KD-트리 기반 방법이 전수 탐색 방법 대비 일정 수준의 계산 효율 개선 효과를 유지함을 확인할 수 있다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F9.jpg
Fig. 9.

Performance comparison of wall-distance computation methods for NASA CRM

종합적으로, NASA CRM 형상과 같이 대규모 비정렬 격자를 포함하는 문제에서는 MPI/OpenMP 하이브리드 병렬 환경에서 수행된 KD-트리 기반 벽면 거리 계산 방법이 가장 안정적이고 효율적인 성능을 보였으며, 격자 규모가 증가하더라도 전수 탐색 방법 대비 계산 시간을 유의미하게 감소시키는 효과를 나타냈다.

Fig. 10은 NASA CRM 형상에 대해 계산된 벽면 거리 분포(contour)를 나타낸다. NASA CRM 형상은 앞선 두 형상에 비해 날개–동체 결합부를 포함하여 표면 굴곡이 상대적으로 복잡한 형상에 해당한다. 그럼에도 불구하고 제안된 벽면 거리 계산 방법을 적용한 결과, 대부분의 표면 격자에서 계산된 벽면 거리는 격자 생성 시 부여한 첫 번째 격자 레이어 높이의 절반에 해당하는 4.23×10⁻⁷ 수준으로 산출되었음을 확인할 수 있다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F10.jpg
Fig. 10.

Distribution of the first-layer height in the NASA CRM

Coarse 격자에 대해 유동 해석한 결과 받음각 2.75°, 마하수 0.85, 평균 공력 시위 기준 레이놀즈수 2.0×10⁷ 조건에서 수행한 유동 해석 결과로, 날개 상부 표면에서의 압력 계수 분포를 나타낸다. 해당 조건에서 계산된 양력 계수 CL 은 약 0.5918로, 받음각 약 2.7°에서 CL=0.58을 만족하는 것으로 DPW-7에서 타 연구자의 결과와 유사하다. 또한 Fig. 11에 나타난 압력 계수 분포로부터 CRM 날개 상부에서 형성되는 충격파 구조가 적절히 포착되었음을 확인할 수 있다.

https://cdn.apub.kr/journalsite/sites/kscfe/2026-031-01/N0500310105/images/jkscfe_2026_311_065_F11.jpg
Fig. 11.

Pressure coefficient distribution on the upper surface of the NASA CRM wing

4. 결 론

본 연구에서는 RANS 기반 압축성 난류 유동 해석에서 필수적으로 요구되는 벽면 거리 계산의 계산 비용을 절감하기 위해 KD-트리 기반의 2단계 벽면 거리 계산 기법을 제안하고, 그 성능과 정확도를 체계적으로 검증하였다. 제안된 방법은 KD-트리 기반 최근접 탐색을 이용해 계산 대상을 효과적으로 축소한 뒤, 선별된 후보에 대해서만 기하학적으로 정확한 거리 계산을 수행함으로써 계산 효율과 벽면 거리 정확도를 동시에 확보하도록 설계되었다. HB-2, ONERA M6, NASA CRM 형상에 대한 수치 실험 결과, 제안된 기법은 전수 탐색 방식 대비 벽면 거리 계산 시간을 유의하게 감소시켰으며, 특히 MPI/OpenMP 하이브리드 병렬 환경에서 격자 규모가 증가하더라도 안정적인 성능 향상을 유지함을 확인하였다.

또한 CRM 형상과 같이 표면 굴곡이 큰 복잡한 형상에서도 격자 중심점으로부터 벽면까지의 거리가 첫 번째 격자 요소 높이의 절반 수준으로 일관되게 계산되어, 제안된 기법이 기하학적 신뢰성을 유지함을 확인할 수 있다. 이러한 결과로부터, KD-트리 기반 2단계 벽면 거리 계산 기법은 대규모 비정렬 격자계를 사용하는 RANS 유동 해석에서 실용적이고 효율적인 대안이 될 수 있음을 보여준다. 향후 연구에서는 MPI 환경에서의 확장성 특성에 대한 추가 분석과 탐색 파라미터 설정에 대한 체계적인 검토를 통해 적용 범위를 확장할 예정이다.

Acknowledgements

본 연구는 (주)LIG넥스원(75776-1)의 지원으로 수행되었으며, 지원에 감사드립니다.

References

1

2000, Spalart, P., “Trends in Turbulence Treatments,” AIAA Paper, No. 2000-2306.

10.2514/6.2000-2306
2

1996, Breuer M., Lakehal D. and Rodi W., “Flow around a Surface-Mounted Cubical Obstacle: Comparison of LES and RANS results,” Comput. Three-Dimens. Complex Flows, pp.22-30.

10.1007/978-3-322-89838-8_4
3

2002, Fares E. and Schroder W., “A Differential Equation to Determine the Wall Distance,” Int. J. for Numerical Methods in Fluids, Vol.39, pp.743-762.

10.1002/fld.348
4

2024, Zhang X., Diskin B. and Nielsen E. J., ”A Wall-Distance Method for Turbulence Modeling,” Proc. 12th Int. Conf. Comput. Fluid Dyn. (ICCFD12).

5

2013, Roget B., “Wall Distance Search Algorithm using Voxelized Marching Spheres,” J. Comput. Phys., pp.76-94.

10.1016/j.jcp.2013.01.035
6

2005, Tucker, P.G., “Computations of Wall Distances Based on Differential Equations,” AIAA J., Vol.43, No.3, pp.539-549.

10.2514/1.8626
7

1992, Spalart, P.R. and Allmaras, S.R., “ A One-Equation Tubulence Model for Aerodynamic Flows,” AIAA Paper 92-0439.

8

1994, Menter, F.R.,, “Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications,” AIAA J., Vol.32, No. 8, pp.1598-1605.

10.2514/3.12149
9

2022, Park, J.S., “pyBaram: Parallel compressible flow solver in high-performance Python for teaching and research,” Software X, Vol.20, 101272.

10.1016/j.softx.2022.101272
10

2012, Park, J.S. and Kim, C., “Multi-dimensional Limiting Process for Finite Volume Method on Unstructured Grid,” Comp. and fluid., Vol.65, pp.8-24.

10.1016/j.compfluid.2012.04.015
11

2024, Choi, S., Kim, D., Park, J. and Park, J.S., “Robust and accurate Roe-type Riemann solver with compact stencil: Rotated RoeM scheme,” J. Comput. Phys., Vol.505, No.15.

10.1016/j.jcp.2024.112913
12

2003, Eberly, D., Geometric Tools for Computer Graphics, Morgan Kaufmann, pp.374-382.

13

2020, Virtanen, P., Gommers, R., Oliphant, T.E., Haberland, M., Reddy, T., Cournapeau, D., et al., “SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python,” Nature Methods, Vol.17, No.3, pp.261-272.

14

1963,Gray, J.D. and Lindsay, E.E., “Force tests of standard hypervelocity ballistic models HB-1 and HB-2 at Mach 1.5 to 10,” AEDC-TDR-63-137.

15

1979, Schmitt, V. and Charpin, F., “Pressure distributions on the ONERA M6 wing at transonic Mach numbers,” AGARD Rep., No. AR-138.

16

2023, Hue, D., Sartor, F., Petropoulos, I. and Fournis, C.,, “DPW-7: Steady and Unsteady Computations of the Common Research Model at Different Reynolds Numbers,” AIAA J., Vol.60, No.6, pp.1857-1871.

10.2514/1.C037231
페이지 상단으로 이동하기