Рассмотрите проблему нахождения решения системы нелинейных уравнений, якобиан которых разрежен. Размерность проблемы в этом примере 1000. Цель состоит в том, чтобы найти x таким образом что F (x) = 0. Принимая n = 1000, нелинейные уравнения
Чтобы решить большую нелинейную систему уравнений, F (x) = 0, можно использовать доверительную область отражающий алгоритм, доступный в fsolve
, крупномасштабный алгоритм (Крупномасштабный по сравнению с Алгоритмами Средней шкалы).
function [F,J] = nlsf1(x) % Evaluate the vector function n = length(x); F = zeros(n,1); i = 2:(n-1); F(i) = (3-2*x(i)).*x(i)-x(i-1)-2*x(i+1) + 1; F(n) = (3-2*x(n)).*x(n)-x(n-1) + 1; F(1) = (3-2*x(1)).*x(1)-2*x(2) + 1; % Evaluate the Jacobian if nargout > 1 if nargout > 1 d = -4*x + 3*ones(n,1); D = sparse(1:n,1:n,d,n,n); c = -2*ones(n-1,1); C = sparse(1:n-1,2:n,c,n,n); e = -ones(n-1,1); E = sparse(2:n,1:n-1,e,n,n); J = C + D + E; end
xstart = -ones(1000,1); fun = @nlsf1; options = optimoptions(@fsolve,'Display','iter',... 'Algorithm','trust-region',... 'SpecifyObjectiveGradient',true,'PrecondBandWidth',0); [x,fval,exitflag,output] = fsolve(fun,xstart,options);
Отправная точка дана, а также имя функции. Метод по умолчанию для fsolve
является доверительным резким искривлением области, таким образом, необходимо задать 'Algorithm'
как 'trust-region'
в аргументе options
в порядке выбрать алгоритм доверительной области. Установка опции Display
к 'iter'
заставляет fsolve
отображать вывод в каждой итерации. Установка 'SpecifyObjectiveGradient'
к true
, fsolve
причин, чтобы использовать якобиевскую информацию, доступную в nlsf1.m
.
Команды отображают этот вывод:
Norm of First-order Iteration Func-count f(x) step optimality 0 1 1011 19 1 2 16.1942 7.91898 2.35 2 3 0.0228027 1.33142 0.291 3 4 0.000103359 0.0433329 0.0201 4 5 7.3792e-07 0.0022606 0.000946 5 6 4.02299e-10 0.000268381 4.12e-05 Equation solved, inaccuracy possible. fsolve stopped because the vector of function values is near zero, as measured by the value of the function tolerance. However, the last step was ineffective.
Линейная система (приблизительно) решена в каждой главной итерации с помощью предобусловленного метода сопряженных градиентов. Установка PrecondBandWidth
к 0 в options
означает, что диагональный предварительный формирователь используется. (PrecondBandWidth
задает пропускную способность матрицы перед созданием условий. Пропускная способность 0 средних значений там является только одной диагональю в матрице.)
От значений оптимальности первого порядка быстро происходит линейная сходимость. Количество итераций метода сопряженных градиентов (CG), требуемых на главную итерацию, является низким, самое большее пять для проблемы 1 000 размерностей, подразумевая, что линейные системы не очень трудно решить в этом случае (хотя больше, работа требуется как сходимость, прогрессирует).
Если вы хотите использовать трехдиагональный предварительный формирователь, т.е. матрицу перед созданием условий с тремя диагоналями (или пропускная способность одной), установите PrecondBandWidth
на значение 1
:
options = optimoptions(@fsolve,'Display','iter','SpecifyObjectiveGradient',true,... 'Algorithm','trust-region','PrecondBandWidth',1); [x,fval,exitflag,output] = fsolve(fun,xstart,options);
В этом случае вывод
Norm of First-order Iteration Func-count f(x) step optimality 0 1 1011 19 1 2 16.0839 7.92496 1.92 2 3 0.0458181 1.3279 0.579 3 4 0.000101184 0.0631898 0.0203 4 5 3.16615e-07 0.00273698 0.00079 5 6 9.72481e-10 0.00018111 5.82e-05 Equation solved, inaccuracy possible. fsolve stopped because the vector of function values is near zero, as measured by the value of the function tolerance. However, the last step was ineffective.
Обратите внимание на то, что несмотря на то, что то же количество итераций происходит, количество итераций PCG понизилось, так меньше, работа делается на итерацию. Смотрите Предобусловленный Метод сопряженных градиентов.
Установка PrecondBandWidth
к Inf
(это - значение по умолчанию) означает, что решатель использует факторизацию Холесского, а не PCG.