MATLABで学ぶ長松昭男著「モード解析入門」 2.5節:減衰系の強制振動

はじめに

振動工学やモード解析をご存じの方であれば、一度は読んだことがあると言われているバイブル的な書籍「モード解析入門」について再復習しております。

[商品価格に関しましては、リンクが作成された時点と現時点で情報が変更されている場合がございます。]

【中古】 モード解析入門 長松昭男 【古本・古書】
価格:9400円(税込、送料別) (2021/10/13時点)

 

再復習で学んだことをMATLABコードを示しながら解説したいと思います。

今回は2.5節の「減衰系の強制振動」のP.50から復習します。

まだ過去の記事を読んでいない方、この章や節以外を復習したい方は下記記事を参照ください。下記記事が目次となっており、各章や節にアクセスできるようにしています。

長松昭男著「モード解析入門」のMATLABコード付き解説
長松昭男先生の「モード解析入門」をMATLABコード付きで解説しています。

応答

減衰系の強制振動の運動方程式は式(2.87)に示します。

$$
m\ddot{x(t)} + c\dot{x(t)} +kx(t) = F e^{jωt} \tag{2.87}
$$

式(2.87)を満たす解は式(2.76)です。

$$
x(t) =X e^{jωt} \tag{2.76}
$$

式(2.76)を式(2.87)に代入して、整理すると式(2.89)を得られる。

$$
(-ω^2 \frac{m}{k}+jω \frac{c}{k} +1) X =\frac{F}{k} \tag{2.89}
$$

振動工学ではおなじみですが、下記に置き換えます。Ω=√(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.47)を式(2.49)に代入します。

$$
(-(\frac{ω}{Ω})^2 +2jζ\frac{ω}{Ω} +1) X =\frac{F}{k} = X_{st} \tag{2.91}
$$

ここで、下記のβというものを導入する。

$$
β=\frac{ω}{Ω}
$$

ここからは、振動工学ではおなじみの式になります。

$$
\frac{X}{X_{st}} = \frac{1}{ 1 – β^2 + 2jζβ } \tag{2.92}
$$

Xstって何?

Xstは静的な変位です

静的な変位って………?

フックの法則はF=kXですよね?
簡単に説明するために、これを静的な変位と考えましょう
この時の変位をXstとすると、

Xst=F/kですよね

式(2.92)は下記を意味します。
(動的な力が作用したときの変位)/(静的な力が作用したときの変位)

つまり、静的な変位に対して、動的な変位がどの程度の倍率なのかを直感的に把握することができます。

さっそく、式(2.92)を計算してみましょう。
結果は下図です。

上記の結果からわかることは下記です。
・角周波数ωが小さいときは、静的な変位と動的な変位は等しい。
・共振周波数Ωよりも少し小さい角周波数ωのときに、応答が最大となる。
・共振周波数Ωと角周波数ωが一致するときが最大応答ではない。
・共振周波数Ωと角周波数ωが一致するときに位相角は-90°になる。
・共振周波数Ωよりも角周波数ωがはるかに大きくなると、応答は0に漸近する。
・共振周波数Ωよりも角周波数ωがはるかに大きくなると、応答が遅れ位相角が-180°になる。

 

簡単な式なので、MATLABプログラムは不要かもしれませんが一応載せておきます。

clear all;clc;close all;


k=1000; %ばね定数
m=1; %質量
% % c=1; %減衰c
% % cc=2*sqrt(m*k); %臨界減衰係数
Omega=sqrt(k/m); %固有値

freq=0.001:0.001:10; %周波数
w=2*pi*freq; %角周波数

beta=w/Omega;
% zeta=c/cc;
zeta=0:0.1:1;
j=sqrt(-1);

RGB=jet(length(zeta));

figure
for ii1=1:length(RGB)
X_Xst=1./(1-beta.^2+2*j*zeta(ii1)*beta);
subplot(121)
plot(beta,abs(X_Xst),'color',RGB(ii1,:))
hold on
grid on
subplot(122)
plot(beta,angle(X_Xst),'color',RGB(ii1,:))
hold on
grid on
end
subplot(121)
legend(split(num2str(zeta)))
ylim([0 4])
xlabel('β=ω/Ω')
ylabel('|X|/X_{st}')
subplot(122)
legend(split(num2str(zeta)))
xlabel('β=ω/Ω')
ylabel('angle( |X|/X_{st} )')

 

今日の解説はここまで。

Good luck

 

コメント