Это учебное пособие содержит несколько примеров использования двух нелинейных решателей оптимизации. fminunc и fminconи как задать параметры. Принципы, описанные в этом учебном пособии, применимы к другим нелинейным решателям, таким как fgoalattain, fminimax, lsqnonlin, lsqcurvefit, и fsolve.
Примеры учебных пособий охватывают следующие задачи:
Минимизация целевой функции
Минимизация одной и той же функции с дополнительными параметрами
Минимизация целевой функции с помощью ограничения
Получение более эффективного или точного решения путем предоставления градиентов или гессена, или путем изменения опций
Рассмотрим проблему нахождения минимума функции
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 дефолт 'quasi-newton' алгоритм. Этот шаг обеспечивает одинаковую работу учебного пособия в каждой версии MATLAB.
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, а затем вложенную функцию.
Рассмотрим целевую функцию из предыдущего примера.
x2 + y2 )/20.
Параметризуйте функцию с помощью (a, b, c) следующим образом:
x-a) 2 + (y-b) 2 )/c.
Эта функция является сдвинутой и масштабированной версией исходной целевой функции.
Рассмотрим целевую функцию файла MATLAB с именем bowlpeakfun определяется следующим образом.
type bowlpeakfunfunction 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
Рассмотрим nestedbowlpeak , которая реализует цель как вложенную функцию.
type nestedbowlpeakfunction [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
Оба подхода дают одинаковые ответы, так что вы можете использовать тот, который вы находите наиболее удобным.
Рассмотрим предыдущую проблему с ограничением:
y2 )/20,
2/2≤2 xy/2 + (x + 2) 2 + (y-2).
Набор зависимостей представляет собой внутреннюю часть наклонного эллипса. Просмотр контуров целевой функции, выводимой на печать вместе с наклонным эллипсом.
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
Можно более эффективно и точно решать задачи оптимизации, предоставляя градиенты. Этот пример, как и предыдущий, решает проблему ограничения неравенства
y2 )/20,
2/2≤2 xy/2 + (x + 2) 2 + (y-2).
Для обеспечения градиента f (x) кfmincon, запишите целевую функцию в виде файла MATLAB.
type onehumpfunction [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 tiltellipsefunction [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');
Задайте параметры для использования данных градиента в функциях цели и ограничения. Примечание: эти параметры должны быть включены, или информация о градиенте будет игнорироваться.
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
В этом примере по-прежнему используются градиенты и решается та же проблема с зависимостями
y2 )/20,
2/2≤2 xy/2 + (x + 2) 2 + (y-2).
В этом случае более точное решение достигается путем переопределения критериев завершения по умолчанию (options.StepTolerance и options.OptimalityTolerance). Значения по умолчанию для fmincon алгоритм внутренней точки: options.StepTolerance = 1e-10 и options.OptimalityTolerance = 1e-6.
Переопределите эти два критерия завершения по умолчанию.
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, лямбда) оценивает гессен лагранжиана; см. Hessian для алгоритма внутренней точки fmincon.
Решатели вычисляют значения lambda.ineqnonlin и lambda.eqlin; функция Гессена сообщает решателям, как использовать эти значения.
Этот пример имеет одно ограничение неравенства, поэтому гессен определяется, как указано в hessfordemo функция.
type hessfordemofunction 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