%図13.3
%
clear all; close all;

%図13.3のプロット

%パラメータの値を与える
K = 1; %K=1
omegan = 10; %ω_n=10
zeta = 0.5; %ζ=0.5

%伝達関数の分子・分母多項式を与える
num = [ 0 0 K*omegan^2 ]; %ζ=0.5の場合の分子多項式
den = [ 1 2*zeta*omegan omegan^2 ]; %ζ=0.5の場合の分母多項式

%積分ゲインのパラメータを与える
Ki1 = 1.0; %K_i=1
Ki2 = 3.0; %K_i=3
Ki3 = 9.0; %K_i=9

%コントローラの分子・分母多項式を与える．
numc1 = [ 0 Ki1 ]; %K_i=1の場合の分子多項式
numc2 = [ 0 Ki2 ]; %K_i=3の場合の分子多項式
numc3 = [ 0 Ki3 ]; %K_i=9の場合の分子多項式
denc = [ 1 0 ]; %分母多項式

%制御対象とコントローラの伝達関数表現を与える
sys = tf( num, den ); %制御対象の伝達関数表現
c1 = tf( numc1, denc ); %K_i=1.0の場合のコントローラの伝達関数表現
c2 = tf( numc2, denc ); %K_i=3の場合のコントローラの伝達関数表現
c3 = tf( numc3, denc ); %K_i=9の場合のコントローラの伝達関数表現

%開ループ伝達関数を求める
sysL1 = c1*sys; %K_i=1の場合の開ループ伝達関数
sysL2 = c2*sys; %K_i=の3の場合の開ループ伝達関数
sysL3 = c3*sys; %K_i=9の場合の開ループ伝達関数

%角周波数の範囲を指定
w = logspace(-1, 3, 1000); %対数的に等間隔なベクトルの生成(10^{-1}から10^{3}で1000点)

%図13.3のプロット

%図13.3(a)
%ゲインと位相の計算
[gainL1 phaseL1 wL1] = bode(sysL1, w); %ゲインと位相
gainL1 = gainL1(:); %K_{i}=1の場合のゲインの配列を1次元ベクトルに変換
gainL1dB = 20*log10(gainL1); %K_{i}=1の場合のゲインをデシベル値に変換
phaseL1 = phaseL1(:); %K_{i}=1の場合の位相の配列を1次元ベクトルに変換

%図13.3(a)のプロット
subplot(2,3,1) %複数の図を並べるためのコマンド．2行3列の1行1列目という意味
semilogx(wL1, gainL1dB); %ゲイン線図をプロット
xlim([10^(-1),10^3]) %横軸（角周波数）の範囲の指定
ylim([-50 50]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1) 10^(2) 10^(3)]) %横軸の目盛りの値の設定
yticks([-40 -20 0 20 40]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
ylabel('Gain [dB]'); %縦軸のラベル表示
title('K_i = 1'); %タイトルの表示

subplot(2,3,4) %複数の図を並べるためのコマンド．2行3列の2行1列目という意味
semilogx(wL1, phaseL1); %位相線図をプロット
xlim([10^(-1),10^3]) %横軸（角周波数）の範囲の指定
ylim([-270 -90]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1) 10^(2) 10^(3)]) %横軸の目盛りの値の設定
yticks([-270 -225 -180 -135 -90]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
ylabel('Phase [deg]'); %縦軸のラベル表示

%図13.3(b)
%ゲインと位相の計算
[gainL2 phaseL2 wL2] = bode(sysL2, w); %ゲインと位相
gainL2 = gainL2(:); %K_i=3の場合のゲインの配列を1次元ベクトルに変換
gain2dB = 20*log10(gainL2); %K_i=3の場合のゲインをデシベル値に変換
phaseL2 = phaseL2(:); %K_i=3の場合の位相の配列を1次元ベクトルに変換

%図13.3(b)のプロット
subplot(2,3,2) %複数の図を並べるためのコマンド．2行3列の1行2列目という意味
semilogx(wL2, gain2dB); %ゲイン線図をプロット
xlim([10^(-1),10^3]) %横軸（角周波数）の範囲の指定
ylim([-50 50]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1) 10^(2) 10^(3)]) %横軸の目盛りの値の設定
yticks([-40 -20 0 20 40]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
%ylabel('Gain [dB]'); %縦軸のラベル表示
title('K_i = 3'); %タイトルの表示

subplot(2,3,5) %複数の図を並べるためのコマンド．2行3列の2行2列目という意味
semilogx(wL2, phaseL2); %位相線図をプロット
xlim([10^(-1),10^3]) %横軸（角周波数）の範囲の指定
ylim([-270 -90]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1) 10^(2) 10^(3)]) %横軸の目盛りの値の設定
yticks([-270 -225 -180 -135 -90]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
%ylabel('Phase [deg]'); %縦軸のラベル表示

%図13.3(c)
%ゲインと位相の計算
[gainL3 phaseL3 wL3] = bode(sysL3, w); %ゲインと位相
gainL3 = gainL3(:); %K_i=9の場合のゲインの配列を1次元ベクトルに変換
gainL3dB = 20*log10(gainL3); %K_i=9の場合のゲインをデシベル値に変換
phaseL3 = phaseL3(:); %K_i=9の場合の位相の配列を1次元ベクトルに変換

%図13.3(c)のプロット
subplot(2,3,3) %複数の図を並べるためのコマンド．2行3列の1行3列目という意味
semilogx(wL3, gainL3dB); %ゲイン線図をプロット
xlim([10^(-1),10^3]) %横軸（角周波数）の範囲の指定
ylim([-50 50]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1) 10^(2) 10^(3)]) %横軸の目盛りの値の設定
yticks([-40 -20 0 20 40]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
%ylabel('Gain [dB]'); %縦軸のラベル表示
title('K_i = 9'); %タイトルの表示

subplot(2,3,6) %複数の図を並べるためのコマンド．2行3列の2行3列目という意味
semilogx(wL3, phaseL3); %位相線図をプロット
xlim([10^(-1),10^3]) %横軸（角周波数）の範囲の指定
ylim([-270 -90]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1) 10^(2) 10^(3)]) %横軸の目盛りの値の設定
yticks([-270 -225 -180 -135 -90]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
%ylabel('Phase [deg]'); %縦軸のラベル表示

Published with MATLAB® R2016b