Understanding Misalignment-Induced Gravity Critical Speed (Subharmonic Resonance) in Rotating Machinery via an Equivalent Base-Excited Pendulum Model

Sang-Jin Park; Huinam Rhee; Daesic Jang; Jeonghan Lee; Jin-Ho Park

doi:10.7734/COSEIK.2026.39.3.143

Preview

Research Paper

Journal of the Computational Structural Engineering Institute of Korea. 30 June 2026. 143-151
https://doi.org/10.7734/COSEIK.2026.39.3.143

Understanding Misalignment-Induced Gravity Critical Speed (Subharmonic Resonance) in Rotating Machinery via an Equivalent Base-Excited Pendulum Model

등가 지지점 가진 진자 모델을 이용한 회전기기의 오정렬 유발 중력 위험속도(저조파 공진) 현상의 이해

Sang-Jin Park¹

Huinam Rhee²^*

Daesic Jang³

Jeonghan Lee⁴

Jin-Ho Park⁵

박 상진¹

이 희남²^*

장 대식³

이 정한⁴

박 진호⁵

¹Senior Researcher, Advanced Fuel Cycle Technology Division, Korea Atomic Energy Research Institute, Daejeon, 31145, Korea

²Professor, School of Mechanical and Aerospace Engineering, Sunchon National University, Suncheon, 57922, Korea

³Senior Engineer, Nuclear System Integrity Sensing & Diagnosis Division, Korea Atomic Energy Research Institute, Daejeon, 31145, Korea

⁴Principal Engineer, Nuclear System Integrity Sensing & Diagnosis Division, Korea Atomic Energy Research Institute, Daejeon, 31145, Korea

⁵Principal Researcher, Nuclear System Integrity Sensing & Diagnosis Division, Korea Atomic Energy Research Institute, Daejeon, 31145, Korea

¹한국원자력연구원 선진핵주기기술개발부 선임연구원

²순천대학교 기계우주항공공학부 명예교수

³한국원자력연구원 기기안전진단연구부 선임기술원

⁴한국원자력연구원 기기안전진단연구부 책임기술원

⁵한국원자력연구원 기기안전진단연구부 책임연구원

^{*Corresponding Author}

License (open-access, https://creativecommons.org/licenses/by-nc/4.0/):

This is an Open-Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (https://creativecommons.org/licenses/by-nc/4.0) which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

ABSTRACT

Machine learning (ML)-based condition-based maintenance (CBM) of rotating machinery requires reliable training datasets, which are difficult to obtain from real operating environments. Because fault data are particularly difficult to obtain from actual operating environments, several studies have focused on generating fault data using computational and numerical methods. Among various mechanical faults in rotating machinery, misalignment faults are caused by the asymmetric stiffness of rotating components, which results in characteristic 2X vibration components and gravity critical speed phenomena. The dynamic behavior caused by misalignment in rotating machinery is fundamentally equivalent to base-excited pendulum (BEP) motion as a type of parametric excitation. Therefore, to artificially generate misalignment fault data and evaluate their physical accuracy and validity, it is necessary to thoroughly understand the fundamental dynamics of BEP motion. Against this backdrop, this study investigates the feasibility and validity of generating time-domain misalignment fault data by simulating BEP motion using a multibody dynamic analysis framework.

Keywords

misalignment

parametric excitation

gravity critical speed

Mathieu equation

multibody dynamics

base-excited pendulum

최근 기계학습(ML) 기반의 회전기기 상태기반정비(CBM) 진단 기술 개발을 위해 학습 데이터 확보가 중요하며, 실제 운전 현장으로부터 확보하기 힘든 결함 데이터를 전산 및 수치적으로 생산하려는 연구들이 수행되고 있다. 회전기기의 다양한 기계적 결함 들 중 오정렬 결함 현상은 회전 부품의 비대칭 강성(asymmetric stiffness)으로부터 기인되어 2X 진동과 gravity critical speed 현상을 발생시킨다. 흥미로운 것은 이러한 회전기기의 오정렬 결함 현상이 파라메트릭 가진 현상(parametric excitation)으로서 본질적으로 지지점 가진 진자(BEP) 운동과 등가적이다는 것이다. 그러므로 회전기기의 오정렬 결함 데이터를 인공적으로 생성하고 그 물리적 정확도와 타당성을 평가하기 위해서는 기초적인 BEP 운동 현상을 충분히 이해할 필요가 있다. 본 연구에서는 이러한 배경을 바탕으로 다물체 동역학(MBD) 해석을 바탕으로 BEP 운동을 시뮬레이션하여 오정렬 결함의 시계열 데이터 생산 가능성과 타당성을 평가하였다

키워드

오정렬

파라메트릭 가진

중력 위험속도

Mathieu 방정식

다물체 동역학

지지점 가진 진자

MAIN

1. 서 론
2. 오정렬 결함과 지지점 가진 진자(base-excited pendulum, BEP) 운동
2.1 비대칭 강성의 회전부품을 갖는 회전기기 시스템
2.2 오정렬에 의한 저조파 공진 현상
2.3 BEP 모델과 안정성 해석
2.4 Floquet 이론 기반 안정성 해석
2.5 안정성 해석 결과(Mathieu 안정성 차트)
3. BEP 운동의 다물체 동역학 시뮬레이션
3.1 시뮬레이션 모델과 조건
3.2 시뮬레이션 결과
3.3 오정렬 유발 중력 위험속도(저조파 공진) 현상에 대한 고찰
4. 결 론

1. 서 론

최근 산업 전반에서 시간기반정비(time-based maintenance, TBM) 체계에 의해 발생되는 경제적 비효율성을 개선시키고, 기계 설비들의 고장 예방과 예측을 통한 운정 효율성을 향상시키기 위한 목적으로 상태기반정비(condition-based maintenance, CBM) 기술이 빠르게 발전하고 있다(Jardine et al., 2006; Mobley, 2002; Park et al., 2019). 특히, 다양한 기계설비들 중 여러 산업 분야에서 널리 쓰이는 펌프, 압축기, 팬과 같은 회전기기 시스템들은 근원적으로 회전 현상에 의해 쉽게 진동 문제가 발생되는 설비들로서 오랫동안 진동 기반의 감시진단 기술이 발전되어 왔으며, 근래에는 상술한 배경과 함께 기계학습(machine-learning, ML) 기술의 폭발적인 발전 추세에 따라 ML 기반의 CBM 기술 개발이 집중적으로 이뤄지고 있다.

하지만, 이러한 ML 기반의 CBM 기술 성능은 학습에 사용되는 데이터의 양과 품질에 크게 의존한다. 실제 운전 환경에서는 설비의 신뢰성 확보를 위해 결함 상태에서의 장시간 운전이 제한되므로, 다양한 결함 데이터를 충분히 확보하는 데 근본적인 한계가 존재하며 이로 인한 학습 데이터 불균형 및 일반화 성능 저하 문제가 지속적으로 제기되고 있다(He and Garcia, 2009). 이러한 한계를 극복하기 위한 대안으로, 최근에는 전산 구조 해석을 통해 결함 데이터를 인공적으로 생성하려는 연구들이 보고되고 있고, 이를 통한 ML 모델의 개발 가능성과 성능 역시 이미 제시되었다(Kahr et al., 2022; Qin et al., 2024; Song et al., 2022). 전통적으로 전산 구조 해석 및 시뮬레이션 기법은 주로 기계 설비의 설계 목적으로 활용되어 왔는데(Kim, 2023), 근래에는 감시진단 기술 개발 목적으로서 활용되는 움직임이 나타나고 있는 것이다.

한편, 회전기기의 오정렬(misalignment) 결함은 회전기기의 대표적인 기계적 결함으로서, 불평형(unbalance) 결함과 함께 가장 빈번히 발생하는 결함들 중 하나이다(Jeong and Kim, 2008). 이러한 오정렬 결함은 구동부(drive part)와 피구동부(driven part)가 커플링(coupling)으로 연결되어 회전력이 전달될 때, 두 회전체의 회전 중심선(shaft centerline)이 설계 의도대로 일치하지 않는 상태를 의미한다. 즉 구동부(예: 모터, 터빈)에서 나온 토크가 피구동부(예: 펌프, 컴프레서)로 전달되는 과정에서 두 축의 중심선이 동축(coaxial)을 이루지 못하거나, 중심선은 평행하더라도 기울어짐(tilt)이 존재하여 각 변형(angular offset)이 생긴 경우이다. 오정렬이 존재하면 커플링과 베어링 등에 주기적인 반력이 발생하여 결과적으로 회전기기의 고진동, 베어링 하중 증가, 커플링/씰(seal) 마모, 발열 등의 문제를 일으키기에 회전기기의 진단 기술에서 반드시 검출해야 하는 대상이다. 그러나 오정렬 현상은 피상적으로 회전축 간의 엇갈림으로서 단순해 보이지만, 그 특성과 메커니즘은 결코 단순하지 않으며 오래된 회전체동역학의 학술 분야에서 관련 연구가 최근까지도 보고되고 있음이 이를 반증한다(Avendano and Childs, 2013; Lees, 2007; Redmond, 2007). 따라서, 상술한 배경을 바탕으로 ML 학습용 결함 데이터 생성과 활용을 위해 오정렬 결함에 대한 전산 구조 해석 기반의 접근은 피상적으로 단순하게 생각되지만, 생성된 데이터의 타당성과 물리적 정합성에 대한 평가와 검토는 오정렬 현상에 대한 충분한 이해 없이는 불가능하다.

본 연구의 주된 목적은 오정렬 유발 동역학과 지지점 가진 진자(BEP) 운동 간의 물리적 등가성을 바탕으로 ML 기반 상태기반정비(CBM)를 위한 물리 기반 데이터 생성 프레임워크를 제시하는 데 있다.

즉, 본 연구에서는

1) 회전체 오정렬 동역학과 BEP 운동 간의 물리적 등가성을 직관적으로 제시하고,

2) 다물체 동역학(MBD) 시뮬레이션을 통해 파라메트릭 공진 거동이 이론과 일관되게 재현됨을 확인하며,

3) 이를 기반으로 시간 영역 결함 데이터 생성 가능성을 제시하고자 한다.

2. 오정렬 결함과 지지점 가진 진자(base-excited pendulum, BEP) 운동

오정렬 결함은 크게 평행형 오정렬(parallel misalignment)와 각도형 오정렬(angular misalignment)로 구분되며, 유형별 세부적 진동 응답은 서로 차이를 보이나 공통적으로 결함이 심화되었을 때 2X 진동(입력 회전수 2배 진동)이 나타나는 것으로 알려져 있다(Bently and Hatch, 2002; Muzynska, 2005). 이러한 진동 특성은 크게 베어링의 비선형성과 회전부품의 비대칭 강성에 의한 것으로 밝혀져 있으며(Avendano and Childs, 2013; Bently and Hatch, 2002; Muzynska, 2005), 본 연구에서는 후자의 발생론만을 집중적으로 다룬다.

2.1 비대칭 강성의 회전부품을 갖는 회전기기 시스템

회전체 시스템에서 오정렬을 포함하여 회전자(disc)의 균열 등이 발생하면 비대칭 강성이 유발되는 것으로 알려져 있다(Avendano and Childs, 2013;, Muzynska, 2005). 비대칭 강성의 회전부품(예: 커플링, 회전축, 회전자 등)을 갖는 회전체 시스템은 Fig. 1과 같다. 질량 $M$ 을 갖는 회전자는 $x y z$ 고정 좌표계의 중심 $O$ 에 중앙지지형(center-hung type)으로서 놓여 지며, 비대칭 강성( $k_{ζ} \neq k_{η}$ )을 갖는 회전축에 대하여 굽힘 모드의 주축(principal axis) 방향으로 정렬된 $ζ η z$ 회전 좌표계는 회전입력 $Ω$ 에 따라 회전한다. 마지막으로 회전축 양단은 등방성(isotropic) 베어링에 의해 지지되며 회전체 시스템에 작용하는 정적 외력으로서 - $y$ 방향의 중력과 비감쇠를 고려하면 운동방정식은 회전 좌표계 상에서 식 (1)과 같이 기술된다(Lee, 1993).

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F1.jpg

Fig. 1.

Coordinate of rotating machinery with rotating part with asymmetric stiffness

(1)

\ddot{ζ} - 2 Ω \dot{η} + (ω_{ζ}^{2} - Ω^{2}) ζ = - g \cos Ω t, \ddot{η} + 2 Ω \dot{ζ} + (ω_{η}^{2} - Ω^{2}) η = g \sin Ω t .

여기서, $ω_{ζ}^{2} = \frac{k_{ζ}}{M}, ω_{η}^{2} = \frac{k_{η}}{M}$ 이고, 식 (1)의 해는

(2)

ζ = ζ_{f} e^{j Ω t} + ζ_{b} e^{- j Ω t}, \bar{ζ} = {\bar{ζ}}_{f} e^{- j Ω t} + \bar{ζ} e^{j Ω t},

식 (2)와 같고, 여기서 $ζ_{f}$ 와 $ζ_{b}$ 는

(3)

ζ_{f} = \frac{- g (ω_{ζ}^{2} - ω_{η}^{2})}{4 Ω^{2} (ω_{ζ}^{2} + ω_{η}^{2}) - 2 ω_{ζ}^{2} ω_{η}^{2}}, ζ_{b} = \frac{- g \{8 Ω Ω^{2} - (ω_{ζ}^{2} + ω_{η}^{2})\}}{4 Ω^{2} (ω_{ζ}^{2} + ω_{η}^{2}) - 2 ω_{ζ}^{2} ω_{η}^{2}},

와 같다. 그러므로 고정 좌표계에서의 진동 응답 $r_{s}$ 는

(4)

r_{s} = ζ e^{j Ω t} = ζ_{f} e^{j 2 Ω t} + ζ_{b}

와 같다. 결국 식 (4)를 통해 비대칭 강성을 갖는 회전체 시스템은 중력 하에서 $ζ_{b}$ 크기를 갖는 정적 처짐과 함께 반경 $ζ_{b}$ 와 $2 Ω$ 회전수를 갖는 원형 whirl의 거동을 한다는 것을 알 수 있다.

2.2 오정렬에 의한 저조파 공진 현상

한편 회전 좌표계 상에서 기술된 운동 방정식 식 (1)은 정지(고정) 좌표계 상에서

(5)

\ddot{x} + \{\frac{ω_{ζ}^{2} + ω_{η}^{2}}{2} + \frac{ω_{ζ}^{2} - ω_{η}^{2}}{2} \cos 2 Ω t\} x + \{\frac{ω_{ζ}^{2} - ω_{η}^{2}}{2} \sin 2 Ω t\} y = 0, \ddot{y} + \{\frac{ω_{ζ}^{2} + ω_{η}^{2}}{2} - \frac{ω_{ζ}^{2} - ω_{η}^{2}}{2} \cos 2 Ω t\} y + \{\frac{ω_{ζ}^{2} - ω_{η}^{2}}{2} \sin 2 Ω t\} x = - g

와 같이 기술된다(Lee, 1993). 이때 수식에 대한 물리적 직관을 위해

(6)

ω_{0} = (ω_{ζ} + ω_{η}) / 2, Δ ω = (ω_{η} - ω_{ζ}) / 2

식 (6)의 조건과

(7)

Δ ω / ω_{0} ≪ 1

식 (7)의 가정을 고려하면 식 (5)는

(8)

\ddot{x} + \{ω_{0}^{2} - 2 ω_{0} Δ ω \cos (2 Ω t)\} x - \{2 ω_{0} Δ ω \sin (2 Ω t)\} y = 0, \ddot{y} - \{2 ω_{0} Δ ω \sin (2 Ω t)\} x + \{ω_{0}^{2} + 2 ω_{0} Δ ω \cos (2 Ω t)\} y = - g

와 같이 정리된다. 식 (8)은 $x y$ 방향의 운동이 서로 연립된(coupled) 2계 선형 상미분 방정식인 것을 볼 수 있으며, 시스템 계수(강성항)가 시간에 따라 $2 Ω$ 의 주파수로 변조되는 비자율(non-autonomous) 형태 인 것을 파악할 수 있다. 만일 $ω_{0}$ 에 비해 주파수 불일치량 $Δ ω$ 가 작다면 약연립(weakly-coupled) 시스템으로서 실질적으로 연립항을 무시할 수도 있다.

식 (1)과 (8)은 동일한 물리적 현상을 좌표계에 따라 각각 다르게 나타내고 있으며 물리적 해석 역시 크게 달라진다. 구체적으로 식 (3)에서 분모가 0이 되는 조건으로부터,

(9)

Ω_{s u b} = (ω_{ζ} + ω_{η}) / 4 = ω_{0} / 2

와 같이 저조파 공진(subharmonic resonance)의 위험속도 $Ω_{s u b}$ 가 닫힌 해(closed-form solution) 형태로서 해석적으로 도출되며 이를 통해 비대칭 강성을 갖는 회전기기는 고유진동수 $ω_{0}$ 의 절반에서 공진이 추가적으로 발생한다는 물리적 사실을 파악할 수 있다(Lee, 1993). 식 (1)을 관찰하면, 이 현상은 중력이 마치 외부 조화 가진력으로서 작용하여 시스템의 고유치(eigenvalue)와 외력의 주파수가 서로 일치하여 나타나는 공진 현상으로 해석될 수 있다. 그러므로 비대칭 강성을 갖는 회전체 시스템에서 나타나는 저조파 위험속도 현상은 중력에 의해 유발되기에 gravity critical speed(혹은 gravitational whirl, secondary critical speed, one-half critical speed) 라고도 불린다(Avendano and Childs, 2013; Ishida et al., 2008; Messal and Bonthron, 1972). 그러나 이러한 저조파 공진을 유발하는 중력항은 정지 좌표계에서 단순한 정적 하중이며, 식 (8)을 바탕으로 외부 조화 가진력은 존재하지 않는다. 따라서, 동일한 현상을 정지 좌표계에서 서술할 때, 회전좌표계에서의 해석과 달리 수학적으로 닫힌 해가 산출되지 않는 동시에 물리적으로는 파라메트릭 공진(parametric resonance)에 의한 것으로 설명될 수 있으며, 실제 중력 효과를 무시할 수 있는 수직 회전체 시스템에서도 비대칭 강성을 갖는 베어링으로부터 기인된 파라메트릭 공진에 의해 저조파 공진 현상이 발생된다는 실험 연구 결과도 보고되어 있다(Messal and Bonthron, 1972).

2.3 BEP 모델과 안정성 해석

지지점 가진 진자(BEP) 모델은 길이 $l$ 과 질량 $m$ 을 갖는 비감쇠 진자 시스템이 수직 방향으로 $Y \cos (Ω t)$ 의 형태의 지지점 가진 조건 하에 놓였다고 가정하면, 미소각 기반의 선형화 가정 하에서 각변위 𝜃에 대한 BEP 운동은

(10)

\ddot{θ} + [ω_{swing}^{2} - q Ω^{2} \cos (Ω t)] θ = 0,

와 같이 나타낼 수 있으며, 여기서

(11)

ω_{swing}^{2} = \frac{g}{l}, q = \frac{Y}{l}

식 (11)과 같고, $Y \cos (Ω t)$ 는 수직 진폭, $g$ 는 중력가속도이다. 식 (10)에 무차원 시간 $τ = Ω t$ 을 적용하면,

(12)

\frac{d^{2} θ}{d τ^{2}} + [{(\frac{ω_{swing}}{Ω})}^{2} - q \cos (τ)] θ = 0,

와 같이 나타낼 수 있다. 한편, 파라메트릭 공진의 수학적 모델로서 감쇠계수 $c$ 를 포함한 표준화된 Mathieu 방정식은

(13)

\ddot{x} + c \dot{x} + (a - ε \cos t) x = 0,

와 같다(Park, 2013). 연계항을 무시한다면, 오정렬 결함에 대한 식 (8)과 BEP 운동에 대한 식 (12)는 모두 식 (13)의 형태를 따르고 있음을 볼 수 있다. 그러므로 두 현상은 비록 서로 다른 물리적 대상에서 나타나지만 수학적으로 서로 등가적 운동이라는 것을 알 수 있다.

2.4 Floquet 이론 기반 안정성 해석

식 (13)과 같이 표준화된 Mathieu 방정식의 해는 일반적으로 섭동법(perturbation method), 다중스케일법(multiple scale method), 그리고 Floquet 이론을 통해 해석되거나 안정성이 분석된다. 이들 중에서 Floquet 이론은 상태 벡터(state vector)를 도입하여 $z_{1} = x, z_{2} = \dot{x}$ 를 통해 식 (13)과 같은 2계 방정식을 식 (14)와 같이 1계 방정식으로 치환한다(Park, 2013).

(14)

\dot{z} (t) = F (t) z (t),

여기서,

(15)

z (t) = [\begin{array}{l} z_{1} \\ z_{2} \end{array}], F (t) = [\begin{matrix} 0 & 1 \\ - (a - ε \cos t) & - c \end{matrix}]

이다. 식 (15)에서 계수 행렬 $F (t)$ 는 주기 $T = 2 π$ 의 주기성을 갖는다는 것을 알 수 있다. 즉 다시 말해, $F (t + T) = F (t)$ 를 만족한다. 식 (13)은 2계 선형 상미분 방정식이므로 일반 해 $z (t)$ 는

(16)

z (t) = A u (t) + B v (t) = P (t) e^{R t} z_{0}, (z_{0} = [\begin{array}{l} A \\ B \end{array}])

와 같이 초기 조건에 의해 결정되는 상수 $A$ , $B$ 와 함께 서로 독립적인 기본 해 $u (t)$ 와 $v (t)$ 의 선형 조합으로서 나타낼 수 있으며, 동시에 주기적 행렬 $P$ , 상수 행렬 $R$ 그리고 초기 상태벡터 $z_{0}$ 의 조합으로도 나타낼 수 있다. Floquet 이론에 의하면 식 (16)는 식 (14)의 해이므로, 주기 $T$ 가 지난 후에도 동일하게 식 (14)의 해이다. 그러므로

(17)

Φ (t + T) = Φ (t) ∙ Φ (T) = Φ (t) M,

가 만족된다. 여기서, 기본 해의 행렬 $Φ (t)$ 는

(18)

Φ (t) = [u (t) v (t)] = P (t) e^{R t},

와 같고 행렬 $M$ 은

(19)

M = [\begin{array}{l} a_{11} & a_{12} \\ a_{21} & a_{22} \end{array}],

식 (19)이다. 특히 행렬 $M$ 은 한 주기 $T$ 이후의 상태를 나타내기에 Floquet 이론에서 모노드로미 행렬(monodromy matrix)로 알려져 있으며, 행렬 내부 성분(entry) $a_{i} j (i, j = 1, 2)$ 는 시스템변수 $a$ $ε$ , $c$ , 에 의해 결정되는 상수이다. 결과적으로 식 (17)과 식 (18)로부터

(20)

M = e^{R T},

식 (20)의 관계식를 얻을 수 있으며, 행렬 $M$ 과 $R$ 의 고유치를 각각 𝜇와 𝜆라고 한다면,

(21)

λ_{i} = \frac{1}{T} \ln (μ_{i}), (i = 1, 2)

식 (21)의 관계식을 도출할 수 있다. 이때, 𝜆는 Floquet 이론에서 Floquet 지수(Floquet exponent)라 불리는 값으로서 아래 식 (22)와 같은 기준으로 시스템의 안정성을 판별한다(Park, 2013).

(22)

|μ_{i}| > 1 \to unstable, |μ_{i}| < 1 \to stable, |μ_{i}| = 1 \to neutral .

2.5 안정성 해석 결과(Mathieu 안정성 차트)

상기한 Floquet 이론을 통해 수치적으로 계산한 안정성 해석 결과(Mathieu 안정성 차트)는 Fig. 2와 같다. Fig. 2(a)는 비감쇠 시스템( $c$ =0)에 대한 계산 결과이며, Fig. 2(b)와 Fig. 2(c)는 각각 $c$ =0.1 및 $c$ =0.2에 대한 결과이다. Fig. 2에서 파란색은 불안정 영역이고, 노란색은 안정 영역을 나타내고 있다. 안정성 해석 결과에 따르면, 시스템의 고유진동수에 해당하는 가로축 변수 $a$ 가 증가할수록 불안정 영역은 낮은 𝜀 값에 대하여 그 영역이 좁아지는 것을 볼 수 있고, 이는 다시 말해 작은 𝜀에 대하여 고차(high-order) 파라메트릭 공진의 발생이 어렵다는 것을 뜻한다. 또한, $a$ = 1/4의 경우는 가장 낮고 상대적으로 넓은 불안정 문턱(instability tongue)을 보이는 특징을 갖는 것을 알 수 있고, 이러한 특징으로 인해 모든 가진 조건들 중에서 파라메트릭 공진이 가장 쉽게 발생할 수 있기 때문에 이 경우는 특별히 파라메트릭 주공진(principal parametric resonance)의 영역으로 불린다. 그러므로 식 (12)에 의해 BEP 운동의 주공진 조건은 저차 근사화를 통해서

(23)

a = {(\frac{ω_{swing}}{Ω})}^{2} = \frac{1}{4} \pm \frac{ε}{2},

와 같으므로(Park, 2013), $Ω = 2 ω_{s w i n g}$ 에서 발생되는 것을 알 수 있다. 불안정 문턱은 감쇠계수 $c$ 가 증가함에 따라 $a$ 의 모든 범위에서 증가하는 것을 알 수 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F2.jpg

Fig. 2.

Mathieu stability charts (blue = unstable, yellow = stable) (a) $c$ =0.0, (b) $c$ =0.1, (c) $c$ =0.2

3. BEP 운동의 다물체 동역학 시뮬레이션

다물체 동역학(multibody dynamics, MBD) 기법은 운동학적 제약조건(kinematic constraints) 하에서 설정된 운동방정식을 바탕으로 질점의 상태(위치, 속도 등)을 수치적으로 계산하는 전산 해석 방법으로서 본 연구에서는 상용 MBD 코드 ADAMS(MSC Software Corporation, 2024)를 통해 BEP 운동의 시뮬레이션이 수행되었다. MBD는 해석적인 해를 구하기 어려운 비선형 거동을 효율적으로 수치해석할 수 있으며, 복잡한 회전체 시스템으로의 확장이 용이하고, ML 학습에 활용 가능한 시간 영역 데이터를 직접 생성할 수 있다는 장점이 있다.

3.1 시뮬레이션 모델과 조건

BEP 시스템의 파라메트릭 가진 현상을 시뮬레이션하기 위해 구축된 MBD 모델은 Fig. 3과 같다. BEP 시스템 줄과 무게추는 모두 강체로 모델링되었고, 시스템 고유진동수 $f_{0}$ 은 비감쇠 및 미소각(small angle) 기반의 선형 시스템 기준으로 1 Hz이다. BEP 시스템의 회전점(pivot point)는 1자유도의 회전연결(revolute joint)를 통해 모델링되었으며, 회전점(지지점)은 $y$ 방향으로 $Y \cos (Ω t)$ 의 형태로서 수직 조화 운동이 적용되었다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F3.jpg

Fig. 3.

Multibody dynamic simulation model for base-excited pendulum

MBD 시뮬레이션 조건으로서 explicit 적분 기법 기반의 GSTIFF I3 적분기가 적용되었으며 적분 시간 간격은 1 × 10^-3 sec이다.

3.2 시뮬레이션 결과

구축된 모델을 통한 시뮬레이션 결과는 Figs. 4, 5, 6, 7, 8과 같다. 우선, Fig. 4는 비감쇠 조건 하에서 비수직 진동( $q$ = 0)과 초기각 $θ_{0}$ = 175°의 대변형(large deformation)을 고려한 자유 진동 시뮬레이션 결과이다. 진자 운동은 대각도(대변형)가 고려되면 아래의 식 (24)와 같이 기하학적 비선형 시스템이 되기 때문에,

(24)

\ddot{θ} + ω_{0}^{2} \sin θ = 0,

감쇠 및 외력이 없는 경우 그 자유 진동의 시간 응답은 Jacobi 타원함수(sn, cn)로 표현될 수 있으며, 구축된 MBD 시뮬레이션의 BEP 모델이 이러한 비선형 시스템의 응답을 매우 정확히 계산하고 있음을 Fig. 4를 통해 알 수 있고, 이를 통해 구축된 MBD 모델 타당성과 정확도를 기초적으로 확인할 수 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F4.jpg

Fig. 4.

Simulation result of pendulum motion under large deformation condition of $θ_{0}$ = 175°

Fig. 5는 비감쇠( $c$ =0) 조건 하에서 초기 각변위 $θ_{0}$ =0, 초기 각속도 ${\dot{θ}}_{0}$ = 10deg/sec, $q = Y / l =$ 4.026 × 10^-3 그리고 $a$ = 1/4( $Ω = 2 ω_{0}$ )의 파라메트릭 주공진에 대한 시뮬레이션 결과이다. Fig. 5(a)의 시간 이력 결과와 Fig. 5(b)의 위상 평면 결과는 모두 시간에 따라 BEP 각도가 최대 20.41°까지 증가한 이후 다시 감소하는 것을 볼 수 있는데, 이는 주공진 영역 내부에서 시간이 증가함에 따라 무한히 발산해야 하는 이론적 결과와 다르다. 이러한 경향성은 구축된 MBD 모델이 Fig. 4와 같이 비선형 시스템이기 때문이다. 다시 말해, 각도가 증가함에 따라 BEP 시스템의 고유진동수가 감소하고 이로 인해 식 (23)의 불안정 영역을 이탈했기 때문에 BEP의 운동은 무한히 증폭되지 않는다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F5.jpg

Fig. 5.

Simulation result of BEP motion under principal parametric resonance condition ( $c$ =0, $ε_{0} = Y / l =$ 4.026 × 10^-3, $a$ = 1/4, $θ_{0}$ =0, ${\dot{θ}}_{0}$ = 10deg/sec), (a) time-angle graph, (b) phase plane

이를 보다 수학적으로 기술하면, 식 (13)은 진자 운동의 등가(유효) 강성을 바탕으로 성립되는 $\sin θ \approx (1 - θ^{2} ⁄ 8) θ$ 와 비감쇠 가정 하에서

(25)

\ddot{θ} + (a^{'} - ε^{'} \cos t) θ = 0 [a^{'} = a κ (θ), ε^{'} = ε κ (θ), κ (θ) = 1 - θ^{2} / 8]

식 (25)와 같이 기술될 수 있다. 이때, 만일 상기 케이스와 같이 BEP 운동이 $a$ = 1/4에서 가진되었다면 각도 $θ_{m a x}$ 는 식 (23)의 주공진 불안정영역의 경계(boundary)에서 발생하므로

(26)

1 - κ (θ_{\max}) = 2 ε ∙ κ (θ_{\max}),

와 같은 관계식이 성립되며, 최종적으로 식 (26)으로부터

(27)

θ_{\max} = \sqrt{\frac{16 ε}{1 + 2 ε}} = \sqrt{\frac{32 q}{1 + 4 q}},

와 같은 관계식이 도출된다. 그러므로 식 (27)과 상기 시뮬레이션 조건을 통해 계산된 BEP의 최대각(bounded angle) 이론치는 20.40°로서 생성된 전산 데이터는 0.04%의 오차를 보인다. 따라서, 이러한 경향성은 구축된 MBD 모델이 물리학적으로 엄밀하고 정확히 시스템의 응답을 계산한다는 것을 반증하고 있다. 또한, 동일한 조건에서 수직 가진 크기 $q$ 를 증가시킨 경우에 대한 시뮬레이션 결과는 Fig. 6에 제시되었다. Fig. 5의 케이스에 대한 수직 가진 크기의 2배값 $2 q$ 부터 2배씩 증가할수록 BEP의 최대 운동 각도가 28.71°, 40.13°, 48.10° 로서 증가하고 식 (27) 기반의 이론치와 각각 0.31%, 0.67%, 1.14%의 오차를 갖는 것으로 확인되어 정량적으로 매우 정확하다는 것을 알 수 있으며, $q$ 값에 따른 오차 증가는 이론적 접근에서 저차 근사화 영향성으로 파악된다. 이러한 경향성은 Fig. 2에서 확인할 수 있듯이 파라메트릭 가진 크기 $q$ 가 커질수록 불안정성 영역이 넓어져 BEP 시스템의 고유진동수가 다소 변화하여도 불안정성 영역 내부에 위치하기에 파라메트릭 공진이 발생하는 것이다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F6.jpg

Fig. 6.

Simulation result of BEP motion under principal parametric resonance condition with increase of $q$

더불어, Fig. 7은 비감쇠( $c$ =0) 조건 하에서 $q = Y / l =$ 4.026 × 10^-3 및 $a$ = 1( $Ω = ω_{0}$ )에 대한 BEP 시뮬레이션 결과이다. Fig. 7(a)에서 볼 수 있듯이 시간이 증가함에 따라 진폭은 변하지 않으며, Fig. 7(b)의 위상평면 결과 역시 완벽히 닫힌 원(closed circle)의 형태를 보여주고 있다. 이는 Fig. 2 (a) 상에서 볼 수 있듯이 매우 작은 𝜀 값에서 $a$ = 1에 대한 불안정성 영역에 매우 좁기 때문이며 매우 작은 고유진동수의 변동에도 불안정성 영역을 벗어나기에 자유 진동이 일어나는 것으로 판단할 수 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F7.jpg

Fig. 7.

Simulation result of BEP motion under 2^nd order parametric resonance condition ( $c$ =0, $q = Y / l =$ 4.026 × 10^-3, $a$ = 1, $θ_{0}$ =0, ${\dot{θ}}_{0}$ = 10deg/sec), (a) time-angle graph, (b) phase plane

마지막으로 Fig. 5의 계산과 동일한 $q$ , $a$ 의 주공진 조건에서 모달 감쇠 0.8%를 고려했을 경우, Fig. 8(a)의 시간이력과 Fig. 8(b)의 위상평면 결과 모두에서 파라메트릭 공진이 발생되지 못 하고 감쇠되는 것을 볼 수 있다. 이러한 이유는 상기한 바와 같이 시스템 감쇠에 따라 불안정성 문턱이 증가하기 때문이다. 따라서, 회전기기에서도 오정렬에 의한 파라메트릭 공진을 회피하기 위해서는 시스템 감쇠를 증가시킬 필요가 있다는 것을 알 수 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-03/N0040390301/images/Figure_jcoseik_39_03_01_F8.jpg

Fig. 8.

Simulation result of BEP motion under principal parametric resonance condition ( $q = Y / l =$ 4.026 × 10^-3, $a$ = 1/4, $θ_{0}$ =0, ${\dot{θ}}_{0}$ = 10deg/sec) with modal damping of 0.8%, (a) time-angle graph, (b) phase plane

결론적으로 BEP 운동에 대한 모든 MBD 시뮬레이션 결과는 Mathieu 안정성 차트에서 확인할 수 있는 파라메트릭 공진 특성을 잘 따르며, 특히 $a$ = 1/4의 주공진 조건에서 파라메트릭 공진이 상대적으로 쉽게 발생됨을 확인하였다. 또한, 구축된 MBD 시뮬레이션 모델은 선 형 이론에서 고려되지 못 한 비선형 시스템의 거동 역시 정확히 계산하고 있음을 알 수 있었으며, 이를 통해 MBD 기반의 전산 해석 기법이 회전기기 시스템의 ML 학습용 결함 데이터 생성을 위한 합리적 대안으로서 고려되어도 충분하다는 것으로 판단된다(Kahr et al., 2022; Song et al., 2022).

3.3 오정렬 유발 중력 위험속도(저조파 공진) 현상에 대한 고찰

수행된 Mathieu 안정성 차트 해석과 BEP 운동의 MBD 시뮬레이션을 통해 파라메트릭 공진은 주공진 조건에서 쉽게 발생될 수 있음을 파악했음으로, 회전기기에서 오정렬 유발 저조파 공진 역시 결과적으로 파라메트릭 주공진 조건에서 발생될 수 있음을 추론할 수 있다. 구체적으로 식 (8)은 $ω_{0}$ 에 비해 주파수 불일치량 $Δ ω$ 가 작은 약연립(weakly-coupled) 시스템으로서 연립항을 무시한 가정 하에서 무차원 시간 $τ = 2 Ω t$ 를 적용하면,

(28)

\frac{d^{2} x}{d τ^{2}} + \{{(\frac{ω_{0}}{2 Ω})}^{2} - \frac{ω_{0}}{Ω} Δ ω \cos (τ)\} x = 0, \frac{d^{2} y}{d τ^{2}} + \{{(\frac{ω_{0}}{2 Ω})}^{2} - \frac{ω_{0}}{Ω} Δ ω \cos (τ)\} y = - g

식 (28)과 같이 변형된다. 그러므로 표준화 Mathieu 방정식을 바탕으로

(29)

\frac{a}{4} = δ = {(\frac{ω_{0}}{2 Ω})}^{2},

와 같은 것을 파악할 수 있다. 따라서, 식 (29)의 관계식으로부터 유도되는 주공진 조건 $δ$ =1/4을 통해 $Ω = ω_{0} / 2$ 에서 발생된다는 것을 파악할 수 있고, 이는 식 (9)와 일치하는 것을 알 수 있다.

결론적으로 오정렬에 의한 중력 위험속도 현상은 관찰 좌표계에 따라 물리적 해석이 달라지지만 완벽히 동일한 운동임을 규명할 수 있으며, 근원적으로 BEP 운동의 파라메트릭 공진 현상과 동일하다는 것을 파악할 수 있다.

4. 결 론

최근 기계학습 (ML) 기반의 기계설비 진단 기술 개발에 요구되는 비정상(결함) 데이터의 확보를 위해 전산 시뮬레이션 모델의 활용이 하나의 대안으로서 주목받고 있으며(Kahr et al., 2022; Qin et al., 2024; Song et al., 2022), 시뮬레이션 모델로부터 생성된 결과 데이터를 ML 학습 목적으로 활용하기 위해서는 결과에 대한 물리적 타당성과 정합성을 반드시 판단해야 하므로 대상 설비와 결함 현상에 대한 도메인 지식(domain knowledge) 기반의 연구와 평가 역시 중요해지고 있다.

또한, 회전체 시스템에서 회전부품의 강성 이방성에 기인된 오정렬 현상은 2X 진동과 함께 주 위험속도(primary critical speed)의 1/2 속도 지점에서 추가적으로 중력 위험속도 (저조파 공진) 현상을 발생시킨다. 흥미로운 것은 이 현상이 회전 좌표계 상에서 중력에 의한 강제 진동 현상으로 해석되지만 정지 좌표계에서는 파라메트릭 공진에 의한 것으로 해석된다는 사실이다.

본 연구에서는 오정렬 현상과 그 진동 특성을 이해하고 이에 대한 합리적 결함 모사와 생성 데이터의 타당성 평가를 위해 근원적으로 등가적 물리적 현상이라고 할 수 있는 지지점 가진 진자(BEP) 운동에 대하여 다물체 동역학(MBD) 기반의 시뮬레이션을 수행하였으며, 이를 통해 오정렬에 의해 발생되는 중력 위험속도 현상이 파라메트릭 공진에 의한 것이며 근원적으로 BEP 운동과 같다는 것을 제시하였다. 결국 이러한 등가적 접근을 활용하여 오정렬 현상을 모사하여 전산적으로 생성한 결함 데이터의 정확성과 물리적 타당성 평가가 가능할 것으로 판단된다.

따라서, 제안한 MBD 기법의 전산 해석을 통한 ML 학습용 결함 데이터셋 생성 접근론은 비정상 데이터 확보 현안에 대한 하나의 대안으로서 충분히 합리적인 것으로 판단되며, 본 연구 결과는 회전기기 시스템의 결함 데이터에 대한 인공적 생성 및 ML 진단 모델 개발 연구와 회전기기의 강연립(strongly-coupled) 시스템에 대한 오정렬 거동 연구의 기초적 자료로서 활용될 수 있을 것으로 사료된다.

Acknowledgements

본 연구는 기후에너지환경부(MCEE)와 한국에너지기술평가원(KETEP)의 지원을 받아 수행한 연구 과제입니다(No. 20224B10100060).

References

Avendano, R.D., Childs, D.W. (2013) One Explanation for Two Times Running Speed Response Due to Misalignment in Rotors Connected by Flexible Couplings, J. Eng. Gas Turbines Power, 135(6), pp.1~10.

10.1115/1.4023232

Bently, D.E., Hatch, C.T. (2002) Fundamentals of Rotating Machinery Diagnostics, Minden, NV: Bently Pressurized Bearing Press, p.726.

He, H., Garcia, E.A. (2009) Learning from Imbalanced Data, IEEE Trans. Knowl. & Data Eng., 21(9), pp.1263~1284.

10.1109/TKDE.2008.239

Ishida, Y., Liu, J., Inoue, T., Suzuki, A. (2008) Vibrations of an Asymmetrical Shaft with Gravity and Nonlinear Spring Characteristics (Isolated Resonances and Internal Resonances), J. Vibr. & Acoust., 130, pp.041004-1~041004-8.

10.1115/1.2889475

Jardine, A.K.S., Lin, D., Banjevic, D. (2006) A Review on Machinery Diagnostics and Prognostics Implementing Condition-based Maintenance, Mech. Syst. & Signal Proc., 20(7), pp.1483~1510.

10.1016/j.ymssp.2005.09.012

Jeong, J.H., Kim, J.W. (2008) Equipment Diagnosis through Vibration Measurement, Fire Prot. Technol., 44, pp.30~43.

Kahr, M., Kovács, G., Loinig, M., Brückl, H. (2022) Condition Monitoring of Ball Bearings Based on Machine Learning with Synthetically Generated Data, Sens., 22(7), p.2490.

10.3390/s2207249035408105PMC9002605

Kim, S.J. (2023) Design and Dynamic Analysis of Rolling Element Bearings for High-Speed Turbo-Machineries, J. Comput. Struct. Eng. Inst. Korea, 36(1), pp.15~23.

10.7734/COSEIK.2023.36.3.185

Lee, C.W. (1993) Vibration Analysis of Rotors, Springer, p.315.

10.1007/978-94-015-8173-8

Lees, A.W. (2007) Misalignment in Rigidly Coupled Rotors, J. Sound & Vibr., 305, pp.261~271.

10.1016/j.jsv.2007.04.008

Messal, E.E., Bonthron, R.J. (1972) Subharmonic Rotor Instability Due to Elastic Asymmetry, J. Eng. Ind., 94(1), pp.185~192.

10.1115/1.3428110

Mobley, R.K. (2002) An Introduction to Predictive Maintenance, 2nd ed., Butterworth-Heinemann, p.437.

10.1016/B978-075067531-4/50006-3

MSC Software Corporation (2024) Adams 2024: Automated Dynamic Analysis of Mechanical Systems, Hexagon AB, Stockholm, Sweden.

Muzynska, A. (2005) Rotordynamics, CRC Press, Boca Raton, FL, p.1128.

Park, C.I. (2013) Stability Analysis of Mathieu Equation by Floquet Theory and Perturbation Method, Trans. KSNVE, 23(8), pp.734~741.

10.5050/KSNVE.2013.23.8.734

Park, J.H., Yoon, D.B., Lee, J.H., Choi, Y.C. (2019) Condition Monitoring and Diagnosis in Mechanical Components of Nuclear Power Plant, J. Korean Soc. Nondestruct. Test., 39(6), pp.369~377.

10.7779/JKSNT.2019.39.6.369

Qin, Y., Fu, X., Li, X., Li, H. (2024) ADAMS Simulation and HHT Feature Extraction Method for Bearing Faults of Coal Shearer, Processes, 12(1), p.164.

10.3390/pr12010164

Redmond, I. (2007) Shaft Misalignment and Vibration: A Model, Saudi Aramco J. Technol., 4, pp.41~51.

Song, M.M., Xiong, Z.C., Zhong, J.H., Xiao, S.G., Tang, Y.H. (2022) Research on Fault Diagnosis Method of Planetary Gearbox based on Dynamic Simulation and Deep Transfer Learning, Sci. Rep., 12, 17023.

10.1038/s41598-022-21339-536220866PMC9553873

Journal of the Computational Structural Engineering Institute of Korea ISSN:1229-3059(Print) 2287-2302(Online) 한국전산구조공학회 논문집

Preview

Understanding Misalignment-Induced Gravity Critical Speed (Subharmonic Resonance) in Rotating Machinery via an Equivalent Base-Excited Pendulum Model

ABSTRACT

MAIN

Fig. 1.

Coordinate of rotating machinery with rotating part with asymmetric stiffness

(1)

(2)

(3)

(4)

(5)

(6)

(7)

(8)

(9)

(10)

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(18)

(19)

(20)

(21)

(22)

(23)

Fig. 2.

Mathieu stability charts (blue = unstable, yellow = stable) (a) c=0.0, (b) c=0.1, (c) c=0.2

Fig. 3.

Multibody dynamic simulation model for base-excited pendulum

(24)

Fig. 4.

Simulation result of pendulum motion under large deformation condition of θ0 = 175°

Fig. 5.

Simulation result of BEP motion under principal parametric resonance condition (c=0, ε0=Y/l= 4.026 × 10-3, a= 1/4, θ0=0, θ˙0= 10deg/sec), (a) time-angle graph, (b) phase plane

(25)

(26)

(27)

Fig. 6.

Simulation result of BEP motion under principal parametric resonance condition with increase of q

Fig. 7.

Simulation result of BEP motion under 2nd order parametric resonance condition (c=0, q=Y/l= 4.026 × 10-3, a= 1, θ0=0, θ˙0= 10deg/sec), (a) time-angle graph, (b) phase plane

Fig. 8.

Simulation result of BEP motion under principal parametric resonance condition (q=Y/l= 4.026 × 10-3, a= 1/4, θ0=0, θ˙0= 10deg/sec) with modal damping of 0.8%, (a) time-angle graph, (b) phase plane

(28)

(29)

Acknowledgements

References

Mathieu stability charts (blue = unstable, yellow = stable) (a) $c$ =0.0, (b) $c$ =0.1, (c) $c$ =0.2

Simulation result of pendulum motion under large deformation condition of $θ_{0}$ = 175°

Simulation result of BEP motion under principal parametric resonance condition ( $c$ =0, $ε_{0} = Y / l =$ 4.026 × 10^-3, $a$ = 1/4, $θ_{0}$ =0, ${\dot{θ}}_{0}$ = 10deg/sec), (a) time-angle graph, (b) phase plane

Simulation result of BEP motion under principal parametric resonance condition with increase of $q$

Simulation result of BEP motion under 2^nd order parametric resonance condition ( $c$ =0, $q = Y / l =$ 4.026 × 10^-3, $a$ = 1, $θ_{0}$ =0, ${\dot{θ}}_{0}$ = 10deg/sec), (a) time-angle graph, (b) phase plane

Simulation result of BEP motion under principal parametric resonance condition ( $q = Y / l =$ 4.026 × 10^-3, $a$ = 1/4, $θ_{0}$ =0, ${\dot{θ}}_{0}$ = 10deg/sec) with modal damping of 0.8%, (a) time-angle graph, (b) phase plane