2017年6月25日日曜日

pythonで一次遅れのボード線図を描く

pythonを使って制御工学の問題を解いていこうと思っています。少しづつ更新していきますので、今続きがなくても、たまに確認してみてくださいね。

ボード線図を描く


今回は一次遅れ系のボード線図を描くことを考えてみます。

伝達関数


まず一次遅れの伝達関数の一例は次の式で表されます。

\begin{equation} G(s)=\frac{1}{s+1} \end{equation}

周波数伝達関数


 ここでこの伝達関数のボード線図を描く前にすることは$s=j\omega$ と置き換えをします。

\begin{equation} G(j\omega)=\frac{1}{1+j\omega} \end{equation}


この置き換えをした後の伝達関数を周波数伝達関数と呼びます。

すると伝達関数は$\omega$を変数とする、複素数と見ることができます。

ところでボード線図は、ゲイン曲線と位相曲線の二つに分けて描きます。

ゲイン線図は周波数伝達関数からゲイン$20\log |G(j\omega)|$を計算したものを横軸が対数目盛の片対数グラフとしたものです。
横軸が角周波数$\omega$で、縦軸はゲインです。

複素数の実部と虚部そして有理化


複素数からゲインを求めるためには複素数の実部と虚部を明らかにしておく必要があります。以下ではその過程を述べてみたいと思います。
一次遅れの伝達関数は分母が複素数の形をした分数の形をした関数でこのような姿の関数を有理関数と呼びます。
複素数の分母を複素数ではないものにすることを有理化と呼びます。以下で有理化の作業を実際にやってみます。

有理化には分母の複素数の共役複素数を分母と分子にかけます。同じものを分母と分子にかけますから。かける前とかけた後も関数自身は同じものです。

\begin{equation} \begin{split} G(j\omega)&=\frac{1}{1+j\omega}\frac{1-j\omega}{1-j\omega}\\ &=\frac{1-j\omega}{1+\omega^2} \end{split} \end{equation}

以上が有理化と言われる作業です。複素数であった分母からは複素数がなくなりました。
さて、話を戻して実部と虚部を求めてみましょう。

複素数を眺めてみたときに虚数単位のjが付いていない項が実部で、jが付いてる項のjの係数を虚部と言います。
したがって、今の場合は、$\frac{1}{1+\omega^2}$が実部で、$-\frac{\omega}{1+\omega^2}$が虚部と言うことができます。

 ここで今後の記述が大変なので実部と虚部を以下のように置き換えておきます。

\begin{equation} \begin{split} x=&\frac{1}{1+\omega^2}\\ y=-&\frac{\omega}{1+\omega^2} \end{split} \end{equation}

これで、ボード線図を描くための準備がほぼ整いました。
ゲインと位相を計算します。

 ゲインの計算


ゲインを計算するためには複素数の大きさを求めなければなりません。大きさは以下の式から計算します。
\begin{equation} |G(j\omega)|=\sqrt{x^2+y^2}=(x^2 + y^2)^{\frac{1}{2}} \end{equation}

ゲインは以上の対数をとって20倍したものと定義されています。

\begin{equation}
20 \log |G(j\omega)|=20 \log (x^2 + y^2)^{\frac{1}{2}}=10 \log (x^2 + y^2)
\end{equation}

$x$と$y$に元の式を代入して整理してみましょう。

\begin{equation}
10 \log (x^2 + y^2)= 10 \log \left\{\left(\frac{1}{1+\omega^2}\right)^2 + \left(-\frac{\omega}{1+\omega^2}\right)^2\right\}
\end{equation}

さて、さらに式の整理を進めましょう

\begin{equation}
\begin{split}
&10 \log \left\{\left(\frac{1}{1+\omega^2}\right)^2 + \left(-\frac{\omega}{1+\omega^2}\right)^2\right\}=10 \log \left\{ \frac{1+\omega^2}{(1+\omega^2)^2}\right\}\\
&=10 \log \frac{1}{1+\omega^2}
\end{split}
\end{equation}

だいぶ簡単になりましたね。ここで対数法則を使って、分数を引き算にします。

\begin{equation}
10 \log \frac{1}{1+\omega^2} = 10 ( \log 1 - \log (1+ \omega^2) )
\end{equation}

$\log 1$は0なので、さらに簡単になって以下のようになりました。

\begin{equation}
10 ( \log 1- \log (1+\omega^2) ) = -10 \log(1+\omega^2)
\end{equation}

位相の計算


ゲインの計算がひとまず終わったので位相の計算をします。

位相は先ほど定義した実部$x$、虚部$y$を使って説明すると、以下のように表されます。

\begin{equation}
\tan^{-1} \frac{y}{x}
\end{equation}

ここでもに元の式を代入して計算を進めてみましょう。

\begin{equation}
\begin{split}
\tan^{-1} \frac{y}{x}&=\tan^{-1} \frac{-\frac{\omega}{1+\omega^2}}{\frac{1}{1+\omega^2}}\\
&=\tan^{-1} (-\omega)\\
&=-\tan^{-1} \omega
\end{split}
\end{equation}

 数値計算!


#ここからソース
#ここまでソース

さて、上のソースが今回ボード線図を描くのに使用したソースコードです。
本記事はpythonでボード線図を描ききるのが主眼です。

import numpy as np
import matplotlib.pyplot as plt

でpythonに数値配列用のライブラリ超有名なnumpyを読み込みます。使うときはnpとして使います。

次に、グラフを描くためのライブラリであるmatplotlib.pyplotを読み込みます。
私は長年、仕事にgnuplotを使ってきましたが、だいぶ慣れてきたので、グラフ描画をpyplotでやってもいいかなと思い始めています。


描きあがったボード線図が以下です。

0 件のコメント:

コメントを投稿