制約付き非線形問題の解法、問題ベース

典型的な最適化問題

この例では、問題ベースのアプローチを使用して制約付き非線形最適化問題を解く方法を示します。例は典型的なワークフローを示します。つまり、目的関数を作成し、制約を作成して、問題を解き、結果を調べます。

メモ:

目的関数または非線形制約が初等関数で構成されていない場合、fcn2optimexpr を使用して、その非線形関数を最適化式に変換しなければなりません。この例の終わりにある fcn2optimexpr を使用した代替定式化、または非線形関数から最適化式への変換を参照してください。

この問題に対するソルバーベースのアプローチについては、[最適化] ライブエディタータスクまたはソルバーを使用した制約付き非線形問題を参照してください。

問題の定式化: Rosenbrock 関数

1 + Rosenbrock 関数の対数を最小化する問題を考えます。

$f (x) = \log (1 + 100 {(x_{2} - x_{1}^{2})}^{2} + (1 - x_{1})^{2}),$

"単位円板" (原点を中心とした半径 1 の円板) で最小化します。つまり、 $x_{1}^{2} + x_{2}^{2} \leq 1$ という条件で関数 $f (x)$ を最小にする $x$ を求めます。この問題は非線形制約付きの非線形関数の最小化です。

Rosenbrock 関数は最適化の標準的なテスト関数であり、点 [1,1] で一意に最小値 0 に到達するため、 $f (x)$ は同じ点で同じ最小値に到達します。この関数は曲線の深い谷に奥行きのない最小値をもつため、最小値の検索がアルゴリズムによっては困難になります。この問題の解は、点 [1,1] ではありません。この点は制約を満たしていないからです。

この図は単位円板内の関数 $f (x)$ の 2 つのビューを示します。等高線は表面プロットの下に描きます。

rosenbrock = @(x)log(1 + 100*(x(:,2) - x(:,1).^2).^2 + (1 - x(:,1)).^2); % Vectorized function

figure1 = figure('Position',[1 200 600 300]);
colormap('gray');
axis square;
R = 0:.002:1;
TH = 2*pi*(0:.002:1); 
X = R'*cos(TH);
Y = R'*sin(TH); 
Z = rosenbrock([X(:),Y(:)]); % Z = f(x)
Z = reshape(Z,size(X));

% Create subplot
subplot1 = subplot(1,2,1,'Parent',figure1);
view([124 34]);
grid('on');
hold on;

% Create surface
surf(X,Y,Z,'Parent',subplot1,'LineStyle','none');

% Create contour
contour(X,Y,Z,'Parent',subplot1);

% Create subplot
subplot2 = subplot(1,2,2,'Parent',figure1);
view([234 34]);
grid('on');
hold on

% Create surface
surf(X,Y,Z,'Parent',subplot2,'LineStyle','none');

% Create contour
contour(X,Y,Z,'Parent',subplot2);

% Create textarrow
annotation(figure1,'textarrow',[0.4 0.31],...
    [0.055 0.16],...
    'String',{'Minimum at (0.7864,0.6177)'});

% Create arrow
annotation(figure1,'arrow',[0.59 0.62],...
    [0.065 0.34]);

title("Rosenbrock's Function: Two Views")

hold off

関数ハンドル rosenbrock は、任意の数の 2-D 点での関数 $f (x)$ を同時に計算します。このベクトル化は関数のプロットを高速化します。また、その他のコンテキストでは複数の点での関数の評価を高速化するのに役立つ場合があります。

関数 $f (x)$ は "目的関数" と呼ばれます。目的関数は、最小化する関数です。不等式 $x_{1}^{2} + x_{2}^{2} \leq 1$ は "制約" と呼ばれます。制約は、ソルバーが最小値を探す $x$ の集合の範囲を限定します。任意の数の制約 (不等式または方程式) をもつことができます。

最適化変数を使用した問題の定義

最適化に対する問題ベースのアプローチでは、最適化変数を使用して目的関数と制約を定義します。これらの変数を使用して式を作成するには、2 つのアプローチがあります。

多項式関数や三角関数などの初等関数の場合、変数で式を直接記述します。
他のタイプの関数の場合、fcn2optimexpr を使用して関数を最適化式に変換します。この例の終わりにある「fcn2optimexpr を使用した代替定式化」を参照してください。

この問題では、目的関数と非線形制約の両方が初等であるため、最適化変数で式を直接記述できます。'x' という名前の 2 次元最適化変数を作成します。

x = optimvar('x',1,2);

最適化変数の式として目的関数を作成します。

obj = log(1 + 100*(x(2) - x(1)^2)^2 + (1 - x(1))^2);

obj を目的関数とする prob という名前の最適化問題を作成します。

prob = optimproblem('Objective',obj);

最適化変数の多項式として非線形制約を作成します。

nlcons = x(1)^2 + x(2)^2 <= 1;

非線形制約を問題に含めます。

prob.Constraints.circlecons = nlcons;

問題を確認します。

show(prob)

  OptimizationProblem : 

	Solve for:
       x

	minimize :
       log(((1 + (100 .* (x(2) - x(1).^2).^2)) + (1 - x(1)).^2))


	subject to circlecons:
       (x(1).^2 + x(2).^2) <= 1

問題を解く

最適化問題を解くには、solve を呼び出します。問題には初期点、つまり最適化変数の初期値を指定する構造体が必要です。[0 0] という $x$ 値を持つ初期点構造体 x0 を作成します。

x0.x = [0 0];
[sol,fval,exitflag,output] = solve(prob,x0)

Solving problem using fmincon.

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

sol = struct with fields:
    x: [0.7864 0.6177]

fval = 0.0447

exitflag = 
    OptimalSolution

output = struct with fields:
              iterations: 23
               funcCount: 44
         constrviolation: 0
                stepsize: 7.7862e-09
               algorithm: 'interior-point'
           firstorderopt: 8.0000e-08
            cgiterations: 8
                 message: 'Local minimum found that satisfies the constraints....'
            bestfeasible: [1x1 struct]
     objectivederivative: "reverse-AD"
    constraintderivative: "closed-form"
                  solver: 'fmincon'

解の検証

解には exitflag = OptimalSolution と表示されます。この終了フラグは、解が局所的最適解であることを示します。より良い解を見つける方法については、ソルバーが成功する場合を参照してください。

終了メッセージは、解が制約を満たしていることを示します。解が実際に実行可能であることをいくつかの方法で確認できます。

output 構造体の constrviolation フィールドで報告された実行不可能性を確認する。

infeas = output.constrviolation

infeas = 0

実行不可能性が 0 の場合、解が実行可能であることを示します。

解での実行不可能性を計算する。

infeas = infeasibility(nlcons,sol)

infeas = 0

ここでも、実行不可能性が 0 の場合、解が実行可能であることを示します。

x のノルムを計算して 1 以下であることを確認する。

nx = norm(sol.x)

nx = 1.0000

output 構造体は、反復の回数 (23)、ソルバー (fmincon)、および関数評価の回数 (44) など、解法プロセスに関する詳細を示します。これらの統計量の詳細については、許容誤差と停止条件を参照してください。

`fcn2optimexpr` を使用した代替定式化

より複雑な式の場合、目的関数または制約関数の関数ファイルを記述し、fcn2optimexpr を使用して最適化式に変換します。たとえば、disk.m ファイルに非線形制約関数の基底があるとします。

type disk

function radsqr = disk(x) 

radsqr = x(1)^2 + x(2)^2;

この関数ファイルを最適化式に変換します。

radsqexpr = fcn2optimexpr(@disk,x);

さらに、プロットルーチンの最初に定義された関数ハンドル rosenbrock を最適化式に変換します。

rosenexpr = fcn2optimexpr(rosenbrock,x);

これらの変換後の最適化式を使用して最適化問題を作成します。

convprob = optimproblem('Objective',rosenexpr,'Constraints',radsqexpr <= 1);

新しい問題を表示します。

show(convprob)

  OptimizationProblem : 

	Solve for:
       x

	minimize :
       log(((1 + (100 .* (x(2) - x(1).^2).^2)) + (1 - x(1)).^2))


	subject to :
       (x(1).^2 + x(2).^2) <= 1

新しい問題を解きます。解は基本的に以前と同じです。

[sol,fval,exitflag,output] = solve(convprob,x0)

Solving problem using fmincon.

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

sol = struct with fields:
    x: [0.7864 0.6177]

fval = 0.0447

exitflag = 
    OptimalSolution

output = struct with fields:
              iterations: 23
               funcCount: 44
         constrviolation: 0
                stepsize: 7.7862e-09
               algorithm: 'interior-point'
           firstorderopt: 8.0000e-08
            cgiterations: 8
                 message: 'Local minimum found that satisfies the constraints....'
            bestfeasible: [1x1 struct]
     objectivederivative: "reverse-AD"
    constraintderivative: "closed-form"
                  solver: 'fmincon'

サポートされている関数の一覧については、最適化変数および式でサポートされる演算を参照してください。