Параметрированное скалярное уравнение

Параметрированное уравнение, чтобы решить

$$\sinh (x)-3x=a$$,

где $$a$$ числовой параметр, который идет от 0 до 5. В $$a=0$$, одно решение этого уравнения $$x=0$$когда $$a$$ не является слишком большим в абсолютном значении, уравнение имеет три решения. Чтобы визуализировать уравнение, создайте левую сторону уравнения как анонимная функция. Постройте функцию.

fun = @(x)sinh(x) - 3*x;
t = linspace(-3.5,3.5);
plot(t,fun(t),t,zeros(size(t)),'k-')
xlabel('x')
ylabel('fun(x)')

Когда $$a$$ является слишком большим или слишком маленьким, существует только одно решение.

Основанный на проблеме Setup

Чтобы создать целевую функцию для подхода, основанного на проблеме, создайте выражение оптимизации expr в переменной x оптимизации.

x = optimvar('x');
expr = sinh(x) - 3*x;

Создайте и постройте решения

Запуск с начального решения $$x=0$$ в $$a=0$$, найдите решения для 100 значений $$a$$ от 0 до 5. Поскольку fun скалярная нелинейная функция, solve вызовы fzero как решатель.

Настройте объект задачи, опции и структуры данных для содержания статистики решения.

prob = eqnproblem;
options = optimset('Display','off');
sols = zeros(100,1);
fevals = sols;
as = linspace(0,5);

Решите уравнение в цикле, начинающем с первого решения sols(1) = 0.

for i = 2:length(as)
    x0.x = sols(i-1); % Start from previous solution
    prob.Equations = expr == as(i);
    [sol,~,~,output] = solve(prob,x0,'Options',options);
    sols(i) = sol.x;
    fevals(i) = output.funcCount;
end

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

subplot(2,1,1)
plot(as,sols,'ko')
xlabel 'a'
ylabel('Solution(x)')
subplot(2,1,2)
plot(fevals,'k*')
xlabel('Iteration Number')
ylabel('Fevals')

Скачок в решении происходит рядом $$a=2.5$$. В той же точке, количестве вычислений функции, чтобы достигнуть решения увеличения от приблизительно 15 до приблизительно 40. Чтобы понять почему, исследуйте более подробный график функции. Постройте функцию и каждую седьмую точку решения.

figure
t = linspace(-3.5,3.5);
plot(t,fun(t));
hold on
plot([-3.5,min(sols)],[2.5,2.5],'k--')
legend('fun','Maximum a','Location','north','autoupdate','off')
for a0 = 7:7:100
    plot(sols(a0),as(a0),'ro')
    if mod(a0,2) == 1
        text(sols(a0) + 0.15,as(a0) + 0.15,num2str(a0/7))
    else
        text(sols(a0) - 0.3,as(a0) + 0.05,num2str(a0/7))
    end
end
plot(t,zeros(size(t)),'k-')
hold off

Как $$a$$ увеличения, сначала решения перемещаются налево. Однако, когда $$a$$ выше 2.5, около предыдущего решения больше нет решения. fzero требует, чтобы дополнительные вычисления функции искали решение и находит решение около x = 3. После этого значения решения перемещаются медленно направо как $$a$$ увеличения далее. Решатель требует только приблизительно 10 вычислений функции для каждого последующего решения.

Выберите Different Solver

fsolve решатель может быть более эффективным, чем fzero. Однако fsolve может стать всунутым локальный минимум и может не решить уравнение.

Настройте объект задачи, опции и структуры данных для содержания статистики решения.

probfsolve = eqnproblem;
sols = zeros(100,1);
fevals = sols;
infeas = sols;
asfsolve = linspace(0,5);

Решите уравнение в цикле, начинающем с первого решения sols(1) = 0.

for i = 2:length(as)
    x0.x = sols(i-1); % Start from previous solution
    probfsolve.Equations = expr == asfsolve(i);
    [sol,fval,~,output] = solve(probfsolve,x0,'Options',options,'Solver','fsolve');
    sols(i) = sol.x;
    fevals(i) = output.funcCount;
    infeas(i) = fval;
end

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

subplot(2,1,1)
plot(asfsolve,sols,'ko',asfsolve,infeas,'r-')
xlabel 'a'
legend('Solution','Error of Solution','Location','best')
subplot(2,1,2)
plot(fevals,'k*')
xlabel('Iteration Number')
ylabel('Fevals')

fsolve несколько более эффективно, чем fzero, требование приблизительно 7 или 8 вычислений функции на итерацию. Снова, когда решатель не находит решения около предыдущего значения, решатель требует, чтобы намного больше вычислений функции искало решение. На этот раз поиск неудачен. Последующие итерации требуют немногих вычислений функции по большей части, но не удаются найти решение. Error of Solution постройте показывает значение функции $$fun(x)-a$$.

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

rng default % For reproducibility
sols = zeros(100,1);
fevals = sols;
asfsolve = linspace(0,5);

Решите уравнение в цикле, начинающем с первого решения sols(1) = 0.

for i = 2:length(as)
    x0.x = sols(i-1); % Start from previous solution
    probfsolve.Equations = expr == asfsolve(i);
    [sol,~,exitflag,output] = solve(probfsolve,x0,'Options',options,'Solver','fsolve');
    
    while exitflag <= 0 % If fsolve fails to find a solution
        x0.x = 5*randn; % Try again from a new start point
        fevals(i) = fevals(i) + output.funcCount;
        [sol,~,exitflag,output] = solve(probfsolve,x0,'Options',options,'Solver','fsolve');
    end
    
    sols(i) = sol.x;
    fevals(i) = fevals(i) + output.funcCount;
end

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

subplot(2,1,1)
plot(asfsolve,sols,'ko')
xlabel 'a'
ylabel('Solution(x)')
subplot(2,1,2)
plot(fevals,'k*')
xlabel('Iteration Number')
ylabel('Fevals')

На этот раз, fsolve восстанавливается с плохой начальной точки рядом $$a=2.5$$ и получает решение, похожее на то, полученное fzero. Количество вычислений функции для каждой итерации обычно равняется 8, увеличиваясь до приблизительно 30 в точке, где решение переходит.

Преобразуйте целевую функцию Используя fcn2optimexpr

Для некоторых целевых функций или версий программного обеспечения, необходимо преобразовать нелинейные функции в выражения оптимизации при помощи fcn2optimexpr. Смотрите Поддерживаемые Операции на Переменных и выражениях Оптимизации. В данном примере преобразуйте исходный функциональный fun используемый в графическом выводе к выражению оптимизации expr:

expr = fcn2optimexpr(fun,x);

Остаток от примера является точно тем же самым после этого изменения в определении expr.