MATLABで学ぶ振動工学 無周期運動/減衰振動

はじめに

振動工学やモード解析をご存じの方であれば、一度は読んだことがあると言われているバイブル的な書籍「モード解析入門」を再復習しております。今回はモード解析入門の2.3節の減衰系の自由振動について復習します。なお、振動工学の参考書はどの本も似たり寄ったりで内容が重複しています。例えば、「実用モード解析入門」は長松先生の息子さんが著者ということもあり、重複している部分が多いです。

各参考書で重複している部分は本記事で大まかに理解できると思いますので、他の参考書から本記事へアクセスした場合も最後まで読んでいただけると幸いです。

参考書解説一覧
MATLABとNVH(Noise, Vibration and Harshness)をこよなく愛する人生の先輩、それが私(MATLABパイセン)です。NVHを極めるにあたり、周辺の知識も網羅的に勉強することになり、その知識を共有すべく、本HPを運営しています。日本の製造業を応援すべく、機械エンジニアの「車輪の再開発」の防止と業務効率化の手助けをモットーに活動。専門である振動騒音工学や音響工学について、プログラムを示しながら解説。大学の授業よりもわかりやすく説明することを目指す。質問頂ければ、回答したいと思います。

運動方程式

粘性減衰項を含む運動方程式は式(2.4)となります。

$$
m\ddot{x(t)} + c\dot{x(t)} +kx(t) = 0 \tag{2.4}
$$

ちなみに、粘性減衰は速度に比例する減衰で、FEM解析でよく使われる構造減衰は変位に比例する減衰です。なので、どの減衰を使うかによって運動方程式はかわります。
初歩的な機械系の書籍では「粘性減衰」を単に「減衰」と記述しているものもあります。

モード解析入門の2.3節では速度に比例する粘性減衰を用いて、「減衰系の自由振動」を説明しているということを頭の片隅にでも覚えておいてください。

式(2.4)を満たす解は式(2.42)です。

$$
x =x_0 e^{λt} \tag{2.42}
$$

式(2.42)を式(2.4)に代入すると式(2.45)が得られます。

\[
\left\{
\begin{array}{l}
x_0&=0\\
mλ^2+cλ+k&=0
\end{array}
\right.
\tag{2.45}
\]

x0=0は系が静止していることを意味するので、振動した状態であるためには式(2.45)の下式が成立している必要があります。

式(2.45)の下式をλについて解くと、中学校で習う「解の公式」から式(2.46)が得られます。

\begin{align*}
x & = {-c \pm \sqrt{c^2-4mk} \over 2m} \\
& = -{c \over 2m}\pm \sqrt{({c \over 2m})^2-{k \over m}} \\
& = \sqrt{{k \over m}}({-c \over {2 \sqrt{mk}}}) \pm \sqrt{{k \over m}} \sqrt{{c^2 \over 4mk}-1}  \\
\tag{2.46}
\end{align*}

振動工学ではおなじみですが、下記に置き換えます。Ω=√(k/m)は減衰のない角共振振動数ですよね。なので、Ωを不減衰固有角振動数と呼びます。

\[
\left\{
\begin{array}{l}
c_c=2 \sqrt{mk}=2m\sqrt{\frac{k}{m}}=2mΩ\\
ζ=c/c_c
\end{array}
\right.
\tag{2.47}
\]

すると、式(2.46)は式(2.48)で書き換えられます。

$$
λ = -Ωζ \pm Ω\sqrt{ζ^2-1} \tag{2.48}
$$

\[
\left\{
\begin{array}{l}
λ_1 = -Ωζ + Ω\sqrt{ζ^2-1} \\
λ_2 = -Ωζ – Ω\sqrt{ζ^2-1}
\end{array}
\right.
\tag{2.48}
\]

式(2.48)の一般的な解は式(2.49)となる。ここで、X1とX2は未定係数です。

$$
x = X_1 e^{λ_1 t} + X_2 e^{λ_2 t} \tag{2.49}
$$

この式(2.49)は式(2.28)の平方根の中身の正負によって、大きく現象が変化するので、場合分けして考えることになります。

無周期運動(式(2.48)の平方根が正のとき)

ζ≥1の場合

λ1とλ2は負になります。
なので、とりあえず下記がどうなるのかだけを考えてみます。
$$
X_1 e^{λ_1 t}
$$

高校生の知識でy=e^(-t)はt=0のときy=1なのは容易にわかりますね。なので、X1はt=0のときの切片に相当すると思われます。試しに、X1を1~10の間で変化させてみます。
(λはとりあえず-10とします。)

予想どおりですね。

clear all;close all;

t=0:0.001:1;
L1=-10;
X1=1:10;

X1_eL1t=zeros(length(X1),length(t));
RGB=jet(length(X1));
figure
for ii1=1:length(X1)
X1_eL1t=X1(ii1)*exp(L1*t);
plot(t,X1_eL1t,'color',RGB(ii1,:))
hold on
end
legend(split(num2str(X1)))

では、λが変化したときどうなるか、0に収束する時間が変わるだけですよね。
高校生の知識でわかりますよね。

やっぱりそうですよね。

clear all;close all;

t=0:0.01:100;
L1=[-10 -1 -0.1 -0.01];
X1=1;

X1_eL1t=zeros(length(L1),length(t));
RGB=jet(length(L1));
figure
for ii1=1:length(L1)
X1_eL1t=X1*exp(L1(ii1)*t);
plot(t,X1_eL1t,'color',RGB(ii1,:))
hold on
end
legend(split(num2str(L1)))

つまり、ζ≥1の場合は振動しません。

今日の解説はここまで。

Good luck

 

コメント