В этом примере показано, как неоднократно решать уравнение, когда параметр изменяется путем запуска последующих решений с предыдущей начальной точки. Часто, этот процесс приводит к эффективным решениям. Однако решение может иногда исчезать, требуя запуска от новой точки или точек.
Параметризованное уравнение, чтобы решить
,
где числовой параметр, который идет от 0 до 5. В , одно решение этого уравнения когда не является слишком большим в абсолютном значении, уравнение имеет три решения. Чтобы визуализировать уравнение, создайте левую сторону уравнения как анонимная функция. Постройте функцию.
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)')
Когда является слишком большим или слишком маленьким, существует только одно решение.
Использовать fun
как целевая функция в подходе, основанном на проблеме, преобразуйте fun
к выражению оптимизации в переменной x
оптимизации.
x = optimvar('x');
expr = fcn2optimexpr(fun,x);
Запуск с начального решения в , найдите решения для 100 значений от 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')
Скачок в решении происходит рядом . В той же точке, количестве функциональных оценок, чтобы достигнуть решения увеличения от приблизительно 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
Как увеличения, сначала решения перемещаются налево. Однако, когда выше 2.5, около предыдущего решения больше нет решения. fzero
требует, чтобы дополнительные функциональные оценки искали решение и находит решение около x = 3
. После этого значения решения перемещаются медленно направо как увеличения далее. Решатель требует только приблизительно 10 функциональных оценок для каждого последующего решения.
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 функциональных оценок на итерацию. Снова, когда решатель не находит решения около предыдущего значения, решатель требует, чтобы намного больше функции evalutions искало решение. На этот раз поиск неудачен. Последующие итерации требуют немногих функциональных оценок по большей части, но не удаются найти решение. Error of Solution
постройте показывает значение функции .
Чтобы попытаться преодолеть локальный минимум не быть решением, ищите снова от различной стартовой точки когда 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
восстанавливается с плохой начальной точки рядом и получает решение, похожее на то, полученное fzero
. Количество функциональных оценок для каждой итерации обычно равняется 8, увеличиваясь до приблизительно 30 в точке, где решение переходит.