고체 로켓 추진기관에서 실리카/페놀릭 열반응 해석 연구

1. 서론

고체 로켓 추진기관에서 노즐은 연소가스를 운동에너지로 효과적으로 전환시켜주는 역할을 한다[1]. 추진기관의 연소 중 노즐의 역할을 유지하기 위해서는 공력형상을 최대한 유지해주어야 한다. 이뿐만 아니라 구조물의 온도에 대해 일정수준 이하로 제한해야 한다. 따라서 노즐 열설계자는 적절한 형상 및 재료를 선택해야 한다[2].

노즐에 사용되는 내열재료는 위치 별 특성을 고려하여 선택된다. 노즐목의 경우는 탄소/탄소 복합재료 또는 그라파이트를 사용한다. 그 외 확대부와 축소부는 목적에 따라 일반적으로 탄소 섬유에 페놀릭 수지를 함침시킨 탄소/페놀릭 또는 실리카 섬유에 페놀릭 수지를 함침시킨 실리카/페놀릭 재료를 사용한다. 제한된 내면형상, 무게 등을 만족시키는 노즐 열설계를 위해서는 내열재료에 대한 적절한 열반응(내부 열반응, 표면반응) 예측이 필요하다. 이러한 이유로 내열재료에 대한 수치해석 연구가 진행되고 있다[3,4]. 또한 재료에 대한 열물성 연구도 수행되었다[5].

본 연구에서는 고체추진기관의 노즐 내열재료에 사용되는 실리카/페놀릭 복합재료의 열반응 예측을 위한 수치해석을 수행하였다. 실리카/페놀릭이 고온, 고압의 연소가스에 노출되면 재료의 온도상승에 따른 내부 열반응이 발생한다. 바인더(페놀릭)의 열분해와 더불어 열분해 완료 후 생성된 숯과 실리카 섬유의 화학반응이 발생하며, 섬유가 용융되어도 계속 진행된다. 내부 열반응 뿐만 아니라 표면삭마(표면열반응, 표면후퇴)가 발생한다. 연소가스의 산화물질과 숯의 산화반응, 용융층에서의 실리카-숯 화학반응 및 부피 수축, 용융층의 증발 및 전단이 발생한다. 본 연구에서는 내부열반응 및 표면삭마를 고려한 열전도 해석을 수행하였으며, 이때 사용된 경계조건은 가장 지배적인 대류열전달을 고려하였다. 코어유동은 1차원 등엔트로피 가정을 통한 이론해를 적용하고, 벽면에서의 대류열전달계수를 위하여 경계층 적분방정식을 해석하였다[6,7]. 복잡한 물리적 현상을 수식화하기 위해서는 적절한 가정이 필요하며, 본 논문에서의 내부 열반응 모델링에 적용된 부분은 Bae[8]의 문헌을 인용하였다.

2. 실리카/페놀릭 내부열반응 모델

2.1 실리카/페놀릭 내부 열반응 과정

실리카/페놀릭 재료는 실리카 섬유와 페놀릭 바인더로 구성된다. 이 재료는 실리카 섬유의 질량분율에 따라 다른 열분해 거동을 나타낸다. 실리카의 질량분율이 낮은 경우에는 고온의 연소가스에 노출 시 페놀릭 바인더의 열분해가 발생하고, 표면 산화반응에 의한 삭마가 발생한다. 이는 탄소/페놀릭 복합재료와 비슷한 거동이다. 하지만 실리카 섬유의 질량분율이 큰 경우(페놀릭 바인더 질량분율이 작은 경우)에는 얇은 실리카 용융층이 숯층을 덮고 연소가스와의 산화반응을 방지하는 역할을 한다. 후자의 경우가 전형적인 고체추진기관에서 사용되는 실리카/페놀릭 복합재료이다. Fig. 1은 실리카/페놀릭의 열반응 과정을 나타낸다. 열반응에 의한 재료층 변화는 처녀층, 열분해층, 숯층, 용융층으로 구분할 수 있다. 우선 재료가 고온, 고압의 연소가스에 노출되면 재료 중에 함유된 수분이 증발한다. 이후 약 600 K부터 페놀릭 바인더의 열분해가 시작되어 열분해층을 형성하며, 열분해 가스가 분출되기 시작한다. 약 1,200 K에서 열분해가 완료되어 숯층을 형성하게 되며, 이후 약 1,300 K부터 탄화된 숯(탄소)과 실리카 섬유의 화학반응이 발생한다. 계속된 온도 상승에 의해 섬유가 연화되고, 약 2,000 K에서 섬유가 완전 용융되어 용융층을 형성한다. 이때 숯과 실리카 섬유의 화학반응은 계속 진행되고, 용융된 섬유는 숯의 공극을 채운된다. 재료의 표면에서는 용융된 섬유의 전단 및 화학반응으로 표면 열반응이 발생하고, 약 2,400 K부터는 용융된 섬유가 증발한다.

Fig. 1

Thermal response of silica/phenolic composite material.

2.4. 용융 실리카-숯 화학 반응 모델

약 2,000 K부터 실리카 섬유는 용융되어 용융층을 형성한다. 온도가 높아짐에 따라 실리카 섬유는 연화되면서 용융되기 때문에 그 용융점을 명확히 규정하는 것이 어렵다. 또한 실제로는 실리카섬유가 연화되는 온도와 실리카-숯 화학반응 시작온도가 비슷하여 숯층과 용융층의 구분은 모호하다. 하지만 본 연구에서는 해석 모델의 편의를 위해 숯층과 용융층을 따로 구분하였다. 용융층에서는 숯층에서 시작된 실리카-숯 화학반응이 계속 진행된다. Eq. 8과 Eq. 9은 각각 용융층에서 반응식을 나타낸다.

$\[ Si{O}_2^S\to Si{O}_2^L \]$

(8)

$\[ C^S+Si{O}_2^L\to {SiO}^G+{CO}^G \]$

(9)

용융층은 두 단계인 가상의 비수축 용융층과 실제 수축 용융층으로 나누어 모델링할 수 있다. Fig. 2는 비수축 용융층과 수축 용융층을 나타낸다. 여기서 q_w는 전면에서의 열유속, q_zB는 용융층 후면에서의 열유속, r_w*는 비수축 용융층에서의 반경, r_w는 수축 용융층에서의 반경, r_rB는 비수축 용융층 후면에서의 반경, ∆η_*는 비수축 용융층에서의 층두께, ∆η는 수축 용융층에서 층두께, A_w는 단면적, ρ_z*는 비수축 용융층 혼합물의 밀도, ρ_{z exh*}는 비수축 용융층에서 소진상태의 밀도, ρ_z는 수축 용융층에서 혼합물의 밀도, ϵ_z*는 비수축 용융층에서 공극율을 의미한다.

Fig. 2

Non contracted melt layer and contracted melt layer.

비수축 용융층과 수축 용융층의 질량이 같다는 성질을 이용하여 용융속도를 구할 수 있다. Eq. 10은 용융층 후면의 용융속도를 나타내며, 용융층의 두께가 얇다는 가정과 에너지 평형식을 이용하여 구할 수 있다. 여기서 $$ y_{{\text{SiO}}_2}$$는 비수축 용융층에서 실리카 섬유 질량분율, $$ ∆Q_{{SiO}_2,melt}$$는 실리카 용융 흡열량 의미한다.

$\[ V_{zB,melt}=\frac{∆{\eta }_{*}}{∆t}=\frac{q_w-q_{zB}}{y_{{SiO}_2}∆Q_{{SiO}_2,melt}} \]$

(10)

실제 수축 용융층은 공극이 존재하지 않음으로 밀도가 아닌 질량변화를 아레니우스식을 이용하여 모델링 하며, 이를 비수축 용융층과의 질량 등치 관계를 이용하여 비수축 용융층에서의 겉보기 밀도 변화로 변환 시킬 수 있다. Eq. 11는 가상의 비수축 용융층에서 용융 실리카 섬유와 숯의 화학반응에 의한 겉보기 밀도 변화를 아레니우스 식으로 표한 식이다. Eq. 12와 Eq. 13은 각각 가상의 비수축 용융층 공극률과 용융층에서 생성된 반응가스 질량 유속을 나타낸다. 표면삭마 측면에서는 용융층에서 실리카-숯의 화학반응으로 인해 줄어든 질량만큼 표면이 후퇴한다.

$\[ -\left(\frac{d{\rho }_{z*}}{dt}\right)=K_c{e}^{-\frac{E_c}{\overline{R}T}}{\rho }_0\left({ϵ}_{c exh}-{ϵ}_{z*}\right) \]$

(11)

$\[ {ϵ}_{z*}={ϵ}_{z*1}+\int K_c{e}^{-\frac{E_c}{\overline{R}T}}\left({ϵ}_{c exh}-{ϵ}_{z*}\right)dt \]$

(12)

$\[ G_z=-\frac{∆{\eta }_{*}}{∆\eta }\int \frac{\partial {\rho }_{z*}}{\partial t}dr \]$

(13)

3. 실리카/페놀릭 복합재료 표면 열반응

4. 유한차분법을 이용한 열전도 해석

4.2 열전도 유한차분법

실리카/페놀릭의 열분해를 고려한 1차원 비정상 열전도 지배방정식인 Eq. 21과 Eq. 22를 유한차분법을 이용하여 해석을 수행하였다. 해석을 위한 이산화(Discretization) 작업에 있어서 공간(반경방향)은 중앙차분, 시간은 후진차분을 이용한 BTCS(Backward-time Centered-cpace)방법 기반으로 수행하였다[11]. Fig. 3은 노즐에서 실리카/페놀릭 열전도 유한차분법 격자를 나타낸다.

Fig. 3

Grids for finite difference method.

유한차분법을 이용하여 온도분포를 구한 다음 시간단계 내에서 발생한 삭마두께를 구해야한다. 시간단계에서 발생한 표면삭마는 Eq. 23를 이용하여 얻을 수 있으며, 시간에 대해 전진차분으로 근사화하였다. 발생한 삭마두께를 고려하여 다음 시간단계에서 격자크기를 재구성하여 계산을 수행하였다.

$\[ \begin{array}{c}{\delta }_{ablation}={\int }_{t_0}^{t_0+∆t}{V}_{w,oxi}+V_{zA,melt}+V_{zA,chem}\hfill \\ +V_{w,evap}+V_{w,shear}dt\hfill \end{array} \]$

(23)

5. 해석결과 및 고찰

본 논문에서 제시된 실리카/페놀릭 열반응 수치해석의 타당성을 검증하기 위해 내열성능평가모터인 TPEM(Thermal Performance Evaluation Motor)을 이용한 결과와 비교하였다[12]. TPEM모터는 평가할 내열재료를 목삽입재로 제작하여 내열성 및 단열성을 평가하고, 해당 내열재의 열반응 모델링의 적합성을 평가할 수 있다. Fig. 4는 TPEM 중 본 시험에서 사용된 10인치급 모터를 나타내며, 그레인은 후미연소형(End burning)이다. 사용된 모터의 제원 및 시험 성능은 Table 1과 같다. Table 2는 사용된 추진제의 연소가스 산화몰분율을 나타내며, CEA 프로그램을 사용하였다[13]. 목삽입재로 사용된 실리카/페놀릭 복합재료는 Flat Laminate 공법으로 제작하였다.

Fig. 4

Schematic of TPEM-10.

Table 1.

Dimension and performance about TPEM-10.

Table 2.

Mole fraction of oxidizing species of combustion gas at chamber and mass fraction of aluminum.

k,시험에서의 열반응 결과(삭마두께, 숯깊이)를 해석결과와 비교분석했다. Fig. 5는 실리카/페놀릭 열반응에 대한 시험결과와 해석결과를 나타낸다. 노즐목에서 해석결과 대비 삭마두께는 7%, 숯깊이는 11% 차이를 나타냈으며, 전반적으로 잘 일치하는 것을 확인 할 수 있었다.

Fig. 5

Results of thermal response of silica/phenolic.

해석결과와 시험결과의 차이는 여러 가지 원인이 있다. 실리카/페놀릭은 탄소/페놀릭의 열반응 대비 용융/고체 실리카-숯 화학 반응이 추가 발생하고 더 다양한 영향에 의한 표면삭마가 나타나므로 그 열반응 현상이 더욱 복잡하다. 용융층 구분의 모호성, 내부 열반응 및 표면삭마에 대한 열반응 상수의 불확실성, 용융층에서의 숯 노출 정도의 불확실성 등 실리카/페놀릭 재료의 열반응 특성치의 불명확성 영향으로 인해 시험과 해석결과의 차이가 나타난 것으로 추정된다. 이뿐만 아니라 경계조건 측면에서의 차이 발생 원인은 노즐과 같은 급가속 유동이 발생하는 경우 내부유동은 층류화 현상이 발생하여, 대류열전달계수에 영향을 미치고 노즐목 반경에 대한 노즐목 곡률반경의 비, 경계층 발달길이 등 다양한 영향 인자들로부터 기인된 것으로 판단된다[14,15].

Fig. 6은 노즐목에서 연소시간에 따른 온도변화를 나타내며, 각 시점에서의 숯 형성 위치는 세로방향의 선으로 표시하였다. 재료의 특성상 단열효과로 인하여 실리카/페놀릭으로 구성된 목삽입재의 배면에서 온도 상승은 발생하지 않았다. Fig. 6의 해석결과로부터 모든 시간에 대해 숯깊이를 기준으로 배면 방향으로 온도 구배가 작아지는 경향을 나타냈다. 또한 연소시간에 따라 숯깊이가 커지며 전반적으로 재료의 온도가 상승하는 것을 확인할 수 있다. 추후 온도해석에 대한 타당성은 목삽입재에 깊이별 온도센서를 장착한 후 온도측정 결과를 비교하여 평가할 수 있을 것으로 판단된다.

Fig. 6

Temperature distribution by thickness at nozzle throat.

6. 결론

고체 로켓 추진기관에 사용되는 내열복합재료인 실리카/페놀릭의 열반응에 대한 수치해석을 수행하였다. 이에 대한 타당성을 확인하기 위하여 내열성능평가모터인 TPEM-10을 이용한 열반응 결과와 비교하였다. 그 결론은 다음과 같다.

• ① 실리카/페놀릭의 내부 열반응과 표면삭마를 고려한 열전도 해석을 유한차분법을 통하여 수행하였고, 경계조건은 적분방정식을 이용한 대류열전달계수를 적용하여 안정적인 해를 얻을 수 있었다.
• ② 본 논문에서 기술한 해석방법을 이용한 결과와 내열성능평가 모터인 TPEM-10의 시험결과를 비교하였으며, 시험결과와 해석결과의 열반응 결과가 전반적으로 잘 일치하였다. 시험결과와 해석결과의 차이는 실리카/페놀릭 재료의 열반응 특성치의 불명확성 영향으로 인해 시험과 해석결과의 차이가 나타난 것으로 추정된다. 이뿐만 아니라 경계조건 측면에서의 차이 발생 원인은 급가속 유동에 의한 층류화, 노즐목 반경에 대한 노즐목 곡률반경의 비, 경계층 발달길이 등 다양한 인자가 있다.

Acknowledgments

이 논문은 한국추진공학회 2017년도 춘계학술대회(2017. 5. 31-6. 2, 라마다프라자 제주호텔) 발표논문을 심사하여 수정 ․ 보완한 것임.

References

• Reydellet, D., "Design Methods in Solid Rocket Motors", Advisory Group for Aerospace Research and Development, AD-A199356, (1988).
• ZieBland, H., and Parkinson, R.C., "Heat Transfer in Rocket Engines", Advisory Group for Aerospace Research and Development, AGARD-AG-148-71, (1971).
• Shi, S., Liang, J., Yi, F., and Fang, G., "Modeling of one-dimensional thermal response of silica-phenolic composites with volume ablation", Journal of Composite Materials, 47(18), p2219-2235, (2013). [https://doi.org/10.1177/0021998312454907]
• Bae, J.Y., Song, J.W., Kim, T.W., Hahm, H.C., and Cho, H.H., "Study on the Ablation of Silica-Phenolic Ablator by Numerical Method", KIMST Annual Conference Proceedings, Gyeongju, Korea, p2010-2013, Jun), (2012.
• Bae, J.Y., Park, S.K., Kim, T.W., Hahm, H. C., Bae, J.C., and Cho, H.H., "Kinetic Constants Measurement of Silica/Phnenolic Materials", 16th Guided Weapons Conference, Daejon, Korea, p326-330, Sep), (2012.
• Yu, M.S., Cho, H.H., Hwang, K.Y., and Bae, J.C., "Hybrid method for jet vane thermal analysis in supersonic nozzle flow", Journal of Thermophysics and Heat Transfer, 20(3), p614-617, (2006). [https://doi.org/10.2514/1.17675]
• Seo, S.K., Hahm, H.C., and Kang, Y.G., "Analysis of Boundary layer in Solid Rocket Nozzle and Numerical Analysis for Thermal response of Carbon/Phenolic using Finite Difference Method", 2016 KSPE Fall Conference, Jungsun, Korea, p748-755, Dec), (2016.
• Bae, J.C., "Analysis Model for the Indepth Thermal Response of Charring Ablators", Agency for Defense Development, ADDR-421-130474, (2013).
• Lapp, P., and Quesada, B., "Analysis of Solid Rocket Motor Nozzle", 28th Joint Propulsion Conference and Exhibit, Nashville, T.N., U.S.A., AIAA 92 3616, Jul), (1992.
• Boyarintsev, V.I., and Zvyagin, Yu. V., "Turbulent Boundary Layer on Reacting Graphite Surface", 5th Int. Heat Transfer Conference, Tokyo, Japan, p264-268, Sep), (1974.
• Levy, D., "Introduction to numerical analysis", Department of Mathematics and Center for Scientific Computation and Mathematical Modeling, University of Maryland, College Park, M.D., U.S.A., (2010).
• Hahm, H.C., "A Study on the Thermal Response Characteristics of Carbon/Carbon Composites for Nozzle Throat Insert", Journal of the Korean Society of Propulsion Engineers, 10(1), p30-37, (2006).
• McBride, B.J., and Gordon, S., "Computer Program for Calculation of Complex Chemical Equilibrium Composition and Applications, II. Users Manual and Program Description", NASA RP-1311, (1996).
• Torii, S., and Yang, W.J., "Numerical study on laminarizing gas flow in strongly heated pipe", International Journal of Heat and Mass Transfer, 40(13), p3105-3117, (1997).
• Hahm, H.C., and Kang, Y.G., "Comparative Studies of Heat Transfer Coefficients for Rocket Nozzle", Journal of the Korean Society of Propulsion Engineers, 16(2), p42-50, (2012). [https://doi.org/10.6108/KSPE.2012.16.2.042]