초공동 수중운동체 주위 공동 특성과 추력 전산 해석

Numerical Analysis of Cavity Characteristics and Thrust for Supercavitating Underwater Vehicle

  • ABSTRACT

    Cavitation is used in various fields. This study examined the drag reduction of an underwater vehicle using cavitation. In this study, the natural partial cavitation analysis results were verified using CFD code with the Navier-Stokes equation based on a mixture model. The momentum and continuity equations in the mixture phase were separately solved in the liquid and vapor phases. The solver employs an implicit preconditioning algorithm in curvilinear coordinates. The results of a computational analysis showed good agreement with the experiment. A computational analysis was also performed on the supercavity. The study investigated the cavity characteristics and drag of an underwater vehicle and studied the speed required to achieve a supercavity. Finally, a 1DOF analysis was carried out to investigate the thrust system for a supercavity. As a result, one of the methods for determining a suitable thrust system for a supercavitating underwater vehicle was presented.


  • KEYWORD

    공동 , 초공동 , 수중운동체 , 추력 , 1자유도

  • 1. 서 론

    공동은 국부적인 압력 강하로 인한 기화 현상을 의미하는데, 일반적인 기계공학의 범주에서는 공동의 생성과 소멸 과정에서부터 기인하는 기계의 부식, 소음 등의 문제로 부정적인 인식이 있다. 하지만 최근에는 이러한 공동 현상 이점을 이용하여 수중운동체 항력저감에 활용하는 연구가 국내외에서 활발하게 진행중이다. 이러한 수중운동체 주위 공동유동장 해석 방법 중 가장 보편적인 방법은 각 상의 에너지 방정식을 하나의 혼합류(Mixture flow)에 대해 계산하는 방법으로 예조건화(Presonditioning)와 이중 시간방법(Dual time stepping method)을 이용한 부분-초공동 수중운동체 연구를 수행하였다(Kunz et al., 1999). 초공동 수중운동체의 경우 연소가스에 의한 추진을 포함한다고 알려져 있어 고온, 고압의 연소가스 해석을 위해서 압축성 효과를 고려한 연구가 필요하다. 하지만 기존 해석 모델의 경우 압축성 효과를 고려하지 못 하였는데 최근에는 압축성과 온도 변화 까지 모두 고려한 완전 압축성 모델을 이용한 전산 해석이 수행 되었다 (Linau et al., 2003; Owis and Nayfeh, 2003). 국내에서는 비점성 공동 해석을 국외 연구 결과 및 다양한 실험과 비교 검증하는 전산 해석이 수행되었으며(Kim et al., 2013b), 점성을 고려한 3차원 내비어-스톡스(Navier-Stokes) 방정식을 이용하여 수중운동체 주위 다상유동장 전산 해석 역시 수행 되었다(Kim et al., 2013a; Jin et al., 2013). 최근에는 완전 압축성 효과를 반영한 다양한 유동 현상에 대한 연구(Ha and Park, 2016)가 수행되었으며, 그 영역이 공동-초공동 다상 유동 해석 영역 까지 확대 되었다(Park, 2016).

    초공동(Supercavitation)은 공동이 수중운동체를 완전히 덮게 되는 현상을 말하며, 물과 공기의 밀도차이에 의해 수중운동체의 항력(Drag)이 크게 감소하여 속도 증가를 가져오게 된다. 초공동 어뢰의 경우 러시아에서는 이미 개발되어 운용하고 있는 것으로 알려져 있으며, 미국, 독일 그리고 최근에는 중국도 다양한 연구를 수준급으로 진행 하고 있는 것으로 알려져 있다. 대부분의 정보가 군사적 목적에 의한 것이므로 현재 단계에서는 제한적이며 일부 공개된 정보에 따라 초공동 수중운동체 유사 모형에 대한 점섬, 비점성 유동장 해석과 초공동 로켓 시스템에 대한 일부 연구가 국내에 알려져 있다(Kim et al., 2013b; Kim et al., 2013c; Ahn et al., 2014).

    초공동 어뢰의 경우 선두부로부터 공동이 발생하게 되는데, 이때 발생하는 공동의 길이 및 직경의 정도에 따라 수중운동체의 초공동화가 결정이 된다. 이러한 공동을 발생 시킬 수 있는 조건은 인위적인공동(Artificial cavitation)을 이용하지 않는 다면 수중운동체의 자체 속도 증가가 유일한 방법이고 이를 달성하기 위해 현재 다양한 방법을 이용한 추력 시스템에 대한 연구가 진행 중이다.

    본 연구에서는 3차원 내비어-스톡스 방정식을 이용한 부분공동 전산해석을 수행하고 검증하여 이후 다양한 속도 조건에 따른 수중운동체 주위 다상 유동 전산 해석을 통해 공동의 길이 및 직경을 계산하여 초공동 달성 조건에 대해 확인 한 후, 이러한 속도 조건을 얻기 위한 추력에 대한 선행 연구를 수행하였다.

    2. 지배 방정식 및 수치해석 방법

    본 연구에서는 수중운동체 항력 저감을 위한 공동 관련 연구에 필요한 무차원 변수를 확인하고 무차원 변수 값 변화에 따른 공동의 형상적 특징에 대한 연구와 이를 확인하기 위해 자체개발 다상유동 전산해석 프로그램을 이용한 전산 해석을 수행하였다. 해당 프로그램은 3차원 내비어-스톡스 방정식을 기반으로 다상(Multiphase)에 대해서는 각각의 연속방정식을 따로 풀고 운동량 방정식의 경우에는 혼합(Mixture)류 방정식을 사용하는 모델을 사용하였다.

    image
    image
    image

    여기서 하첨자 ν, l 그리고 m은 기상, 액상 및 혼합을 의미하고 p, u, ρ, α는 압력, 속도, 밀도, 체적분율(Volume fraction)을 의미하며 r는 물리적 시간, t에 대응하는 의사 시간(Pesudo-time)을 의미한다. 혼합류에서의 밀도의 정의와 체적분율은 다음과 같다.

    image
    image

    위 식을 기반으로 전산 해석을 위한 일반 곡선 좌표계(Generalized curvilinear coordinate system)로 변경하는 식은 다음과 같다.

    image

    여기서, 이고 , , 는 각 ξ , η, ζ방향에 대한 대류 플럭스 항이며, , , 는 그에 대응하는 점성 플럭스 항이고, Γe는 자코비안(Jacobian) 행렬이다. Γ는 예조건화(Pre-conditioning) 행렬을 의미하고 이 기법을 통해서 전체 계산과정에서 안정성을 높여주었다(Owis and Nayfeh, 2004).

    공동 모델의 경우에는 기화(Evaporation)와 응축(Condensation)의 비율을 나타내는 과 을 사용하여 기상과 액상 사이의 물질 전달(Mass transfer)에 대해서 수학적 모델링 하게 되는데, Merkle에 의해 종속 변수에 단점 개선이 증명된 모델을 사용하였다(Park et al., 2009; Ha and Park., 2012).

    image
    image

    여기서 하첨자 ν는 기상, l은 액상을 의미하며 kVkl은 스케일링 상수이다.

    3. 전산 해석 결과

       3.1 자연 부분공동(Natural partial-cavitation) 전산 해석

    공동에 의한 항력저감으로 얻을 수 있는 고속 수중운동체의 경우 그 주위를 감싸는 공동은 공동수(Cavitation number)라는 무차원 수로 정의 될 수 있으며, 공동 수(σ)에 따라 공동에 대한 물리적인 특성을 평가 할 수 있으며 다음은 그 정의 이다.

    image

    여기서, P는 수중운동체 주위의 압력, PV는 유체의 증기압 이고 ρ는 유체의 밀도, V는 유동 속도이다

    Fig. 1은 반구형 선두부(Hemispherical forebody)를 가지는 실린더 형상 수중운동체에서 서로 다른 공동 수(σ)에 대한 전산 해석 결과이다. 표면 압력 계수(Surface pressure coefficient)에 대한 전산해석 결과가 실험(Rouse and McNown, 1948)결과와 실린더 길이와 직경을 무차원화 한 s/d 값에 따라 유사함을 확인 할 수 있는데 특히 공동 발생 구간에서는 거의 오차 없이 일치함을 보이고 있다. s/d 값이 3이상인 경우 낮은 공동수(σ)에서 일부 차이를 보이고 이러한 차이는 Fig. 2에서 두르러 진다. Fig. 2는 선두부 형상이 다른(0-caliber) 경우에 대한 전산 해석 결과이고 공동 수에 따라 대체로 잘 일치 함을 보이고 있는데, 일부 공동수 영역에 대해서 공동 끝단 구간에서 Owis & Neyfeh의 해석과 차이를 보이고 그 이유는 Owis & Neyfeh의 해석에 사용된 공동 모델(Cavitation model)이 공동 끝단에서 압력 회복에 의한 상변화를 정확하게 모델링 하지 못한 결과로 판단된다. 본 연구에서는 개선된 모델을 사용함으로써 이러한 문제를 해결하여 모든 구간에서 향상된 결과를 보이고 있다.

    Fig. 3는 반구형 실런더에 대해 공동의 체적분율(Volume fraction)을 Kunz et al.(1999)의 전산해석 결과와 비교한 것으로 공동 시작점 및 전체 공동 크기가 상당히 잘 일치함을 보이고 있으며 공동 끝점의 차이는 공동 후미부에서 발생하는 재유입류(Re-entrant jet)등의 영향으로 판단된다. Figs. 1-3을 통한 자연 부분 공동 관련 실험과 다른 연구자 해석 결과와의 비교 결과, 공동의 크기 및 길이 뿐 아니라 표면 압력 계수를 포함한 공동수(σ)를 포함한 유체력 예측에 있어서도 본 연구의 전산 해석 프로그램 사용이 적합하다는 것을 확인 하였다.

       3.2 초공동 전산 해석

    부분공동 전산해석을 통한 해석 프로그램의 적합성을 바탕으로 실제 수중운동체와 유사한 모형에 대해서 공동 특성을 살펴보고자 하였다. 속도 변화에 따라 공동의 직경 및 길이가 변화하고 일정 속도 이상의 구간에서 공동이 수중운동체 전체를 덮는 초공동에 대한 전산해석을 수행하였다. Fig. 4은 해석 대상 모델이며, 그 형상은 러시아에서 이미 개발한 시크발(Shkval) 초공동 어뢰 형상을 단순화 한 것이다. 해당 수중운동체의 선두부는 판형(Diks-type)으로 그 지름(DC)이 5mm이고, 몸체의 지름은 3×DC, 전체 길이는 30×DC이다. 공동의 직경(DS)은 그 값이 최대가 되는 값으로 정의 하였으며, 길이(LS)는 선두부로부터 공동의 끝단까지의 값이며, 각 값들은 선두부(Cavitator) 지름(DC)에 기반한 무차원 수로 나타 내었다. 해석은 수심 1.0m인 수온 15.0℃물에 대해서 수행하였고 자유수면에 의한 공동 및 유동 변화는 없으며 수중운동체는 충분히 잠긴 상태로 가정하였고 전체 압력장은 중력가속도(g)의 영향을 받는다.

    Fig. 5는 서로 다른 속도 조건에 따른 항력(Drag) 해석 결과이다. 수중운동체 전체 항력은 속도가 증가 함에 따라 큰 선두부 저항에 의해 선형적으로 증가하지만 몇 차례 큰 변곡점을 지나게 되는 것을 확인 할 수 있고, 80m/s 이상의 속도 조건에서는 완전한 초공동이 형성되어 다시 선형적으로 증가 함을 알 수 있다. 속도가 60m/s와 80m/s 구간에서의 항력은 여러차례 급격히 변하되는데 그 이유는 공동이 수중운동체 몸체에서 차지하는 비율이 높아지면서 부분 공동에서 초공동으로 변하기 때문이라 판단된다. Fig. 6은 항력의 총합을 압력(Pressure)과 점성(Viscous)에 대한 값으로 나눈 것이고 초공동이 확연하게 구분되는 60m/s이전과 80m/s이후 영역의 점성 저항이 크게 감소하는 것을 확인 할 수 있다. Fig. 7에서는 전 속도 구간에서의 공동의 길이 및 직경에 대해서 체적분율로 나타내었고 그에 해당하는 속도와 항력 값을 나타내었다. 이를 통해서 공동이 해당 수중운동체를 충분히 덮을 수 있는 속도 조건은 75m/s 이상 조건이며 75m/s 조건에서 공동이 모델의 선두부 직경 대비 3.83배의 직경과 37.17배의 길이를 가질 때 초공동 수중운동체가 된다는 것을 알 수 있고, 이때 공동의 길이를 수중운동체의 몸체 길이 대비로 나타내면 1.15배의 길이에 해당한다. Fig. 8은 선두부 직경 대비 공동의 직경을 속도에 따라 나타낸 것으로 속도가 증가함에 따라 공동의 직경이 커지는 것을 확인 할 수 있었다. Fig. 9는 선두부 직경 대비 공동의 길이에 해당하는 결과인데, 여기서 공동의 길이는 선두부로 부터 이어지는 공동의 길이로 Fig. 6에서 보여지는 70m/s 이하의 속도 조건에서 후미부 공동의 길이는 제외된다. 70m/s와 75m/s 속도 구간에서 공동의 길이가 갑자기 증가하는 것은 초공동으로 인한 후미부 공동과 길이가 연결되었기 때문이다. Figs. 5-9를 종합하면, 해당 수중운동체 모델의 경우에는 75m/s 이상의 속도를 가질때 초공동이라 할 수 있지만, 75m/s 구간에서, 후미부 일부 영역에서의 공동 직경 값이 상대적으로 작아 외란 등의 외부 요인에 의해 파괴 될 수 있을 수 있는 등의 불안요소를 가지고 있다. 따라서 안정적으로 완전히 발달한 초공동 효과를 가질 수 있는 80m/s이상의 속도 구간이 필요 할 것으로 판단할 수 있다.

       3.3 추력을 이용한 수중운동체의 초공동화 연구

    인위 공동(분사공동)을 제외한 자연 초공동 수중운동체 달성을 위해서는 후미부 추력에 의한 수중운동체의 속도 증가가 요구되는데, 수중운동체의 경우 단순 프로펠러 구동 만으로는 초공동에 필요한 속도 조건을 얻을 수 없다. 따라서 별도의 추력 시스템에 대한 요구가 있으며 이때 실제 수중운동체 크기 및 무게를 고려한 제한 된 추력 시스템이 필요하다. 여기서는 수중 운동체가 65kg의 질량을 가지는 모델에 대해서 추력 시간에 대한 1자유도를 가지는 전산해석을 수행 하였으며, 그 지배 방정식은 다음과 같다.

    image

    위 식 (10)은 6자유도를 고려할 수 있는 지배 방정식이며, 여기서 A는 질량 행렬(Mass matrix), s는 직선 및 회전 속도를 의미하고 b는 관성 커플링을 의미한다. 그리고 f는 힘과 모멘트를 의미하고 해당 연구에서는 수중운동체에 작용하는 힘은 추력으로 1자유도 만을 고려한 해석을 수행한다.

    1.5초와 2초의 추력 시간 조건에 대한 해석을 수행하였으며 2초 후에는 추진체 질량 6kg가 감소한 59kg가 수중운동체의 최종 질량이 된다. 연구 목적에 따라 초공동을 달성하기에 적합한 추력은 일정 값 이상에서 선형적으로 증가하게 되는데, 이는 실제 초공동 어뢰의 추력 시스템과 유사하게 모델링 되었다(Kim et al., 2013c).

    Fig. 10은 시간에 따라 변하는 1.5초와 2.0초의 추력에 대한 조건을 나타낸 것이다. 각 조건에 해당하는 시간 이후에는 급격히 차단하여 추력이 없는 상태에서 수중운동체는 감속 운동을 하게 된다.

    해당 수중운동체 모델의 경우 자연 초공동 전산 해석의 결과로 미루어 볼 때, 75m/s에서 초고동이 달성되고 80m/s 이상의 구간에서 공동이 안정적으로 완전히 발달하여 초공동 이 되는 것을 확인 하였는데, Fig. 11의 결과에서 보면 추력이 1.5초 동안 지속되는 경우 해당 수중운동체의 최고 속도는 76.1m/s, 추력이 2.0초 동안 지속되는 경우는 최고 속도가 88.1m/s 이므로 완전한 초공동 달성을 위해서는 1.5초 보다는 2.0초 동안 추력이 지속된 시간 조건이 타당하다고 할 수 있다. 두 가지 경우 모두 속도가 추력의 지속시간 동안에는 증가하고 지속시간 이후에서는 급격히 감소하는 것을 확인 할 수 있었다.

    4. 결 론

    본 연구에서는 자연 부분, 초월 공동 해석을 3차원 네비어-스톡스 방정식 기반의 다상유동 전산해석 프로그램을 이용하여 수행하였다. 우선적으로 자연 부분공동에 대해서 표면 압력 계수와 공동의 형상을 이용한 비교를 통해 신뢰성을 검증 하였으며, 그 적합성을 바탕으로 정지상태부터 100m/s 이상 까지 다양한 속도 조건에 대한 초공동 유동 전산해석을 수행하였다. 초공동 어뢰 형상의 수중운동체 주위 공동 유동 특성을 예측 하였는데, 이를 통해 안정적인 초공동 달성을 위한 속도를 계산 하였으며, 해당 수중운동체 모델의 경우 75m/s 이상 조건에서 초공동화가 달성되었고 이때 공동은 선두부 직경 대비 3.83배 직경과 37.17배 길이의 특성을 보였다. 이러한 속도를 달성하기 위해 1자유도를 고려한 추력의 지속시간 연구를 수행하여 적합한 시간 조건 역시 제시 되었다. 이러한 연구는 초공동을 달성하기 위한 임의 수중운동체 모델에 대한 추력 시스템에 대한 연구이며, 이는 앞서 제시한 속도에 따른 공동 특성과 함께 고려되어 초공동 수중운동체의 제원과 추력 시스템을 설계 및 결정에 적용이 가능 한지에 대한 향후 연구가 필요할 것이다.

  • 1. Ahn B.K., Kim J.H., Choi J.K., Kim H.T., Nah Y.I., Lee D.H. 2014 Numerical Analysis of Supercavitating Flow Based on Viscous/Inviscid Method [Journal of the Korea Institute of Military Science and Technology] Vol.17 P.25-32 google doi
  • 2. Ha C.T., Park W.G. 2012 Sensitivity Evaluation of Empirical Coefficients in Cavitation Models [Korea-Japan CFD Workshop] google
  • 3. Ha C.T., Park W.G. 2016 Evaluation of a New Scaling Term in Preconditioning Schemes for Computations of Compressible Cavitating and Ventilated Flows [Ocean Engineering] Vol.126 P.432-466 google doi
  • 4. Jin M.S., Park W.G., Jung C.M. 2013 Numerical Analysis of Cavitating Flow Past an Axisymmetric Cylinder by Comparison with Experiments [Journal of Mechanical Science and Technology] Vol.27 P.3673-3681 google doi
  • 5. Kim D.H., Park W.G., Jung C.M. 2013a Numerical Simulation of Cavitating Flow Past Axisymmetric Body [International Journal of Naval Architecture and Ocean Engineering] Vol.4 P.256-266 google
  • 6. Kim J.H., Jang H.G., Ahn B.G., Lee S.C. 2013b A Numerical Analysis of the Supercavitating Flow around Three-dimensional Axisymmetric Cavitators [Journal of the Society of Naval Architects of Korea] Vol.50 P.160-166 google doi
  • 7. Kim K.M., Lee H.J., Khil T.O. 2013c Supercavitation Rocket System [Journal of the Korea Institute of Military Science and Technology] Vol.16 P.867-880 google doi
  • 8. Kunz R., Boger D., Stinebring D., Chyczewski T. 1999 A Preconditioned Navier-Stokes Method for Two-Phase Flows with Application to Cavitation Prediction P.676-688 google
  • 9. Lindau J.W., Venkateswaran S., Kunz R.F., Merkle C.L. 2003 Multi-phase Computations for Underwater Propulsive Flows [Proceedings of the 16th AIAA Computational Fluid Dynamics Conference] P.4105 google
  • 10. Owis F.M., Nayfeh A.H. 2003 Computations of the Compressible Multi-phase Flow over the Cavitating High-speed Torpedo [Journal of Fluids Engineering] Vol.125 P.459-468 google doi
  • 11. Owis F.M., Nayfeh A.H. 2004 Numerical Simulation of 3-D Incompressible, Multi-phase Flows over Cavitatiing Projectiles [European Journal of Mechanics B/Fluids] P.339-351 google
  • 12. Park W.G. 2016 Compressible Multiphase Flow Analysis for 6DOF Water Entry Behavior and Ventilated Cavitation [11th Asian Computational Fluid Dynamics Conference] google
  • 13. Park W.G., Ha C.T., Merkle C.K. 2009 Multiphase Flow Analysis of Cylinder using a New Cavitation Model [7th International Symposium on Cavitation, CAV2009] P.99 google
  • 14. Rouse H., McNown J.S. 1948 Cavitation and Pressure Distribution: Head Forms at Zero Angle of Yaw google
  • [] 
  • [] 
  • [] 
  • [] 
  • [] 
  • [] 
  • [] 
  • [] 
  • [] 
  • [Fig. 1] Surface pressure coefficient on the hemispherical cylinder
    Surface pressure coefficient on the hemispherical cylinder
  • [Fig. 2] Surface pressure coefficient on the 0-caliber(blunt fore) cylinder
    Surface pressure coefficient on the 0-caliber(blunt fore) cylinder
  • [Fig. 3] Comparison with other numerical result at σ = 0.3
    Comparison with other numerical result at σ = 0.3
  • [Fig. 4] Modeling of underwater vehicle and generated supercavity
    Modeling of underwater vehicle and generated supercavity
  • [Fig. 5] Predicted total drag with various speed conditions
    Predicted total drag with various speed conditions
  • [Fig. 6] Predicted presssure and viscous drag with various speed conditions
    Predicted presssure and viscous drag with various speed conditions
  • [Fig. 7] Predicted supercavites with various speed conditions
    Predicted supercavites with various speed conditions
  • [Fig. 8] Predicted cavity diameter with various speed conditions
    Predicted cavity diameter with various speed conditions
  • [Fig. 9] Predicted cavity length with various speed conditions
    Predicted cavity length with various speed conditions
  • [] 
  • [Fig. 10] Limited thrust system for underwater vehicle
    Limited thrust system for underwater vehicle
  • [Fig. 11] Predicted velocity of underwater vehicle with diffrent thrust conditions
    Predicted velocity of underwater vehicle with diffrent thrust conditions