Решение нелинейной системы без и включая якобиан

Этот пример показывает сокращение вычислений функции, когда вы предоставляете производные для системы нелинейных уравнений. Как объяснено в Writing Vector and Matrix Objective Functions, якобиан J(x) системы уравнений F(x) является Jij(x)=Fi(x)xj. Предоставьте эту производную как второй выход вашей целевой функции.

Для примера, multirosenbrock функция является n-мерное обобщение Розенбрка функции (см. Решение ограниченной нелинейной задачи, основанной на проблеме) для любого положительного, четного значения n:

F(1)=1-x1F(2)=10(x2-x12)F(3)=1-x3F(4)=10(x4-x32)F(n-1)=1-xn-1F(n)=10(xn-xn-12).

Решение системы уравнений F(x)=0 есть точка xi=1, i=1n.

Для этой целевой функции все якобианские термины Jij(x) являются нулем, кроме членов, где i и j отличаются самое большее одним. Для нечетных значений i<n, ненулевые условия

Jii(x)=-1J(i+1)i=-20xiJ(i+1)(i+1)=10.

The multirosenbrock helper функция в конце этого примера создает целевую функцию F(x) и его якобиан J(x).

Решите систему уравнений, начиная с точки xi=-1.9 для нечетных значений i<n, и xi=2 для четных значений i. Определить n=64.

n = 64;  
x0(1:n,1) = -1.9; 
x0(2:2:n,1) = 2;
[x,F,exitflag,output,JAC] = fsolve(@multirosenbrock,x0);
Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.

Исследуйте расстояние вычисляемого решения x от истинного решения и количества вычислений функции, которые fsolve требуется вычислить решение.

disp(norm(x-ones(size(x))))
     0
disp(output.funcCount)
        1043

fsolve находит решение и принимает для этого более 1000 вычислений функции.

Снова решить систему уравнений, на этот раз используя якобиан. Для этого установите 'SpecifyObjectiveGradient' опция для true.

opts = optimoptions('fsolve','SpecifyObjectiveGradient',true);
[x2,F2,exitflag2,output2,JAC2] = fsolve(@multirosenbrock,x0,opts);
Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.

Снова исследуйте расстояние вычисляемого решения x2 от истинного решения и количества вычислений функции, которые fsolve требуется вычислить решение.

disp(norm(x2-ones(size(x2))))
     0
disp(output2.funcCount)
    21

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

Функция помощника

Этот код создает multirosenbrock вспомогательная функция.

function [F,J] = multirosenbrock(x)
% Get the problem size
n = length(x);  
if n == 0, error('Input vector, x, is empty.'); end
if mod(n,2) ~= 0
   error('Input vector, x ,must have an even number of components.');
end
% Evaluate the vector function
odds  = 1:2:n;
evens = 2:2:n;
F = zeros(n,1);
F(odds,1)  = 1-x(odds);
F(evens,1) = 10.*(x(evens)-x(odds).^2); 
% Evaluate the Jacobian matrix if nargout > 1
if nargout > 1
   c = -ones(n/2,1);    C = sparse(odds,odds,c,n,n);
   d = 10*ones(n/2,1);  D = sparse(evens,evens,d,n,n);
   e = -20.*x(odds);    E = sparse(evens,odds,e,n,n);
   J = C + D + E;
end
end

См. также

Похожие темы