EEIC後期実験演習例 Python Control版

前にjupyterで書いてみたは良いが重くてGithubで公開をためらった代物です。 Markdown変換してみました。

ossyaritoori.hatenablog.com

はじめに：Environment settings
Initialization
Saving

はじめに：Environment settings

Using anaconda you just need to activate following command.

pip install python-control

For further information see here. folllowing command is additional.

Textbook for this experiment is here.

# !python -m pip install --upgrade pip
# !pip install slycot 

# I failed to install these command

Initialization

Import matlab like package from control.matlab

from control.matlab import *    # MATLAB-like functions

basic feature

At first you can use tf function like matlab.

下記は演習2の(1)の問題を例に取り上げます。

# Q2-1
num1 = [1,2]
den1 = [3,1]
sys = tf(num1,den1)
print(sys) # this is the system output

 s + 2
-------
3 s + 1

Load matplotlib for drawing figure

import matplotlib.pyplot as plt
%matplotlib nbagg

# import scipy
from scipy import arange

# step responce 
# Generate step responses
plt.figure(1)
(y1a, T1a)=step(sys,T = arange(0, 10, 0.01))
plt.plot(T1a,y1a)

<IPython.core.display.Javascript object>

[<matplotlib.lines.Line2D at 0x15ae2ea4438>]

Error メモ

sysの関数を変更すると？ときどき'tuple' object is not callableというエラーが起きます。イマイチ自信がないですがどこか変更してはいけない部分を変更しようとしているようです。

一度カーネルをRestartすると治るのでとりあえずはこれで勘弁してください。

# bode plot
bode1 = bode(sys)

<IPython.core.display.Javascript object>

演習4：2次系の応答

数式も打ち込みたいけど時間がないので割愛。減衰比（ダンピング）を変えたらどうなる？という話です。

共振周波数を変えるVerもありましたが，ダンピングの例だけ示します。

Gs = []
K = 1
w_n = 10
Ys = []
Ts =[]
legend = []
# changing zeta

plt.figure(3)
for z in range(0,8):
    zeta = z * 0.2
    Gs.append(tf([K*w_n*w_n],[1,2*zeta*w_n,w_n*w_n]))
    (y1, T1)=step(Gs[-1],T = arange(0, 1, 0.001))
    Ys.append(y1)
    Ts.append(T1)
    # plot
    plt.plot(Ts[-1],Ys[-1])
    legend.append('\zeta is '+str(zeta))
    
plt.legend(legend)
plt.show()

<IPython.core.display.Javascript object>

どうせならbode線図も見てみましょう

plt.figure(4)
for G in Gs:
    bode(G)

<IPython.core.display.Javascript object>

演習5：ゼロ点の性質

最終値の定理から行き着く値は1/15で0.06くらいになることはわかりきっている。ではゼロ点はどこに寄与するのか？

よしんば過渡応答であったとしてどのように？

という設問であるが，私が思うに分子がs+1になるような極端な例がないとこの設問の意図は組みにくいであろう。早急に実験レポートに反映されるべきである。

# 
sys5 = [tf([1],[1,8,15]), tf([1,4],[4,32,60]) , tf([1,50],[50,400,750]), tf([1,1],[1,8,15])]

plt.figure(5)
for sys in sys5:
    (y1, T1)=step(sys,T = arange(0, 1, 0.001))
    plt.plot(T1,y1)

plt.legend(['G1','G2','G3','Additional'])
plt.show()

<IPython.core.display.Javascript object>

レポートの解答で無難に言うならば，「遅いゼロ点はオーバーシュートを引き起こす。早い零点はまぁ無視しても良い。」ということになるが，これはどうしてか。せっかくなのでBode線図も見てみよう

plt.figure(6)
for sys in sys5:
    print(sys)
    bode(sys)

plt.legend(['G1','G2','G3','Additional'])

<IPython.core.display.Javascript object>

      1
--------------
s^2 + 8 s + 15


      s + 4
-----------------
4 s^2 + 32 s + 60


       s + 50
--------------------
50 s^2 + 400 s + 750


    s + 1
--------------
s^2 + 8 s + 15






<matplotlib.legend.Legend at 0x15ae3a986a0>

遅いゼロ点ではゲインが不用意に盛り上がっていることが確認できる。不要な周波数成分を増幅することにより変な振動が起きてしまう → オーバーシュートが出てくるといった感じである。

直感的に説明すると，「分子が1の場合，安定極であれば分母は単調増加して高周波ほどゲインが下がるのに対し，分子にsの項が加わると分子も増加する要素が加わるためsの寄与が相対的に高い遅い零点ではある周波数周辺でゲインが盛り上がる → オーバーシュート」という感じであろうか。

演習８

結構難しい気がする。問題の意図が汲めない… 一応低域の位相を遅らせることと引き換えに高域の位相を少し早めているって言うことなんだろうが。誰か解説手伝って。

余談：zpk記述

因数分解されたシステムは記述が面倒くさい。 zpk2tfというscipyの関数を用いることでゼロ点，極，DCゲインの３つからtfを生成可能だが，生成結果はtupleであり２変数を受け取るtfに直接食わせられない。

そんなときは*を用いた展開表現を用いることで，tupleを複数argumentとして認識させることが出来る。

参考

Gp =tf(*zpk2tf([],[0,-2,-10],40))
Gc =tf(*zpk2tf([-0.2,-1],[-0.02,-5],5))

plt.figure(7)
bode(Gp)
bode(Gp*Gc)
bode(feedback(Gp*Gc))

plt.legend(['Plant','Plant+C', 'Feedback'])
plt.show()

<IPython.core.display.Javascript object>

根軌跡

rlocusやってrlockfindがない！ってOctaveで毎年なるんですがpythonにもなかった…

tau = 0.01
plant = tf([1],[1,3,2,0])
gc = tf([-1.5],[-1.0/tau])

x1,y1 = rlocus(plant*gc)

<IPython.core.display.Javascript object>

Saving

データとか保存用になにかする（予定）

粗大メモ置き場

個人用，たまーに来訪者を意識する雑記メモ