ディジタル制御（1987年）その３v2 pythonで記述してみる（５.単一ループディジタル制御）

2024.3.29 ディジタル制御（1987年）その３v2 pythonで記述してみる（５.単一ループディジタル制御） 下 条 誠 電気通信大学名誉教授 https://researchmap.jp/read0072509/ “高橋安人，ディジタル制御 ，岩波書店， 1987 The University of Electro-Communications https://www.docswell.com/user/m_shimojo Department of Mechanical Engineering and Intelligent System

• https://researchmap.jp/read0072509/
• https://www.docswell.com/user/m_shimojo

目 次 1. 緒 言 1.1 ディジタル制御工学 4.4 4.5 1.2 1.3 1.4 1.5 5. 連続時間プラントの離散時間形 5.1 フィードバック制御系 5.2 １次プラントの制御 ディジタル制御系の構成 量子化、サンプリングおよびホールド 例－ロボットアームのデッドビート応答 プログラム作成入門 2. ｚ変換と反復計算 2.1 ｚ領域における時系列 2.2 入力と出力 2.3 極とゼロ 2.4 2.5 調和振動 リカーシブ・アルゴリズム 5.3 5.4 5.5 時系列からの伝達関数 観測器による状態推定 追従系の設計 PID制御則 PIDゲイン調整 6. 線形2乗最適制御 6.1 行列リカチ(Riccati)式 6.2 慣性体のLQ制御 6.3 I動作を含むLQ制御 6.4 観測器を含むLQI制御 6.5 ０型プラントのLQI制御 ３. 線形離散時間系 3.1 3.2 3.3 状態空間系 伝達関数 同伴形 3.4 応答の算定 7.1 LQ制御系の固有値 3.5 安定に関する定理 7.2 7.3 7.4 7.5 LQ制御系の根軌跡 固有値を指定の制御系および状態推定系設計 有限整定系 多変数制御系 ４. 連続時間プラントの離散時間形 4.1 プラント微分方程式と等価の差分式 4.2 パルス伝達関数と拡張ｚ変換 4.3 解析的手法によるｚ変換と拡張ｚ変換 7. 8. 固有値を指定した制御系 同定および適応制御 2 8.1 8.2 簡単な適応サーボ系 模型規範適応制御 8.3 8.4 8.5 MRASによる同定と適応制御への応用 模型追従適応制御 最小2乗同定アルゴリズム 9. 周波数領域における諸関係 9.1 振動状定常状態 9.2 周波数応答 9.3 制御系の周波数応答 9.4 観測器を含む状態ベクトルフィードバック系 9.5 離散形式のフーリエ変換 10. ランダム信号系 10.1 乱数から造成する正規性白色ノイズ 10.2 10.3 10.4 有色ノイズと自己相関 共分散行列 連続時間ノイズを含む離散時間形 10.5 最適フィルタ

コメント 3 本書の面白いところは、全てのアルゴリズムが、HP-Basi言語と付属の基本行列 演算で作られている点です。正準形変換の方法やリカチ行列の解法などの記述が 多々あり大変興味深いものです。これらは現在ではpython-controlなどを使えば 簡単に解けます。しかし、その解法の基本原理まで踏み込み、意外と簡単なアル ゴリズムで計算できることを本書を通じて実感できることは楽しいです。なお本 解説は本人の再学習のためであり、主な内容のメモとプログラムの記述だけです。 デジタル制御の解説ではありません。あとpython-controlなどでの検証を少し加 えました。 また、作成したpythonプログラムにはミスなどあると思います。書籍の結果とは 比べチェックはしましたが完全ではありません。その点ご注意ください。 （Pythonのprogramは下記のGitHubに格納してあります） https://github.com/m4881shimojo/YTakahashi_Digital-control/tree/main 高橋 安人氏（1912-1996)は、日本の機械工学者。工学博士（東京帝国大学）、 カリフォルニア大学バークレー校名誉教授、豊橋技術科学大学名誉教授。当初は 伝熱工学を専攻していたが、太平洋戦争後は日本の制御工学の草分けとして活躍 した。https://ja.wikipedia.org/wiki/高橋 安人

• https://github.com/m4881shimojo/YTakahashi_Digital-control/tree/main
• https://ja.wikipedia.org/wiki/

5. 5.1 連続時間プラントの離散時間形 フィードバック制御系 5.2 １次プラントの制御 5.3 追従系の設計 5.4 PID制御則 5.5 PIDゲイン調整 5章は主に古典的手法を、 ６～７章では現代制御的 手法を対象とする 4

5.1 フィードバック制御系 V(z) 外乱 応答 設定点 R(z) + E(z) ディジタル制御 5 U(z) + Gc(z) + サンプラホールド 付きプラントGp(z) Y(z) 𝐺c 𝑧 制御法則 𝐺𝑝 𝑧 プラント伝達関数 𝐺 𝑧 = 𝐺𝑝 𝑧 ∙ 𝐺𝑐 𝑧 図5-1 目的 単一ループディジタルフィードバック制御系 開ループ伝達関数 1. 追跡 (Tracking) Type1 2. 自動調節 (Regulation) Type0 𝑌 𝑧 = 𝐺𝑝 𝑧 𝐺(𝑧) 𝑅 𝑧 + 𝑉(𝑧) 1 + 𝐺(𝑧) 1 + 𝐺𝑝 𝑧 (5-1) 𝐸 𝑧 = 𝐺𝑝 𝑧 1 𝑅 𝑧 + 𝑉(𝑧) 1 + 𝐺(𝑧) 1+𝐺 𝑧 (5-2)

5.1 フィードバック制御系 ⚫ ランプ入力 ⚫ 単位ステップ入力 𝑅 𝑧 = 𝐸 𝑧 = 𝑧 𝑧−1 𝑧 1 𝑧 − 1 1 + 𝐺(𝑧) 6 𝑉 𝑧 =0 (5-4) 𝑇𝑧 𝑅 𝑧 = 𝑧−1 2 𝑇𝑧 1 𝐸 𝑧 = 𝑧 − 1 2 1 + 𝐺(𝑧) 𝑉 𝑧 =0 (5-6) lim 𝑒 𝑘 = lim 𝑧 − 1 𝐸(𝑧) 定常偏差は： 𝑒0 = 𝑘→∞ 𝑧→1 𝑧 1 𝑒0 = lim 𝑧 − 1 𝑧→1 𝑧 − 1 1 + 𝐺(𝑧) 1 𝑒0 = 1 + 𝐺(1) 𝑧→1 2 1 1 + 𝐺(𝑧) 𝑇 𝑒1 = lim 𝑧→1 (𝑧 − 1)𝐺(𝑧) (5-5) G(z)有限ならばe0は有限 （G(z)がz=1に極を持たない） Type0 𝑒1 = lim 𝑧 − 1 𝑇𝑧 𝑧−1 (5-7) G(z)がz=1に単極を持てばe1は有限 Type1

5.2 １次プラントの制御 5.1 フィードバック制御系 5.2 １次プラントの制御 5.3 追従系の設計 5.4 PID制御則 5.5 PIDゲイン調整 古典的手法 での解法 7

5.2 １次プラントの制御（例） 対象システム例： 𝐺 𝑠 = 𝑝 U(z)+ E(z) 1 𝑠+𝑎 （一次プラント） U(z) Gc(z) Gp(z) Y(z) 𝐺c 𝑧 制御法則 𝐺𝑝 𝑧 プラント伝達関数 𝑌 𝑧 = 単一ループディジタルフィードバック制御系 Step1: 離散時間形へ変換 𝐺𝑝 𝑠 = 1 𝑠+𝑎 𝐺 𝑧 𝑈(𝑧) 1+𝐺 𝑧 𝐺 𝑧 = 𝐺𝑐 𝑧 𝐺𝑝 𝑧 Step2: 制御法則の決定 𝐺𝑐 𝑧 ディジタル制御 その２v2 4.3節 𝐺𝑝 𝑧 8 type1の制御系 とする 極 z=1 𝐺𝑐 z = 𝐾𝑐 を含む 𝑧−𝑑 𝑧−1 d:設計の自由度を持たすために導入

5.2 １次プラントの制御（例） 𝓛 変換形 Step1: 9 𝓩 変換形 1 𝐺𝑝 𝑠 = 𝑠+𝑎 1 1 − 𝑒 −𝑎𝑇 𝐺𝑝 z = 𝑎 𝑧 − 𝑒 −𝑎𝑇 ディジタル制御 その２ (4-3) 𝐺𝑝 z = 𝐾𝑝 = 1−𝑝 𝑎 𝐾𝑝 𝑧 − 𝑝 (5-8) 𝑝 = 𝑒 −𝑎𝑇 Step2: Type 0 𝐺𝑐 z = 𝐾𝑐 𝐺 𝑧 = 比例動作のみ 仕様 要求 とす る （極 z=1 を含む） Type 1 𝐺𝑐 z = 𝐾𝑐 ランプ入力での 定常偏差 𝑧−𝑑 𝑧−1 𝐾𝑝 𝐾𝑐 𝑧−𝑝 𝐺 z = 𝐺𝑐 z 𝐺𝑝 (𝑧) 𝑧−𝑑 𝐺 z = 𝐾𝑐 𝐾𝑝 𝑧−1 𝑧−𝑝 𝑇 𝑇 1−𝑝 = (𝐾 = 𝐾𝑐 𝐾𝑝 ) 𝑧→1 (𝑧 − 1)𝐺(𝑧) 𝐾 1−𝑑 有限値になる (z-1)は相殺 𝑒1 = lim (5-10) (5-7)

5.2 １次プラントの制御（例） 閉ループ特性式： 𝐹 𝑧 = 𝑧 2 + 𝐾 − 1 − 𝑑 𝑧 + 𝑝 − 𝐾𝑑 10 (5-12) (5-10)より 1+𝐺 z =1+K 𝑧−𝑑 𝒛−𝟏 𝒛−𝒑 +𝑲 𝒛−𝒅 = 𝑧−1 𝑧−𝑝 𝑧−1 𝑧−𝑝 ①閉ループの極をz1, z2とすると 𝐹 𝑧 = 𝑧 − 𝑧1 𝑧 − 𝑧2 = 𝑧 2 − 𝑧1 + 𝑧2 𝑧 + 𝑧1 𝑧2 ②これをKとdについて解くと 𝐾 = 1 + p − 𝑧1 + 𝑧2 , 𝑑 = 𝑝 − 𝑧1 𝑧2 Τ𝐾 (5-14) 極を指定し て制御系を 作る 次頁にこの1次系の挙動を示す。 a=1, T=0.1 ① ② ③ ④ case1: case2: case3: case4: K=1.6 , d=0.6 K=1.905, d=0.475 K=2.5, d=0.3 K=3.0, d=0.3 書籍にあり 次頁に結果を示す。

5.2 １次プラントの制御（例） 閉ループ特性式： 1 + 𝐺 z = 1 + K 11 𝑧−𝑑 𝒛−𝟏 𝒛−𝒑 +𝑲 𝒛−𝒅 = 𝑧−1 𝑧−𝑝 𝑧−1 𝑧−𝑝 𝐹 𝑧 = 𝑧 2 + 𝐾 − 1 − 𝑑 𝑧 + 𝑝 − 𝐾𝑑 (5-12) 閉ループの極をz1, z2とすると 𝐹 𝑧 = 𝑧 − 𝑧1 𝑧 − 𝑧2 = 𝑧 2 − 𝑧1 + 𝑧2 𝑧 + 𝑧1 𝑧2 これをKとdについて解くと 𝐾 = 1 + p − 𝑧1 + 𝑧2 , 𝑑 = 𝑝 − 𝑧1 𝑧2 Τ𝐾 (5-14) 次頁にこの1次系の挙動を示す。（a=1, T=0.1） ① case1: K=1.6 , d=0.6 ② case2: K=1.905, d=0.475 ③ case3: K=2.5, d=0.3 ④ case4: K=3.0, d=0.3 書籍にあり 次頁に結果を示す。

case1: F-pole:特性方程式の根 5.2 １次プラントの制御 imag d=0.6;kgain=1.6 real case2: imag d=0.475;kgain=1.905 real DGC_no3/ p83v4 12

case3: F-pole:特性方程式の根 5.2 １次プラントの制御 imag d=0.3;kgain=2.5 real case4: imag d=0.3;kgain=3 （書籍に無し） 単位円外 real DGC_no3/ p83v4 13

5.2 １次プラントの制御（例） 14 書籍の図 page84 “高橋安人,ディジタル制御 ,岩波書店,1987

Gs→Gz変換 step応答での比較 応答は同じ（あたりまえか） K=2 Gp(s) Gc(s) G=Gc*Gp Gp(z) Gc(z) 極 z=1 を含む G/(1+G) 上式の離散化 15

5.3 追従系の設計 5.1 フィードバック制御系 5.2 １次プラントの制御 5.3 追従系の設計 5.4 PID制御則 5.5 PIDゲイン調整 例）2次系、サーボ機構 16

5.3 追従系の設計1/7 17 page85 物体の位置を目標値R(z)に追従させる制御系（サーボ機構）の設計 例：次のシステムを考える 𝑦′ u b m 𝑚 𝑑2 𝑦 ′ 𝑑𝑡 ′ 2 𝑑𝑦 ′ +𝑏 ′ = 𝑢 𝑑𝑡 𝑡 = 𝑏Τ𝑚 𝑡 ′ 𝑦 = 𝑏 2 Τ𝑚 𝑦 ′ 運動方程式： ラプラス変換式： ｚ変換式： 𝑑 2 𝑦 𝑑𝑦 + =𝑢 𝑑𝑡 2 𝑑𝑡 𝑌(𝑠) 1 𝐺𝑝 𝑠 = = 𝑈(𝑠) 𝑠(𝑠 + 1) 𝐺𝑝 𝑧 = 𝐾𝑝 𝑧+𝑞 𝑧−1 𝑧−𝑝 m: ｂ: y‘: u: 質量 減衰係数 直線変位 力

5.3 追従系の設計2/7 18 Z変換式の求め方（部分分数展開して変換） 𝑑2 𝑦 𝑡 𝑑𝑦 𝑡 + = 𝑢(𝑡) 𝑑𝑡 2 𝑑𝑡 微分式 𝐺𝑝 𝑠 = Laplace変換式 1 1 =− + 𝑠 1+𝑠 𝑠 1 𝑧 → 𝑠 𝑧−1 注 1 𝐺 𝑧 = 1 − 𝑧 −1 𝒵 ℒ −1 𝐺 𝑠 𝑠 1− 𝑧 −1 −1 展開 1 𝑧 → 1+𝑠 𝑧−𝑝 2 1−𝑧 𝑇 𝑧−1 1 1+𝑠 変換 −1 𝑧−1 𝑧−𝑝 𝑇 𝑧−1 𝑧+𝑞 𝐺𝑝 𝑧 = −1 + + = 𝐾𝑝 𝑧−1 𝑧−𝑝 𝑧−1 𝑧−𝑝 𝑝 = 𝑒 −𝑇 𝐾𝑝 = 𝑝 + 𝑇 − 1 1 + 𝑠2 1 𝑇𝑧 → 𝑠2 𝑧−1 1−𝑧 1 ｚ変換式 1 𝑌 𝑠 = 𝑈 𝑠 𝑠 1+𝑠 𝑞= 1 − 𝑝 − 𝑝𝑇 𝑝+𝑇−1 注）ディジタル制御 その２ (4-3)

5.3 追従系の設計3/7 𝑧+𝑞 𝑧−1 𝑧−𝑝 𝑧+𝑏 𝐺𝑐 𝑧 = 𝐾𝑐 𝑎 <1 𝑧+𝑎 制御対象： 𝐺𝑝 𝑧 = 𝐾𝑝 制御方式： 19 (5-16) (5-17) 特性式： 1 + 𝐺𝑝 𝑧 𝐺𝑐 𝑧 → 𝐹 𝑧 = 𝑧 3 + 𝑎 1 𝑧 2 + 𝑎 2 𝑧 + 𝑎(3) (5-18) 特性式の固有値を指定することで、制御式のa, b, Kc が以下の手順で決まる 𝑎 1 =a−p−1+k 𝑎 2 = −ap − a + p + K b + q 𝑎 2 𝑎 3 = ap + Kbq a(1),(2),a(3) 決定 制御式のa, b, Kcが 決まる a(1), a(2), a(3)は 特性値の固有値を指 定すると決まる 𝐹 1 = 1 + 𝑎 1 + 𝑎 2 + 𝑎(3) (K = 𝐾𝑝 𝐾c ) K 1 + b = 𝐹(1)/(1 + 𝑞) 1 a= 𝑎 1 𝑞 + 𝑎(3) + 𝑞 1 + 𝑝 − 𝐾 1 + 𝑏 𝑞 𝑝+𝑞 K=a 1 +1+p−a

5.3 追従系の設計4/7 R(z) + + 20 Limitter E(z) 設定点 ディジタル制御 Gc(z) 𝐺𝑐 𝑧 = 𝐾𝑐 𝑧+𝑏 , 𝑧+𝑎 U(z) U’(z) out in 𝑎 <1 サンプラホールド 付きプラントGp(z) 𝐺𝑝 𝑧 = 𝐾𝑝 Y(z) 応答 𝑧+𝑞 𝑧−1 𝑧−𝑝 固有値を指定 特性式： 𝐹 𝑧 = 𝑧 − 𝑧1 𝑧 − 𝑧2 𝑧 − 𝑧3 = 𝑧 3 +𝑎 1 𝑧 2 + 𝑎 2 𝑧 + 𝑎(3) ① case (a): z1=z2=z3=0 ② case (b): z1=0.5, z2=0.5+0.5j, z3=0.5-0.5j 次頁に結果を示す。但し、本システムのシミュレーションは書籍では未掲載のため自 作した（←ご注意ください！）。書籍の結果も併せて記載する。ステップ応答はほぼ 同じだが、ランプ応答は書籍の方が良いように思える。

5.3 追従系の設計5/7 21 特性式の固有値 ① case (a): z1=z2=z3=0 DGC_no3/ p87a

5.3 追従系の設計6/7 22 特性式の固有値 ② case (b): z1=0.5, z2=0.5+0.5j, z3=0.5-0.5j DGC_no3/ p87b

5.3 追従系の設計7/7 23 書籍の図 page87 “高橋安人,ディジタル制御 ,岩波書店,1987

今回の制御方式Gc(z)の特性 24 Gc(z)の役割がイメージできないので、その応答を計算してみた 𝐺𝑐 𝑧 = 𝐾𝑐 case(a) 𝑧+𝑏 𝑧+𝑎 𝑎 <1 計算に用いたa, bは, 図5-3で用いた値である 図の右はstep応答、左は sine波に対する応答である。 これをみる限りでは、 ① 減衰するbang-bang 振動？ ② 位相を進める効果？ case(b) があるのか？？ DGC_no3/p87Gc.ipynb

制御方式Gc(z)のBode線図（付録参照のこと） Gc(z)がどのような特性になるのかBode線図を描いてみた zをexp(jωT)に置き換え、周波数特性を求めた 𝑧 = 𝑒 𝑗𝜔𝑇 1 + 𝑏𝑒 −𝑗𝜔𝑇 𝐺𝑐 𝑧 = 𝐾𝑐 1 + 𝑎𝑒 −𝑗𝜔𝑇 25 意味があるのか？ 𝑎 <1 𝑧 −1 = 𝑒 −𝑗𝜔𝑇 𝑇 = 𝑑𝑡 p87a p87b 計算に用いた、a, b, dtはcase(a), case(b)で用いた値である DGC_no3/ GcBode.ipynb

5.4 PID制御則 5.1 フィードバック制御系 5.2 １次プラントの制御 5.3 追従系の設計 5.4 PID制御則 5.5 PIDゲイン調整 26

5.4 PID制御則 偏差E(s)から制御入力U(s) R(s) + E(s) Gc(s) U(s) 27 Gp(s) Y(s) を算出する方式の１つ PID制御則： 𝐺𝑐 (𝑠) = kc= ゲイン Ti= I 動作時間 Td= D動作時間 𝑈(𝑠) 1 = 𝑘𝑐 1 + + T𝑑 𝑠 𝐸(𝑠) 𝑇𝑖 𝑠 𝐸 𝑠 = 𝑅 𝑠 − 𝑌(𝑠) 刻み差d(k)で考える： 𝑑 𝑘 = 𝑘𝑐 𝑑 𝑘 = 𝑢 𝑘 − 𝑢(𝑘 − 1) 𝑒 𝑘 − 𝑒(𝑘 − 1) + 𝑇 𝑒 𝑘 + 𝑒(𝑘 − 1) 𝑇𝑑 + 𝑒 𝑘 − 2𝑒(𝑘 − 1) + 𝑒(𝑘 − 2 𝑇𝑖 2 𝑇 𝒊𝒇: 𝑟 𝑘 − 2 = 𝑟 𝑘 − 1 = 𝑟(𝑘) 𝑑 𝑘 = 𝑘𝑐 𝑦 𝑘 − 1 − 𝑦(𝑘) + 𝑇: sampling interval 仮定：3刻み期間不変 𝑇 𝑦 𝑘 − 1 + 𝑦(𝑘) 𝑇𝑑 𝑟 𝑘 − + 2𝑦 𝑘 − 1 − 𝑒 𝑘 − 2 + 𝑒(𝑘 𝑇𝑖 2 𝑇 𝑑 𝑘 = 𝐾𝑝 𝑦 𝑘 − 𝑦(𝑘 − 1) + 𝐾𝑖 𝑟 𝑘 − 𝑦(𝑘) + 𝐾𝑑 2𝑦 𝑘 − 1 − 𝑦(𝑘 − 2) − 𝑦(𝑘) 𝑇 1 𝑘 = 𝑘 𝑖 𝑐 𝑘𝑝 = 𝑘𝑐 − 𝑘𝑖 𝑇𝑖 2 𝑘𝑑 = 𝑘𝑐 𝑇𝑑 𝑇

5.4 PID制御則 PID制御則： 𝐺𝑐 (𝑠) = 𝑈(𝑠) 1 = 𝑘𝑐 1 + + T𝑑 𝑠 𝐸(𝑠) 𝑇𝑖 𝑠 28 kc= ゲイン Ti= I 動作時間 Td= D動作時間 (5-19) 𝐸 𝑠 = 𝑅 𝑠 − 𝑌(𝑠) 「速度アルゴリズム」と「トータルアルゴリズム」の２方式がある （１）速度アルゴリズム：制御出力の1刻み差d(k)を与える方法 次節に例 𝑑 𝑘 = 𝑢 𝑘 − 𝑢(𝑘 − 1) 𝑑 𝑘 = 𝐾𝑝 𝑦 𝑘 − 𝑦(𝑘 − 1) + 𝐾𝑖 𝑟 𝑘 − 𝑦(𝑘) + 𝐾𝑑 2𝑦 𝑘 − 1 − 𝑦(𝑘 − 2) − 𝑦(𝑘) (5-21) （２）トータルアルゴリズム：操作量u(k)そのものを与える方法 𝑑 𝑘 → 1 − 𝑧 −1 𝑈 𝑧 , 𝑦 𝑘 →𝑌 𝑧 , 𝑧 𝑈 𝑧 = −𝐾𝑝 𝑌 𝑧 + 𝐾𝑖 𝑧−1 𝑅 𝑧 − 𝑌 𝑧 𝑟(𝑧) → 𝑅(𝑧) より − 𝐾𝑑 1 − 𝑧 −1 𝑌(𝑧) (5-23)

5.5 PIDゲイン調整 5.1 フィードバック制御系 5.2 １次プラントの制御 5.3 追従系の設計 5.4 PID制御則 5.5 PIDゲイン調整 29

5.5 PIDゲイン調整 P、I、D動作の調節方法にはアナログ制御系についてはジーグラ・ニコルス (Ziegler Nichols)公式が知られている。 高橋安人らは、ジーグラ・ニコルス(Ziegler Nichols)公式をデジタル制御系へ拡 張した。これは、計算機シミュレーションによるものである。但し、この公式 は一つの目安に過ぎないとしている。真の好ましい応答を得るためには、Ki, Kp, Kdのいろいろな値を試す必要があるとした。 次では、一例としてむだ時間つき2時定数のプロセスをPID制御した結果を示す 制御対象： 𝐾𝑒 −𝑠𝐿 𝐺 𝑠 = 1 + 𝑇1 𝑠 1 + 𝑇2 𝑠 30

5.5 PIDゲイン調整（例） 𝐾𝑒 −𝑠𝐿 𝐺 𝑠 = 1 + 𝑇1 𝑠 1 + 𝑇2 𝑠 r(k) y(k) step状外乱 (step:25で) V=0.5を入力 （書籍と同じ） 速度アルゴリズム例 この例は、page70で用いた モデルを利用した。プロセス 特性は表4-5で求めた時系列 形を用いている。すなわち、 g1,g2,g3,g4,g5,g6,pを利用 した。 d(k) u(k) ディジタル制御 その２ (4-4) DGC_no3/page95v4 31

DGC_no3/page95v4 書籍の例 Kp=1.5; Ki=0.5; Kd=0.5 Kp=1.5; Ki=0.35; Kd=0.5 パラメータを 変化させた Kp=1.8; Ki=0.6; Kd=0.5 Kp=2.38; Ki=1.06; Kd=2.19 32

5.5 PIDゲイン調整 33 書籍の図 ステップ外乱応答 “高橋安人,ディジタル制御 ,岩波書店,1987 page94

5.5 PIDゲイン調整 34 𝐾𝑒 −𝑠𝐿 𝐺 𝑠 = 1 + 𝑇1 𝑠 1 + 𝑇2 𝑠 r(k) y(k) ここでは指令値r(k)変化 への追従性を検査した この例は、page70で用いた モデルを利用した。 但し、 プロセス特性は表4-5で求め た時系列形を用いている。 すなわち、 g1,g2,g3,g4,g5,g6,pを利用 した。 d(k) u(k) ディジタル制御 その２ (4-4) DGC_no3/page95v4

5.5 PIDゲイン調整 35 書籍のリスト page95 “高橋安人,ディジタル制御 ,岩波書店,1987 130以前は記載無し

5.5 PIDゲイン調整（制御モデル） むだ時間つき2時定数のプロセス（例） page235 5.5節での制御モデル 加熱室とその室内に加熱物体を置き、 物体を加熱する。室内の温度をx1 (t) 、 物体温度をx2 (t)とする 𝐾𝑒 −𝑠𝐿 𝐺 𝑠 = 1 + 𝑇1 𝑠 1 + 𝑇2 𝑠 36 𝑢 𝑡 𝑅1 𝑣 𝑡 ～ ～ ～ 𝐶1 𝑅2 𝐶2 1. 熱収支 𝑥2 𝑡 − 𝑥1 𝑡 𝑑𝑥1 𝑡 −𝑥1 𝑡 + 𝑣(𝑡) 𝑐 = + 加熱室： 1 𝑑𝑡 𝑅1 𝑅2 𝑥1 𝑡 − 𝑥2 𝑡 𝑑𝑥 𝑡 物体 ： 𝑐2 2 = 𝑑𝑡 𝑅2 𝑥1 𝑡 + 𝑢(𝑡) 2. C1>>C2ならば、加熱室の熱収支で物体の存在を無視できる 𝑋1 (𝑠) 𝑅1 = , 𝑈(𝑠) 𝑇1 𝑠 + 1 𝑋2 (𝑠) 𝑅2 = 𝑋1 (𝑠) 𝑇2 𝑠 + 1 𝑥2 𝑡 𝑇1 = 𝐶1 𝑅1 C1：加熱室の熱容量 R1：加熱室壁の熱伝達抵抗 C2：加熱物体の熱容量 R2：物体表面の熱伝達抵抗 v(t)：外気温度（外乱） 𝑇2 = 𝐶2 𝑅2 3. 加熱物体の温度を出力y(t)とすると、全体の伝達関数は下記となる 𝐺 𝑠 = 𝑌(𝑠) 𝑅1 𝑅2 ≃ 𝑈(𝑠) 𝑇1 𝑠 + 1 𝑇2 𝑠 + 1 𝑌 𝑠 = 𝑋2 (𝑠) この様な関係は、圧力、拡散、混合、 反応過程にも適用されるそうです

5.5 PIDゲイン調整（制御モデル） 1 + 𝑇1′ 𝑠 1 + 𝑇2′𝑠 ⋯ 1 + 𝑇𝑚′ 𝑠 𝐺 𝑠 =K 1 + 𝑇1 𝑠 1 + 𝑇2 𝑠 ⋯ 1 + 𝑇𝑛 𝑠 G(s)の単位ステップ応答は0から立上り、 ゲインKへ暫定して頭を打つ（０型） 37 page236 化学工業プラントの特性はしばしば この形の伝達関数で近似される (n>m) 左図に単位ステップ応答の例を示す。こ の4次系で短い時定数2と1の和3をむだ時 間に置き換え、長い時定数6と10はその ままにすると、次の単位ステップ応答は 4次系のものとほとんど重なる 𝐾𝑒 −3𝑠 𝐺 𝑠 = 1 + 10𝑠 1 + 6𝑠 工業用プロセスではこの原理を使って、 ステップ応答の実験結果に合うように K,L,T1,T2を決めることが良くあるとの こと 𝐾𝑒 −𝑠𝐿 𝐺 𝑠 = 1 + 𝑇1 𝑠 1 + 𝑇2 𝑠

付録 z伝達関数のBode線図 この形の伝達関数の周波数特性が 知りたい として、5.3節では自作programを 用いて応答＆BODE線図を示した。 この付録では low-pass & high-pass filterを例に取り、 自作programと既存programの解析結果 を比較し、その信頼性を検証する 𝐺𝑐 𝑧 = 𝐾𝑐 𝑧+𝑏 𝑧+𝑎 1 + 𝑏𝑒 −𝑗𝜔𝑇 𝐺𝑐 𝑧 = 𝐾𝑐 1 + 𝑎𝑒 −𝑗𝜔𝑇 𝑧 −1 = 𝑒 −𝑗𝜔𝑇 (𝑇: 𝑠𝑎𝑚𝑝𝑙𝑖𝑛𝑔 𝑖𝑛𝑡𝑒𝑟𝑣𝑎𝑙) 𝐺𝑐 𝑧 = 𝐾𝑐 𝑧+𝑏 , 𝑧+𝑎 38 𝑎 <1

位相遅れ＆位相進み(s空間) tsample=0.1 s空間での伝達関数 (python control) T=0.1; alfa=10 Gd=(matlab.tf([T,1],[alfa*T,1]))*alfa; Gd#位相遅れ T=1;bata=0.1 Gl=matlab.tf([T,1],[beta*T,1]); Gl #位相進み matlab.bode(Gd); 位相遅れ 位相遅れ(Gd) 位相進み(Gl) matlab.bode(Gl); 位相進み 39

位相遅れ＆位相進み(Z空間) 40 G(s).ipynb

自作Bode関数との比較 41 「 python-controlで計算したBODE線図」と「自作プログラムのBODE線図」 との比較 位相遅れ： 𝐺𝑑 (𝑠) = 𝛼 𝑇𝑠 + 1 𝛼𝑇𝑠 + 1 matlab.bode これがpython-control で可能ならば自作する 必要はないです 𝑇 = 0.1 𝛼 = 10 𝐺𝑑 (𝑧) = 𝐺𝑐 𝑧 = 𝑧+𝑏 𝑧+𝑎 𝑇𝑠 + 1 𝐺 (𝑠) = 𝑙 位相進み： 𝛽𝑇𝑠 + 1 𝑧−𝑏 𝑧−𝑎 𝑇=1 𝛽 = 0.1 自作Bode関数 ←何らかの方法があるかと思いますが、 自作のほうが簡単なので作ってみました 結論：両者（ｚ空間）はほぼ一致した。但し、連続形（ｓ空間） との比較では離散化するサンプリング間隔が長くなるにつれ精度 は低下する（当たり前か） DGC_no3/ GcBode.ipynb

自作Bode関数との比較：位相遅れ 1/3 matlab.bode 自作Bode関数 42 DGC_no3/GcBode.ipynb 位相遅れ(dt＝0.1）

自作Bode関数との比較：位相遅れ ２/3 matlab.bode 自作Bode関数 43 DGC_no3/GcBode.ipynb 位相遅れ(dt＝0.01）

自作Bode関数との比較：位相遅れ 3/3 matlab.bode 自作Bode関数 44 DGC_no3/GcBode.ipynb 位相遅れ(dt＝0.001）

自作Bode関数との比較：位相進み matlab.bode 自作Bode関数 45 DGC_no3/GcBode.ipynb 位相進み(dt＝0.001）

電気回路：位相遅れ＆位相進み 電気回路で、「位相遅れ」、「位相進み」の 回路を作ってみました。 Ltspiceによるシミュレーション結果も併せて 記します 46

伝達関数： 1 + 𝑗𝜔𝐶𝑅2 1 + 𝑗𝜔𝑇 𝐺= = 1 + 𝑗𝜔𝐶 𝑅1 + 𝑅2 1 + 𝑗𝜔𝛼𝑇 𝛼= Hz 位相遅れ回路 𝐺= 𝑇𝑠 + 1 𝛼𝑇𝑠 + 1 𝐶 𝑅1 + 𝑅2 𝑅1 9 = + 1 = + 1 = 10 𝐶𝑅2 𝑅2 1 𝑇 = 𝐶𝑅2 = 10 × 10−6 × 1 × 103 = 1 × 10−2 𝑓1 = 1 1 = = 1.5915𝐻𝑧 = 10 𝑟𝑎𝑑/𝑠 2𝜋𝐶 𝑅1 + 𝑅2 2𝜋 × 10 × 10−6 × 9 + 1 × 103 𝑓2 = 1 1 = = 15.915𝐻𝑧 = 100 𝑟𝑎𝑑/𝑠 2𝜋𝐶 𝑅2 2𝜋 × 10 × 10−6 × 1 × 103 47

1 + 𝑗𝜔𝐶𝑅2 𝑅1 ∙ 伝達関数： 𝐺 = 𝑅1 𝑅2 𝑅1 + 𝑅2 1 + 𝑗𝜔𝐶 𝑅1 + 𝑅2 𝛽= 𝑅 𝑅2 𝐶 𝑅 1+ 𝑅 1 𝐶𝑅2 2 = 位相進み回路 Hz 𝐺= 𝑅1 𝑇𝑠 + 1 ∙ 𝑅1 + 𝑅2 𝛽𝑇𝑠 + 1 𝑅1 1 1 = = 𝑅1 + 𝑅2 1 + 9 10 𝑇 = 𝐶𝑅2 = 10 × 10−6 × 1 × 103 = 1 × 10−2 𝑓1 = 1 1 = = 1.768𝐻𝑧 = 100/9 𝑟𝑎𝑑/𝑠 2𝜋𝐶 𝑅2 2𝜋 × 10 × 10−6 × 9 × 103 𝑓2 = 1 1 = = 17.684𝐻𝑧 = 1000/9𝑟𝑎𝑑/𝑠 9∗1 𝑅1 𝑅2 −6 3 2𝜋 × 10 × 10 × 9 + 1 × 10 2𝜋𝐶 𝑅 + 𝑅 1 2 48

49 つづく ここまで5章です。次回は6章以降に進みます。 今回PID制御の勉強で「 Pythonによる制御工学入門 著者：南裕樹」を参考にさせて • いただきました。分かり易く、大変参考になりました。感謝いたします そこで下記の記述にびっくり致しました。 • https://y373.sakura.ne.jp/minami/pyctrl 渡すのが、アドレス or 内容なのかの違いで しょうか？明示的でないのが不安です ココ • それ以来、計算結果に齟齬があった場合、私のアルゴリズムの間違いなのか、python言語 へ理解不足の結果なのか、疑心暗鬼となりました。今回GitHubに置いたプログラムリスト には、このような重大なミスがない事を願っていますが.. • 今回ｚ空間でのGc(z)の振る舞いが気になりました。 しかし、Z空間からｓ空間の変換が Python-controlではできず？ Gc(z)のままでBODE線図を作ってみることにしました。 ただ、いろいろと奥が深そうな気がしています。 • あと、時定数とradのことを書こうと思ったのですが、それはまた後日機会があったら。

• https://y373.sakura.ne.jp/minami/pyctrl

50