Матлаб как найти максимум функции

Main Content

syms x
num = 3*x^2 + 6*x -1;
denom = x^2 + x - 3;
f = num/denom

f1 = 
6 x+6x2+x-3-2 x+1 3 x2+6 x-1x2+x-32

f1 = 
-3 x2+16 x+17x2+x-32

crit_pts = 
(-133-83133-83)

fplot(f)
hold on
plot(double(crit_pts), double(subs(f,crit_pts)),'ro')
title('Maximum and Minimum of f')
text(-4.8,5.5,'Local minimum')
text(-2,4,'Local maximum')
hold off

f2 = diff(f1);
inflec_pt = solve(f2,'MaxDegree',3);
double(inflec_pt)

ans = 3×1 complex

  -5.2635 + 0.0000i
  -1.3682 - 0.8511i
  -1.3682 + 0.8511i

idx = imag(double(inflec_pt)) == 0;
inflec_pt = inflec_pt(idx)

inflec_pt = 
-139 16954-2197181/3-16954-2197181/3-83

fplot(f,[-9 6])
hold on
plot(double(inflec_pt), double(subs(f,inflec_pt)),'ro')
title('Inflection Point of f')
text(-7,1,'Inflection point')
hold off

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

Общие сведения

Итак, в этом блоке ничего про Matlab не будет, лишь информация о понятии оптимизации. Это понятие сводится к терминам минимума и максимума функции, или, если коротко — экстремумам.

Под минимумом понимают такое значение функции, которое в некоторой окрестности этой функции, принимает наименьшее значение из всех возможных значений в этой окрестности. Соответственно максимум — это наибольшее значение функции в какой-либо окрестности.

Если не понятно — вот простой пример с всеми известной параболой:

У этой функции есть один минимум, и он находится в точке x = 0. Эта точка называется точкой минимума, а само значение этой функции есть минимум (он тоже равен 0). Максимумов у этой функции нет, но если бы функцию перевернули вверх ногами, то он бы появился.

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

И теперь, когда вы хоть как то познакомились с понятиями оптимизации и нахождения экстремумов функции, можем перейти к нахождению максимумов и минимумов в Matlab. Далее, как обычно мы будем разбирать примеры различной сложности. Часть примеров будет содержать в себе стандартные команды Matlab для нахождения минимума и максимума функции, а другая часть — это реализация метода с нуля.

Стандартные методы Matlab

Разберем 2 задачи нахождения минимума в Matlab:

1 пример. Вычислить минимум функции f(x)=-x^1/x, определив графически интервал его локализации. Вычисления провести с минимальным шагом по аргументу 1*10^-5

Для начала создадим скрипт, который отобразит эту функцию. Вот код для этого:

x = 0.00001:0.00001:10;
y = -x.^(x.^(-1));
plot (x,y);
hold on;
grid on;

Запускаем скрипт и получаем:
По графику функции делаем вывод, что имеется один минимум, и его координаты находятся в интервале 2.5 — 3, то есть мы сократим интервал поиска минимума.

Теперь создадим еще один скрипт, дадим ему название first.m и пропишем в него функцию:

function fun=first(x)
fun = -x.^(x.^(-1));
end

Таким образом в этом m-файле мы определили функцию. Теперь в командном окне мы пропишем следующий код:

>> [x,y] = fminbnd(@first,2.5,3)

И получаем такие значения:
x = 2.7183 — координата точки минимума
y = -1.4447 — значение минимума

В этой части кода мы использовали стандартный метод Matlab для нахождения минимума функции — fminbnd. мы передаем 3 параметра — саму функцию и интервалы для поиска минимума. Стоит отметить, что этот метод подходит только для функций, зависящих от одной переменной.

Итак, для этой задачи мы создали 2 скрипт-файла, которые вы можете скачать в конце статьи.

2 пример. Вычислить минимум функции двух переменных x⁴+y⁴-2x²+4xy-2y²+1 с точность 1*10-5.
Координаты начальной точки поиска [1.0,-1.0].

Для начала построим график функции от двух переменных — для этого создадим новый скрипт и пропишем там этот код:

[x y] = meshgrid(-2:0.1:2, -2:0.1:2);
z = x.^4 + y.^4 - 2*x.^2 + 4 * x.*y - 2*y.^2 + 1;
surf(x,y,z);

Функция surf позволяет строить трехмерные графики и отображать глубину значений функции для лучшего понимания. Запускаем скрипт — в итоге получился такой график:
Как видно из графика, имеется два участка, где присутствует локальный минимум (темно-синие участки), и наша задача найти координаты и значения двух этих точек. Воспользуемся стандартными инструментами Matlab и создадим новый скрипт с именем second.m, в котором и пропишем код:

function fun=second(x)
fun = x(1)^4 + x(2)^4 - 2*x(1)^2 + 4 * x(1)*x(2) - 2*x(2)^2 + 1;
end

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

>> [z,f,exitflag,output] = fminsearch(@second, [1.0,-1.0], optimset('TolX',1e-5))

Получаем такой вывод:

z = 1.4142   -1.4142
f = -7.0000
exitflag = 1
output = 
    iterations: 40
    funcCount: 74
    algorithm: 'Nelder-Mead simplex direct search'
    message: [1x196 char]

Для нахождение минимумов в Matlab на этот раз мы использовали функцию fminsearch. Эта функция реализует симплекс — метод Нелдера-Мида. В выводе мы получили несколько переменных: в z записались значения координат точек минимума, в f само значение этого минимума. А в переменных exitflag и output помещены условия прерывания процесса поиска и информация об оптимизации соответственно.

В итоге у нас опять получилось 2 m-файла.

Метод Ньютона Matlab

А теперь попробуем сами реализовать метод Ньютона для оптимизации функции.

3 пример. Методом Ньютона найти точку минимума x* и минимальное значение f* функции f(x)=(x-2)⁴-lnx на отрезке xє[2;3] c точностью 10^-7

Начнем с того, что создадим новый скрипт и назовем его Newton.m. Затем пропишем в нем код:

function [Xk, Yk] = Newton(f,diap)
    a = diap(1); % границы
    b = diap(2);
    df = char(diff(sym(f))); % символьно ищем первую 
    ddf = char(diff(sym(df))); % и вторую производные
    F = inline(f); % преобразуем в функции
    F1 = inline(df);
    F2 = inline(ddf);
    eps = 0.0000001; % задаем точность
    if F(a)*F2(a) > 0 % проверка с какой границы начинать искать
        Xk = b;
    else
        Xk = a;
    end
    while abs(F1(Xk)) > eps
        X0 = Xk; % X0 - значение предыдущего шага
        Xk =  X0 - (F1(X0)/(F2(X0))); % расчет нового значения
        Yk = F(Xk);
    end
end

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

>> fun = '(x-2)^4 - log(x)';
>> diap = [2,3];
>> [Xk, Yk] = Newton(fun, diap);

В итоге получилось:
Xk = 2.4663
Yk = -0.8554

Не будем приводить график, но вы сами можете проверить, что значения найдены правильно. Важно сказать, что этот метод позволяет найти только локальный экстремум, и если на выбранном диапазоне есть несколько экстремумов, то метод может найти не тот, который нужен вам.

Также, очень важно задавать как можно узкий диапазон поиска, иначе метод может работать некорректно, особенно это проявляется с периодическими функциями по типу cos(x) и т.п.

Заключение

Ну что ж, в этой статье мы рассмотрели некоторые методы для нахождения экстремумов в Matlab. Мы использовали как стандартные методы, так и реализовали метод Ньютона в среде Matlab. Их исходники чуть ниже.
Скачать исходники

На этом сегодня все, оставляйте ваши комментарии и задавайте вопросы.

МАТЛАБ.
ОПТИМИЗАЦИЯ, ЛИНЕЙНОЕ ПРОГРАММИРОВАНИЕ.

Поиск
экстремума функции одной переменной.
Поиск локального (в пределах определенного
интервала) минимума осуществляет функция
[x,y]=fminbnd(имяФункции,границыИнтервала).
Если же надо найти локальный минимум,
то исследуемую функцию берут с другим
знаком и ищут ее минимум, который
фактически является ее максимумом. При
этом значение функции в максимуме это
найденное значение с другим знаком.

Поиск
экстремума функции нескольких переменных.
Если функция вычисляется от нескольких
переменных, то для поиска ее минимума
используется (после приблизительной
оценки начальных приближений переменных
в предполагаемом минимуме) команда
[x,y]=fminsearch(имяФункции,
векторНачальныхПриближенийПеременных).

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

Пусть
функция L(x₁,x₂,x₃,x₄)=x₁+x₂+x₃+x₄
и надо найти ее максимум и соответствующие
ее максимуму значения переменных.

Пусть
даны ограничения в виде неравенств.
x₁≥0,
x₂≥0,
x₃≥0,
x₄≥0.

Пусть
также даны еще ограничения в виде
неравенств.

3x₁-x₂≤7,
x₂-2x₃≤-1,
4x₃-x₄≤3,
5x₁+2x₄≥14.

Для
решения задачи линейного программирования
(нахождения минимума) используется
функция [x,
L,
f]=linprog(c,
A,
b,
A1,
b1,Lx,
Rx)
где

x
вектор значений переменных, полученный
в качестве ответа; L
значение функции в минимуме;

f
параметр, характеризующий вычислительный
процесс (если он ноль то решение
приостановлено после достижения
максимального числа итераций, если
положителен то все нормально решено,
если отрицателен то решения не найдено);

c
функция цели представленная в виде
вектора коэффициентов (в нашем случае
[1 1 1 -1] но так как нам нужен максимум, а
функция [x,
L,
f]=linprog(c,
A,
b,
A1,
b1,
Lx,
Rx)
ищет минимум, то в выражении для функции
поменяем знак, поэтому вектор коэффициентов
будет [-1 -1 -1 1]);

A,
b
система ограничений, заданная в матричном
виде Ax≤b
(это в нашей задаче ограничения 3x₁-x₂≤7,
x₂-2x₃≤-1,
4x₃-x₄≤3,
5x₁+2x₄≥14,
но так как 5x₁+2x₄≥14
не подходит, то надо поменять знак в
этом неравенстве и тогда оно будет
-5x₁-2x₄≤-14,
в таком случае матрица А состоит из
коэффициентов (при переменных) в этих
неравенствах, а столбец b
состоит из правых (не содержащих
переменных) частей неравенств;

A1,
b1
система равенств вида A1x=b
(в нашей задаче такой системы ограничений
нету, но могла бы быть);

Lx,
Rx
относятся к ограничениям в виде Lx≤x≤Rx,
Lx≤x,
x≤Rx
(в нашей задаче есть ограничения вида
Lx≤x,
поэтому вектор Lx
будет 0 0 0 0).

При
использовании функции linprog
в списке аргументов вместо тех, которые
не указаны, ставятся пустые квадратные
скобки.

Рассмотрим
еще пример для тренировки. Пусть дана
функция W=x₁+x₂+3x₃-x₄
и надо найти ее максимум. Пусть даны
ограничения в виде неравенств.

x₁-5x₂+4x₃≤5,
x₁-2x₂-3x₃≤4,
x₁+6x₂+5x₃≤4,
x₂+x₃≤1.
И пусть есть еще ограничения в виде
неравенств. x₁≥0,
x₂≥0,
x₃≥0,
x₄≥0.

Решим
задачу с помощью функции [x,
L,
f]=linprog(c,
A,
b,
A1,
b1,
Lx,
Rx).
Для ее применения нам надо подготовить
аргументы функции. Вектор с коэффициентов
функции равен (1 1 3 -1), но так как функция
linprog
ищет минимум, а нам нужен максимум, то
поменяем знак у коэффициентов функции
и тогда вектор с будет (-1 -1 -3 1). Матрица
А состоит из коэффициентов при переменных
в системе неравенств x₁-5x₂+4x₃≤5,
x₁-2x₂-3x₃≤4,
x₁+6x₂+5x₃≤4,
x₂+x₃≤1.
Так как все они со знаком ≤ то ничего
менять не требуется. Вектор b
состоит из правых частей этих же
неравенств. Матрица A1
из системы A1х=b1
для нас не актуальна (нет таких условий).
Lx,
Rx
относятся к ограничениям в виде Lx≤x≤Rx,
Lx≤x,
x≤Rx.
У нас есть ограничения в виде неравенств.
x₁≥0,
x₂≥0,
x₃≥0,
x₄≥0.
Тогда вектор Lx=[0;
0; 0; 0]. Что касается остальных типов
ограничений, то так как у нас нет таких
ограничений, то для нас они не важны.

Задачи
нелинейного программирования.
Нелинейная задача отличается от линейной.
В ней есть система линейных неравенств
Ax≤b,
система линейных равенств A1x=b1,
система ограничений вида lx≤x≤rx,
lx≤x,
x≤rx.
Однако помимо этих уже знакомых (из
линейного программирования) условий,
есть еще система нелинейных ограничений
вида g1(x)≤0,
и система нелинейных равенств g2(x)=0.
Для решения задачи нелинейного
программирования (после выбора вектора
начальных приближений х0) используется
функция [x,
y,
f]=fmincon(F,
x0,
A,
b,
A1,
b1,
Lx,
Rx,
G),
где G
это М-функция, вычисляющая левые части
нелинейных ограничений g1(x)≤0,
g2(x)=0.
Если какие-то аргументы не определены,
то вместо них в вызове функции ставятся
квадратные скобки.

Рассмотрим
для примера следующую задачу. пусть
дана функция x₁²+x₂²=1
при условии, что x₁x₂=4,
x₁≥0,
x₂≥0.
Вектор начальных приближений примем
равным х0=(0.1, 0.1).