2次遅れ系は、振動的な応答を示すシステムの基本モデルです。ばね-ダンパ系、RLC回路、モータ制御系など多くの工学システムが2次系で記述されます。減衰比 $\zeta$ と固有振動数 $\omega_n$ の2つのパラメータで応答特性が完全に決まります。
本記事では、2次遅れ系の伝達関数、減衰比による応答の分類、過渡応答特性の計算公式、Pythonでの可視化までを解説します。
本記事の内容
- 2次遅れ系の標準形と2つのパラメータ
- 減衰比 $\zeta$ による応答の分類
- オーバーシュート・整定時間・立ち上がり時間の公式
- Pythonでの可視化
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
2次遅れ系の伝達関数
2次遅れ系の標準形は
$$ \begin{equation} G(s) = \frac{\omega_n^2}{s^2 + 2\zeta\omega_n s + \omega_n^2} \end{equation} $$
ここで
- $\omega_n$ : 固有角振動数（自然振動数）[rad/s]
- $\zeta$ : 減衰比（ダンピング比）[無次元]
定常ゲインは $G(0) = 1$ なので、単位ステップ入力に対する定常値は $1$ です。
極の位置
特性方程式 $s^2 + 2\zeta\omega_n s + \omega_n^2 = 0$ の解は
$$ \begin{equation} s = -\zeta\omega_n \pm \omega_n\sqrt{\zeta^2 – 1} \end{equation} $$
減衰比による応答の分類
不足減衰（$0 < \zeta < 1$）
$\zeta^2 – 1 < 0$ のため、極は複素共役になります。
$$ s = -\zeta\omega_n \pm j\omega_n\sqrt{1 – \zeta^2} = -\sigma \pm j\omega_d $$
ここで $\sigma = \zeta\omega_n$（減衰係数）、$\omega_d = \omega_n\sqrt{1 – \zeta^2}$（減衰固有振動数）です。
ステップ応答は
$$ \begin{equation} y(t) = 1 – \frac{e^{-\sigma t}}{\sqrt{1-\zeta^2}} \sin(\omega_d t + \phi) \quad (t \geq 0) \end{equation} $$
ここで $\phi = \arctan\frac{\sqrt{1-\zeta^2}}{\zeta}$ です。振動しながら定常値に収束します。
臨界減衰（$\zeta = 1$）
極は重根 $s = -\omega_n$（2重極）です。
$$ y(t) = 1 – (1 + \omega_n t) e^{-\omega_n t} \quad (t \geq 0) $$
振動せずに最も速く定常値に収束する減衰条件です。
過減衰（$\zeta > 1$）
極は2つの異なる負の実数になります。
$$ s_1 = -\zeta\omega_n + \omega_n\sqrt{\zeta^2 – 1}, \quad s_2 = -\zeta\omega_n – \omega_n\sqrt{\zeta^2 – 1} $$
$$ y(t) = 1 – \frac{s_2 e^{s_1 t} – s_1 e^{s_2 t}}{s_2 – s_1} \quad (t \geq 0) $$
振動せずにゆっくりと定常値に収束します。臨界減衰より遅い応答です。
過渡応答特性
不足減衰（$0 < \zeta < 1$）の場合の主要な過渡応答特性を求めます。
オーバーシュート $M_p$
ステップ応答の最大値が定常値をどれだけ超えるかを百分率で表します。
$$ \frac{dy}{dt} = 0 \text{ の最初の解は } t_p = \frac{\pi}{\omega_d} = \frac{\pi}{\omega_n\sqrt{1-\zeta^2}} $$
このとき
$$ \begin{equation} M_p = e^{-\frac{\zeta\pi}{\sqrt{1-\zeta^2}}} \times 100 \quad [\%] \end{equation} $$
オーバーシュートは減衰比 $\zeta$ のみで決まり、$\omega_n$ には依存しません。
| $\zeta$ | $M_p$ |
|---|---|
| 0.1 | 72.9% |
| 0.3 | 37.2% |
| 0.5 | 16.3% |
| 0.7 | 4.6% |
| 0.9 | 0.2% |
ピーク時間 $t_p$
応答が最大値に達する時間です。
$$ \begin{equation} t_p = \frac{\pi}{\omega_d} = \frac{\pi}{\omega_n\sqrt{1-\zeta^2}} \end{equation} $$
立ち上がり時間 $t_r$
応答が定常値の 10% から 90%（または 0% から 100%）に達するまでの時間です。0% から 100% の定義では
$$ \begin{equation} t_r = \frac{\pi – \phi}{\omega_d} = \frac{\pi – \arctan\frac{\sqrt{1-\zeta^2}}{\zeta}}{\omega_n\sqrt{1-\zeta^2}} \end{equation} $$
整定時間 $t_s$
応答が定常値の $\pm 2\%$（または $\pm 5\%$）以内に収まる時間です。包絡線 $\frac{e^{-\sigma t}}{\sqrt{1-\zeta^2}}$ が許容範囲内に入る条件から
$$ \begin{equation} t_s \approx \frac{4}{\zeta\omega_n} \quad (2\% \text{ 基準}) \end{equation} $$
$$ t_s \approx \frac{3}{\zeta\omega_n} \quad (5\% \text{ 基準}) $$
Pythonでの実装
減衰比による応答の変化
import numpy as np
import matplotlib.pyplot as plt
import control as ctrl
wn = 2.0 # 固有角振動数
zetas = [0.1, 0.3, 0.5, 0.7, 1.0, 1.5, 2.0]
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))
t = np.linspace(0, 10, 1000)
for zeta in zetas:
G = ctrl.tf([wn**2], [1, 2*zeta*wn, wn**2])
t_out, y_out = ctrl.step_response(G, T=t)
label = f'ζ={zeta}'
if zeta < 1:
label += ' (underdamped)'
elif zeta == 1:
label += ' (critical)'
else:
label += ' (overdamped)'
ax1.plot(t_out, y_out, linewidth=2, label=label)
ax1.axhline(y=1.0, color='k', linestyle='--', alpha=0.3)
ax1.set_xlabel('Time [s]')
ax1.set_ylabel('Output y(t)')
ax1.set_title(f'Step Response (ωn={wn} rad/s)')
ax1.legend(fontsize=8)
ax1.grid(True, alpha=0.3)
ax1.set_xlim([0, 10])
# 極の配置
colors = plt.cm.viridis(np.linspace(0, 0.9, len(zetas)))
for zeta, color in zip(zetas, colors):
G = ctrl.tf([wn**2], [1, 2*zeta*wn, wn**2])
poles = ctrl.poles(G)
ax2.plot(poles.real, poles.imag, 'x', markersize=10,
markeredgewidth=2, color=color, label=f'ζ={zeta}')
# 固有振動数の円
theta = np.linspace(0, 2*np.pi, 100)
ax2.plot(wn*np.cos(theta), wn*np.sin(theta), 'k:', alpha=0.3, label=f'|s|=ωn={wn}')
ax2.axvline(x=0, color='k', linewidth=0.5)
ax2.axhline(y=0, color='k', linewidth=0.5)
ax2.set_xlabel('Real')
ax2.set_ylabel('Imaginary')
ax2.set_title('Pole Locations')
ax2.legend(fontsize=7, loc='upper left')
ax2.grid(True, alpha=0.3)
ax2.set_aspect('equal')
ax2.set_xlim([-5, 1])
ax2.set_ylim([-3, 3])
plt.tight_layout()
plt.show()
過渡応答特性の計算と可視化
import numpy as np
import matplotlib.pyplot as plt
import control as ctrl
wn = 3.0
zeta = 0.3
G = ctrl.tf([wn**2], [1, 2*zeta*wn, wn**2])
# 過渡応答特性の計算
wd = wn * np.sqrt(1 - zeta**2)
sigma = zeta * wn
Mp = np.exp(-zeta * np.pi / np.sqrt(1 - zeta**2)) * 100
tp = np.pi / wd
tr = (np.pi - np.arctan2(np.sqrt(1 - zeta**2), zeta)) / wd
ts = 4 / (zeta * wn)
print(f"減衰固有振動数 ωd = {wd:.3f} rad/s")
print(f"オーバーシュート Mp = {Mp:.1f}%")
print(f"ピーク時間 tp = {tp:.3f} s")
print(f"立ち上がり時間 tr = {tr:.3f} s")
print(f"整定時間 ts(2%) = {ts:.3f} s")
# 可視化
t = np.linspace(0, 6, 1000)
t_out, y_out = ctrl.step_response(G, T=t)
plt.figure(figsize=(10, 6))
plt.plot(t_out, y_out, 'b', linewidth=2, label='Step response')
# 定常値
plt.axhline(y=1.0, color='k', linestyle='--', alpha=0.3)
# ±2%帯
plt.axhline(y=1.02, color='gray', linestyle=':', alpha=0.5)
plt.axhline(y=0.98, color='gray', linestyle=':', alpha=0.5)
plt.fill_between(t, 0.98, 1.02, alpha=0.1, color='green', label='±2% band')
# オーバーシュート
y_peak = 1 + Mp / 100
plt.plot(tp, y_peak, 'ro', markersize=8, zorder=5)
plt.annotate(f'Mp={Mp:.1f}%\ntp={tp:.2f}s',
xy=(tp, y_peak), xytext=(tp + 0.3, y_peak + 0.05),
fontsize=10, arrowprops=dict(arrowstyle='->', color='red'))
# 立ち上がり時間
plt.axvline(x=tr, color='orange', linestyle='--', alpha=0.7, label=f'tr={tr:.2f}s')
# 整定時間
plt.axvline(x=ts, color='purple', linestyle='--', alpha=0.7, label=f'ts(2%)={ts:.2f}s')
# 包絡線
t_env = np.linspace(0, 6, 500)
env_upper = 1 + np.exp(-sigma * t_env) / np.sqrt(1 - zeta**2)
env_lower = 1 - np.exp(-sigma * t_env) / np.sqrt(1 - zeta**2)
plt.plot(t_env, env_upper, 'r:', alpha=0.4, linewidth=1)
plt.plot(t_env, env_lower, 'r:', alpha=0.4, linewidth=1)
plt.xlabel('Time [s]')
plt.ylabel('Output y(t)')
plt.title(f'2nd Order System: ωn={wn}, ζ={zeta}')
plt.legend(loc='upper right', fontsize=9)
plt.grid(True, alpha=0.3)
plt.xlim([0, 6])
plt.tight_layout()
plt.show()
オーバーシュートと減衰比の関係
import numpy as np
import matplotlib.pyplot as plt
zeta = np.linspace(0.01, 0.99, 200)
Mp = np.exp(-zeta * np.pi / np.sqrt(1 - zeta**2)) * 100
plt.figure(figsize=(8, 5))
plt.plot(zeta, Mp, 'b', linewidth=2)
plt.xlabel('Damping ratio ζ')
plt.ylabel('Overshoot Mp [%]')
plt.title('Overshoot vs Damping Ratio')
plt.grid(True, alpha=0.3)
# 代表的な値をマーク
for z in [0.1, 0.3, 0.5, 0.7, 0.9]:
mp = np.exp(-z * np.pi / np.sqrt(1 - z**2)) * 100
plt.plot(z, mp, 'ro', markersize=6)
plt.annotate(f'({z}, {mp:.1f}%)', xy=(z, mp),
xytext=(z + 0.05, mp + 3), fontsize=9)
plt.xlim([0, 1])
plt.ylim([0, 80])
plt.tight_layout()
plt.show()
まとめ
本記事では、2次遅れ系の応答特性について解説しました。
- 2次遅れ系は $G(s) = \frac{\omega_n^2}{s^2 + 2\zeta\omega_n s + \omega_n^2}$ で、$\zeta$ と $\omega_n$ の2パラメータで特性が決まる
- $0 < \zeta < 1$ は不足減衰（振動）、$\zeta = 1$ は臨界減衰、$\zeta > 1$ は過減衰
- オーバーシュートは $M_p = e^{-\zeta\pi / \sqrt{1-\zeta^2}} \times 100\%$ で $\zeta$ のみに依存
- 整定時間は $t_s \approx 4 / (\zeta\omega_n)$ で、$\zeta$ と $\omega_n$ の積に反比例
次のステップとして、以下の記事も参考にしてください。