Modeling Basic - 2 (damped pendulum)
1. Introduction
Simple pendulum model 에서는 감쇠가 없는 Simple Pendulum에 대해서 운동방정식을 도출하고 해를 계산해 봤다. 이번 페이지에서는 Damped Pendulum에 대해서 운동방정식을 도출하고, 해를 구해보자. 이전 페이지를 참고 없이, 본 페이지만 보고도 알 수 있도록 식을 다시 작성하였다.
2. Modeling
Damped Pendulum에 작용하는 힘은 Fig 1. 처럼나타낼 수 있다.
이 때, 진자의 길이는 $l[\mathrm{m}]$, 질량은 $m[\mathrm{kg}]$, 중력 가속도는 $g[\mathrm{m/s^2}]$,
수직선과 진자가 이루는 각도는 $\theta[\mathrm{rad}]$, 저항력은 $F_d[\mathrm{N}]$이다.
2.1 Derivation of the Differential Equation
미분방정식을 구하기 위해서는 뉴턴의 제 2법칙 $F=ma$을 적용해 보면, 식 (2.1)을 만족해야 한다.
\[ml\ddot{\theta} = -mg\sin{\theta} - F_d\tag{2.1}\]$F_d$는 속도에 따라 변하므로 식 (2.2)처럼 나타낼 수 있다.
\[F_d = b\dot{\theta} \tag{2.2}\]여기에서 $b[\mathrm{Ns/m}]$는 damping coefficient를 의미한다. Damping coefficient는 공기저항, 마찰 등을 의미한다.
식을 정리하기 위해 식 (2.1)의 양변을 $ml$로 나누고 식 (2.2)를 대입하면, 식 (2.3) 과 같이 표현할 수 있다.
\[\ddot{\theta} + \frac{b}{ml}\dot{\theta} + \frac{g}{l} \sin\theta = 0 \tag{2.3}\]2.2 Linearization of the pendulum model
각도가 작을 때 ($\theta \approx 0$), Taylor Expression를 이용하면 식 (2.4)과 같이 나타낼 수 있다;.
\[\sin\theta \approx \theta \tag{2.4}\]따라서, 식 (2.5)과 같이 선형 미분방정식 형태로 나타낼 수 있다.
\[\ddot{\theta} + \frac{b}{ml}\dot{\theta} + \frac{g}{l} \theta = 0 \tag{2.5}\]식 (2.5), 를 보면 미분방정식에 질량 $m$이 $\dot{\theta}$ 항에 포함되어 있는 것을 알 수 있다. 이러한 damping coefficient가 있기 때문에 시간이 지나면 $\theta$는 0에 수렴한다.
3. Solution of the Odinary Differential Equation
위에서 도출한 선형 미분방정식 식 (3.1)의 해를 구해보자. 이 방정식은 감쇠가 포함된 선형 2차 미분방정식이다.
\[\ddot{\theta} + \frac{b}{ml} \dot{\theta} + \frac{g}{l} \theta = 0 \tag{3.1}\]특성 방정식의 해를 구하면, $r$은 식 (3.2)처럼 산출할 수 있다.
\[r^2 + \frac{b}{ml} r + \frac{g}{l} = 0 \tag{3.2}\]근을 구해보면, 식 (3.3)과 같다.
\[r_{1,2} = \frac{-\frac{b}{ml} \pm \sqrt{\left(\frac{b}{ml}\right)^2 - 4 \frac{g}{l}}}{2} \tag{3.3}\]$r_{1,2}$가 허근, 중근, 실근을 갖는 3가지 경우를 생각해 볼 수 있다.
3.1 Underdamped System ($\frac{b^2}{m^2l^2} < 4\frac{g}{l}$)
Underdamped system은 진동하면서 감쇠하는 시스템이며, 가장 일반적인 형태이다. $r$이 허수일 경우를 의미한다. $\theta$의 시간에 대한 변화는 식 (3.4)와 같다.
\[\theta(t) = e^{-\frac{b}{2ml} t} (C_1 \cos\omega' t + C_2 \sin\omega' t) \tag{3.4}\]이 때, $\omega’$는 감쇠 주파수를 의미하며 식 (3.5)와 같다.
\[\omega' = \sqrt{\frac{g}{l} - \left(\frac{b}{2ml}\right)^2} \tag{3.5}\]3.2 Critical Damped System ($\frac{b^2}{m^2l^2} = 4\frac{g}{l}$)
Critical system은 진동하지 않고 감쇠하는 시스템이며 $r$이 중근일 경우를 의미한다. 가장 빠르게 멈추는 특징이 있다. $\theta$의 시간에 대한 변화는 식 (3.6)와 같다.
\[\theta(t) = (C_1 + C_2 t) e^{-\frac{b}{2ml} t} \tag{3.6}\]3.3 Overdamped System ($\frac{b^2}{m^2l^2} > 4\frac{g}{l}$)
Overdamped system은 진동하지 않고 감쇠하는 시스템이며 $r$이 실근일 경우를 의미한다. $\theta$의 시간에 대한 변화는 식 (3.7)와 같다.
\[\theta(t) = C_1 e^{\displaystyle{\left(-\frac{b}{2ml} + \sqrt{\left(\frac{b}{2ml}\right)^2 - \frac{g}{l}} \right)t}} + C_2 e^{\displaystyle{\left(-\frac{b}{2ml} - \sqrt{\left(\frac{b}{2ml}\right)^2 - \frac{g}{l}} \right)t}} \tag{3.7}\]4. Simulation result
Fig 2. ~ Fig 4.는 matlab을 이용하여 식 (3.1)을 시뮬레이션한 결과이다. Damping coefficient는 각각 Fig 2.는 $b=m\sqrt{4gl}-50[\mathrm{Ns/m}]$, Fig 3.은 $b=m\sqrt{4gl}[\mathrm{Ns/m}]$, Fig 4.는 $b=m\sqrt{4gl}+50[\mathrm{Ns/m}]$로 설정하였다. 공통적으로는 $g=9.81[\mathrm{m/s^2}]$, $l=5[\mathrm{m}]$, 초기 값으로 $\theta = 1 [\mathrm{rad}]$, $\dot{\theta} = 0[\mathrm{rad/s}]$ 로 설정하였다.
