Simulink モデルでの fminimax の使用

R2026b

この例では、Simulink® モデルのパラメーターを調整する方法を示します。この例を実行すると、モデル optsim が含まれます。このモデルは Simulink ブロック線図としてモデル化された非線形プロセス プラントを含みます。

飽和要素を含むアクチュエータをもつプラント

このプラントは、アクチュエータに制限をもつ低減衰の 3 次のモデルです。アクチュエータは、飽和制限とスルーレート制限をもちます。アクチュエータの飽和制限により、入力値が最大 2 単位、最小 -2 単位の範囲に制限されます。そしてその変化率も 0.8 単位 / 秒以内としています。ステップ入力のシステムの閉ループ応答を閉ループ応答に示します。(コマンド ラインで optsim と入力するか、モデル名をクリックして) モデルを開き、[シミュレーション] メニューから [実行] を選択すると、この応答を確認できます。応答は Scope ブロックにプロット表示されます。

閉ループ応答

問題は、システムへの単位ステップ入力をトラックするフィードバック制御ループの設計です。閉ループ プラントは、プラントとアクチュエータをまとめて 1 つの階層構造の Subsystem ブロックとして表します。Scope ブロックは、設計プロセスの間、出力軌跡を表示します。

閉ループ モデル

このシステムを最適化するには、0 と 100 の間の任意の時刻 $$t$$ での出力の最大値を最小にします。(対照的に、Simulink モデルでの lsqnonlinの例では、出力と入力信号の間の誤差を最小化して解を求めます。)

この例のコードは、補助関数 runtrackmm に含まれています (の例の終わりに掲載)。目的関数は、単に sim コマンドによって返される出力 yout です。しかし、全タイム ステップで最大出力を最小化すると、いくつかのタイム ステップに対して出力が 1 よりはるかに小さくなる可能性があります。はじめの 20 秒後、0.95 を超える出力を保つために、制約関数 trackmmcon は t = 20 から t = 100 まで制約 yout >= 0.95 を含みます。制約は $$g\le 0$$ の形式でなければならないので、関数における制約は g = -yout(20:100) + 0.95 です。

trackmmobj と trackmmcon は共に現在の PID 値から計算された sim からの結果 yout を使用します。シミュレーションを 2 回呼び出さないようにするために、runtrackmm では入れ子関数になっており、yout の値が目的関数と制約関数で共有されるようになっています。シミュレーションは現在の点が変化した場合にのみ呼び出されます。

runtrackmm を呼び出します。

[Kp,Ki,Kd] = runtrackmm

                  Objective        Max     Line search     Directional 
 Iter F-count         value    constraint   steplength      derivative   Procedure 
    0      5              0       1.11982                                            
    1     11          1.184       0.07978            1           0.482     
    2     17          1.012       0.04285            1          -0.236     
    3     23         0.9995      0.007058            1         -0.0186    Hessian modified twice  
    4     29         0.9997     9.708e-07            1         0.00716    Hessian modified  
    5     29         0.9997     9.708e-07            1          -0.162    Hessian modified twice  

Local minimum possible. Constraints satisfied.

fminimax stopped because the size of the current search direction is less than
twice the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.

<stopping criteria details>

Kp = 
0.5910

Ki = 
0.0606

Kd = 
5.5383

出力の Objective value 列の最後の値は、全タイム ステップに対する最大値が 1 よりわずかに少ない値であることを表しています。

この解は問題の定式化を異なる方法で行った Simulink モデルでの lsqnonlin の解とは異なります。

補助関数

次のコードは補助関数 runtrackmm を作成します。

function [Kp, Ki, Kd] = runtrackmm
mdl = "optsim";
open_system(mdl)                             % Load the model
in = Simulink.SimulationInput(mdl);          % Create simulation input object
in = in.setModelParameter(StopTime="100"); % Stop time 100
pid0 = [0.63 0.0504 1.9688];                 % Initial gain values
a1 = 3; a2 = 43;                             % Initialize model plant variables
pold = [];                                   % Initialize pold
yout = [];                                   % Initialize yout
options = optimoptions(@fminimax,Display="iter",...
    StepTolerance=0.001,FunctionTolerance=0.001);
% Optimize the gains
set_param(mdl,FastRestart="on");             % Fast restart
pid = fminimax(@trackmmobj,pid0,[],[],[],[],[],[],...
    @trackmmcon,options);
set_param(mdl,FastRestart="off");
% Return the gains
Kp = pid(1); Ki = pid(2); Kd = pid(3); 

     function F = trackmmobj(pid)
        % Track the output of optsim to a signal of 1.
        % Variables a1 and a2 are shared with runtrackmm.
        % Variable yout is shared with runtrackmm and 
        % trackmmcon.
        updateIfNeeded(pid)
        
        F = yout.yout;
    end

    function [ineqnonlin,eqnonlin] = trackmmcon(pid)
        % Track the output of optsim to a signal of 1.
        % Variable yout is shared with runtrackmm and
        % trackmmobj.
        updateIfNeeded(pid)
        
        ineqnonlin = -yout.yout(20:100)+.95;
        eqnonlin=[];
    end

    function updateIfNeeded(pid)
         if ~isequal(pid,pold) % Compute only if needed.
              in = in.setVariable(Kp=pid(1),Workspace=mdl);
              in = in.setVariable(Ki=pid(2),Workspace=mdl);
              in = in.setVariable(Kd=pid(3),Workspace=mdl);
              
              % Simulate
              yout = sim(in);
              pold = pid;
         end
    end
end

参考

fminimax

トピック

• 多目的最適化