Пример для Optimization Toolbox™

Этот пример показывает, как использовать два нелинейных решателя оптимизации и как установить опции. Нелинейными решателями, которые мы используем в этом примере, является fminunc и fmincon.

Все принципы, обрисованные в общих чертах в этом примере, применяются к другим нелинейным решателям, таким как fgoalattain, fminimax, lsqnonlin, lsqcurvefit и fsolve.

Пример запускается с минимизации целевой функции, затем продолжает минимизировать ту же функцию с дополнительными параметрами. После этого пример показывает, как минимизировать целевую функцию, когда существует ограничение, и наконец показывает, как получить более эффективное и/или точное решение путем обеспечения градиентов или Гессиана, или путем изменения некоторых опций.

Неограниченный пример оптимизации

Рассмотрите проблему нахождения минимума функции:

xexp(-(x2+y2))+(x2+y2)/20.

Постройте функцию, чтобы понять то, где она минимизирована

f = @(x,y) x.*exp(-x.^2-y.^2)+(x.^2+y.^2)/20;
fsurf(f,[-2,2],'ShowContours','on')

График показывает, что минимум около точки (-1/2,0).

Обычно вы задаете целевую функцию как файл MATLAB. На данный момент эта функция достаточно проста задать как анонимная функция:

fun = @(x) f(x(1),x(2));

Возьмите предположение в решении:

x0 = [-.5; 0];

Установите опции оптимизации не использовать крупномасштабный алгоритм fminunc по умолчанию, поскольку тот алгоритм требует, чтобы градиент целевой функции был обеспечен:

options = optimoptions('fminunc','Algorithm','quasi-newton');

Просмотрите итерации, когда решатель вычисляет:

options.Display = 'iter';

Вызовите fminunc, неограниченный нелинейный минимизатор:

[x, fval, exitflag, output] = fminunc(fun,x0,options);
                                                        First-order 
 Iteration  Func-count       f(x)        Step-size       optimality
     0           3          -0.3769                         0.339
     1           6        -0.379694              1          0.286  
     2           9        -0.405023              1         0.0284  
     3          12        -0.405233              1        0.00386  
     4          15        -0.405237              1       3.17e-05  
     5          18        -0.405237              1       3.35e-08  

Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.

Решатель нашел решение в:

uncx = x
uncx = 2×1

   -0.6691
    0.0000

Значение функции в решении:

uncf = fval
uncf = -0.4052

Мы будем использовать количество функциональных оценок как мера эффективности в этом примере. Общее количество функциональных оценок:

output.funcCount
ans = 18

Неограниченная оптимизация с дополнительными параметрами

Мы теперь передадим дополнительные параметры в качестве дополнительных аргументов к целевой функции. Мы показываем два различных способа сделать это - использование файла MATLAB или использование вложенной функции.

Рассмотрите целевую функцию от предыдущего раздела:

f(x,y)=xexp(-(x2+y2))+(x2+y2)/20.

Мы параметризуем функцию с (a, b, c) следующим образом:

f(x,y,a,b,c)=(x-a)exp(-((x-a)2+(y-b)2))+((x-a)2+(y-b)2)/c.

Эта функция является переключенной и масштабированной версией исходной целевой функции.

Метод 1: Функция файла MATLAB

Предположим, что у нас есть целевая функция файла MATLAB под названием bowlpeakfun, заданный как:

type bowlpeakfun
function y = bowlpeakfun(x, a, b, c)
%BOWLPEAKFUN Objective function for parameter passing in TUTDEMO.

%   Copyright 2008 The MathWorks, Inc.

y = (x(1)-a).*exp(-((x(1)-a).^2+(x(2)-b).^2))+((x(1)-a).^2+(x(2)-b).^2)/c;

Задайте параметры:

a = 2;
b = 3;
c = 10;

Создайте указатель анонимной функции на файл MATLAB:

f = @(x)bowlpeakfun(x,a,b,c)
f = function_handle with value:
    @(x)bowlpeakfun(x,a,b,c)

Вызовите fminunc, чтобы найти минимум:

x0 = [-.5; 0];
options = optimoptions('fminunc','Algorithm','quasi-newton');
[x, fval] = fminunc(f,x0,options)
Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.
x = 2×1

    1.3639
    3.0000

fval = -0.3840

Метод 2: вложенная функция

Рассмотрите следующую функцию, которая реализует цель как вложенную функцию

type nestedbowlpeak
function [x,fval] =  nestedbowlpeak(a,b,c,x0,options)
%NESTEDBOWLPEAK Nested function for parameter passing in TUTDEMO.

%   Copyright 2008 The MathWorks, Inc.

[x,fval] = fminunc(@nestedfun,x0,options);  
    function y = nestedfun(x)
      y = (x(1)-a).*exp(-((x(1)-a).^2+(x(2)-b).^2))+((x(1)-a).^2+(x(2)-b).^2)/c;    
    end
end

В этом методе параметры (a, b, c) видимы к вложенной целевой функции под названием nestedfun. Внешняя функция, nestedbowlpeak, вызывает fminunc и передает целевую функцию, nestedfun.

Задайте параметры, исходное предположение и опции:

a = 2;
b = 3;
c = 10;
x0 = [-.5; 0];
options = optimoptions('fminunc','Algorithm','quasi-newton');

Запустите оптимизацию:

[x,fval] =  nestedbowlpeak(a,b,c,x0,options)
Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.
x = 2×1

    1.3639
    3.0000

fval = -0.3840

Вы видите, что оба метода произвели идентичные ответы, таким образом используйте, какой бы ни один вы находите самыми удобными.

Ограниченный пример оптимизации: Неравенства

Рассмотрите вышеупомянутую проблему с ограничением:

минимизировать xexp(-(x2+y2))+(x2+y2)/20,

подвергающийся  xy/2+(x+2)2+(y-2)2/22.

Ограничительный набор является внутренней частью наклоненного эллипса. Посмотрите на контуры целевой функции, построенной вместе с наклоненным эллипсом

f = @(x,y) x.*exp(-x.^2-y.^2)+(x.^2+y.^2)/20;
g = @(x,y) x.*y/2+(x+2).^2+(y-2).^2/2-2;
fimplicit(g)
axis([-6 0 -1 7])
hold on
fcontour(f)
plot(-.9727,.4685,'ro');
legend('constraint','f contours','minimum');
hold off

График показывает, что самое низкое значение целевой функции в эллипсе происходит около нижней правой части эллипса. Мы собираемся вычислить минимум, который был только построен. Возьмите предположение в решении:

x0 = [-2 1];

Установите опции оптимизации: используйте алгоритм внутренней точки и включите отображение результатов в каждой итерации:

options = optimoptions('fmincon','Algorithm','interior-point','Display','iter');

Решатели требуют, чтобы нелинейные ограничительные функции дали два выходных параметров: один для нелинейных неравенств, второго для нелинейных равенств. Таким образом, мы пишем ограничение с помощью функции deal, чтобы дать оба выходных параметров:

gfun = @(x) deal(g(x(1),x(2)),[]);

Вызовите нелинейный ограниченный решатель. Нет никаких линейных равенств или неравенств или границ, так передача [] для тех аргументов:

[x,fval,exitflag,output] = fmincon(fun,x0,[],[],[],[],[],[],gfun,options);
                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       3    2.365241e-01    0.000e+00    1.972e-01
    1       6    1.748504e-01    0.000e+00    1.734e-01    2.260e-01
    2      10   -1.570560e-01    0.000e+00    2.608e-01    9.347e-01
    3      14   -6.629160e-02    0.000e+00    1.241e-01    3.103e-01
    4      17   -1.584082e-01    0.000e+00    7.934e-02    1.826e-01
    5      20   -2.349124e-01    0.000e+00    1.912e-02    1.571e-01
    6      23   -2.255299e-01    0.000e+00    1.955e-02    1.993e-02
    7      26   -2.444225e-01    0.000e+00    4.293e-03    3.821e-02
    8      29   -2.446931e-01    0.000e+00    8.100e-04    4.035e-03
    9      32   -2.446933e-01    0.000e+00    1.999e-04    8.126e-04
   10      35   -2.448531e-01    0.000e+00    4.004e-05    3.289e-04
   11      38   -2.448927e-01    0.000e+00    4.036e-07    8.156e-05

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.

Решение этой проблемы было найдено в:

x
x = 1×2

   -0.9727    0.4686

Значение функции в решении:

fval
fval = -0.2449

Общее количество функциональных оценок было:

Fevals = output.funcCount
Fevals = 38

Ограничение неравенства удовлетворено в решении.

[c, ceq] = gfun(x)
c = -2.4608e-06
ceq =

     []

С тех пор c (x) близко к 0, ограничение "активно", означая, что ограничение влияет на решение. Вспомните, что неограниченное решение было найдено в

uncx
uncx = 2×1

   -0.6691
    0.0000

и неограниченная целевая функция, как находили, была

uncf
uncf = -0.4052

Ограничение переместило решение и увеличило цель

fval-uncf
ans = 0.1603

Ограниченный пример оптимизации: предоставленные пользователями градиенты

Задачи оптимизации могут быть решены более эффективно и точно если градиенты предоставляются пользователем. Этот пример показывает, как это может быть выполнено. Мы снова решаем ограниченную неравенством проблему

минимизировать xexp(-(x2+y2))+(x2+y2)/20,

подвергающийся  xy/2+(x+2)2+(y-2)2/22.

Чтобы обеспечить градиент f (x) к fmincon, мы пишем целевую функцию в форме файла MATLAB:

type onehump
function [f,gf] = onehump(x)
% ONEHUMP Helper function for Tutorial for the Optimization Toolbox demo

%   Copyright 2008-2009 The MathWorks, Inc.

r = x(1)^2 + x(2)^2;
s = exp(-r);
f = x(1)*s+r/20;

if nargout > 1
   gf = [(1-2*x(1)^2)*s+x(1)/10;
       -2*x(1)*x(2)*s+x(2)/10];
end

Ограничение и его градиент содержатся в файле MATLAB tiltellipse:

type tiltellipse
function [c,ceq,gc,gceq] = tiltellipse(x)
% TILTELLIPSE Helper function for Tutorial for the Optimization Toolbox demo

%   Copyright 2008-2009 The MathWorks, Inc.

c = x(1)*x(2)/2 + (x(1)+2)^2 + (x(2)-2)^2/2 - 2;
ceq = [];

if nargout > 2
   gc = [x(2)/2+2*(x(1)+2);
       x(1)/2+x(2)-2];
   gceq = [];
end

Выскажите предположение в решении:

x0 = [-2; 1];

Установите опции оптимизации: мы продолжаем использовать тот же алгоритм в целях сравнения.

options = optimoptions('fmincon','Algorithm','interior-point');

Мы также устанавливаем опции использовать информацию о градиенте в ограничительных функциях и цели. Примечание: эти опции MUST быть включенными или информация о градиенте будут проигнорированы.

options = optimoptions(options,...
    'SpecifyObjectiveGradient',true,...
    'SpecifyConstraintGradient',true);

На этот раз должно быть меньше функциональных количеств, поскольку fmincon не должен оценивать градиенты с помощью конечных разностей.

options.Display = 'iter';

Вызовите решатель:

[x,fval,exitflag,output] = fmincon(@onehump,x0,[],[],[],[],[],[], ...
                                   @tiltellipse,options);
                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       1    2.365241e-01    0.000e+00    1.972e-01
    1       2    1.748504e-01    0.000e+00    1.734e-01    2.260e-01
    2       4   -1.570560e-01    0.000e+00    2.608e-01    9.347e-01
    3       6   -6.629161e-02    0.000e+00    1.241e-01    3.103e-01
    4       7   -1.584082e-01    0.000e+00    7.934e-02    1.826e-01
    5       8   -2.349124e-01    0.000e+00    1.912e-02    1.571e-01
    6       9   -2.255299e-01    0.000e+00    1.955e-02    1.993e-02
    7      10   -2.444225e-01    0.000e+00    4.293e-03    3.821e-02
    8      11   -2.446931e-01    0.000e+00    8.100e-04    4.035e-03
    9      12   -2.446933e-01    0.000e+00    1.999e-04    8.126e-04
   10      13   -2.448531e-01    0.000e+00    4.004e-05    3.289e-04
   11      14   -2.448927e-01    0.000e+00    4.036e-07    8.156e-05

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.

fmincon оценил градиенты хорошо в предыдущем примере, таким образом, итерации в текущем примере подобны.

Решение этой проблемы было найдено в:

xold = x
xold = 2×1

   -0.9727
    0.4686

Значение функции в решении:

minfval = fval
minfval = -0.2449

Общее количество функциональных оценок было:

Fgradevals = output.funcCount
Fgradevals = 14

Сравните это с количеством функциональных оценок без градиентов:

Fevals
Fevals = 38

Изменение допусков завершения по умолчанию

На этот раз мы решаем ту же ограниченную проблему

минимизировать xexp(-(x2+y2))+(x2+y2)/20,

подвергающийся  xy/2+(x+2)2+(y-2)2/22,

более точно путем переопределения критериев завершения по умолчанию (опции. StepTolerance и опции. OptimalityTolerance). Мы продолжаем использовать градиенты. Значения по умолчанию для алгоритма внутренней точки fmincon являются опциями. StepTolerance = 1e-10, опции. OptimalityTolerance = 1e-6.

Замените два критерия завершения по умолчанию: погрешности на x для завершения и fval.

options = optimoptions(options,...
    'StepTolerance',1e-15,...
    'OptimalityTolerance',1e-8);

Вызовите решатель:

[x,fval,exitflag,output] = fmincon(@onehump,x0,[],[],[],[],[],[], ...
                                   @tiltellipse,options);
                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       1    2.365241e-01    0.000e+00    1.972e-01
    1       2    1.748504e-01    0.000e+00    1.734e-01    2.260e-01
    2       4   -1.570560e-01    0.000e+00    2.608e-01    9.347e-01
    3       6   -6.629161e-02    0.000e+00    1.241e-01    3.103e-01
    4       7   -1.584082e-01    0.000e+00    7.934e-02    1.826e-01
    5       8   -2.349124e-01    0.000e+00    1.912e-02    1.571e-01
    6       9   -2.255299e-01    0.000e+00    1.955e-02    1.993e-02
    7      10   -2.444225e-01    0.000e+00    4.293e-03    3.821e-02
    8      11   -2.446931e-01    0.000e+00    8.100e-04    4.035e-03
    9      12   -2.446933e-01    0.000e+00    1.999e-04    8.126e-04
   10      13   -2.448531e-01    0.000e+00    4.004e-05    3.289e-04
   11      14   -2.448927e-01    0.000e+00    4.036e-07    8.156e-05
   12      15   -2.448931e-01    0.000e+00    4.000e-09    8.230e-07

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.

Мы теперь принимаем решение видеть больше десятичных чисел в решении в порядке видеть более точно значение, которое имеют новые допуски.

format long

Оптимизатор нашел решение в:

x
x = 2×1

  -0.972742227363546
   0.468569289098342

Сравните это с предыдущим значением:

xold
xold = 2×1

  -0.972742694488360
   0.468569966693330

Изменение

x - xold
ans = 2×1
10-6 ×

   0.467124813385844
  -0.677594988729435

Значение функции в решении:

fval
fval = 
  -0.244893137879894

Решение, улучшенное

fval - minfval
ans = 
    -3.996450220755676e-07

(это отрицательно, поскольку новое решение меньше),

Общее количество функциональных оценок было:

output.funcCount
ans = 
    15

Сравните это с количеством функциональных оценок, когда проблема будет решена с обеспеченными пользователями градиентами, но с допусками по умолчанию:

Fgradevals
Fgradevals = 
    14

Ограниченный пример оптимизации с предоставленным пользователями гессианом

Если вы даете не только градиент, но также и Гессиан, решатели еще более точны и эффективны.

Решатель внутренней точки fmincon берет матрицу Гессиана в качестве отдельной функции (не часть целевой функции). Функция Гессиана H (x, lambda) должна оценить Гессиан функции Лагранжа; см. Руководство пользователя для определения этого термина.

Решатели вычисляют значения lambda.ineqnonlin и lambda.eqlin; ваша функция Гессиана говорит решатели, как использовать эти значения.

В этой проблеме мы имеем, но одно ограничение неравенства, таким образом, Гессиан:

type hessfordemo
function H = hessfordemo(x,lambda)
% HESSFORDEMO Helper function for Tutorial for the Optimization Toolbox demo

%   Copyright 2008-2009 The MathWorks, Inc.

s = exp(-(x(1)^2+x(2)^2));
H = [2*x(1)*(2*x(1)^2-3)*s+1/10, 2*x(2)*(2*x(1)^2-1)*s;
        2*x(2)*(2*x(1)^2-1)*s, 2*x(1)*(2*x(2)^2-1)*s+1/10];
hessc = [2,1/2;1/2,1];
H = H + lambda.ineqnonlin(1)*hessc;

В порядке использовать Гессиан, необходимо установить опции соответственно:

options = optimoptions('fmincon',...
    'Algorithm','interior-point',...
    'SpecifyConstraintGradient',true,...
    'SpecifyObjectiveGradient',true,...
    'HessianFcn',@hessfordemo);

Допуски были задержаны к значениям по умолчанию. На этот раз должно быть меньше функциональных количеств.

options.Display = 'iter';

Вызовите решатель:

[x,fval,exitflag,output] = fmincon(@onehump,x0,[],[],[],[],[],[], ...
                                   @tiltellipse,options);
                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       1    2.365241e-01    0.000e+00    1.972e-01
    1       3    5.821325e-02    0.000e+00    1.443e-01    8.728e-01
    2       5   -1.218829e-01    0.000e+00    1.007e-01    4.927e-01
    3       6   -1.421167e-01    0.000e+00    8.486e-02    5.165e-02
    4       7   -2.261916e-01    0.000e+00    1.989e-02    1.667e-01
    5       8   -2.433609e-01    0.000e+00    1.537e-03    3.486e-02
    6       9   -2.446875e-01    0.000e+00    2.057e-04    2.727e-03
    7      10   -2.448911e-01    0.000e+00    2.068e-06    4.191e-04
    8      11   -2.448931e-01    0.000e+00    2.001e-08    4.218e-06

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.

Были меньше и различные итерации на этот раз.

Решение этой проблемы было найдено в:

x
x = 2×1

  -0.972742246093537
   0.468569316215571

Значение функции в решении:

fval
fval = 
  -0.244893121872758

Общее количество функциональных оценок было:

output.funcCount
ans = 
    11

Сравните это с номером с помощью только оценки градиента с теми же допусками по умолчанию:

Fgradevals
Fgradevals = 
    14

Похожие темы

Для просмотра документации необходимо авторизоваться на сайте