このページは機械翻訳を使用して翻訳されました。最新版の英語を参照するには、ここをクリックします。

パターン探索を用いた制約付き最小化、問題ベース

この例では、問題ベースのアプローチでパターン探索を使用して、非線形不等式制約と境界に従って目的関数を最小化する方法を示します。この問題のソルバーベースのバージョンについては、パターン探索を用いた制約付き最小化、ソルバーベースを参照してください。

制約付き最小化問題

この問題では、最小化する目的関数は 2 次元変数 X と Y の単純な関数です。

camxy = @(X,Y)(4 - 2.1.*X.^2 + X.^4./3).*X.^2 + X.*Y + (-4 + 4.*Y.^2).*Y.^2;

この機能は、L.C.W. Dixon と G.P. Szego [1] で説明されているように、「カム」として知られています。

さらに、この問題には非線形の制約と境界があります。

   x(1)*x(2) + x(1) - x(2) + 1.5 <= 0  (nonlinear constraint)
   10 - x(1)*x(2) <= 0                 (nonlinear constraint)
   0 <= x(1) <= 1                      (bound)
   0 <= x(2) <= 13                     (bound)

目的関数のサーフェスプロット上に非線形制約領域をプロットします。制約により、解は両方の赤い曲線の上にある小さな領域に制限されます。

x1 = linspace(0,1);
y1 = (-x1 - 1.5)./(x1 - 1);
y2 = 10./x1;
[X,Y] = meshgrid(x1,linspace(0,13));
Z = camxy(X,Y);
surf(X,Y,Z,"LineStyle","none")
hold on
z1 = camxy(x1,y1);
z2 = camxy(x1,y2);
plot3(x1,y1,z1,'r-',x1,y2,z2,'r-')
xlim([0 1])
ylim([0 13])
zlim([0,max(Z,[],"all")])
hold off

Figure contains an axes object. The axes object contains 3 objects of type surface, line.

最適化変数、問題、制約を作成する

この問題を設定するには、最適化変数 x と y を作成します。変数を作成するときに境界を設定します。

x = optimvar("x","LowerBound",0,"UpperBound",1);
y = optimvar("y","LowerBound",0,"UpperBound",13);

最適化式として目標を作成します。

cam = camxy(x,y);

この目的関数を使用して最適化問題を作成します。

prob = optimproblem("Objective",cam);

2 つの非線形不等式制約を作成し、問題に含めます。

prob.Constraints.cons1 = x*y + x - y + 1.5 <= 0;
prob.Constraints.cons2 = 10 - x*y <= 0;

問題を確認します。

show(prob)

  OptimizationProblem : 

	Solve for:
       x, y

	minimize :
       (((((4 - (2.1 .* x.^2)) + (x.^4 ./ 3)) .* x.^2) + (x .* y)) + (((-4) + (4 .* y.^2)) .* y.^2))


	subject to cons1:
       ((((x .* y) + x) - y) + 1.5) <= 0

	subject to cons2:
       (10 - (x .* y)) <= 0

	variable bounds:
       0 <= x <= 1

       0 <= y <= 13

初期点を設定して解く

初期点を、フィールド x が 0.5 に等しく、y が 0.5 に等しい構造体として設定します。

x0.x = 0.5;
x0.y = 0.5;

patternsearch ソルバーを指定して問題を解決します。

[sol,fval] = solve(prob,x0,"Solver","patternsearch")

Solving problem using patternsearch.
Optimization finished: mesh size less than options.MeshTolerance 
and constraint violation is less than options.ConstraintTolerance.

sol = struct with fields:
    x: 0.8122
    y: 12.3122

fval = 
9.1324e+04

patternsearch は、目的関数値が 9.1324e4 である解点 x = 0.8122、y = 12.3122 を見つけます。

視覚化を追加する

ソルバーの進行状況を観察するには、2 つのプロット関数を選択するオプションを指定します。プロット関数 psplotbestf は、各反復における最適な目的関数の値をプロットし、プロット関数 psplotmaxconstr は、各反復における最大制約違反をプロットします。これら 2 つのプロット関数をセル配列に設定します。また、Display オプションを 'iter' に設定して、コマンドウィンドウにソルバーの進行状況に関する情報を表示します。

options = optimoptions(@patternsearch,...
    "PlotFcn",{@psplotbestf,@psplotmaxconstr},...
    "Display","iter");

options 引数を含めてソルバーを実行します。

[sol,fval] = solve(prob,x0,"Solver","patternsearch","Options",options)

Solving problem using patternsearch.

                                      Max
  Iter   Func-count       f(x)      Constraint   MeshSize      Method
    0         1     0.373958         9.75       0.9086    
    1        18       113581    1.617e-10        0.001   Increase penalty
    2       147        92267            0        1e-05   Increase penalty
    3       373      91333.2            0        1e-07   Increase penalty
    4       638        91324            0        1e-09   Increase penalty
Optimization finished: mesh size less than options.MeshTolerance 
and constraint violation is less than options.ConstraintTolerance.

Figure Pattern Search contains 2 axes objects. Axes object 1 with title Best Function Value: 91324, xlabel Iteration, ylabel Function value contains an object of type scatter. Axes object 2 with title Maximum Constraint Violation: 0, xlabel Iteration, ylabel Constraint violation contains an object of type scatter.

sol = struct with fields:
    x: 0.8122
    y: 12.3122

fval = 
9.1324e+04

非線形制約により、patternsearch は各反復で多くのサブ問題を解決します。プロットと反復表示の両方に示されているように、解決プロセスには反復がほとんどありません。ただし、反復表示の Func-count 列には、反復ごとに多数の関数評価が表示されます。プロットと反復表示の両方から、初期点が実行不可能であり、目的関数が初期点で低いことがわかります。解決プロセス中、目的関数の値は最初は増加し、その後最終値まで減少します。

サポートされていない機能

目的関数または非線形制約関数が最適化変数および式でサポートされる演算ではない場合は、fcn2optimexpr を使用して、問題ベースのアプローチに適した形式に変換します。たとえば、制約 $x y \geq 10$ の代わりに制約 $I_{1} (x) + I_{1} (y) \geq 10$ があり、ここで $I_{1} (x)$ は修正ベッセル関数 besseli(1,x) であるとします。(ベッセル関数はサポートされていない関数です。)次のように fcn2optimexpr を使用してこの制約を作成します。まず、 $I_{1} (x) + I_{1} (y)$ の最適化式を作成します。

bfun = fcn2optimexpr(@(t,u)besseli(1,t) + besseli(1,u),x,y);

次に、制約 cons2 を制約 bfun >= 10 に置き換えます。

prob.Constraints.cons2 = bfun >= 10;

問題を解きます。制約領域が異なるため、解も異なります。

[sol2,fval2] = solve(prob,x0,"Solver","patternsearch","Options",options)

Solving problem using patternsearch.

                                      Max
  Iter   Func-count       f(x)      Constraint   MeshSize      Method
    0         1     0.373958        9.484       0.9307    
    1        18       113581            0        0.001   Increase penalty
    2       183      962.841            0        1e-05   Increase penalty
    3       499      960.942            0        1e-07   Increase penalty
    4       636       960.94            0    8.511e-15   Update multipliers
Optimization finished: mesh size less than options.MeshTolerance 
and constraint violation is less than options.ConstraintTolerance.

Figure Pattern Search contains 2 axes objects. Axes object 1 with title Best Function Value: 960.94, xlabel Iteration, ylabel Function value contains an object of type scatter. Axes object 2 with title Maximum Constraint Violation: 0, xlabel Iteration, ylabel Constraint violation contains an object of type scatter.

sol2 = struct with fields:
    x: 0.4998
    y: 3.9981

fval2 = 
960.9401

参考文献

[1] Dixon, L. C. W., and G .P. Szego (eds.).Towards Global Optimisation 2.North-Holland:Elsevier Science Ltd., Amsterdam, 1978.

参考

patternsearch | solve