Parallel Computing Toolbox を使用した、時間のかかる最適化問題の最小化

この例では次を使用します。

この例では、Optimization Toolbox™ と Global Optimization Toolbox の関数を使用して、時間のかかる最適化問題の最小化を高速で行う方法を説明します。例の前半では、関数を連続的に評価することで最適化問題を解きます。例の後半では、並列 for ループ (parfor) 機能を使用して関数を並列で評価することにより、同じ問題を解きます。最適化関数の処理にかかる時間を、両方のケースで比較します。

時間のかかる最適化問題

この例では、4 つの変数をもつ問題を解きます。一時停止することで、意図的に目的関数と制約関数の処理に時間がかかるようにします。

function f = expensive_objfun(x)
%EXPENSIVE_OBJFUN An expensive objective function used in optimparfor example.

%   Copyright 2007-2013 The MathWorks, Inc.

% Simulate an expensive function by pausing
pause(0.1)
% Evaluate objective function
f = exp(x(1)) * (4*x(3)^2 + 2*x(4)^2 + 4*x(1)*x(2) + 2*x(2) + 1);

function [c,ceq] = expensive_confun(x)
%EXPENSIVE_CONFUN An expensive constraint function used in optimparfor example.

%   Copyright 2007-2013 The MathWorks, Inc.

% Simulate an expensive function by pausing
pause(0.1);
% Evaluate constraints
c = [1.5 + x(1)*x(2)*x(3) - x(1) - x(2) - x(4); 
     -x(1)*x(2) + x(4) - 10];
% No nonlinear equality constraints:
ceq = [];

`fmincon` を使用した最小化

逐次処理で fmincon にかかる時間を測定し、並列実行の時間と比較できるようにします。

startPoint = [-1 1 1 -1];
options = optimoptions('fmincon','Display','iter','Algorithm','interior-point');
startTime = tic;
xsol = fmincon(@expensive_objfun,startPoint,[],[],[],[],[],[],@expensive_confun,options);
time_fmincon_sequential = toc(startTime);
fprintf('Serial FMINCON optimization takes %g seconds.\n',time_fmincon_sequential);

                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       5    1.839397e+00    1.500e+00    3.211e+00
    1      11   -9.760099e-01    3.708e+00    7.902e-01    2.362e+00
    2      16   -1.480976e+00    0.000e+00    8.344e-01    1.069e+00
    3      21   -2.601599e+00    0.000e+00    8.390e-01    1.218e+00
    4      29   -2.823630e+00    0.000e+00    2.598e+00    1.118e+00
    5      34   -3.905338e+00    0.000e+00    1.210e+00    7.302e-01
    6      39   -6.212992e+00    3.934e-01    7.372e-01    2.405e+00
    7      44   -5.948761e+00    0.000e+00    1.784e+00    1.905e+00
    8      49   -6.940062e+00    1.233e-02    7.668e-01    7.553e-01
    9      54   -6.973887e+00    0.000e+00    2.549e-01    3.920e-01
   10      59   -7.142993e+00    0.000e+00    1.903e-01    4.735e-01
   11      64   -7.155325e+00    0.000e+00    1.365e-01    2.626e-01
   12      69   -7.179122e+00    0.000e+00    6.336e-02    9.115e-02
   13      74   -7.180116e+00    0.000e+00    1.069e-03    4.670e-02
   14      79   -7.180409e+00    0.000e+00    7.797e-04    2.815e-03
   15      84   -7.180410e+00    0.000e+00    6.368e-06    3.120e-04

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.

Serial FMINCON optimization takes 18.2876 seconds.

遺伝的アルゴリズムを使用した最小化

遺伝的アルゴリズム (ga) では fmincon よりもずっと多くの関数評価が行われるため、時間のかかる制約をこの問題から取り除き、代わりに制約なしの最適化を実行します。そのために、制約として空行列 [] を渡します。また、妥当な時間内に ga が終了できるように、ga の最大世代数を 15 に制限します。ga の処理にかかる時間を測定して、並列での ga の評価と比較できるようにします。ga の実行には、Global Optimization Toolbox が必要なことに注意してください。

rng default % for reproducibility
try
    gaAvailable = false;
    nvar = 4;
    gaoptions = optimoptions('ga','MaxGenerations',15,'Display','iter');
    startTime = tic;
    gasol = ga(@expensive_objfun,nvar,[],[],[],[],[],[],[],gaoptions);
    time_ga_sequential = toc(startTime);
    fprintf('Serial GA optimization takes %g seconds.\n',time_ga_sequential);
    gaAvailable = true;
catch ME
    warning(message('optim:optimdemos:optimparfor:gaNotFound'));
end

Single objective optimization:
4 Variable(s)

Options:
CreationFcn:       @gacreationuniform
CrossoverFcn:      @crossoverscattered
SelectionFcn:      @selectionstochunif
MutationFcn:       @mutationgaussian

                                  Best           Mean      Stall
Generation      Func-count        f(x)           f(x)    Generations
    1              100      -5.546e+05       1.483e+15        0
    2              147      -5.581e+17      -1.116e+16        0
    3              194      -7.556e+17       6.679e+22        0
    4              241      -7.556e+17      -7.195e+16        1
    5              288      -9.381e+27      -1.876e+26        0
    6              335      -9.673e+27      -7.497e+26        0
    7              382      -4.511e+36      -9.403e+34        0
    8              429      -5.111e+36      -3.011e+35        0
    9              476      -7.671e+36       9.346e+37        0
   10              523       -1.52e+43      -3.113e+41        0
   11              570      -2.273e+45       -4.67e+43        0
   12              617      -2.589e+47      -6.281e+45        0
   13              664      -2.589e+47      -1.015e+46        1
   14              711      -8.149e+47      -5.855e+46        0
   15              758      -9.503e+47       -1.29e+47        0
Optimization terminated: maximum number of generations exceeded.
Serial GA optimization takes 81.5878 seconds.

Parallel Computing Toolbox の設定

Parallel Computing Toolbox™ を使用でき、かつワーカーの並列プールが存在する場合、Optimization Toolbox の関数によって導関数を近似する有限差分化は、parfor 機能を使用して並列的に行われます。同様に、Global Optimization Toolbox の ga、gamultiobj および patternsearch の各ソルバーは、関数を並列に評価します。parfor 機能を使用するために、関数 parpool を使用して並列環境を設定します。この例のパブリッシュ先となるコンピューターには 4 つのコアがあるため、parpool は 4 つの MATLAB® ワーカーを起動します。この例を実行する際に並列プールが既に存在する場合は、そのプールを使用します。詳細については、parpool のドキュメンテーションを参照してください。

if max(size(gcp)) == 0 % parallel pool needed
    parpool % create the parallel pool
end

並列な `fmincon` を使用した最小化

関数 fmincon を並列で使用して時間のかかる最適化問題を最小化するには、目的関数と制約関数が並列で評価されうること、そして fmincon にその並列な機能性を可能な限り使用させることを明示的に指定する必要があります。現在、有限差分を並列で実行できます。fmincon の処理にかかる時間を測定して、fmincon の逐次的実行と比較できるようにします。

options = optimoptions(options,'UseParallel',true);
startTime = tic;
xsol = fmincon(@expensive_objfun,startPoint,[],[],[],[],[],[],@expensive_confun,options);
time_fmincon_parallel = toc(startTime);
fprintf('Parallel FMINCON optimization takes %g seconds.\n',time_fmincon_parallel);

                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       5    1.839397e+00    1.500e+00    3.211e+00
    1      11   -9.760099e-01    3.708e+00    7.902e-01    2.362e+00
    2      16   -1.480976e+00    0.000e+00    8.344e-01    1.069e+00
    3      21   -2.601599e+00    0.000e+00    8.390e-01    1.218e+00
    4      29   -2.823630e+00    0.000e+00    2.598e+00    1.118e+00
    5      34   -3.905338e+00    0.000e+00    1.210e+00    7.302e-01
    6      39   -6.212992e+00    3.934e-01    7.372e-01    2.405e+00
    7      44   -5.948761e+00    0.000e+00    1.784e+00    1.905e+00
    8      49   -6.940062e+00    1.233e-02    7.668e-01    7.553e-01
    9      54   -6.973887e+00    0.000e+00    2.549e-01    3.920e-01
   10      59   -7.142993e+00    0.000e+00    1.903e-01    4.735e-01
   11      64   -7.155325e+00    0.000e+00    1.365e-01    2.626e-01
   12      69   -7.179122e+00    0.000e+00    6.336e-02    9.115e-02
   13      74   -7.180116e+00    0.000e+00    1.069e-03    4.670e-02
   14      79   -7.180409e+00    0.000e+00    7.797e-04    2.815e-03
   15      84   -7.180410e+00    0.000e+00    6.368e-06    3.120e-04

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.

Parallel FMINCON optimization takes 8.79291 seconds.

並列な遺伝的アルゴリズムを使用した最小化

関数 ga を並列で使用して時間のかかる最適化問題を最小化するには、目的関数が並列で評価されうること、そして ga にその並列な機能性を可能な限り使用させることを明示的に指定する必要があります。並列 ga を使用するには、[ベクトル化] オプションを既定値 (つまり [off]) に設定することも必要です。ここでも、ga の処理にかかる時間を測定して、ga の逐次的実行と比較できるようにします。ga では乱数発生器が使用されるため、この実行は逐次的実行とは異なることがありますが、時間のかかる関数評価の回数は両方の実行で同じです。ga の実行には、Global Optimization Toolbox が必要なことに注意してください。

rng default % to get the same evaluations as the previous run
if gaAvailable
    gaoptions = optimoptions(gaoptions,'UseParallel',true);
    startTime = tic;
    gasol = ga(@expensive_objfun,nvar,[],[],[],[],[],[],[],gaoptions);
    time_ga_parallel = toc(startTime);
    fprintf('Parallel GA optimization takes %g seconds.\n',time_ga_parallel);
end

Single objective optimization:
4 Variable(s)

Options:
CreationFcn:       @gacreationuniform
CrossoverFcn:      @crossoverscattered
SelectionFcn:      @selectionstochunif
MutationFcn:       @mutationgaussian

                                  Best           Mean      Stall
Generation      Func-count        f(x)           f(x)    Generations
    1              100      -5.546e+05       1.483e+15        0
    2              147      -5.581e+17      -1.116e+16        0
    3              194      -7.556e+17       6.679e+22        0
    4              241      -7.556e+17      -7.195e+16        1
    5              288      -9.381e+27      -1.876e+26        0
    6              335      -9.673e+27      -7.497e+26        0
    7              382      -4.511e+36      -9.403e+34        0
    8              429      -5.111e+36      -3.011e+35        0
    9              476      -7.671e+36       9.346e+37        0
   10              523       -1.52e+43      -3.113e+41        0
   11              570      -2.273e+45       -4.67e+43        0
   12              617      -2.589e+47      -6.281e+45        0
   13              664      -2.589e+47      -1.015e+46        1
   14              711      -8.149e+47      -5.855e+46        0
   15              758      -9.503e+47       -1.29e+47        0
Optimization terminated: maximum number of generations exceeded.
Parallel GA optimization takes 15.2253 seconds.

逐次実行と並列実行の時間の比較

X = [time_fmincon_sequential time_fmincon_parallel];
Y = [time_ga_sequential time_ga_parallel];
t = [0 1];
plot(t,X,'r--',t,Y,'k-')
ylabel('Time in seconds')
legend('fmincon','ga')
ax = gca;
ax.XTick = [0 1];
ax.XTickLabel = {'Serial' 'Parallel'};
axis([0 1 0 ceil(max([X Y]))])
title('Serial Vs. Parallel Times')

parfor によって関数評価を並列利用することで、fmincon と ga の効率はどちらも向上しました。こうした向上は多くの場合、時間のかかる目的関数と制約関数により適しています。