В примере Нелинейные уравнения с Аналитическим якобианом, функциональным nlsf1 вычисляет якобиевский J, разреженная матрица, наряду с оценкой F. Что, если код, чтобы вычислить якобиан не доступен? По умолчанию, если вы не указываете, что якобиан может быть вычислен в nlsf1 (путем установки 'SpecifyObjectiveGradient' опция в options к true), fsolve, lsqnonlin, и lsqcurvefit вместо этого использует конечное дифференцирование, чтобы аппроксимировать якобиан.
Для этого конечного дифференцирования, чтобы быть максимально эффективными, необходимо предоставить шаблон разреженности якобиана установкой JacobPattern к разреженной матрице Jstr в options. Таким образом, предоставьте разреженную матрицу Jstr чьи ненулевые записи соответствуют ненулям якобиана для всего x. Действительно, ненули Jstr может соответствовать надмножеству ненулевых местоположений J; однако, в целом вычислительная стоимость разреженной процедуры конечной разности увеличится с количеством ненулей Jstr.
Обеспечение шаблона разреженности может решительно уменьшать время, должен был вычислить конечное дифференцирование на больших проблемах. Если шаблон разреженности не обеспечивается (и якобиан не вычисляется в целевой функции ни один), затем, в этой проблеме с 1 000 переменных, код конечного дифференцирования пытается вычислить все записи 1000 на 1000 в якобиане. Но в этом случае существует только 2 998 ненулей, существенно меньше, чем 1 000 000 возможных ненулей, которые код конечного дифференцирования пытается вычислить. Другими словами, эта проблема разрешима, если вы обеспечиваете шаблон разреженности. В противном случае большинство компьютеров исчерпывает память, когда полное плотное конечное дифференцирование предпринято. На самых небольших проблемах не важно обеспечить структуру разреженности.
Предположим разреженная матрица Jstr, вычисленный ранее, было сохранено в файле nlsdat1.mat. Следующий драйвер вызывает fsolve примененный nlsf1a, который является nlsf1 без якобиана. Разреженное конечное дифференцирование используется, чтобы оценить разреженную якобиевскую матрицу по мере необходимости.
function F = nlsf1a(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;
xstart = -ones(1000,1);
fun = @nlsf1a;
load nlsdat1 % Get Jstr
options = optimoptions(@fsolve,'Display','iter','JacobPattern',Jstr,...
'Algorithm','trust-region','SubproblemAlgorithm','cg');
[x,fval,exitflag,output] = fsolve(fun,xstart,options);В этом случае отображенный вывод
Norm of First-order
Iteration Func-count f(x) step optimality
0 5 1011 19
1 10 16.1942 7.91898 2.35
2 15 0.0228025 1.33142 0.291
3 20 0.00010336 0.0433327 0.0201
4 25 7.37924e-07 0.0022606 0.000946
5 30 4.02301e-10 0.000268382 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.В качестве альтернативы возможно выбрать разреженный прямой линейный решатель (т.е. разреженную QR-факторизацию) путем указания на “полный” предварительный формирователь. Например, если вы устанавливаете PrecondBandWidth к Inf, затем разреженный прямой линейный решатель используется вместо предобусловленной итерации метода сопряженных градиентов:
xstart = -ones(1000,1); fun = @nlsf1a; load nlsdat1 % Get Jstr options = optimoptions(@fsolve,'Display','iter','JacobPattern',Jstr,... 'Algorithm','trust-region','SubproblemAlgorithm','factorization'); [x,fval,exitflag,output] = fsolve(fun,xstart,options);
и получившееся отображение
Norm of First-order
Iteration Func-count f(x) step optimality
0 5 1011 19
1 10 15.9018 7.92421 1.89
2 15 0.0128161 1.32542 0.0746
3 20 1.73502e-08 0.0397923 0.000196
4 25 1.10732e-18 4.55495e-05 2.74e-09
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.При использовании разреженного прямого решателя нет никаких итераций CG. Заметьте что итоговая оптимальность и f (x) значение (который для fsolve, f (x), сумма квадратов значений функции), ближе, чтобы обнулить, чем использование метода PCG, который часто имеет место.