Compressible-incompressible Two-phase Research of the Cryogenic Cavitation Flow about Thermal effects and several Cavitation models 초 록 일반적으로, 극저온유체에서의캐비테이션현상은상당한온도효과와유체물성치의강한변화를발생시킨다. 이러한온도효과를이해하고정량화하는것은액체수소와액체산소를가압하는액체로켓터보시스템을디자인하는데중요하다. 이현상을정확하고강건하게시뮬레이션하기위해서는, 수치기법이압축성 - 비압축성이상유동장의공존을다룰수있어야하고, 극저온유체의상태방정식이정확하게모델링되어야한다. 본연구에서는, 본연구그룹에서개발한강건한압축성 - 비압축성이상유동수치기법이이용되었고, regression analysis 를이용하여상태방정식을모델링하였다. 이를바탕으로, 액체질소에서의캐비테이션현상해석결과를수록하였다 ABSTRACT Generally, cavitation in cryogenic fluids generates substantial thermal effects and strong variations in fluid properties. To simulate this phenomena accurately and robustly, numerical method should be able to treat the co-existence of compressible and incompressible two-phase flow fields, and equation of state of cryogenic fluid should be accurately modeled. In this paper, a robust compressible-incompressible two-phase scheme which was developed by our research group is used and equation of state is modeled using a regression analysis. And then, numerical simulations of cavitation in liquid nitrogen are presented. Key Words : Cryogenic Fluid( 극저온유체 ), Cavitation( 공동 ), Thermal Effects( 온도효과 ), Compressible-incompressible Flow( 압축성 - 비압축성유동 ) 1. 서론일반적으로, 극저온유체의캐비테이션현상은온도효과와유체물성치의강한변화로인하여캐비테이션특성이작동유체가물일때와는다른양상을보인다. 이러한극저온캐비테이션현상은주로액체로켓의산화제와연료를가압하는터보펌프인듀서주변에서발생하게되는데, 캐비테이션이발생하면터보펌프가엔진연소부에필요한압력을상승시켜주지못하게되어문제가발생하기도하고, 주기적인캐비테이션진동으로인한구조적문제점을발생시키기도한다. 따라 서이러한극저온캐비테이션현상을정확하게예측하는것은터보펌프성능을향상시키는데있어아주중요하다. 극저온유체는유체의임계점근처에서작동영역을가지게되므로액체와기체의부피비가물에비해서작아지게되고, 이에따라캐비테이션발생시일정량의기체를생성하는데물보다더많은양의액체가필요하게된다. 따라서 evaporation cooling effect가현저하게발생하여국부적으로유체의온도를하강시키게된다. 한편, 캐비테이션과관련된유체의특성중하나인 101
vaporization pressure 는온도의함수로서유체가 임계온도근처에서는온도변화에매우민감하게반응한다. 따라서전술한국부적온도하강은 vaporization pressure을낮추게되고따라서캐비테이션발생이억제되게되는데이를극저온유체캐비테이션온도효과라고한다. 이러한극저온캐비테이션에대한연구는여러연구그룹에의해연구되어왔다. Ruggeri and Moore 와이후에 Hord 는다양한형상에대한극저온실험을통하여경계에서의열소산효과가이러한캐비테이션억제의주된원인이라고분석하였다. Cooper 는에너지방정식을고려하지않고비압축성유동으로만간주하여해석하였다. Hosangadi et al. 연구그룹은압축성지배방정식을이용하여압축성-비압축성동시해석이가능한알고리즘을통해해석을수행하여비교적정확한결과를보여주고있지만에너지방정식에있어압력에의한일을무시하였다. 본연구에서는완전한형태의지배방정식으로부터압축성-비압축성의넓은마하수를해석할수있도록예조건화를도입하였고, 본연구그룹에서개발한정확하고안전한이상유동 RoeM 수 치기법 을사용하였다. 또한극저온유체의특징을정확하게모사하기위하여 NIST database 를이용한 regression analysis를하여상태방정식을근사하였다. 이를바탕으로다양한형상에서의캐비테이션현상을해석하였고, 온도효과를분석하였다. 2. 극저온이상유동시스템 2.1 지배방정식균질혼합류모델 (Homogeneous Equilibrium Model) 을이용한압축성 / 비압축성이상유동지배방정식은혼합류의질량, 운동량, 에너지보존식과기체상의질량보존식으로구성된다. 각각의상을구분하는함수로는여러가지방법이있지만본연구에서는질량비율 (mass fraction) 을이용하여액체상과기체상을구분하였다. 예조건화된 2차원이상유동지배방정식은다음과같이표현된다. (1) 보존량벡터 와원시변수벡터 는다음과 같다., (2), (3) 위식에서 은혼합유체의밀도이고, 은 기체상의질량비율이다. 저마하수에서의시스템의 stiffness 문제를해결하기위하여시스템을예조건화하였다. 한편식 (1) 에서비점성및점성플럭스벡터는다음과같다. (4) (5) (6) (7) 또한상변화를묘사하는 source 항으로서 는 다음과같이나타난다. (8) sensible energy와 latent energy를포함하는 Total internal energy에대해에너지방정식을구성하여상변화시발생하는에너지변화를고려할수있었다. 여기서사용한상변화모델은 Merkle et al., Kunz et al., Singhal et al. 의모델을사용하였다. Merkle et al. Kunz et al. 102
Singhal et al. 성치를계산한결과를비교한예로써, MEOS를비선형해법을통해계산한결과와비교해높은수준의정확도를확보할수있음을보여준다. 한편혼합류의밀도는다음과같이제어체적 내에서각상이차지하고있는영역에각각의상태방정식을적용하여구해진다.. (9) 이제각계산격자내에서는각상의압력과온도가같다는동역학적, 열역학적평형을가정함으로써, 전체 2상유동의지배방정식시스템이닫히게된다.,. (10) 2.2 상태방정식극저온유체는잠열이상온의다른유체에비해훨씬작기때문에작은온도변화에도물성치의변화가민감하게나타난다. 이러한특성때문에극저온유체유동을해석하는데정확한상태방정식을사용해야하며기존의극저온유체의상태방정식으로넓은온도및압력범위에서높은정확도를보장하는 MEOS(Multiparameter EOS) 가널리사용된다. 본연구에서는지배방정식을예조건화하는과정에서보존변수중밀도를압력으로대체하였으나위의열거한상태방정식들은압력을밀도와온도의함수로표현하는형태로되어있다. 그러므로음속을계산하는과정에서밀도를계산하기위해상태방정식에대한비선형해법을필요로한다. 이러한과정은수치해석과정의효율성을크게감소하는원인으로작용하기때문에본연구에서는상태방정식의물성치에대한 regression analysis를수행하여상태방정식으로대체하였다. Fig 1은 Nitrogen에대해 MEOS와본연구에서적용한상태방정식을통해음속, 엔탈피등의물 Fig 1. Comparison of MEOS and Regression model (Example : Nitrogen, 70K) 3. 수치기법 3.1 이상유동 RoeM 수치기법기체상에서 RoeM수치기법 을유도하는과정과같이 2상유동 Roe 수치기법으로부터 RoeM 수치기법을유도하였다. 이때, 제어함수 f와 g에는다음과같이새롭게고안된 2상유동충격파포착항이사용된다., (11), (12). (13) 최종적으로접촉불연속면을정확히포착하면서도팽창영역의불안전성을제거하기위해 signal velocity b1과 b2를도입하고, 최종형태에서음속과마하수를재조정하는방식에따라시스템고유값을속도의차수로재조정한예조건화된 2상유동 RoeM 수치기법은다음과같이정리된다., (14), (15) 103
(16),, (17) 석결과를관찰한다. 격자계는 Fig 2와같고, 유 동조건은 Table 1과같다. Table 1. Cryogenic fluid flow condition ( ) Run Number 290C 83.06 23.9 1.70 1.13x10 여기서, (18) 3.2 이상유동 LU-SGS 시간적분법 LU-SGS 는 Block diagonal solver 가필요없는내재적시간적분기법으로, 효율성이좋기 때문에널리사용되고있다. 본연구에서는예조 건화행렬을포함하는비정상 2 차원점성문제에 대해, 비점성항과점성항모두를내재적으로처 리한 LU-SGS 를사용하였고, 다음과같다., (19). (20) Fig 2. Grid system for Hord hydrofoil, (21), (22), (23) 위식은 block diagonal solver가필요하지않으며, 주어진계산격자에서단순한 형태에대한 solver만을필요로한다. 그런데 는결정된형태이므로, 본연구에서는 x에대한일반해를미리구해이용하였다. 따라서시간적분과정에서실질적인 matrix solver는사용되지않았다. 본연구에서는이러한 LU-SGS를바탕으로해석을수행하였다. 4. 수치해석결과 4.1 Hord hydrofoil cryogenic cavitation problem Hord 연구그룹에서수행한 hydrofoil 실험결과와, 같은 problem에대한 Hosangadi 연구그룹의해석결과와비교하여 cryogenic cavitation solver의성능을검증하고, cavitation model에따른수치해 Fig 3. 1-phase pressure contour Fig 4을보면 contour에서나타나는액체와기체의경계가표면에수직하게분포하는것을알수있다. 이것은작동유체를물로하였을때와는다른결과로, 극저온유체를이용했을때나타나는캐비테이션억제현상이다. 극저온유체는물에비해기체상과액체상의부피비가작기때문에같은부피의기체상으로변하기위해더많은액체가필요하게되고, 이에따른잠열로인해국부적으로온도가낮아진다. 그에따라 vaporization pressure가낮아져캐비테이션현상이억제되는것으로나타난다. 이는온도를나타낸 Fig 5을통해알수있다. 또한 Fig 6을보면극저온유체에서의캐비테이션현상은 Mach number변화가 1.5에달하는상당한압축성효과를동반하고있는것을알수있다. 또기체상과액체상의경계면이유선과일 104
치하지않아경계면을통과하여유동이진행되는특성을확인할수있다. Fig 8과 Fig 9을보면실험치에비해 pressure와 temperature 모두 expansion을더크게예측하고있음을알수있다. Hosangadi 그룹의결과와비교하면세모델모두 cavitation길이를더길게예측하고있고, Merkle model이액체와기체의경계면을다른 model에비해날카롭게잡아내고있다. Singhal model의경우는 closure region에서의급격한압력회복을예측하지못하고있다. 전체적으로, Hosangadi 그룹에비해실험치에근접한결과를나타내고있음을확인할수있었다. Fig 7. Pressure contour for cryogenic cavitation Fig 8. Pressure depression comparison Fig 4. VOF contour for cryogenic cavitation Fig 5. Temperature contour for cryogenic cavitation Fig 6. Mach contour for cryogenic cavitation Fig 9. Temperature depression comparison 5. 결론액체로켓터보펌프등에서발생하는극저온유체의캐비테이션현상을해석하기위하여온도효과를모사할수있는지배방정식및상태방정식을구성하였고, 본연구그룹에서개발된비압축성에서압축성에이르는넓은범위의마하수에서정확하고안정적으로작동하는수치기법을사용하여물과극저온유체의캐비테이션현상을해 105
석하였다. 작동유체의차이로인하여극저온유체에서는물과는다른캐비테이션특성, 즉온도효과에의한캐비테이션억제효과가발생하였고, 이로인하여압력계수의형태도현저한차이가나타났다. 또한물과는달리압축성영역에해당하는유동이나타남으로인하여압축성-비압축성동시해석의중요성을확인할수있었다. 후기본연구는한국연구재단을통해교육과학기술부의우주기초원천기술개발사업 (NSL, National Space Lab, 과제번호 20090091724), 2011년도 2단계두뇌한국21사업및항공우주신기술연구소의지원으로수행되었습니다. 참고문헌 (1) Ruggeri, R. S., and Moore, R. D., 1969, "Method for Prediction of Pump Cavitation Performance for Various Liquids, Liquid Temperature, and Rotation Speeds", NASA, TN D-5292 (2) Hord, J., 1973, "Cavitation in Liquids Cryogenics II- Hydrofoil", NASA CR-2156 (3) Cooper, P., 1967, "Analysis of Single and Two-Phase Flows in Turbopump Inducers", J. Eng. Power, 89, pp. 577-588 (4) Ashivin Hosangadi and Vineet Ahuja, 2005, "Numerical Study of Cavitation in Cryogenic Fluids", J. Fluids Engineering, 127, pp. 267-281 (5) Ihm, S. and C. Kim, 2008, "Computations of Homogeneous Equilibrium Two-Phase Flows with Accurate and Efficient Shock-Stable Schemes", AIAA journal, 46(12), pp. 3012-3037 (6) NIST Reference Fluid Thermodynamic and Transport Properties Database(REFPROP): version 8.0," NIST standard Reference Database 23[online database], http://www.nist.gov/srd/nist23.htm (7) Merkle, C. L., Feng, J. Z., and Buelow, P.E.O, 1998, "Computational Modeling of the Dynamics of Sheet Cavitation", Proceedings of the 3rd International Symposium on Cavitation, Grenoble, France, April 7-9 (8) Robert F. Kunz, David A. Boger, David R. Stinebring, Thomas S. Chyczewski, Howard J. Gibeling, 1999, "A preconditioned Navier-Stokes method for two-phase flows with application to cavitation prediction" AIAA journal, pp. 3329-3341 (9) Ashok K. Singhal, Mahesh M. Athavale, Huiying Li, Yu Jiang, 2002, "Mathematical basis and validation of the fill cavitation model", Journal of fluids engineering, Vol. 124, pp. 617-624 (10) Roland Span, Eric W. Lemmon, Richard T Jacobsen, 2000, "A reference equation of state for the thermodynamic properties of nitrogen for temperature from 63.151 to 1000K and pressures to 2200Mpa", J. Phys. Chem. Ref. Data, Vol 29, pp.1361-1433 (11) J. W. Leachman, R. T Jacobsen, S. G. Penoncello, E. W. Lemmon, 2009, "Fundamental equation of state for parahydrogen, normal hydrogen, and ortho hydrogen", J. Phys. Chem. Ref. Data Vol 38, pp. 721-748 (12) Kim, S.-s., Chongam Kim, 2003, "Cures for the shock instability: Development of a shock-stable Roe scheme" Journal of Computational Physics, 185(2): pp. 342-374 (13) Yoon, S. and A. Jameson, 1988, Lower-upper symmetric-gauss-seidel method for the Euler and Navier-Stokes equations. AIAA journal, 26(9): pp. 1025-1026 106