Main Content

2 次元の半無限制約

以下を最小化する x の値を探します。

f(x) = (x1 – 0.2)2 + (x2– 0.2)2 + (x3– 0.2)2,

ここで、

K1(x,w)=sin(w1x1)cos(w2x2)11000(w150)2sin(w1x3)x3+...                 sin(w2x2)cos(w1x1)11000(w250)2sin(w2x3)x31.5,

すべての値に対して、w1 と w2 の範囲は以下になります。

1 ≤ w1 ≤ 100,
1 ≤ w2 ≤ 100,

x = [0.25,0.25,0.25] を開始値とします。

半無限制約が 2 次元、すなわち、行列であることに注意してください。

まずはじめに目的関数を計算するファイルを記述します。

function f = myfun(x,s)
% Objective function
f = sum((x-0.2).^2);

次に mycon.m と呼ばれる制約のファイルを記述してください。mycon が呼び出される度に半無限制約の表面プロットを描くコードが含まれています。X が最小化されていくとき、どのように制約が変化しているのかを見ることができます。

function [c,ceq,K1,s] = mycon(X,s)
% Initial sampling interval
if isnan(s(1,1)),
   s = [2 2];
end

% Sampling set
w1x = 1:s(1,1):100;
w1y = 1:s(1,2):100;
[wx,wy] = meshgrid(w1x,w1y);

% Semi-infinite constraint 
K1 = sin(wx*X(1)).*cos(wx*X(2))-1/1000*(wx-50).^2 -...
       sin(wx*X(3))-X(3)+sin(wy*X(2)).*cos(wx*X(1))-...
       1/1000*(wy-50).^2-sin(wy*X(3))-X(3)-1.5;

% No finite nonlinear constraints
c = []; ceq=[];

% Mesh plot
m = surf(wx,wy,K1,'edgecolor','none','facecolor','interp');
camlight headlight
title('Semi-infinite constraint')
drawnow

次に最適化ルーチンを呼び出します。

x0 = [0.25, 0.25, 0.25];    % Starting guess
[x,fval] = fseminf(@myfun,x0,1,@mycon)

9 回の反復後、解は以下になります。

x
x =
    0.2523    0.1715    0.1938

解 x での関数値および半無限制約の最大値は、次のとおりです。

fval
fval =
    0.0036

ゴールは半無限制約 K1(x,w) ≤ 1.5 を満たすように目的関数 f(x) を最小化することでした。解 xmycon を評価し、行列 K1 の最大要素を見つけると、制約が簡単に満たされることが示されます。

[c,ceq,K1] = mycon(x,[0.5,0.5]);  % Sampling interval 0.5
max(max(K1))
ans =
   -0.0370

この mycon の呼び出しによって、解 x の半無限制約を示した次の表面プロットが作成されます。

Wavy surface less than or equal to zero, with some points attaining zero

参考

関連するトピック