Matlab_Simulink：PIDゲイン最適化プログラム③_入力飽和考慮（無料公開）

無料教材：matlab

実現したいこと

  • Simulinkとmファイルを用いたPID制御シミュレータの構築
  • 上記プログラムのPIDゲインの最適化
  • fminsearchを使っての最適化
  • 過大な入力が印加されないような最適化
  • 入力飽和にペナルティを与えて最適化
    ※厳密にはサチュレーションブロックを用いない最適化

Matlabのバージョン

Matlab2021a
※ダウンロード形式として過去のバージョンも用意

必要なtoolbox

ダウンロードURL

【ダウンロードリンク】

PID_Control_fmin_satu - Google ドライブ 

※ご利用中のMatlabバージョンのフォルダをダウンロードしてください。
　例：2021aバージョンであれば，上記リンクの「2021a」を選択
※上記プログラムの利用で生じたトラブルは一切の責任を負いかねます。

フォルダ構成

上記ダウンロードし展開。
（下記のようなフォルダ構成になっていることを確認）

実行方法

「Main_PID_fmin_input.m」を実行
※PCスペックに依って最適化計算に5分程度時間がかかります。

Simulinkファイル（PID_fmin.slx）

mファイルソースコード（Main_PID_fmin_input_satu.m）

※下記赤文字が過去の記事からの変更点

%% 初期化
clc % コマンドウィンドウの初期化
clear % ワークスペースの初期化
close all % グラフを全部閉じる

%% 変数宣言
% ===== 制御対象 =====
K = 1;
A = 1;
B = 1;
C = 1;
% ===== 目標値 =====
r_val = 10; % ステップ状の目標値
% ===== Simulink関係 =====
N_File = 'PID_fmin.slx'; % simulinkファイル名
FinalTime = 30; % シミュレーション終了時刻[s]
SamplingTime = 0.1; % サンプリング時間[s]
% ===== 推定システムゲイン =====
K_hat = K; % fminsearchで偏差と入力の値のオーダーを一致させることに利用
% ===== 入力の上限値 =====
th_u = 11; % u(t) > th_u の場合，最適化評価規範fを2倍悪化させる 
 
%% fminsearch
ini_x = [0.5 0.5, 0.5]; % PIDゲインの初期値
[x, fval] = fminsearch(@(x) func_fmin_input_satu(x,N_File,K_hat,th_u),ini_x);

%% Simulinkの実行
open(N_File); % Simulinkを起動
sim(N_File); % Simulinkの実行

%% グラフ化
% ===== Simulinkデータ格納 =====
t = ScopeData.time; % 時刻情報
y = ScopeData.signals(1).values(:, 1); % 出力
r = ScopeData.signals(1).values(:, 2); % 目標値
u = ScopeData.signals(2).values(:, 1); % 入力
% ===== グラフ描画 =====
% ----- y -----
figure;
subplot(211);
plot(t,y);
hold on;
plot(t,r,'--');
grid on;
xlabel('$ t {\rm [s]} $', 'interpreter', 'latex');
ylabel('$ y(t) $', 'interpreter', 'latex');
legend('$ y(t) $', '$ r(t) $', 'interpreter', 'latex');
% ----- u -----
subplot(212);
plot(t,u);
hold on;
grid on;
xlabel('$ t {\rm [s]} $', 'interpreter', 'latex');
ylabel('$ u(t) $', 'interpreter', 'latex');

mファイルソースコード（func_fmin_input_satu.m）

function f = func_fmin_input_satu(x,N_File,K_hat,th_u)

%% 変数宣言
Kp = x(1);
Ki = x(2);
Kd = x(3);
assignin('base','Kp',Kp);
assignin('base','Ki',Ki);
assignin('base','Kd',Kd);

% Simulink実行
sim(N_File);

%% Scopeデータの取得
t = ScopeData.time; % 時刻情報
y = ScopeData.signals(1).values(:, 1); % 出力
r = ScopeData.signals(1).values(:, 2); % 目標値
u = ScopeData.signals(2).values(:, 1); % 入力

 

%% 評価規範
% 入力の整定値u(end), K_hatを導入することで，偏差(r-y)と同じオーダーにする
 f = sum( ( r-y ).^2 ) + ( K_hat.*(u(end)-u)).^2 ); 
 
% u(t) > th_u の場合，最適化評価規範fを2倍悪化させる
if max(u) > th_u
    f = 2*f;
end

 

%% コマンドウィンドウ表示
fprintf('[Kp, Ki, Kd, fval] = %.3g, %.3g, %.3g, %.3g\n',Kp, Ki, Kd, f);

実行結果

上記結果の通り，設定した飽和値th_u=11を超えない出力結果を得ることができた。

まとめ・考察

上記赤文字の箇所の通り，if文を用いて，指定した飽和th_uを超えた場合に評価値 fを2倍悪化させている。

もちろん f=\inftyとして最大の悪化値を与えても良いのだが，その場合，PIDゲインの最適化が進まないため，今回は上記の対応を実施した。
（ f=\inftyの場合，例えば K_Pを小さくすると評価値が良くなるか悪化するか分からず各ゲイン修正の勾配情報が失われる。一方， f=2fとすることで， K_Pを小さくした場合でも fの値に変化があるため，修正のための勾配情報は維持される）

本記事は，全く学術的ではなく力業のアプローチであるが，取り急ぎの最適化などには活用できるのではないかと考える。

次回の記事では，その他評価規範による出力結果の違いについてついて執筆予定である。

併せて確認推奨の過去記事

forfree.hatenablog.jp

forfree.hatenablog.jp

forfree.hatenablog.jp

forfree.hatenablog.jp