「周波数応答」の版間の差分

編集の要約なし
編集の要約なし
16行目: 16行目:




図\ref{sin_res_2nd}は,伝達関数<math>\displaystyle \frac{1}{s^2 + 0.5s + 1}</math>に<math>\omega=</math>0.1, 1, 10 [rad/s]の3種類の正弦波を入力した時の応答をMATLAB(Octave)の\texttt{lsim}関数で計算した結果である(リストsample2_1).
図\ref{sin_res_2nd}は,伝達関数<math>\displaystyle \frac{1}{s^2 + 0.5s + 1}</math>に<math>\omega=</math>0.1, 1, 10 [rad/s]の3種類の正弦波を入力した時の応答をMATLAB(Octave)の<tt>lsim</tt>関数で計算した結果である(リストsample2_1).


グラフから,過渡状態はおおよそ開始から10秒強で,その後は定常状態となっていることがわかる.定常状態でのゲインと位相は周波数によって変化していることがわかる.
グラフから,過渡状態はおおよそ開始から10秒強で,その後は定常状態となっていることがわかる.定常状態でのゲインと位相は周波数によって変化していることがわかる.
93行目: 93行目:
=== MATLAB(Octave)での複素数の取扱い ===
=== MATLAB(Octave)での複素数の取扱い ===


MATLAB(Octave)では複素数は以下のように入力し,そのまま実数と同じように扱える.なお,入力する際は,虚数単位は<tt>i</tt>でもjでもよい.
MATLAB(Octave)では複素数は以下のように入力し,そのまま実数と同じように扱える.なお,入力する際は,虚数単位は<tt>i</tt>でも<tt>j</tt>でもよい.


  1 + 2i
  1 + 2i
103行目: 103行目:


と出力される.
と出力される.
しかし,変数として\texttt{i}\texttt{j}を使っていると,虚数単位と混同するので,複素数を表す関数\texttt{complex()}を使用した方が間違いがない.また,絶対値は関数\texttt{abs()},偏角は関数\texttt{angle()}で求められる(単位はラジアン).以下のこれらの使用例を示す.
しかし,変数として<tt>i</tt><tt>j</tt>を使っていると,虚数単位と混同するので,複素数を表す関数<tt>complex()</tt>を使用した方が間違いがない.また,絶対値は関数<tt>abs()</tt>,偏角は関数<tt>angle()</tt>で求められる(単位はラジアン).以下のこれらの使用例を示す.


  > z = complex(1, -1)
  > z = complex(1, -1)
114行目: 114行目:
=== MATLAB(Octave)での有理式の取扱い ===
=== MATLAB(Octave)での有理式の取扱い ===


伝達関数は一般に有理式(分子・分母が多項式)で表される(むだ時間など一部除く).MATLAB(Octave)では多項式を係数を降べきの順に並べたベクトルで表し,有理式は分子・分母それぞれの係数ベクトルのペアで表現する.例えば,伝達関数<math>\displaystyle G(s) = \frac{s + 2}{s^2 + 3s + 2} </math>は,次のように二つのベクトル\texttt{num}\texttt{den}のペアで表す.
伝達関数は一般に有理式(分子・分母が多項式)で表される(むだ時間など一部除く).MATLAB(Octave)では多項式を係数を降べきの順に並べたベクトルで表し,有理式は分子・分母それぞれの係数ベクトルのペアで表現する.例えば,伝達関数<math>\displaystyle G(s) = \frac{s + 2}{s^2 + 3s + 2} </math>は,次のように二つのベクトル<tt>num</tt><tt>den</tt>のペアで表す.


  > num = [1 2];
  > num = [1 2];
  > den = [1 3 2];
  > den = [1 3 2];


また,MATLAB(Octave)では,多項式の変数に値を代入して計算した結果を返してくれる関数\texttt{polyval()}が用意されている.使用例を以下に示す.このように複素数にも対応するので,周波数応答の計算に利用できる.
また,MATLAB(Octave)では,多項式の変数に値を代入して計算した結果を返してくれる関数<tt>polyval()</tt>が用意されている.使用例を以下に示す.このように複素数にも対応するので,周波数応答の計算に利用できる.


  > f = [1 2 1];
  > f = [1 2 1];
132行目: 132行目:
ここではMATLAB(Octave)によるボード線図のプロットについて,以下の二つの方法を紹介する.
ここではMATLAB(Octave)によるボード線図のプロットについて,以下の二つの方法を紹介する.


* \texttt{polyval}関数を利用して原理に基づいて計算する方法
* <tt>polyval</tt>関数を利用して原理に基づいて計算する方法
* Control System Toolboxの\texttt{bode}関数を利用する方法
* Control System Toolboxの<tt>bode</tt>関数を利用する方法




==== \texttt{polyval}関数を利用した方法 ====
==== <tt>polyval</tt>関数を利用した方法 ====


リストsample2_2に\texttt{polyval}関数を利用してボード線図をプロットするスクリプト例を示す.
リストsample2_2に<tt>polyval</tt>関数を利用してボード線図をプロットするスクリプト例を示す.


変数\texttt{minp}\texttt{maxp}にそれぞれ周波数軸の下限と上限を10のベキ数で指定する.この例では<math>10^{-2}=0.01</math>[rad/s]から<math>10^2=100</math>[rad/s]までを指定している.
変数<tt>minp</tt><tt>maxp</tt>にそれぞれ周波数軸の下限と上限を10のベキ数で指定する.この例では<math>10^{-2}=0.01</math>[rad/s]から<math>10^2=100</math>[rad/s]までを指定している.


グラフ表示のウィンドウを\texttt{subplot}命令で2行1列に分割し,上側にゲイン線図,下側に位相線図をそれぞれ横軸対数表示である\texttt{semilogx}関数を使ってプロットしている.図Bode_2_2に実行結果のグラフを示す.
グラフ表示のウィンドウを<tt>subplot</tt>命令で2行1列に分割し,上側にゲイン線図,下側に位相線図をそれぞれ横軸対数表示である<tt>semilogx</tt>関数を使ってプロットしている.図Bode_2_2に実行結果のグラフを示す.


sample2_2 \texttt{polyval}関数を利用したボード線図の描画スクリプト}
sample2_2 <tt>polyval</tt>関数を利用したボード線図の描画スクリプト}
  num = [1]; %伝達関数の分子
  num = [1]; %伝達関数の分子
  den = [1 0.5 1]; %伝達関数の分母
  den = [1 0.5 1]; %伝達関数の分母
175行目: 175行目:


Bode_2_2.pdf
Bode_2_2.pdf
<math>\displaystyle \frac{1}{s^2 + 0.5s+ 1}</math>のBode線図(\texttt{polyval}関数利用)
<math>\displaystyle \frac{1}{s^2 + 0.5s+ 1}</math>のBode線図(<tt>polyval</tt>関数利用)






リストsample2_3は,MATLAB(Octave)固有の要素ごとの演算を利用して\texttt{for}ループを用いないスクリプト例である.ここで,\texttt{logspace}関数は等比数列を出力する関数である.
リストsample2_3は,MATLAB(Octave)固有の要素ごとの演算を利用して<tt>for</tt>ループを用いないスクリプト例である.ここで,<tt>logspace</tt>関数は等比数列を出力する関数である.


sample2_3 \texttt{polyval}関数を利用したボード線図の描画スクリプト(forループ無しバージョン)
sample2_3 <tt>polyval</tt>関数を利用したボード線図の描画スクリプト(forループ無しバージョン)
  num = [1]; %伝達関数の分子
  num = [1]; %伝達関数の分子
  den = [1 0.5 1]; %伝達関数の分母
  den = [1 0.5 1]; %伝達関数の分母
206行目: 206行目:
  grid('on');
  grid('on');


==== \texttt{bode}関数を利用した方法 ====
==== <tt>bode</tt>関数を利用した方法 ====


MATLABのControl System Toolboxに含まれる\texttt{bode}関数を利用すると,簡単にボード線図を描画できる.リストsample2_4にスクリプト例を,図Bode_2_4に実行結果のグラフを示す.なお,この例では変数\texttt{omega}で周波数ベクトルを指定しているが,省略することもできる.その場合はbode関数が自動的に周波数範囲を設定して描画する.
MATLABのControl System Toolboxに含まれる<tt>bode</tt>関数を利用すると,簡単にボード線図を描画できる.リストsample2_4にスクリプト例を,図Bode_2_4に実行結果のグラフを示す.なお,この例では変数<tt>omega</tt>で周波数ベクトルを指定しているが,省略することもできる.その場合は<tt>bode</tt>関数が自動的に周波数範囲を設定して描画する.


sample2_4 \texttt{bode}関数を利用したボード線図の描画スクリプト
sample2_4 <tt>bode</tt>関数を利用したボード線図の描画スクリプト
  num = [1]; %伝達関数の分子
  num = [1]; %伝達関数の分子
  den = [1 1 1]; %伝達関数の分母
  den = [1 1 1]; %伝達関数の分母
224行目: 224行目:


Bode_2_4.pdf
Bode_2_4.pdf
<math>\displaystyle \frac{1}{s^2 + 0.5s+ 1}</math>のBode線図(\texttt{bode}関数利用)
<math>\displaystyle \frac{1}{s^2 + 0.5s+ 1}</math>のBode線図(<tt>bode</tt>関数利用)


リストsample2_4の例のように\texttt{bode}関数単体でグラフまで描画するが,ゲイン線図だけを描画したい時や,実験結果と重ねてグラフ化したい時などは,次のように,\texttt{bode}関数の左辺に変数をおくと,グラフ表示は行わず周波数応答の計算結果だけが得られる(ゲインは絶対値,位相は度の単位).
リストsample2_4の例のように<tt>bode</tt>関数単体でグラフまで描画するが,ゲイン線図だけを描画したい時や,実験結果と重ねてグラフ化したい時などは,次のように,<tt>bode</tt>関数の左辺に変数をおくと,グラフ表示は行わず周波数応答の計算結果だけが得られる(ゲインは絶対値,位相は度の単位).


  [gain phase omega] = bode(sys, omega);
  [gain phase omega] = bode(sys, omega);
233行目: 233行目:
これを利用したスクリプト例をリストsample2_5に示す.これはボード線図のゲイン線図だけをプロットする例である.実行結果を図Bode_2_5に示す.
これを利用したスクリプト例をリストsample2_5に示す.これはボード線図のゲイン線図だけをプロットする例である.実行結果を図Bode_2_5に示す.


なお,\texttt{bode}関数が多入出力に対応しているため,MATLABの場合左辺の変数\texttt{gain}\texttt{phase}は3次元の行列として出力される(対象システム1入力1出力であっても).よって,\texttt{plot}関数で使用する前にこのスクリプト例のように\texttt{squeeze}関数で1次元ベクトルに変換する必要があるので注意する.
なお,<tt>bode</tt>関数が多入出力に対応しているため,MATLABの場合左辺の変数<tt>gain</tt><tt>phase</tt>は3次元の行列として出力される(対象システム1入力1出力であっても).よって,<tt>plot</tt>関数で使用する前にこのスクリプト例のように<tt>squeeze</tt>関数で1次元ベクトルに変換する必要があるので注意する.


sample2_5 \texttt{bode}関数を利用したゲイン線図のスクリプト
sample2_5 <tt>bode</tt>関数を利用したゲイン線図のスクリプト
  num = [1]; %伝達関数の分子
  num = [1]; %伝達関数の分子
  den = [1 0.5 1]; %伝達関数の分母
  den = [1 0.5 1]; %伝達関数の分母
256行目: 256行目:


Bode_2_5.pdf
Bode_2_5.pdf
<math>\displaystyle \frac{1}{s^2 + 0.5s+ 1}</math>のBode線図(\texttt{bode}関数利用,ゲイン線図のみ)
<math>\displaystyle \frac{1}{s^2 + 0.5s+ 1}</math>のBode線図(<tt>bode</tt>関数利用,ゲイン線図のみ)




262行目: 262行目:
=== 周波数応答(ベクトル軌跡)のプロット ===
=== 周波数応答(ベクトル軌跡)のプロット ===


ベクトル軌跡は周波数応答関数を<math>\omega = 0 \sim \infty</math>の範囲で複素平面上にそのままプロットしたものである.MATLAB(Octave)では複素数ベクトルを\texttt{plot}関数でプロットすると複素平面上で表されるので,それをそのまま利用できる.リストsample2_6がスクリプトの例であり,実行結果は図vector_loci_2_6である.なお,MATLAB(Octave)には\texttt{nyquist}関数が用意されており,それを利用してもベクトル軌跡を描くことができる.
ベクトル軌跡は周波数応答関数を<math>\omega = 0 \sim \infty</math>の範囲で複素平面上にそのままプロットしたものである.MATLAB(Octave)では複素数ベクトルを<tt>plot</tt>関数でプロットすると複素平面上で表されるので,それをそのまま利用できる.リストsample2_6がスクリプトの例であり,実行結果は図vector_loci_2_6である.なお,MATLAB(Octave)には<tt>nyquist</tt>関数が用意されており,それを利用してもベクトル軌跡を描くことができる.


sample2_6 二次遅れ系のベクトル軌跡を描くスクリプト
sample2_6 二次遅れ系のベクトル軌跡を描くスクリプト
305行目: 305行目:


よって,むだ時間要素のボード線図を描くスクリプト例はリストsample2_7となる.実行結果は図Bode_2_7である.
よって,むだ時間要素のボード線図を描くスクリプト例はリストsample2_7となる.実行結果は図Bode_2_7である.
ここで,関数\texttt{ones}は指定したサイズで要素がすべて1の行列を生成する関数であり,関数\texttt{size}で周波数ベクトルの大きさを取り出して利用している.
ここで,関数<tt>ones</tt>は指定したサイズで要素がすべて1の行列を生成する関数であり,関数<tt>size</tt>で周波数ベクトルの大きさを取り出して利用している.


sample2_7 むだ時間要素のボード線図を描くスクリプト
sample2_7 むだ時間要素のボード線図を描くスクリプト
342行目: 342行目:
実験科目や卒業研究などで,実験装置での実測結果と理論式から求めた理論曲線を比較する場面がよく出てくる.リストsample2_8は実験で求めた周波数応答と理論曲線を重ねて描くスクリプトの例である.
実験科目や卒業研究などで,実験装置での実測結果と理論式から求めた理論曲線を比較する場面がよく出てくる.リストsample2_8は実験で求めた周波数応答と理論曲線を重ねて描くスクリプトの例である.


実験による測定結果の周波数,ゲイン,位相をそれぞれ\texttt{omega\_e}, \texttt{gaindB\_e}, \texttt{phase\_e}という変数に代入しておき,理論曲線と重ねている.この際,実験データの方はマーカーで示すべきであるので,\texttt{semilog}関数内で線類を指定している.また,関数\texttt{legend}でグラフに凡例を付け加えることができる.実行結果を図Bode_2_8に示す.
実験による測定結果の周波数,ゲイン,位相をそれぞれ<tt>omega_e</tt>, <tt>gaindB_e</tt>, <tt>phase_e</tt>という変数に代入しておき,理論曲線と重ねている.この際,実験データの方はマーカーで示すべきであるので,<tt>semilog</tt>関数内で線類を指定している.また,関数<tt>legend</tt>でグラフに凡例を付け加えることができる.実行結果を図Bode_2_8に示す.