college/ws2012/CE/uebungen/P3/P3Matlab/quasi_newton.m
Michael Scholz 0fe686245c CE update
2013-02-20 16:20:56 +01:00

49 lines
1.1 KiB
Matlab

function [N,x] = quasi_newton(x0)
% Function quasi_newton
% In: Start value
% Out: number of iteration steps, vector of solution
N=1;
xCur = x0;
xOld = x0;
epsilon = 10^(-8);
A = [3,4;1,6];
function [y] = f(x)
yTemp = (A - x(3)*eye(2))*[x(1);x(2)];
lambda = x(1)^2+x(2)^2 - 1;
y = [ yTemp ; lambda];
end
function [y] = jacobi(x)
y1 = [A(1,1)-x(3) , A(1,2) , -x(1)];
y2 = [A(2,1) , A(2,2)-x(3) , -x(2)];
y3 = [2*x(1) , 2*x(2) , 0];
y = [y1;y2;y3];
end
function [y] = newton_next(x)
y = x + getDeltaX(x);
end
function [y] = getDeltaX(x)
y = jac\(-f(x));
end
jac = jacobi(x0);
xCur = newton_next(xOld);
x(:,N) = xCur;
while(abs(norm(xCur)-norm(xOld)) > epsilon && N < 60)
N = N + 1;
xOld = xCur;
% anpassen der angenäherten jacobi matrix
% ---------------------------------------
% Schritt 1 - delta X ausrechnen
delta = getDeltaX(xOld);
% Schritt 2 - x^k ausrechnen
xCur = newton_next(xOld);
% Schritt 3 & 4 - f(x^k) = f(xCur) und anpassen der Jacobi Matrix
jac = jac + (1/norm(delta)^2)*(f(xCur) - f(xOld) - jac*delta) * delta';
x(:,N) = xCur;
end
end