MATLABで学ぶ振動工学 基礎加振による応答

はじめに

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

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

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

 

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

今回は2.5.6項の「基礎加振による応答」のP.63から復習します。

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

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

基礎加振とは

モード解析入門の解説をする前に「基礎加振」について説明します。

基礎加振とは下図の基礎(地面もしくは床)の部分が振動することで、質量mが加振されることを意味します。

図 1自由度系の振動モデル(斜線部分が基礎や床と呼ばれている部分)

基礎加振による応答

下図のように基礎が振動する場合を考えましょう。

ここで、絶対空間内の質量mの変位をx、振動する基礎から見た質量の変位をy、基礎の振動をxbとします。

$$
y = x-x_b \tag{2.110}
$$

質量の慣性力をfm、ばねの復元力をfk、粘性抵抗力をfcとすると式(2.111)の関係式で表せます。

$$
f_m = -m \ddot{x} , f_k = -ky , f_c = -c \dot{y}  \tag{2.111}
$$

この1自由度のバネマス系には、外力は作用していないので、力は釣り合っている状態になります。従って、式(2.112)が得られます。

$$
m \ddot{x} +c \dot{y} +ky =0 \tag{2.112}
$$

これまでの式で式(2.112)を整理すると式(2.113)が得られます。

$$
m \ddot{x} +c \dot{x} +kx =c \dot{x_b} +kx_b = B(k+jcω)e{jωt} \tag{2.113}
$$

以下省略

最終的に、xとyは下式で表せます。

$$
x = |X|e^{j(ωt+atan(cω/k)+atan(-2ζβ/(1-β^2)))} \tag{2.117}
$$

$$
y = \frac{Bβ^2}{(1-β^2+2jζβ)}e{jωt} \tag{2.120}
$$

参考文献のモード解析入門には、不減衰(c=0、ζ=0)のときなどの場合分けをして、考察をしています。

参考文献のような場合分けは、頭のいい人なら理解できるのかもしれません。
しかし、より簡単に理解するには、MATLABで式をコード化して、とりあえずパラメータスタディするのが良いと思っています。

なので、とりあえずMATLABでパラメータスタディしてみましょう。
まず、ζをいじってみましょう。

図 ζ=0~1 @31Hz

位相が少し変化したぐらいですかね。。。。

clear all;clc;close all;


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

% zeta=c/cc;
zeta=0:0.1:1;
j=sqrt(-1);
B=1;
RGB=jet(length(zeta));
t=0:0.0001:1;
freq=31; %周波数


figure
for ii1=1:length(RGB)
    w=2*pi*freq; %角周波数
    beta=w/Omega;
    y=B*beta.^2./(1-beta.^2 + 2*j.*zeta(ii1)*beta).*exp(j*w.*t);
    plot(t,real(y),'color',RGB(ii1,:))
hold on
grid on
end
legend(split(num2str(zeta)))
xlim([0 0.3])
xlabel('Time [s]')
ylabel('Ampli. [m]')

そもそもyは基礎から見た質量の変位なので、あまり変化が見られないのかもしれませんね。
振動特性を把握するには、周波数次元での応答倍率を確認するのが良いですよね。

ただ、この応答倍率を導出するのは難しいですよね。
導出過程は理解できていませんが、基礎の変位加振振幅に対する変位応答倍率λbは下式で表せるようです。

$$
λ_b = \frac{|X|}{B} = \frac{B\sqrt{1+(2ζβ)^2}}{B\sqrt{(1-β^2)^2+(2jζβ)^2}} \tag{2.123}
$$

この式は前回紹介した下記式と同じですよね。

$$
λ = sqrt{\frac{1+(2ζβ)^2}{ 1 – β^2 + 2jζβ }} \tag{2.107}
$$

従って、基礎の変位加振振幅に対する変位応答倍率λbは、質量mに作用する力が基礎に伝搬する応答倍率λと等しいことがわかります。

では、簡単にζを変化させた時のλbを図示してみましょう。結果は下図になります。

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:20; %周波数
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)
Lamda=sqrt( (1+(2*zeta(ii1).*beta).^2) ./ ((1-beta.^2).^2+(2*zeta(ii1).*beta).^2) );
plot(beta,Lamda,'color',RGB(ii1,:))
hold on
grid on
end
plot([0 sqrt(2)],[1 1],'k--','linewidth',3)
plot([sqrt(2) sqrt(2)],[0 1],'k--','linewidth',3)
legend(split(num2str(zeta)))
text(1.5,1,'(√2 , 1)')
ylim([0 3])
xlim([0 3])
xlabel('β=ω/Ω')
ylabel('λ')

 

今日の解説はここまで。

Good luck

 

コメント