Modal Analysis for the Rotating Cantilever Beam with a Tip Mass Considering the Geometric Nonlinearity

기하학적 비선형성을 고려한 종단 질량을 갖는 회전하는 외팔보의 모달 분석

Abstract

In this paper, a new dynamic model for modal analysis of a rotating cantilever beam with a tip-mass is developed. The nonlinear strain such as von Karman type and the corresponding linearized stress are used to consider the geometric nonlinearity, and Euler-Bernoulli beam theory is applied in the present model. The nonlinear equations of motion and the associated boundary conditions which include the inertia of the tip-mass are derived through Hamilton's principle. In order to investigate modal characteristics of the present model, the linearized equations of motion in the neighborhood of the equilibrium position are obtained by using perturbation technique to the nonlinear equations. Since the effect of the tip-mass is considered to the boundary condition of the flexible beam, weak forms are used to discretize the linearized equations. Compared with equations related to stiffening effect due to centrifugal force of the present and the previous model, the present model predicts the dynamic characteristic more precisely than the another model. As a result, the difference of natural frequencies loci between two models become larger as the rotating speed increases. In addition, we observed that the mode veering phenomenon occurs at the certain rotating speed.

1. 서 론

회전하는 외팔보는 전통적인 산업분야인 가스터빈이나 풍력발전기, 헬리콥터 등에 사용되는 블레이드 뿐만 아니라 비교적 근래에 등장하기 시작한 드론이나 곤충형 로봇의 날개와 같이 다양한 분야에 응용되기 때문에 현재에도 많은 연구자들의 관심을 받아 오고 있다. 설계자들은 에너지 효율과 출력을 높이기 위해서 회전 구조물을 경량화, 고속화 하고 동시에 정확히 제어하기를 원한다. 이를 위해서 회전운동의 영향을 받는 구조물의 정확한 동적특성을 예측할 수 있는 모델링 방법이 필요하며 몇 가지 모델링 방법이 오래도록 사용되어 오고 있다.

선형 변형률 이론에 기반한 모델링 방법(1~5)에서는 회전운동에 의한 외팔보의 강성변화 효과를 고려하기 위해서 구심력에 의해 발생하는 인장력을 잠재에너지에 직접적으로 내포시켰다. 이러한 방법을 사용하여 Hoa(1)와 Huang 등(2)은 회전속도에 따라 외팔보의 고유진동수가 변하는 현상을 관찰하였다. 그리고 Banerjee와 Kennedy(3)는 더 높은 회전속도에 대한 고유진동수 변화를 관찰하여 코리올리 효과가 회전하는 외팔보의 모달 특성에 미치는 영향을 분석하였다. Kim과 Kang(4)은 축방향 운동을 무시하고 회전축 방향으로 운동하는 면외방향(out-of-plane) 변형만 고려하여 종단 질량이 있는 회전하는 외팔보의 진동을 분석 하였다. 진동 해석 교재(5)에도 사용되어온 이 방법은 회전하는 외팔보의 동적 모델을 쉽고 직관적으로 이해할 수 있는 장점이 있다.

축변형 변수 대신 인장변형 변수를 사용하여 회전하는 외팔보의 동적모델을 구축하는 방법(6~10)은 앞선 방법보다 체계적으로 선형화된 운동방정식을 구할 수 있다. Kane 등(6)에 의해서 소개된 이 방법은 Yoo 등(7)에 의해서 명확하고 간결하게 정식화 하는 방법이 만들어졌다. 미소공간에서 인장변형 변수는 고전적인 축변위와 횡변위 제곱의 합으로 정의되고 축방향으로 임의의 위치 x까지 적분을 하여 인장변형 변수의 기하학적 제한조건을 구한다. 따라서, 운동에너지가 미적분 형태로 나타나지만, 이들 연구들에서는 최종적으로 미적분 운동방정식을 선형화하여 동적 특성을 연구하였다. 결과적으로 선형 변형률 이론에 기반한 모델링 방법(1~5)과 운동방정식은 형태가 비슷하지만, 회전에 의한 강성변화 항을 인위적이지 않고 좀 더 자연스럽게 고려할 수 있다.

비선형 변형률 이론을 사용하여 회전하는 외팔보의 좀 더 정확한 동적 응답을 예측하기 위한 시도가 계속되어 오고 있다(11~15). Simo와 Vu-quoc(11)는 회전에 의한 강성변화가 빔의 기하학적 비선형성에 의거 한다는 관점으로 비선형 지배방정식을 유도하였다. Pesheck(12)는 본-카르만 변형률(von-Karman strain)을 사용하여 축변위와 회전 면외방향(out-of-plane) 변위만을 고려하여 비선형 운동방정식을 유도하고 섭동법으로 정적 평형상태에서 선형화된 운동방정식을 구하였다. 이후, Sharf(13)와 Wang(14)은 축변위와 면내방향(in-plane) 굽힘 변위가 연성되어 있는 비선형 방정식을 유도하여 선형모델의 시간 응답과 비교하여 차이점을 관찰하였다. 비교적 최근에 Kim 등(15)은 선형모델과 비선형 모델의 동적 특성은 회전속도가 높아질수록 차이가 커진다는 것을 해석적으로 증명하였다. 이와 같이, 비선형 모델이 기존의 선형 모델보다 더 정확하지만, 아직까지 기하학적 비선형성을 고려한 종단 질량이 있는 회전하는 외팔보의 동적 특성을 분석한 연구는 찾기 어렵다.

이 논문에서는 종단질량을 갖는 회전하는 외팔보의 모달 해석을 위한 새로운 모델링 방법을 제안한다. 제안하는 모델의 운동 방정식을 유도하기 위해서 유연체 외팔보의 자유단에 집중질량에 대한 운동에너지를 고려하였다. 유연체 외팔보는 비선형 본-카르만 변형률을 사용하고 오일러-베르누이 빔 이론(Euler-Bernoulli beam theory)으로 모델링 하였다. 유도된 비선형 방정식은 섭동법(perturbation method)을 사용하여 정상 상태에서 정적 평형 방정식과 평형 위치 근방에서의 동적 운동방정식을 구하였다. 그리고 현재 제안하는 모델과 기존 선형 모델의 회전속도 및 집중질량에 따른 고유진동수 변화를 비교하고 모드 근접(mode veering) 현상을 관찰하였다.

 

2. 운동방정식

실제 회전하는 블레이드는 외팔보와 같이 단순한 형상은 아니지만, 이 연구에서는 회전운동이 유연체 외팔보에 미치는 영향성에 집중하기 위해 블레이드를 단순한 외팔보로 가정하였다. 수학적인 모델을 수립하기 위한 회전하는 외팔보의 개략도를 Fig. 1에 나타내었다. 반지름 𝑎를 갖는 강체 허브는 원점 O'를 중심으로 회전하며, 유연체 외팔보는 강체허브의 점 O에 부착되어 있다. 외팔보는 길이 L, 단면적 A, 면적 관성 모멘트는 Iz이며, 밀도 ρ, 영률 E를 갖는 정사각 단면 보이다. 외팔보의 변위 벡터는 공간 변수 x와 시간 변수 t의 함수이며, 축변위는 u(x,t), 회전운동의 접선방향으로 변형하는 면내(in-plane)굽힘 변위는 v(x,t)로 나타냈다. 이 논문에서는 논의를 단순하게 하기 위해서 면외 굽힘 변위 w는 고려하지 않았다.

Fig. 1Rotating cantilever beam considering tip-mass with rotating speed Ω

유연체 외팔보의 변형 전 임의의 위치를 P라고 하면, 변형 후 위치는 P* 로 표현 할 수 있다. 이때, 임의의 점 P*에 대한 속도 벡터는 다음과 같이 표현 할 수 있다.

그리고 유연체 외팔보의 자유단에 위치한 집중 질량 mt의 속도 벡터는,

여기서 하첨자 L은 자유단 x = L을 나타내며, uL = u(L,t) 과 vL = (L,t)이다

기하학적인 비선형성을 고려한 유연체 외팔보를 모델링하기 위해 이 논문에서는 본-카르만 변형률 이론을 사용하고 외팔보의 길이가 두께에 비해서 충분히 길다고 가정하여 전단 변형과 회전 관성을 무시한 오일러-베르누이 빔 이론을 적용한다. 본-카르만 변형률 εx 와 이에 대한 선형화된 응력 σx 는 다음과 같이 표현 할 수 있다.

비선형 변형률과 선형화된 응력을 사용하는 이유는 둘 다 비선형을 사용하여 유도한 운동방정식보다 훨씬 단순하기 때문이다. 또한, 모달 특성 분석을 위해 정상상태에서 구한 선형화된 운동방정식은 두 경우 모두 동일하다. 한편, 운동방정식을 구하기 위해서 다음과 같은 해밀턴의 원리를 사용하였다.

여기서 δK 와 δP, δW는 각각 운동에너지와 잠재에너지, 비-보존력의 변분을 나타내며,

와 같다. 여기서, pv는 횡방향으로 가해지는 분포력이다. 변분된 에너지 식 (6)들을 식 (5)에 대입하여 운동방정식을 유도하면,

여기서 (˙)은 시간 t에 대한 미분 연산자를 나타낸다. 기하학적 비선형성을 고려한 효과는 식 (8)에서 좌측의 마지막 항 −∂/∂x[ EA(∂u/∂x(∂v/∂x) ]으로 나타난다. 한편, 이에 대한 경계조건은 다음과 같다.

종단 질량은 자유단 x = L 에서 식 (10)과 같이 축방향 힘으로 작용하며, 면내 방향 경계조건에서는 식 (12)와 같이 전단력을 발생시키는 것을 알 수 있다.

유도된 비선형 운동방정식을 고유치 문제로 정식화하기 위해서 회전속도가 일정하고 외력이 없을 때, 섭동법을 사용하여 평형 위치 근방에서 선형화된 운동방정식을 구할 수 있다. 섭동법에 따르면, 축방향과 횡방향의 해는 정적 평형 위치 근방에서 다음과 같이 정적해와 동적해로 나눌 수 있다(12).

여기서 u0 와 v0 는 시간 t에 대해 무관하며, 등속회전 시 정적 평형 위치에서의 해다. 그리고 ∆u와 ∆v는 평형 위치 근방에서의 섭동 변수이다. 식 (13)과 식 (14)를 비선형방정식 (7)과 (8)에 대입하면, 정적 평형방정식은 다음과 같이 유도된다.

그리고 식 (15)에 대한 축방향 경계조건은 식 (17), (18)과 같다. 횡방향의 정적해 v0는 식 (16)과 같이 0이므로 그에 대한 경계조건은 생략한다.

정적 평형 방정식 (15)를 경계조건 식 (17), (18)을 이용해서 축방향 변위 u0의 엄밀해를 구하면,

여기서, 이다. 그리고 이에 대한 평형 위치근방에서 축방향 운동 ∆u 와 면내방향 운동 ∆v 로 자유진동하는 운동방정식은 다음과 같다.

여기서, 기호를 간단하게 표현하기 위해서 섭동변수의 기호 ∆은 생략 하였다. 식 (21)로 주어지는 면내 방향 운동방정식의 Fx는 식 (8)의 비선형 항이 정적 평형 위치 근방에서 선형화된 항이다. 이 항의 물리적인 의미는 회전에 의한 원심력 때문에 발생하는 축방향 힘이다. 이 힘은 식 (19)를 x에 대해서 미분하여 다음과 같이 유도 할 수 있다.

 

3. 이산화 방정식 유도

갤러킨 방법을 통해 근사해를 구하기 위해서 변분된 운동방정식이 필요하다. 축방향과 면내방향에 대한 가중함수는 각각 와 로 표현하였다. 종단 질량의 관성 효과가 내포된 경계조건을 고려하기 위해서 약한 형태(weak form)를 사용한다. 선형화된 운동방정식 (20), (21)에 가중함수를 곱하고 L부터 0까지 적분 한 후, 탄성항에 대해 부분적분을 수행한다. 그리고 경계조건 식 (9) ~ (12)을 적용하여 약한 형태의 변분된 운동방정식을 유도하면,

갤러킨 방법을 사용하기 위한 근사해는 공간변수와 시간변수가 분리된 형태의 선형 조합으로 이루어져 있다. 이들 근사 함수들은 시도함수라고 불리며, 다음과 같이 표현 할 수 있다.

여기서 N은 기저함수의 개수이며, 시간에 대한 함수 와 은 계산 되어야 할 임의의 함수이다. Uj와 Vn는 정지된 외팔보의 모드 함수를 사용하였다. 시도함수에 대한 가중함수는 다음과 같다.

여기서 와 는 임의의 시간에 대한 함수이다. 식 (25)와 (26)을 식 (23)과 (24)에 대입하면, 유한 차원에서 이산화된 운동방정식을 구할 수 있다.

여기서,

식 (27)과 (28)은 서로 연성되어 있으므로 해를 구하기 위해서는 하나의 행렬로 구성해야 한다. 이산화된 운동방정식 (27)과 (28)을 행렬-벡터 형태로 표현하면 다음과 같다.

여기서, M는 질량 행렬, G는 자이로스코픽 행렬이며, K는 강성행렬이다. 마지막으로 T는 계산되어야 할 미지수 벡터로 다음과 같다.

 

4. 수치 결과

4.1 모델링 비교

이 절에서는 기존의 선형 모델링 방법(1~10)을 통해 구한 운동방정식과 제안하는 모델의 운동방정식을 비교하고 고유진동수의 차이점을 조사한다. 회전속도가 등속인 상태에서 기존의 선형 모델링 방법을 통해서 유도한 운동방정식은 형태상 현재 모델의 운동방정식 (20), (21)과 비슷하다. 하지만, 식 (22)로주어지는 종단 질량이 고려된 원심력에 의한 축방향 힘 Fx는 다르다. 기존의 선형 모델링 방법으로 이 힘을 유도하면,

여기서 기호 Fs는 이 연구에서 사용한 Fx와 구별하기 위해서 도입하였다. 식 (32)는 참고문헌(4)에서 사용한 ‘회전에 의한 원심력’ 항과 동일하다. 집중 질량 mt가 없다면, 기존문헌들(1~3,6~9)에서 사용한 원심력에 의한 축방향 힘과 동일하다.

이 연구에서 제안하는 방법과 기존 연구에서 사용한 방법의 차이는 원심력에 의한 축방향 힘을 계산할 때 나타난다. 기존 연구들에서는 식 (15)와 (18)로 주어지는 정적 평형 방정식과 경계조건에서 u0가 매우 작다고 가정하여 무시하고 미분방정식을 풀면 식 (32)를 유도할 수 있다. 현재 모델링 방법은 축방향 힘 Fx를 계산할 때, 축방향 변위 u0를 무시하지 않았기 때문에 현재 모델링 방법이 기존의 선형 모델링 방법보다 좀 더 정확하게 고유진동수를 예측 할 수 있다.

현재 모델과 기존 모델의 고유진동수의 차이를 확인하기 전에 현재 모델의 정확한 수치해를 구하기 위한 수렴성 조사를 수행한다. 논의를 간편하게 하기 위해서 다음의 무차원화된 변수를 도입하였다.

여기서 제어할 수 있는 무차원화 변수는 종단질량 μt, 회전속도 γ, 강체허브의 반경 δ, 세장비 α이다. 수렴성 조사에 사용된 무차원화 변수의 값은 μt = 1, γ = 10, δ = 1, α = 70이다. 기저 함수의 개수 N을 증가시키면서 구한 무차원화된 진동수 T𝜔n의 수렴성을 조사한 자료는 Table 1에 나타냈다. 기저함수를 9개 사용했을 때와 10개 사용했을 때의 고유진동수의 상대적인 오차는 0.02 % 근방으로써, 이 연구에서는 이를 충분히 수렴하였다고 판단하여 앞으로의 분석에서 총 10개의 기저 함수를 사용하였다. 참고로 이러한 방법은 기존문헌(10)에서 사용되었다.

Table 1Convergence tests of the dimensionless natural frequency of the axial and chordwise motion, T𝜔n, when μt = 1, γ = 10, δ = 1, α = 70

무차원화된 회전속도( 0 ≤ γ ≤ 50 )에 따라 현재 모델과 기존 모델의 동적 특성 차이를 조사하기 위해서 μt = 1, δ = 1, α = 70일 때, 축방향과 면내방향의 1차부터 8차까지 고유진동수에 대한 변화를 Fig. 2에 그렸다. 회전속도가 증가 할수록 현재 모델과 기존 모델의 고유진동수 차이가 점점 커지는 것을 관찰할 수 있다. 중요한 점은 현재 모델의 고유진동수는 기존 모델 보다 좀 더 높게 평가된다. 그 이유는 앞에서 설명한 현재 모델과 기존 모델의 원심력에 의한 축방향 힘을 계산하는 방법의 차이 때문이다. 현재 모델의 Fx는 축방향 변형 u0를 고려하기 때문에 이를 무시한 기존 모델의 Fs보다 좀 더 크게 계산된다. 따라서, Fx에 영향을 받는 현재 모델의 굽힘 강성은 Fs에 영향을 받는 기존 모델의 굽힘 강성보다 높기 때문에 회전 속도가 높아질 수록 현재 모델의 고유진동수가 더 높게 계산 된다.

Fig. 2Comparison of eight lowest dimensionless natural frequencies T𝜔n for the variation of the dimensionless rotating speed g between present model (solid line) and previous model (dotted line), when μt =1, δ = 1, α = 70)

두 모델의 차이점은 회전속도 증가뿐만 아니라 종단 질량의 크기 변화에 따라서도 나타난다. 무차원화된 종단 질량의 크기 변화( 0 ≤ μt ≤ 1 )에 따른 두 모델의 고유진동수 차이를 조사하기 위해서 γ = 50, δ = 1, α = 70일 때, 1차부터 4차까지 현재 모델과 기존 모델의 고유진동수의 궤적 변화는 Fig. 3에 나타냈다. Fig. 3에서 기호 A는 축방향 고유진동수의 궤적을 표시하며, B는 면내방향 굽힘 진동의 고유진동수 궤적을, 그리고 하첨자는 그에 대한 차수를 나타낸다. 이 조사에서도 현재 모델의 고유진동수는 원심력에 의한 축방향 힘을 계산 하는 방법 차이에 의거하여 기존 모델의 고유진동수 보다 높게 예측된다. 한 가지 흥미로운 점은 축방향 고유진동수 A1과 A2에서 무차원화된 종단 질량의 크기 0.2부근을 기점으로 두 모델의 차이가 눈에 띄게 점점 커지지만, 면내방향의 굽힘 고유진동수 B1과 B2에 대한 두 모델의 차이는 그렇지 않다. 이것은 종단 질량이 커질수록 현재 모델의 축방향 고유진동수에 미치는 영향이 더 크다는 것을 시사한다.

Fig. 3Comparison of four lowest dimensionless natural frequencies T𝜔n for the variation of the dimensionless tip-mass ratio μt between present model (solid line) and previous model (dotted line), when γ = 50, δ = 1, α = 70)

4.2 모드 근접 현상

축방향과 횡방향 변위가 연성된 회전하는 외팔보는 회전속도 증가에 따라 고유진동수의 변화를 조사하면 특정 차수의 모드가 교환되는 흥미로운 현상을 관찰할 수 있다. 이러한 현상은 일반적으로 모드 교환 현상이나 모드 근접 현상으로 부른다. 앞 절에서 현재 모델과 기존 모델의 고유진동수 차이를 관찰한 Fig. 2를 보면, γ = 0일 때, 회전속도가 증가하면서 4차 고유진동수가 3차와 2차 고유진동수와 근접하는 지점이 존재하며 이어서 7차 고유진동수는 6차, 5차, 3차 고유진동수와 모드 근접 현상이 발생한다.

모드 근접 현상을 면밀히 검토하기 위해서 무차원화된 회전속도 γ를 0부터 20까지 증가시킬 때, 현재 모델의 1차부터 4차까지의 축방향과 횡방향의 고유진동수 변화를 Fig. 4에 나타냈다. P1-Q1-R1으로 이어지는 1차 고유진동수는 다른 차수의 고유진동수와 근접하지 않지만, P2-Q2-R2로 이어지는 2차 고유진동수는 γ가 10.5 근방일 때 3차 고유진동수와 근접하며, P3-Q3-R3로 이어지는 3차 고유진동수는 γ가 3.05 근방에서 4차 고유진동수와 근접하는 현상이 관찰된다.

Fig. 4Variation of lowest fourth natural frequencies of the rotating cantilever beam with tip-mass, when μt =1, δ = 1, α = 70)

특정 회전 속도에서 각 차수의 모드를 연속적으로 그린 Fig. 5를 참고하면, 모드가 교환되는 현상에 대해서 명확하게 확인할 수 있다. Fig. 4에서 무차원화된 회전속도가 0, 5, 15일 때, 각각의 고유진동수를 붉은색 점으로 표시하였다. 회전속도가 0일 때, 정지된 외팔보의 고유진동수 P1, P2, P3는 굽힘 모드이고 P4는 축방향 모드로써, 모드 형상은 Fig. 5에서 확인할 수 있다. 먼저, Fig. 5에서 2차 고유진동수 P2-Q2-R2의 모드를 확인하면 P2와 Q2의 모드 형상은 상당히 유사하지만 R2의 모드형상과는 다르다. Q2와 Q3모드 형상은 모드 근접 현상이 나타나는 회전속도 이후에 각각 R3와 R2와 같이 뒤바뀌게 된다. 따라서, 정지된 외팔보의 2차 굽힘 모드는 회전속도가 증가하면서 모드 근접 회전속도(γ =10.5) 이후에 3차 모드가 되는 반면 3차 모드는 그 회전속도를 지나면서 2차 모드가 된다. 마찬가지로 3차와 4차 모드도 무차원화 회전속도 γ가 3.05일 때를 기점으로 교환된다. 한편, 4차 모드는 외팔보가 회전하지 않을 때 축 방향 모드로써 Fig. 4에서 P4-Q3-R2로 모드가 교환되면서 이어진다. Fig. 5에서 모드 형상을 확인하면 P4-Q3-R2의 모드 형상이 동일함을 확인할 수 있다.

Fig. 5Mode shapes of the axial and chordwise motions at dimensionless rotating speeds γ = 0 and 5, 10 (when μt = 1, δ = 1, α = 70)

 

5. 결 론

이 논문에서는 종단 질량을 갖는 회전하는 외팔보의 비선형 편미분 방정식을 유도하였다. 현재 모델의 모달 특성을 조사하기 위해서 섭동법을 사용하여 회전속도가 일정한 정상상태에서 선형화된 운동방정식을 구하였다. 현재 모델의 선형화된 운동방정 식은 기존 선형 모델의 운동방정식과 형태가 비슷하지만 원심력에 의한 굽힘 강성 변화를 계산하는 방식에 차이가 있다. 원심력을 계산할 때, 현재 모델은 축 변형을 포함하여 계산하지만, 기존 선형 모델에서는 이를 무시하였다. 따라서, 회전속도에 따라 두 모델의 고유진동수를 비교하면 저속에서는 큰 차이 보이지 않지만, 회전속도가 높아질수록 차이가 커진다. 또한, 종단 질량의 크기에 따라서 현재 모델의 각 차수별 고유진동수는 기존 모델보다 높게 평가된다. 한편, 추가적으로 현재 모델을 사용해서 특정 회전속도를 기점으로 모드가 교환되는 모드 근접 현상을 분석하였다. 이러한 수치 결과를 통해서 우리는 고속으로 회전하는 가스터빈과 같은 시스템의 동적 특성을 정확하게 예측하기 위해서 이 논문의 모델링 방법을 제안한다.