Pythonで学ぶ制御工学　第3弾：伝達関数

# #Pythonで学ぶ制御工学< 伝達関数 >

##はじめに
基本的な制御工学をPythonで実装し，復習も兼ねて制御工学への理解をより深めることが目的である．
その第3弾として伝達関数を扱う．

##伝達関数
**入力と出力の関係を周波数領域で表現したもの**

$P(s)=\frac{Y(s)}{U(s)}$

具体的には，次のような形で伝達関数は入力と出力の関係をむすぶ．

$Y(s)=P(s)U(s)$


##実装方法
ここで，まずPythonでの伝達関数モデルの実装方法を学び，後に前回の制御モデルの伝達関数モデルを実際にPython上で作成してみることとする．なお，伝達関数モデルは上記の$P(s)$のことである．

学習のための簡単な例として，次に示す伝達関数モデルを扱う．

①　$P(s) = \frac{1}{s^2+2s+3}$

②　$P(s) = \frac{s+2}{s^3+5s^2+3s+4}$

③　$P(s) = \frac{s+3}{(s+1)(s+2)^2}$

ここでは，**``control``ライブラリの``tf()``関数で伝達関数の作成**を行う．

以下にソースコードとそのときの出力を示す．
#####ソースコード
```python:trans_function.py
"""
2021/02/14
@Yuya Shimizu

伝達関数モデル
"""
from control import tf, tfdata

##伝達関数を作成
#①式
Np = [0, 1]     #伝達関数の分子多項式の係数（0*s + 1）
Dp = [1, 2, 3]     #伝達関数の分母多項式の係数（1*s^2 + 2*s + 3）

P = tf(Np, Dp)

print(f"①式<伝達関数>\n{P}\n")


#②式
Np = [1, 2]     #伝達関数の分子多項式の係数（1*s + 2）
Dp = [1, 5, 3, 4]     #伝達関数の分母多項式の係数（1*s^3 + 5*s^2 + 3*s + 4）

P = tf(Np, Dp)

print(f"②式<伝達関数>\n{P}\n")


#③式
Np = [1, 3]     #伝達関数の分子多項式の係数（1*s + 3）
Dp = [1, 5, 8, 4]     #伝達関数の分母多項式の係数（(1*s + 1)(1*s + 2)^2  =  1*s^3 + 5*s^2 + 8*s + 4）


P = tf(Np, Dp)

print(f"③式<伝達関数>\n{P}\n")


#③式の別方法
Np1 = [1, 3]     #伝達関数の分子多項式の係数（1*s + 3）
Dp1 = [0, 1]     #伝達関数の分母多項式の係数（0*s + 1）
Np2 = [0, 1]     #伝達関数の分子多項式の係数（0*s + 1）
Dp2 = [1, 1]     #伝達関数の分母多項式の係数（1*s + 1）
Np3 = [0, 1]     #伝達関数の分子多項式の係数（0*s + 1）
Dp3 = [1, 2]     #伝達関数の分母多項式の係数（1*s + 2）

P1 = tf(Np1, Dp1)
P2 = tf(Np2, Dp2)
P3 = tf(Np3, Dp3)

P = P1 * P2 * P3**2

print(f"③式<伝達関数>（分割法）\n{P}\n")


##伝達関数の分子と分母の係数を抽出する
#シンプルに抽出する
numP = P.num
denP = P.den

print(f"<係数抽出>\n\n分子係数：{numP}, 分母係数：{denP}\n\n")


#入れ子を避ける抽出
[[numP]], [[denP]] = tfdata(P)
print(f"<係数抽出>（入れ子を避ける抽出）\n\n分子係数：{numP}, 分母係数：{denP}\n")
```

#####出力
```
①式<伝達関数>

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


②式<伝達関数>

        s + 2
---------------------
s^3 + 5 s^2 + 3 s + 4


③式<伝達関数>

        s + 3
---------------------
s^3 + 5 s^2 + 8 s + 4


③式<伝達関数>（分割法）

        s + 3
---------------------
s^3 + 5 s^2 + 8 s + 4


<係数抽出>

分子係数：[[array([1., 3.])]], 分母係数：[[array([1, 5, 8, 4])]]


<係数抽出>（入れ子を避ける抽出）

分子係数：[1. 3.], 分母係数：[1 5 8 4]

```

分子と分母の多項式係数を引数に伝達関数を作成することができる．また，③式の伝達関数モデルについては，展開してからその係数を与える方法もあるが，簡単に因数分解することで，小さな伝達関数の積で求めることも可能である．
ここでは，さらに係数の抽出についての実装方法も示した．``.num``で分子，``.den``で分母の多項式係数を抽出することができる．しかし，出力からも分かるように，入れ子構造で出力されており，非常に扱いづらい．よって，``control``ライブラリの``tfdata()``関数を使って，うまく取り出せば，単なるリストとして抽出することが可能となる．


##実践
前回の制御モデルから，伝達関数モデルを求める．その導出過程を示した後，Pythonでの実装を行う．
ただし，``tf()``関数は数値でないと受け付けないため，関数を定義することで，擬似的にPythonでの実装を表現し，さらに引数に係数を指定することで，数値の変更を容易にしている．なお，ここでは，重力加速度以外を1として与えている．


####制御モデル①：台車
#####導出過程
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/951300/b8841b0e-7323-80dd-463e-06f5230da18b.png)

#####ソースコード
```python:tf_truck.py
"""
2021/02/14
@Yuya Shimizu

伝達関数モデル（台車）
"""
from control import tf


def truck_tf(M, mu):
    ##伝達関数を作成
    Np = [0, 1]       #分子の多項式係数（0*s + 1）
    Dp = [M, mu]    #分母の多項式係数（M*s + μ）

    P = tf(Np, Dp)

    return P

if __name__ == '__main__':
    M   =    1       #質量 
    mu =    1       #粘性係数
    P = truck_tf(M, mu)

    print(f"伝達関数モデル（台車）\n{P}\n質量: {M}\n粘性係数: {mu}")

```

#####出力
```
伝達関数モデル（台車）

  1
-----
s + 1

質量: 1
粘性係数: 1
```

####制御モデル②：アーム
#####導出過程
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/951300/96e4a1c8-9d8f-4e33-d818-31f397d3ef80.png)


#####ソースコード
```python:tf_arm.py
"""
2021/02/14
@Yuya Shimizu

伝達関数モデル（アーム）
"""
from control import tf


def arm_tf(J, mu, M, g, l):
    ##伝達関数を作成
    Np = [0, 1]       #分子の多項式係数（0*s + 1）
    Dp = [J, mu, M*g*l]    #分母の多項式係数（M*s^2 + μ*s + M*g*l）

    P = tf(Np, Dp)

    return P

if __name__ == '__main__':
    J    =   1       #慣性モーメント
    mu =   1       #粘性係数
    M   =   1       #質量
    g    =   9.8    #重力加速度
    l     =   1      #回転軸から重心までの距離
    P = arm_tf(J, mu, M, g, l)

    print(f"伝達関数モデル（アーム）\n{P}\n慣性モーメント: {J}\n粘性係数: {mu}\n質量: {M}\n重力加速度: {g}\n回転軸から重心までの距離: {l}")

```

#####出力
```
伝達関数モデル（アーム）

      1
-------------
s^2 + s + 9.8

慣性モーメント: 1
粘性係数: 1
質量: 1
重力加速度: 9.8
回転軸から重心までの距離: 1
```


####制御モデル③：RLC回路
#####導出過程
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/951300/d3e51d05-4be8-bae8-15c6-971731adb948.png)


#####ソースコード
```python:tf_RLC.py
"""
2021/02/14
@Yuya Shimizu

伝達関数モデル（RLC回路）
"""
from control import tf


def RLC_tf(R, L, C):
    ##伝達関数を作成
    Np = [0, 1]       #分子の多項式係数（0*s + 1）
    Dp = [L*C, R*C, 1]    #分母の多項式係数（L*C*s^2 + R*C*s + 1）

    P = tf(Np, Dp)

    return P

if __name__ == '__main__':
    R =  1   #抵抗器
    L =  1   #コイル
    C =  1   #コンデンサ
    
    P = RLC_tf(R, L, C)

    print(f"伝達関数モデル（アーム）\n{P}\n抵抗器: {R}\nコイル: {L}\nコンデンサ: {C}")

```

#####出力
```
伝達関数モデル（アーム）

     1
-----------
s^2 + s + 1

抵抗器: 1
コイル: 1
コンデンサ: 1
```

####制御モデル④：増幅回路
#####導出過程
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/951300/75a266e7-9e79-8b58-59ed-988d3d29ef08.png)

#####ソースコード
```python:tf_amplifier.py
"""
2021/02/14
@Yuya Shimizu

伝達関数モデル（増幅回路）
"""
from control import tf


def amplifier_tf(R1, R2):
    ##伝達関数を作成
    Np = [0, -R2]       #分子の多項式係数（0*s - R2）
    Dp = [R1*R2, R1]    #分母の多項式係数（R1*R2*s + R1）

    P = tf(Np, Dp)

    return P

if __name__ == '__main__':
    R1 =  1   #抵抗器1
    R2 =  1   #抵抗器2
    
    P = amplifier_tf(R1, R2)

    print(f"伝達関数モデル（アーム）\n{P}\n抵抗器1: {R1}\n抵抗器2: {R2}")

```

#####出力
```
伝達関数モデル（アーム）

  -1
-----
s + 1

抵抗器1: 1
抵抗器2: 1
```

##感想
``tf()``関数が数値のみを受け付けるということで，文字を含むような一般的な表現はできなかったが，とりあえず伝達関数モデルの作成は身に着けられたと思う．次回は状態空間モデルということで，こちらも楽しみである．

##参考文献
Pyhtonによる制御工学入門　　南 祐樹　著　　オーム社