В этом примере показано, как предоставить якобиевский шаблон разреженности, чтобы решить большую разреженную систему уравнений эффективно. Пример использует целевую функцию, заданную для системы n уравнений,
Уравнения, чтобы решить , . Использование в качестве примера . nlsf1a
функция помощника в конце этого примера реализует целевую функцию .
В примере Большая Разреженная Система Нелинейных уравнений с якобианом, который решает ту же систему, целевая функция, имеет явный якобиан. Однако иногда вы не можете вычислить якобиан явным образом, но можно определить, какие элементы якобиана являются ненулевыми. В этом примере якобиан трехдиагонален, потому что единственные переменные, которые появляются в определении , , и . Таким образом, единственные ненулевые части якобиана находятся на основной диагонали и ее двух соседних диагоналях. Якобиевский шаблон разреженности является матрицей, ненулевые элементы которой соответствуют (потенциально) ненулевым элементам в якобиане.
Создайте разреженное - трехдиагональная матрица из единиц, представляющих якобиевский шаблон разреженности. (Для объяснения этого кода смотрите spdiags
.)
n = 1000; e = ones(n,1); Jstr = spdiags([e e e],-1:1,n,n);
Просмотрите структуру Jstr
.
spy(Jstr)
Установите fsolve
опции, чтобы использовать 'trust-region'
алгоритм, который является единственным fsolve
алгоритм, который может использовать якобиевский шаблон разреженности.
options = optimoptions('fsolve','Algorithm','trust-region','JacobPattern',Jstr);
Установите начальную точку x0
к вектору из –1 записи.
x0 = -1*ones(n,1);
Решите задачу.
[x,fval,exitflag,output] = fsolve(@nlsf1a,x0,options);
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.
Исследуйте получившуюся норму невязки и количество вычислений функции что fsolve
берет.
disp(norm(fval))
1.0522e-09
disp(output.funcCount)
25
Норма невязки является близким нулем, указывая на тот fsolve
решает задачу правильно. Количество вычислений функции справедливо мало, чуть приблизительно две дюжины. Сравните это количество вычислений функции к номеру без предоставленного якобиевского шаблона разреженности.
options = optimoptions('fsolve','Algorithm','trust-region'); [x,fval,exitflag,output] = fsolve(@nlsf1a,x0,options);
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.
disp(norm(fval))
1.0659e-09
disp(output.funcCount)
5005
Решатель достигает по существу того же решения как прежде, но берет тысячи вычислений функции, чтобы сделать так.
Этот код создает nlsf1a
функция.
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; end