MATLAB,Python,Scilab,Julia比較 その79【PID制御⑧】

バックナンバーはこちら。
https://www.simulationroom999.com/blog/compare-matlabpythonscilabjulia-backnumber/

はじめに
登場人物
【再掲】PID制御器のブロック線図と全体構成
Scilabコード
シミュレーション結果
まとめ

はじめに

前回は、Python(Numpy)によるDCモータ状態空間モデルをPID制御のシミュレーションを実施。
MATLABと同様の結果は取得でき、問題無く動作確認はできた。

今回は、これのScilab版

登場人物

博識フクロウのフクさん

イラストACにて公開の「kino_k」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=iKciwKA9&area=1

エンジニア歴8年の太郎くん

イラストACにて公開の「しのみ」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=uCKphAW2&area=1

【再掲】PID制御器のブロック線図と全体構成

太郎くん

今回はScilab版のシミュレーションだね。

フクさん

おなじみのPID制御器のブロック線図と全体構成の再掲。

PID制御器(変形式)ブロック線図、目標値、実値、1/zは前回値保持の意、Ki、Kp、Kd、T1、ΔT、e(t)、u(t)

PID制御器を加味した構成図、指令器0~1の範囲で出力、指令器、制御器、PID制御器、DCモータ、状態空間モデル、target(t)、actual(t)、E(t)、u(t)、ω(t)

太郎くん

ScilabはだいたいMATLABと一緒っぽいからあんまり心配はないかな。

Scilabコード

フクさん

というわけでScilabコードが以下になる。

function [x,y] = statespacemodel(A,B,C,D,u,dt,x) 
    // 状態方程式
    x = x + (A*x + B*u) * dt;
    
    // 出力方程式
    y = C*x + D*u;
endfunction


function [state, u]=PIDController(state, target, actual, Kp, Ki, Kd, dt, t1)
    e_i = target - actual;
    e_p = (t1/dt)*(e_i - state.pzi);
    e_d = (t1/dt)*(e_p - state.dzi);
    u = (dt/t1)*(e_p*Kp + e_i*Ki + e_d*Kd) + state.izi;

    state.pzi = e_i;
    state.dzi = e_p;
    state.izi = u;
endfunction


function [] = statespacemodel_pid()
    K=0.016;
    J=0.000000919;
    R=1.34;
    L=0.00012;

    A=[0,1,0 ; 0,0,K/J ; 0,-K/L,-R/L];
    B=[0 ; 0; 1/L];
    C=[1,0,0;0,1,0;0,0,1];
    D=[0;0;0];
    
    dt = 0.0001;
    t = linspace(0, 1, 10000); // 時間(横)軸
    u = zeros(1,10000);         // 入力信号生成
    u(1,5000:10000)=1;          // 5秒後に0から1へ
    y = zeros(3,length(t));
    x = zeros(3,1);
    
    state.pzi = 0;state.dzi = 0;state.izi = 0;
    ratio = 1/60;
    Kp = 0.80;
    Ki = 0.45;
    Kd = 0.0;
    t1 = 0.005;
    
    omega = 0;
    uPID = zeros(1,10000);
    
    for i = 1:length(t)
        [state,uPID(1,i)] = PIDController( state, u(:,i), omega*ratio, Kp, Ki, Kd, dt, t1 );
        [x,y(:,i)] = statespacemodel(A,B,C,D,uPID(:,i),dt,x);
        omega = y(2,i);
    end
    
    subplot(3,1,1);
    plot(t,uPID,'-r', t, u, '--b');
    p=gca();
    p.tight_limits(:)="on";p.data_bounds(:,2)=[-0.1;1.1];p.data_bounds(:,1)=[0;1];
    subplot(3,1,2);
    plot(t',y(1:2,:)');
    p=gca();
    p.tight_limits(:)="on";p.data_bounds(:,2)=[-5;65];p.data_bounds(:,1)=[0;1];
    subplot(3,1,3);
    plot(t,y(3,:));
    p=gca();p.tight_limits(:)="on";p.data_bounds(:,2)=[-0.1;1];p.data_bounds(:,1)=[0;1];
endfunction

太郎くん

予想通りMATLABと一緒か。
構造体も同じ書き方ができるんだね。

フクさん

構造体については書き方は一緒だが、
概念は大きく異なるようだな。
Scilabの構造体は連想配列として実装されているようだ。

太郎くん

連想配列？
なんだっけ？

連想配列（れんそうはいれつ、英語: associative array）とは、コンピュータプログラミングにおいて、添え字にスカラー数値以外のデータ型(文字列型等)も使用できる配列である。抽象データ型のひとつ。連想リスト、連想コンテナ、辞書（あるいはカタカナでディクショナリ英語: dictionary）、ハッシュ（英語: hash）、マップ（英語: map）とも呼ばれる。
Wikipediaより(https://ja.wikipedia.org/wiki/%E9%80%A3%E6%83%B3%E9%85%8D%E5%88%97)

フクさん

具体的にはこんな書き方だな。
添え字が文字列になる感じ。

state('pzi')=1

太郎くん

あ、こんな書き方はみたことあるな。
配列っぽいんだけど、なんかなんでも入れられちゃう感じ。

フクさん

Scilabで以下にようにすると、連想配列である事実はわかるな。

-->state.pzi = 1;
 
-->state.dzi = 2;
 
-->state.izi = 3;
 
-->state('pzi')
 ans  =
 
    1.  
 
-->state('dzi')
 ans  =
 
    2.  
 
-->state('izi')
 ans  =
 
    3.