1. 서 론

로켓노즐은 고온의 연소가스를 초음속으로 팽창시켜 추력을 발생시키며 이 과정에서 고온/고속의 연소가스 유동에 노출된다. 따라서 고온에서도 형상을 유지하며 노즐의 중요부품이 보호 될 수 있도록 노즐의 각 부위에 내열재를 적절히 적용하는 것이 필요하다. 이러한 노즐의 내열설계를 위해서는 노즐 내의 초음속 유동이 표면에 발생시키는 열전달 양을 예측해야 한다. 그러나 노즐 유동은 속도가 빠르고 경계층이 매우 얇아 유동해석이 어려우며 이에 노즐 내의 열전달을 예측하기 위한 많은 연구가 진행되어 왔다.

노즐의 열전달을 예측하는 방법 중 현재 노즐의 설계에 가장 널리 사용되고 있는 것은 경계층 적분법을 이용한 해석방법이다. 이는 Prandtl의 운동량 적분에 기반을 두고 있으며(Anderson, 2005) 경계층 내의 속도분포를 가정하고 적분함으로써 경계층 내의 보존방정식을 만족시키는 경계층 두께와 벽면에 발생하는 열전달을 역으로 도출하는 방법이다. 이를 초음속 노즐에 적용 한 연구 중 가장 대표적인 것은 Bartz(1965)의 연구로 이 연구에서는 Dittus 등(1930)의 관내유동의 열전달 해석을 위한 수식을 변형하여 면적이 변하는 초음속 노즐 내부의 열전달을 도출하고자 하였다. 또한 Kwong 등(1978)은 경계층 적분법의 마찰계수(C_f ) 및 스탄튼수(St)를 도출함에 있어 역압력구배, 경계층 내 물성변화 및 표면조도의 영향까지 고려한 수식을 사용하고 이를 로켓노즐 실험결과와 비교하였다.

그러나 경계층 적분법은 주어진 경계층 외부유동에 따라 경계층 내부만을 계산하므로 경계층 외부유동에 대해 별도의 모델이 필요하다는 한계가 있다. 일반적으로 초음속 노즐 유동의 코어 유동은 1차원 등엔트로피 모델로 모사할 수 있으며 경계층 적분법을 사용하는 경우 대부분 이 모델을 통해 코어 유동을 계산한다. 이 경우 코어 유동의 압력, 온도 및 속도 등은 노즐 내의 위치에 따라 달라지는 노즐 목 대비 면적비에 따라 계산된다. 그러나 열전달양이 높으며 정확도가 가장 중요한 노즐 목의 경우 곡률이 존재하여 2차원적인 유동의 회전이 발생하는 곳으로 경계층 외부유동의 거동이 1차원 등엔트로피 모델에서 벗어나게 된다. 또한 노즐 벽면에서 경계층이 발달하면 코어 유동이 흐르는 유효 단면적이 좁아지며 이는 코어 유동에 영향을 주게 된다. 두 가지 효과는 경계층 외부유동의 조건이 주어져야하는 경계층 적분법 계산의 정확도를 제한하는 요인 이며 외부유동 계산에서 이를 고려할 경우 경계층 적분법 계산의 정확도를 향상시킬 수 있다.

따라서 본 연구에서는 이와 같은 경계층 적분법의 한계를 극복하고 정확도를 향상시키기 위해 경계층 외부유동의 도출에 Computational Fluid Dynamics(CFD)를 적용하고 이를 경계층 적분법과 통합하였다. 오일러 모델과 SST k-ω 난류 모델을 적용한 두 가지 경우를 CFD로 해석하고 이를 경계층 적분 계산에 적용함으로서 노즐 목에서의 2차원 유동 회전과 경계층 성장에 따른 유효단면적 변화를 고려하고자 하였다. 그리고 선행문헌의 실험결과에 대해 기존의 1차원 등엔트로피 모델을 사용한 해석결과와 비교하여 해석기법의 타당성을 검증하였다.

2. 본 론

2.1 경계층 외부유동의 해석

본 해석의 경계층 외부유동 및 경계층 적분법의 해석에 사용되는 좌표계 및 주요 변수는 Fig. 1에 나타난 것과 같다. 먼저 노즐 내부의 특정 위치의 단면적을 A , 노즐 유동의 비열비를 γ라고 하면 1차원 등엔트로피 가정 하에서 노즐의 경계층 외부유동의 마하수 M_e 은 다음과 같이 표현할 수 있다.

Fig. 1

Coordinate system of nozzle core flow and boundary layer in this research

(1)

(AAth)2=1Me2[2γ+1(1+γ−12Me2)]γ+1γ−1

수식에서 하첨자 th는 노즐 목을, 하첨자 e 는 경계층 외부를 나타낸다. 그리고 이때의 경계층 외부 온도 T_e , 밀도 ρ_e 및 압력 p_e 는 아래와 같은 수식으로 표현된다.

(2a)

T0eTe=1+γ−12Me2

(2b)

ρ0eρe=(1+γ−12Me2)1γ−1

(2c)

p0epe=(1+γ−12Me2)γγ−1

여기서, 하첨자 0은 유동의 정체점에서의 물성을 나타낸다. 노즐 특정 위치에서의 반지름을 통해 면적비를 계산하면 식 (1)에서 마하수를 역산할 수 있다. 그리고 이를 식 (2)에 대입 하여 경계층 외부유동에 대해 경계층 적분법에서 필요한 변수 들을 모두 도출할 수 있게 된다. 실제로 식 (1)의 해를 도출할 때는 면적비와 마하수의 관계가 내재적으로 표현되어있으므로 해를 구하기가 어렵다. 따라서 본 연구에서는 Newton-Rapson 법을 사용하여 주어진 면적비 A/A_th 에 해당하는 M_e 의 해를 도출하였다(Todd, 2016).

CFD를 이용하여 경계층 외부유동을 해석하는 경우는 경계층 외부유동의 도출에 두 가지 방법을 사용하였다. 먼저 오일러 모델을 통해 경계층 외부조건을 도출하는 경우는 축대칭의 비점성 Navier-Stokes 방정식을 해석하며 점성의 영향이 없으므로 벽면에서 no-slip 조건이 적용되지 않는다. 따라서 CFD 해석을 통해 도출된 벽면에서의 마하수 및 온도, 압력, 밀도를 바로 경계층 적분법에 적용하였다.

그러나 경계층 발달에 따른 경계층 외부유동의 변화를 고려하기 위해 SST k-ω 난류모델을 적용한 경우 벽면에서 no-slip 조건이 만족되므로 벽면에서의 속도가 0이 된다. 따라서 이때는 CFD에서 도출된 벽면에서의 값을 경계층 적분에 직접 적용하는 것이 불가능하다. 본 연구에서는 경계층 내에서 높이방향으로는 압력구배가 거의 존재하지 않는 점을 이용하여 경계층 외부의 마하수를 CFD에서 도출된 벽면압력 p_w를 통해 아래와 같이 역산하였다.

(3)

Me=2γ−1{(p0pw)γ−1/γ−1}

2.3 검증을 위한 초음속 노즐해석

본 연구에서 제시한 CFD와 경계층 적분법을 통합한 해석 기법을 검증하기 위해 Back 등(1964)의 실험결과를 이용 하였다. 실험에 사용된 노즐은 30도의 수축각과 15도의 팽창각도를 가지고 있으며 출구 마하수는 2.5이다. 노즐의 형상은 Fig. 2에 나타나있으며 노즐 목의 반지름은 22.9mm 이다. 본 연구에서는 논문에서 수행된 여러 실험케이스 중 TEST 262의 조건을 사용하였으며 이 조건에서는 노즐이 457.2mm의 발달구간을 가지며 노즐의 전압력은 518.5 kPa, 전온도는 843.3K이다. 실험에서는 전온도를 상승시키기 위해 압축공기에 에탄올의 연소가스가 일부 포함되어 있으나 그 영향이 크지 않은 것으로 Back 등(1964)이 보고하였 으므로 본 해석에는 공기의 물성을 사용하였다. 특히 해석에 중요한 비열비 γ는 1.35로 가정하였다. 또한 코드 검증의 측면에서 Nichols 등(2004)의 접근방법과 동일하게 벽면의 온도는 전온도의 1/2인 421.7K으로 고정하고 열유속을 도출 하였다.

Fig. 2

Geometry of supersonic nozzle used in this research

경계층 외부유동의 도출을 위해 Fig. 3과 같이 축대칭으로 작성된 정렬계산격자를 사용하였다. 격자는 각각 길이방향과 높이방향으로 300×100개가 사용되었다. 본 연구에서는 CFD 해석을 통해 직접 경계층을 해석하여 벽면 열전달을 도출하지 않고 CFD 결과를 경계층 적분법에 대입하여 벽면열전달을 도출하므로 격자를 균일하게 배치하는데 초점을 두어 벽면에 많은 격자를 사용하지는 않았다. 해석에는 ANSYS Fluent v.15를 사용하였으며 압축성 유동의 정확도 향상을 위해 density-based solver가 사용되었다. 또한 경계층 적분법은 ANSYS Fluent의 User Defined Function 기능을 이용하여 통합되었다.

Fig. 3

Computational grid for CFD to deduce boundary layer edge values

3. 해석결과 및 고찰

Fig. 4는 각각 1차원 등엔트로피 모델과 CFD를 통해 오일러 모델 및 SST k-ω모델을 해석하여 도출한 노즐표면의 압력을 노즐의 입구 전압력으로 나누어 무차원화 한 것으로 노즐의 후류로 이동함에 따라 유동이 가속되고 팽창하여 압력이 감소하는 것을 확인할 수 있다. Back 등(1964)의 실험에서 측정된 결과와 비교해 보면 계산방법과 관계없이 유동회전이 일어나지 않는 수축부/팽창부의 직선구간에서는 벽면압력의 예측값이 실험과 일치하는 것을 확인할 수 있다. 그러나 유동의 2차원 회전이 일어나는 노즐 목 근처 구간에서는 1차원 등엔 트로피 모델을 이용한 벽면압력은 실험값보다 최대 32.4% 과다 예측되었다. 오일러 모델과 SST k-ω 모델의 경우 유동의 2차원 해석이 가능하므로 노즐 목에서도 실험값과 일치하는 결과를 나타냈다. 오일러 모델과 SST k-ω모델의 벽면압력의 차이는 2% 내외로 SST k-ω로 도출한 압력이 노즐 목 이후로 지속적으로 낮게 나타났다.

Fig. 4

Comparison of pressure distribution along the nozzle wall between experiment of Back et al.(1964) and current calculation

Fig. 5는 도출된 벽면압력으로 역산한 노즐 경계층 외부 유동의 마하수를 나타낸 것이다. 먼저 1차원 등엔트로피 모델은 노즐의 면적비로 마하수를 도출하므로 노즐 목에서의 마하수가 1이 되는 것을 알 수 있다. 그러나 2차원 효과를 고려할 경우 노즐의 유효 노즐 목(effective throat)은 벽면에서의 유동 회전으로 인해 노즐 목 전방에서 발생한다. 이로 인해 유효 노즐 목은 실제 노즐 목보다 작으며 노즐 목에서의 마하수는 1보다 높아진다. 따라서 면적비가 커져 노즐 목의 차이가 상쇄되기 이전에는 CFD를 통해 도출된 마하수가 1차원 등엔트로피 모델에 비해 높게 예측되었다. 또한 SST k-ω모델을 적용한 경우가 오일러 모델에 비해 마하수가 높게 예측되는데 이는 SST k-ω 모델은 CFD에서 경계층 성장에 따른 주유동의 변화까지도 함께 해석하기 때문이다. Fig. 5에 나타난 마하수의 차이는 크지 않지만 이러한 마하수 차이로 인해 경계층 적분법에 경계조건으로 사용되는 압력, 온도 및 밀도가 모두 낮아지게 되므로 경계층 적분 계산결과에는 더 큰 차이를 나타낼 수 있다.

Fig. 5

Computed Mach number distribution along the nozzle

Fig. 6은 각각의 경계층 외부유동 해석방법을 적용하여 경계층 적분을 수행한 결과로 도출된 energy thickness θ_T 이다. 모든 해석에서 θ_T 는 노즐이 수축함에 따라 감소하여 노즐 목 부근에서 최소가 되었다가 노즐이 팽창하며 다시 증가 하는 경향을 나타내고 있다. 그러나 각각의 외부유동 해석방법에 따라 값의 차이가 나타났으며 노즐의 수축부에서 차이가 가장 크게 나타났다. 이는 노즐 입구부분도 곡률을 가지고 있으며 이로 인해 유동의 2차원 효과가 노즐 입구부분에서도 발생할 수 있기 때문으로 판단된다. 입구에서의 오일러 모델과 SST k-ω모델의 차이는 SST k-ω 모델의 경우 노즐 앞부분에 존재하는 발달구간에서 경계층의 발달로 인한 유동의 유효면적 변화가 고려되므로 오일러 모델에서 도출되는 것 보다 마하수가 높을 것이며 이로 인해 경계층이 얇게 예측되어 발생한 것으로 보인다. 또한 앞서 언급한바와 같이 2차원 효과로 인한 유효 노즐 목의 생성이 예측되므로 CFD로 계산된 모델에서는 θ_T 가 최소가 되는 점이 실제 노즐 목보다 전방에 발생하였다. 팽창부 에서는 1차원 등엔트로피 모델로는 예측할 수 없는 노즐 변곡점 에서의 팽창률 변화에 따른 경계층 두께의 변화가 예측되었으며 유동의 마하수가 SST k-ω 모델에서 더 높게 예측되므로 초음속 구간에서 경계층의 두께가 더 빨리 증가하는 것을 알 수 있다.

Fig. 6

Calculated energy thickness θ_T distribution on the nozzle surface

Fig. 7은 본 연구에서 최종적으로 도출한 노즐 벽면의 열유속을 나타낸 것이다. 먼저 Back 등(1964)의 실험에서 도출된 열유속을 살펴보면 앞서 설명한 바와 같이 2차원 효과로 인해 실제 노즐 목 전방에서 생기는 유효 노즐 목에서 열유속이 최대가 되는 형태로 실험결과가 도출되었음을 알 수 있다. 이런 현상은 유동이 흐르는 단면적이 좁아짐에 따라 노즐 목에서 단위 면적을 지나는 질량유속이 최대가 되며 경계층도 가장 얇아지기 때문이다. 이는 식 (17)로 도출되는 열유속에서 질량유속 ρ_eu_e 와 θ_T 의 함수로 정의되는 St의 관계에서 유추할 수 있다.

Fig. 7

Comparison of surface heat flux distribution along the nozzle wall between experiment of Back et al.(1964) and current calculation

반면 1차원 등엔트로피 모델을 적용한 해석결과를 살펴보면 팽창부에서는 실험결과와 유사한 값이 도출되었으나 실제 노즐 목에서 열유속이 최대가 되는 것으로 나타났다. 이는 2차원 효과가 고려되지 못했기 때문이다. 또한 이렇게 노즐 목 계산에서의 오차가 발생하면 노즐 내 유동의 유효단면적 계산 에서의 오류로 누적되어 경계층 외부 마하수와 같은 경계 조건의 해석에 오차가 발생한다. 이로 인하여 노즐의 팽창부 에서의 열전달의 오차가 커지게 되었으며 특히 실제 노즐 목과 그 후방의 2차원 유동 회전구간에서 마하수가 낮게 예측됨에 따라 열유속이 과다 예측되는 결과를 가져왔다.

반면 오일러 모델과 SST k-ω 모델을 적용한 경우 2차원 유동회전과 이로 인한 유효 노즐 목의 형성을 예측할 수 있으므로 팽창부에서의 오차가 감소하는 효과가 나타났다. 두 해석에서 모두 열전달이 최대가 되는 지점은 실제 노즐 목보다 전방에 나타났으며 최대점을 지난 후 더 빨리 감소하였다. SST k-ω 모델로 도출된 열유속은 노즐의 수축부와 목에서 오일러 모델로 계산한 열유속 결과보다 미세하게 높았으나 노즐 팽창부의 변곡점을 지난 이후에는 더 낮게 나타났다. 또한 θ_T 와 마찬가지로 열유속도 노즐 목에서 거리가 먼 하류에서는 다시 1차원 등엔트로피 모델을 사용한 결과에 수렴하는 것으로 나타났다. 그리고 θ_T 에서도 나타났던 노즐 변곡점에서 열전달도 변화하는 것이 예측되었다. 팽창부에서 실험과 해석의 최대오차는 노즐 팽창부의 변곡점 근처에서 관찰되었으며 1차원 등엔트로피 모델은 36.1%, 오일러 모델은 25.1%, SST k-ω 모델은 20.1%로 나타났다.

4. 결 론

본 연구에서는 경계층 적분법의 계산 정확도를 향상시키기 위해 적분 계산의 경계조건으로 사용되는 경계층 외부유동 조 건의 정확성을 향상시키고자 하였다. 노즐 목에서 일어나는 2차원 효과와 경계층이 노즐 코어유동에 미치는 영향을 반영 하기 위해 CFD를 사용하여 경계층 외부유동을 해석하였으며 오일러 모델과 SST k-ω 모델을 적용하였다. 그리고 CFD로 도출된 경계층 외부유동을 경계층 적분법에 적용할 수 있도록 UDF로 연동하는 계산기법을 개발하였다. 개발된 해석기법을 검증하기 위해 Back 등(1964)의 실험에서 사용된 노즐에 대해 노즐 벽면의 열전달을 해석하였다.

해석결과 기존에 사용되던 1차원 등엔트로피 모델에 비해 노즐 벽면의 압력 및 열유속 예측에 있어서 정확도가 향상되는 것을 확인할 수 있었다. 특히 2차원 효과가 크게 발생하는 노즐 목 주위에서 정확도가 향상되었으며 열전달 최대점이 실제 노즐 목보다 전방에 발생하도록 하는 유효 노즐 목의 위치도 예측이 가능한 것이 확인되었다. 결과적으로 기존의 1차원 등엔트로피 모델과 비교하여 열유속 예측의 오차를 최대 16% 감소시킬 수 있는 것으로 나타났다.

본 연구에서 개발된 해석기법은 노즐 목에서의 정확도를 향상시키면서도 경계층 적분법을 이용하여 계산용량을 줄인 것이 특징이다. 따라서 추후 로켓노즐의 내열설계 및 내열재 설계에 유용하게 사용될 것으로 판단된다.