En el ejemplo, la función calcula el jacobiano, una matriz dispersa, junto con la evaluación de.Ecuaciones no lineales con analítica jacobiananlsf1
J
F
¿Qué pasa si el código para calcular el jacobiano no está disponible? De forma predeterminada, si no indica que el jacobiano se puede calcular en (estableciendo la opción en),,, y en su lugar utiliza la diferenciación finita para aproximar el jacobiano.nlsf1
'SpecifyObjectiveGradient'
Opciones
true
fsolve
lsqnonlin
lsqcurvefit
Para que esta diferenciación finita sea lo más eficiente posible, debe suministrar el patrón de la dispersión del jacobiano, estableciendo una matriz dispersa en.JacobPattern
Jstr
Opciones
Es decir, suministrar una matriz dispersa cuyas entradas de distinto cero corresponden a los nonceros del jacobiano para todos.Jstr
x De hecho, los ceros de puede corresponder a un superconjunto de las ubicaciones distinto de cero de; sin embargo, en general, el costo computacional del procedimiento disperso de diferencias finitas aumentará con el número de ceros de.Jstr
J
Jstr
Proporcionar el patrón de dispersión puede reducir drásticamente el tiempo necesario para calcular la diferencia finita en problemas grandes. Si no se proporciona el patrón de dispersión (y el jacobiano no se calcula en la función objetiva), entonces, en este problema con las variables 1000, el código de diferenciación finita intenta computar todas las entradas 1000-por-1000 en el jacobiano. Pero en este caso hay sólo 2998 nonceros, sustancialmente menos que el 1 millón posibles ceros el código de diferenciación finita intenta computar. En otras palabras, este problema es solucionable si usted proporciona el patrón de la dispersión. Si no es así, la mayoría de los equipos se queda sin memoria cuando se intenta realizar una diferenciación finita completa y densa. En la mayoría de los problemas pequeños, no es esencial proporcionar la estructura de la dispersión.
Supongamos que la matriz dispersa, calculada anteriormente, se ha guardado en el archivo.Jstr
nlsdat1.mat
Las siguientes llamadas de conductor aplicadas a, que es sin el jacobiano.fsolve
nlsf1a
nlsf1
Las diferencias finitas dispersas se utilizan para estimar la matriz jacobiana dispersa según sea necesario.
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);
En este caso, la salida mostrada es
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.
Alternativamente, es posible elegir un solucionador lineal directo disperso (es decir, una factorización de QR dispersa) indicando un preacondicionador "completo". Por ejemplo, si se establece en, se utiliza un solucionador lineal directo disperso en lugar de una iteración de degradado conjugada precondicionada: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);
y la visualización resultante es
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.
Cuando se utiliza el solucionador directo disperso, no hay iteraciones CG. Tenga en cuenta que el valor final optimalidad y () (que para, (), es la suma de los cuadrados de los valores de función) están más cerca de cero que el uso de la PCG método, que es a menudo el caso.fxfsolve
fx