The Korean Society for Power System Engineering

Journal Archive

Journal of the Korean Society for Power System Engineering - Vol. 22 , No. 6

[ Research Paper ]
Journal of the Korean Society for Power System Engineering - Vol. 22, No. 6, pp. 11-19
Abbreviation: The J. of the Korean Soc. for Power Syst. Eng.
ISSN: 1226-7813 (Print) 2384-1354 (Online)
Print publication date 31 Dec 2018
Received 08 Jun 2018 Revised 22 Nov 2018 Accepted 26 Nov 2018
DOI: https://doi.org/10.9726/kspse.2018.22.6.011

실베스터-전달강성계수법에 의한 실습선 새동백호 추진축계의 비틀림 자유진동 해석
김명준* ; 왕우경** ; 여동준** ; 최명수**,
*대학원생, 전남대학교 기관시스템공학과
**전남대학교 해양기술학부 교수

Torsional Free Vibration Analysis of Propulsion Shafting of Training Ship SAEDONGBAEK by Sylvester-Transfer Stiffness Coefficient Mehtod
Myung-Jun Kim* ; Woo-Gyeong Wang** ; Dong-Jun Yeo** ; Myung-Soo Choi**,
*Graduate student, Department of Power System Engineering, Chonnam National University.
**Professor, Faculty of Marine Technology, Chonnam National University.
Correspondence to : Myung-Soo Choi : Professor, Faculty of Marine Technology, Chonnam National University. E-mail : engine@jnu.ac.kr, Tel : 061-659-7183


Abstract

In this study, the authors examine the propulsion shafting of the training ship SAEDONGBAEK and perform modeling to analyze the torsional free vibration of the shafting. In this paper, the computational algorithm for analyzing the torsional free vibration of the shafting with a reduction gear is formulated by the sylvester-transfer stiffness coefficient method (S-TSCM) that is a recently developed and a powerful method in free vibration analysis. According to the state of the controllable pitch propeller of the shafting and the temperature of the elastic coupling, the torsional free vibration of the shafting is performed by the S-TSCM. The authors examine the changes of the natural frequencies and natural modes which are the results of the torsional free vibration analysis of the shafting.


Keywords: Torsional vibration, Free vibration, Propulsion shafting, Natural frequency, Sylvester-transfer stiffness coefficient method

1. 서 론

2018년 4월 5일 부산(한진중공업 영도조선소)에서 전남대학교의 새 실습선인 ‘새동백호’의 진수식이 있었다. 승무원과 실습생 110명이 승선할 수 있는 새동백호는 총톤수 2,996 ton, 전장 96.45 m, 선폭 15 m, 깊이 7.6 m이고, 최대 속력은 16 노트(knot)이며 최대 항속거리는 8,400 마일(mile)로 예상되며 현재 건조 중이다.

이 연구에서는 새동백호의 추진축계의 구성을 살펴보고, 선박의 안전 운항을 위해 실베스터-전달강성계수법1-2)으로 가변 피치 프로펠러의 상태 및 탄성 커플링의 온도 변화에 따라 축계의 비틀림 진동 고유진동수와 고유모드를 해석한다.

축계의 비틀림 자유 진동을 해석하는 방법으로는 유한요소법3), 전달행렬법4), 전달강성계수법5) 등 다양한 방법들이 존재한다. 이 연구에서는 최근에 개발된 강력한 자유진동 해석 기법인 실베스터-전달강성계수법으로 감속기(reduction gear)를 갖는 추진축계의 비틀림 자유진동을 해석할 수 있는 알고리즘을 기술한다.


2. 새동백호 추진축계
2.1 추진축계의 구성

Fig. 1은 새동백호 추진축계의 개략도이다. 추진축계는 주기관, 탄성 커플링, 감속기, 중간 축, OD(Oil Distributer) 축, 프로펠러 축 그리고 프로펠러로 구성된다. 주기관은 현대중공업의 7H32/40P, 탄성 커플링은 Centa사의 CX-181, 감속기는 Reintjes사의 LAF 5645이고, 축과 프로펠러의 제작사는 WÄRTSILÄ이다. 새동백호 추진축계의 구체적인 스펙을 Table 1에 제시한다.


Fig. 1 
Propulsion shafting of training ship AEDONGBAEK

Table 1 
Specification of propulsion shafting
Main Engine
Type 4 stroke, single acting,
trunk piston
Max. continuous output 3,500 kW
Max. continuous speed 750 rpm
No. of cylinder 7
Cylinder bore 320 mm
Stroke 400 mm
Ratio of connecting rod 0.23
Mean indicated pressure 27.7 bar
Reciprocating mass 167.42 kg/cyl.
Firing order 1-3-5-7-6-4-2
Flexible coupling
Stiffness 0.2856 MNm/rad (50℃)
Permissible vibratory torque 12.5 kNm
Permissible thermal load 0.6375 kW (50℃)
Reduction gear
Type Vertical, off-set type with hyd. clutch
Gear ratio 3.440:1
Length of shaft
Intermediate shaft 6.450 m
Oil distributer shaft 1.300 m
Propeller shaft 6.230 m
Propeller
Type CPP
Number of blades 4
Diameter 3.3 m
Weight 5,000 kg

2.2 추진축계의 모델링

새동백호 추진축계의 비틀림 자유진동을 해석하기 위하여 추진축계를 다수의 회전 질량과 비틀림 스프링 그리고 기어쌍과 기어축으로 모델링한 것이 Fig. 2이다. Fig. 2에서 절점 1은 캠 구동부이고 절점 2는 플랜지이다. 주기관은 7기통이므로 절점 3부터 절점 9는 순차적으로 크랭크 축과 커넥팅로드와의 연결점이다. 절점 10은 기관의 플라이 휠과 탄성 커플링의 결합점이고, 절점 11은 탄성 커플링과 감속기의 회전 질량과의 결합점이다. 절점 12는 클러치이고, 절점 13은 상부 기어와 하부 기어로 구성된 기어쌍이고, 절점 17은 프로펠러에 상당한다. 절점 14와 절점 15 사이에 중간 축이 있고, 절점 15와 절점 16 사이에 OD 축, 절점 16과 절점 17 사이에 프로펠러 축이 있다.


Fig. 2 
Computational model


3. 비틀림 자유진동 해석 알고리즘

유한요소법, 전달강성계수법, 실베스터-전달강성계수법으로 Fig. 3과 같이 감속기가 있는 축계의 비틀림 자유진동을 해석하는 알고리즘을 기술한다.


Fig. 3 
Analytical model

해석대상 계는 회전 질량, 비틀림 스프링, 2개의 회전 질량J^up,J^down으로 구성된 기어 쌍 및 비틀림 강성 K^3를 갖는 기어 축으로 구성된다.

3.1 유한요소법3)

절점 1, 2, 4, 5에 각각 하나의 회전 질량이 존재하므로, 이들의 질량관성모멘트(Ji)를 가지고 크기 1×1의 요소 질량관성행렬을 구한다.

Jie=Ji i=1,2,4,5(1) 

절점 3의 기어 쌍은 상부와 하부에 각각 질량관성모멘트 J^upJ^down인 2개의 회전 질량으로 구성되므로, 상부 회전 질량을 중심으로 등가화하여 질량관성모멘트를 구하면

J3=J^up+r2×J^down(2) 

여기서 r은 아래 식과 같이 상부 기어와 하부 기어의 회전 속력비이다.

r=ΩdownΩup(3) 

따라서 식 (3)의 등가 질량관성모멘트(J3)를 식(1)에 적용하면 절점 3의 요소 질량관성행렬을 구할 수 있다.

인접한 두 절점을 연결하는 비틀림 스프링 1, 2, 4는 각각의 스프링상수(Ki)를 이용하여 식 (4)와 같은 요소 강성행렬을 유도할 수 있다.

Kie=Ki-Ki-KiKi i=1,2,4(4) 

유한요소해석3)으로부터 하부 기어(J^down)와 회전 질량(J4)을 연결하는 기어 축의 요소 강성행렬은 기어 축의 비틀림 강성이 K^i인 경우 다음과 같다.

Kie=r2K^irK^irK^iK^i i=3(5) 

요소행렬들을 모두 조립하면, Fig. 3의 계를 다음 식과 같은 전체 질량관성행렬(J)과 전체 강성행렬(K)로 나타낼 수 있다.

J=J100000J200000J300000J400000J5(6) 
K=K1-K1000-K1K1+K2-K2000-K2K2+r2K^3rK^3000rK^3K^3+K4-K4000-K4K4(7) 

전체 계의 운동방정식을 행렬로 나타내면

Jθ¨+Kθ=0(8) 

여기서 ϴ는 각 절점의 각변위(ϴi)로 구성되는 전체 각변위벡터(ϴ = {ϴ1, ϴ2, ϴ3, ϴ4, ϴ5}T)이고, θ¨는 전체 각가속도벡터이다.

식 (8)을 일반 고유치 문제로 나타내면

Kψ=λJψ(9) 

가 되고, 여기서 λ는 고유치로서 각고유진동수 ωn의 제곱이고, ψ는 고유벡터이다. 식 (9)의 우변 항을 좌변으로 넘기면

K-λJψ=0(10) 

따라서 진동수 방정식은

detK-λJ=0(11) 

가 되어, 식 (11)에서 고유치(λi)를 구하고, 식 (10)에서 고유벡터(ψi)를 구할 수 있다.

만약 계의 전체 자유도가 커지면 자코비 법(Jacobi method), 역방향 반복법(inverse iteration method), 축소공간 반복법(subspace iteration method) 등의 고유치 문제 해법6)으로 고유치와 고유벡터를 구하는 것이 효과적이다.

3.2 전달강성계수법5)
3.2.1 동강성계수의 정의

절점 i 좌측 및 우측의 토크(T¯i, Ti)와 절점 i의 각변위(ϴi) 사이의 관계를 다음 식과 같이 동강성계수(S¯i, Si)를 사용하여 정의한다.

T¯i=S¯i×θi(12) 
Ti=Si×θi(13) 

기호 TS 위에 “―”를 붙인 것은 절점의 좌측을 뜻하고, “―”가 없는 것은 절점의 우측을 뜻한다. 그리고 하첨자 “i”는 i번째 절점을 의미한다.

3.2.2 동강성계수의 전달

Fig. 4와 같이 절점 i에 질량관성모멘트가 Ji인 회전 질량이 존재하는 경우, 절점 i의 좌측에서 절점 i의 우측으로 동강성계수의 전달식은

Si=S¯i+-ω2Ji i=1,2,4,5(14) 

된다.


Fig. 4 
Rotating masses and torsional spring

Fig. 4에서 절점 i와 절점 i+1 사이에 스프링상수가 Ki인 비틀림 스프링이 존재하는 경우, 절점 i와 절점 i+1의 각변위 사이의 관계는

θi=Vi×θi+1 i=1,2,4(15) 

가 되고, 여기서

Vi=Ki/Gi, Gi=Si+Ki i=1,2,4(16) 

이고, 절점 i의 우측에서 절점 i+1의 좌측으로 동강성계수의 전달식은

S¯i+1=Ki×1-Vi i=1,2,4(17) 

이 된다.

Fig. 5에 나타낸 절점 i의 기어 쌍은 상부 회전 질량을 중심으로 한 등가 질량관성모멘트(Ji=J^up+r2+J^down)를 식 (14)에 적용하면 기어 쌍에 대해 동강성계수를 전달할 수 있다.


Fig. 5 
Gear pair and gear shaft

Fig. 5에서 하부 기어(J^down)와 다음 회전원판(Ji+1) 사이에 연결된 기어 축의 비틀림 스프링상수가 K^i라면, 이 스프링 요소의 좌측과 우측의 토크와 각변위 사이의 관계는 다음 식과 같다.

TiT¯i+1=-r2K^i-rK^irK^iK^iθiθi+1 i=3(18) 

식 (18)Ti식 (13)을 대입하고 정리하면

θi=Vi×θi+1 i=3(19) 

가 되고, 여기서

Vi=-rKi/Gi,Gi=Si+r2Ki i=3(20) 

식 (18)T¯i+1식 (12)의 “i” 대신에 “i+1”을 대입한 식과 식 (19)을 대입하여 정리하면, 기어 축의 좌측에서 우측으로 동강성계수를 전달하는 아래 식을 유도할 수 있다.

S¯i+1=K^i×1+rVi i=3(21) 

Fig. 3의 계산 모델에 상기 식 (14), 식 (17), 식 (21)을 요소 별로 순차적으로 적용하면 최종적으로 마지막 절점인 절점 5의 우측의 동강성계수(S5)를 구할 수 있다.

3.2.3 고유진동수

회전축계는 구속되지 않으므로 마지막 절점의 경계조건은 자유가 된다. 따라서 절점 5 우측의 토크 및 절점 5의 각변위는

T5=0, θ50(22) 

가 되고, 식 (13)으로부터 아래와 같은 전달강성계수법의 진동수 방정식을 유도할 수 있다.

S5ω=0(23) 

전달강성계수법은 식 (23)에 임의의 각진동수(ω~)를 고유진동수 ω 대신에 대입하고 S5의 부호 변화를 통해 고유진동수를 찾아낸다.

한편 전달강성계수법은 식 (23)에서 반대칭 극이 발생할 때, 이것을 고유진동수로 오인하는 문제가 발생한다. 이 문제를 해소하기 위하여 식(24)를 이용하여 고유진동수를 구하지만, 초기 검색 진동수의 간격(Δω0)이 인접한 두 고유진동수 값의 차이보다도 커지면 고유진동수를 찾지 못하는 일이 발생한다. 우리는 인접한 두 고유진동수의 차이를 알지 못하므로, 상기 Δω0를 아주 작은 값으로 계산해야 하며, 이 경우에도 고유진동수를 빼 먹었는지 확인이 되지 않을 뿐더러, Δω0가 작아질수록 계산 시간도 크게 증가하는 문제가 전달강성계수법의 결정적인 취약점이다.

i=14Giω×S5ω=0(24) 
3.2.4 고유모드

마지막 절점의 각변위를 1로 두고(ϴ5 = 1), 나머지 절점의 각변위는 요소에 따라 식 (15) 또는 식 (19)를 사용하면 모든 절점의 각변위를 순차적으로 구할 수 있다.

상부 기어와 하부 기어의 각 변위를 각각 ϴupϴdown이라 하면, 이들 사이의 관계는

θdown=-r×θup=-r×θi i=3(25) 

이므로, 하부 기어의 각변위(ϴdown)는 상부 기어의 각변위(ϴup = ϴ3)를 구한 후에 계산한다.

계산한 각변위 중에서 절대치가 가장 큰 값을 찾고 이를 기준으로 정규화하여 최종적으로 고유모드를 구한다.

3.3 실베스터-전달강성계수법1-2)

식 (9)의 일반 고유치 문제는 아래와 같은 직교성을 만족한다.

ΦTJΦ=I=diag1,1,1,1,1ΦTKΦ=Λ=diagλ1,λ2,λ3,λ4,λ5(26) 

여기서 Φ = [ψ1, ψ1, ψ2, ψ3, ψ4, ψ5]는 1차부터 5차까지의 고유벡터로 구성된 모드행렬이다.

행렬 K - λJ의 앞과 뒤에 모드행렬을 식 (27)과 같이 곱하면

ΦTk-λJΦ=Λ-λI(27) 

대각행렬이 되고 대각항의 요소들은 (λ1 - λ), (λ2 - λ), …, (λ5 - λ)가 된다. 따라서 (K - λJ) 행렬의 고유치는 (λi - λ)가 된다. 그리고 이들 고유치 중에서 음수의 개수는 JK 행렬의 고유치(λi) 중에서 λ보다 작은 고유치의 개수와 같다.

(K - λJ) 행렬의 앞뒤에 아래에 나타낸 행렬 L과 전치행렬인 LT를 곱하면 다음과 같은 대각행렬이 된다.

LK-λJ LT=diagG1,G2,G3,G4,G5(28) 

여기서

L=1V1V1V2V1V2V3V1V2V3V401V2V2V3V2V3V4001V3V3V40001V400001(29) 

실베스터의 관성이론으로부터 (K - λJ)행렬의 고유치 중에서 음의 고유치의 개수는 행렬 diag[G1, G2, G3, G4, G5]의 대각요소의 음수의 개수와 동일하다. 그리고 식 (28)식 (29)Gi, Vi, Si는 3.2절에 소개한 전달강성계수법의 Gi, Vi, Si와 동일하다.

따라서 Fig. 3의 해석 모델에 대해 실베스터-전달강성계수법을 이용하여 고유진동수를 구하는 전산 알고리즘은 다음과 같다.

1) 임의의 각진동수 ω~를 정한다.

2) 요소별 전달강성계수법의 전달식을 적용하여 모든 절점의 GiSi를 계산한다.

3) G1부터 G4 그리고 S5 중에서 0보다 작은 값을 가지는 변수의 개수(m)를 센다.

4) 상기 3)의 개수(m)는 검색 진동수 ω~ 보다 작은 고유각진동수의 개수이다. 진동수 검색 범위 안의 여러 ω~에 대해 반복적으로 m을 구하면 빠뜨리지 않고 모든 고유각진동수를 구할 수 있다.


4. 수치 계산 결과 및 검토
4.1 가변 피치 프로펠러 상태에 따른 추진축계의 비틀림 자유진동 해석

새동백호는 가변 피치 프로펠러가 장착되어 있으므로 프로펠러의 상태에 따라 프로펠러의 질량관성모멘트의 값이 달라진다.

탄성 커플링의 온도를 50℃로 설정하고, 프로펠러가 zero-pitch와 full-pitch일 때, 실베스터-전달강성계수법으로 추진축계의 고유진동수를 각각 계산하였다. 강체모드는 생략하고, 1차부터 5차까지의 고유진동수를 계산한 결과를 Table 2에 제시하였다.

Table 2 
Natural frequencies of computational model (temperature of coupling : 50℃) [rpm]
Order Zero-Pitch Full-Pitch
1 260 236
2 826 820
3 4,071 4,071
4 5,105 5,101
5 9,557 9,557

Table 2로부터 프로펠러의 상태에 따른 고유진동수를 비교해 보면, 프로펠러가 zero-pitch에서 full-pitch로 변경되면 1차 고유진동수는 약 9% 정도 감소되고, 2차 고유진동수는 약 1% 정도 감소되었다. 그러나 3차 이상의 고유진동수에는 거의 변화가 없었다.

zero-pitch 및 full-pitch 2가지 경우에 대해 고유모드를 각각 계산하였다. Fig. 6에는 zero-pitch (o 표시)와 full-pitch(x 표시) 상태의 1차 고유모드를 제시하였다. 1차 고유모드를 살펴보면 탄성 커플링(절점 10과 절점 11사이)에서 진동을 하지 않는 절이 발생하고 있음을 알 수 있었다. 프로펠러의 상태에 따른 고유모드의 값을 비교해 보면 감속기 상부 기어의 좌측(절점 1에서 절점 13)에서 차이가 있음을 알 수 있었다.


Fig. 6 
1st natural mode

Fig. 7은 2차 고유모드를 계산한 결과이다. 2차 고유모드는 탄성 커플링(절점 10과 절점 11 사이)과 프로펠러 축(절점 16과 절점 17 사이)에서 절이 발생하고 있음을 알 수 있었다. 그리고 프로펠러의 피치 값에 무관하게 2차 고유모드는 거의 동일함을 알 수 있었다.


Fig. 7 
2nd natural mode

그리고 저자들은 프로펠러의 상태에 따른 3차부터 5차까지의 고유모드를 구한 후, 비교해 본 결과 거의 차이가 없음을 확인하였다.

4.2 탄성 커플링의 온도 변화에 따른 추진축계의 비틀림 자유진동 해석

새동백호 추진축계의 탄성 커플링은 온도에 따라 비틀림 강성이 변화한다. 커플링 제작사로부터 제공 받은 자료에 의하면 절점 10과 11 사이에 있는 커플링의 비틀림 스프링상수(K10)의 값은 30℃에서 0.3162 MNm/rad이고, 70℃에서 0.2652 MNm/rad이었다.

프로펠러가 full-pitch이고, 커플링의 온도가 30℃와 70℃일 때, 실베스터-전달강성계수법으로 추진축계의 고유진동수를 각각 계산하였다. 강체모드는 생략하고, 1차부터 5차까지의 고유진동수를 계산한 결과를 Table 3에 제시하였다. 계산된 고유진동수를 비교해 보면, 1차 고유진동수는 약 4%의 변화가 있었고, 2차 고유진동수는 약 5%의 변화가 있었다. 그리고 탄성 커플링의 온도가 상승할 때 고유진동수가 감소함을 알 수 있었다. 그러나 3차, 4차, 5차 고유진동수는 거의 변화가 없었다.

Table 3 
Natural frequencies of computational model (state of propeller : full-pitch) [rpm]
Order 30℃ 70℃
1 241 232
2 844 803
3 4,071 4,070
4 5,101 5,101
5 9,557 9,557

30℃와 70℃, 2가지 경우에 대해 고유모드를 각각 계산하였다. Fig. 8에는 탄성 커플링의 온도인 30℃(o 표시)와 70℃(x 표시)의 1차 고유모드를 제시하였다. 1차 고유모드를 살펴보면 탄성 커플링(절점 10과 절점 11사이)에서 절이 발생하고 있음을 알 수 있었다. 탄성 커플링의 온도에 따른 고유모드의 값을 비교해 보면 기어 박스에 해당하는 절점 11과 절점 14 사이에서 약간 차이가 있음을 알 수 있었다.


Fig. 8 
1st natural mode

Fig. 9은 2차 고유모드를 계산한 결과이다. 2차 고유모드는 탄성 커플링(절점 10과 절점 11 사이)과 프로펠러 축(절점 16과 절점 17 사이)에서 절이 발생하고 있음을 알 수 있었다. 그리고 탄성 커플링의 온도와 무관하게 2차 고유모드는 거의 동일함을 알 수 있었다.


Fig. 9 
2nd natural mode

그리고 저자들은 탄성 커플링의 온도에 따른 3차부터 5차까지의 고유모드를 구한 후, 비교해 본 결과 거의 차이가 없음을 확인하였다.


5. 결 론

이 연구에서는 전남대학교 실습선 새동백호의 추진축계의 구성을 살펴보고, 최근에 개발된 강력한 자유진동 해석 기법인 실베스터-전달강성계수법으로 감속기를 갖는 비틀림 축계의 자유진동을 해석할 수 있는 전산 알고리즘을 기술하였다. 그리고 실베스터-전달강성계수법으로 가변 피치 프로펠러의 상태 및 탄성 커플링의 온도 변화에 따라 축계의 비틀림 진동 고유진동수와 고유모드를 계산하였다.

프로펠러가 zero-pitch에서 full-pitch로 변경되면 1차 고유진동수는 약 9% 정도 감소되고, 2차 고유진동수는 약 1% 정도 감소되었다. 그러나 3차 이상의 고유진동수에는 거의 변화가 없었다. 그리고 1차 고유모드는 프로펠러의 상태에 따라 약간 달랐지만, 2차 이상의 고유모드는 거의 차이가 없었다.

탄성 커플링의 온도가 30℃에서 70℃로 증가할 때, 1차 및 2차 고유진동수는 약 4~5%의 감소하였다. 그러나 3차 이상의 고유진동수에는 거의 변화가 없었다. 고유모드는 1차 고유모드는 프로펠러의 상태에 따라 약간 달랐지만, 2차 이상의 고유모드는 거의 동일하였다.


References
1. M. S. Choi, T. Kondou, and Y. Bonkobara, (2012), "Development of Free Vibration Analysis Algorithm for Beam Structures by Combining Sylvester’s Inertia Theorem and Transfer Stiffness Coefficient Method", Journal of Mechanical Science and Technology, 16(1), p11-19.
2. M. S. Choi, T. Kondou, J. H. Byun, and D. J. Yeo, (2015), "Flexural Free Vibration Analysis of Axisymmetric Annular Plates Using Sylvester-Transfer Stiffness Coefficient Method", Journal of The Korean Society for Power System Engineering, 19(6), p60-67.
3. M. Lalanne, and G. Ferraris, (1998), "Rotor Dynamics Prediction in Engineering (2nd edition)", John Wiley and Sons, Chichester, p203-227.
4. J. S. Rao, (1998), "Rotor Dynamics (3rd edition)", New Age International, New Delhi, p9-24.
5. D. H. Moon, M. S. Choi, and J. M. Sim, (1997), "Torsional Vibration Analysis of Shaft System Using Transfer Dynamic Stiffness Coefficient", Journal of The Korean Society for Power System Engineering, 1(1), p91-97.
6. K. J. Bathe, (1996), "Finite Element Procedures", Prentice Hall, New Jersey, p887-978.