Найти безусловный минимум функции одной переменной

clear, clc
close all
 
% задаем функцию и производную
f = @(x) 2*x.^2 + 16./x; % создаем функцию
%f = @(x) 2*(x-3).^2;
h = 1e-6; % приращение (нужно для расчета производной) f' = (f(x+h)-f(x-h))/(2h)
df = @(x) (f(x+h)-f(x-h))/(2*h); % производная (для любой функции)
 
% Шаг 1 (задаем параметры)
x0 = 1;
delta = 1;
eps1 = 0.01;
eps2 = 0.03;
 
x1 = x0;
df1 = df(x1); % Шаг 2
 
% определяем начальный интервал поиска
k = 0;
flag = 1; % если 0, то минимум найден
if df1 < 0 % Шаг 3 - проверяем знак производной
    x2 = x1 + 2^k*delta;
    while df(x1)*df(x2) > 0
        k = k+1;
        x1 = x2;
        x2 = x1 + 2^k*delta;
    end
elseif df1 > 0
    x2 = x1 - 2^k*delta;
    while df(x1)*df(x2) > 0
        k = k+1;
        x1 = x2;
        x2 = x1 - 2^k*delta;
    end
else % df = 0 - попали в точку минимума!
    if df(x1) == 0
        xm = x1;
    else
        xm = x2;
    end
    flag = 0; % минимум нашли
end
 
while flag  % делаем бесконечный цикл, так как мы не знаем, когда найдем минимум
    % Шаг 4
    % интервал поиска
    %[x1, x2]
    f1 = f(x1);
    f2 = f(x2);
    df1 = df(x1);
    df2 = df(x2);
    % Шаг 5
    % считаем параметры
    z = 3*(f1-f2)/(x2-x1) + df1 + df2;
    if x1 < x2
        w = sqrt(z^2 - df1*df2);
    else
        w = -sqrt(z^2 - df1*df2);
    end
    mu = (df2 + w - z)/(df2 - df1 + 2*w);
    % находим минимум куб. сплайна:
    if mu < 0
        xm = x2;
    elseif (mu>=0)&(mu<=1)
        xm = x2 - mu*(x2-x1);
    else
        xm = x1;
    end
    % Шаг 6 - проверяем условие убывания
    while f(xm) > f(x1)
        xm = xm - 1/2*(xm - x1);
    end
    % Шаг 7 - проверяем условия окончания
    if (abs(df(xm)) < eps1) & ( abs((xm-x1)/x1) < eps2 ) % если оба условия выполнены
        flag = 0; % опускаем флаг, чтобы выйти из цикла
    else % если что-то не сошлось - пересчитываем по новой
        if df(xm)*df(x1)<0
            x2 = x1;
            x1 = xm;
        elseif df(xm)*df(x2)<0
            x2 = x2;
            x1 = xm;
        elseif df(x1) == 0
            xm = x1; flag = 0; % минимум нашли
        else
            xm = x2; flag = 0; % минимум нашли
        end
    end
    x_min = xm; % сохраняем результат
end
 
x_min % результат
 
%точное решение
Xmin = fminsearch(f, x0) % ищем минимум начиная от точки х=0
err = abs(x_min - Xmin) % разница между найденным и точным решением
 
% рисуем
fplot(f,[x_min-delta, x_min+delta]) % рисуем функцию на интервале (-1,5)
hold on
plot(x_min,f(x_min),'*r') % рисуем точку минимума
grid on