function []=blatt13()

%Startwerte - z0 muss nicht im zulaessigen Bereich liegen
z0=[1;1];
mu0=1;

[fz,z,mu]=Lagrange_Newton(@f ,@h ,@grad_h ,@grad_L ,@Hess_L, z0, mu0);
g=sprintf('%d ', z(:,end));
Minimum=z(:,end)
Funktionswert=fz(end)

hold on
x=[-1:0.01:1];
y1=sqrt(1-x.^2);
y2=-sqrt(1-x.^2);
l(1)=plot(x,y1,'b');
l(2)=plot(x,y2,'b');
l(3)=plot(z(1,:),z(2,:),'r--o');
legend( l([1 3]),'zulaessiger Bereich','Iterierten','Location','SouthEast');
xlabel('x')
ylabel('y')
title('A 13.5')

end


%Lagrange-N. Verfahren
function [fz,z,mu] = Lagrange_Newton(f,h,grad_h,grad_L,Hess_L,z0,mu0)


n=size(z0);
p=size(mu0);
z(:,1) = z0;
fz(:,1) = f(z0);
mu(:,1) = mu0;

i=1;

while norm( grad_L(z(:,i),mu(:,i)) ) > 10^(-6) || norm( h(z(:,i)) ) > 10^(-6)

   d = -inv( [ Hess_L(z(:,i),mu(:,i)) grad_h(z(:,i));(grad_h(z(:,i)))' 0])*[ grad_L(z(:,i),mu(:,i));h(z(:,i))];
   z(:,i+1) = z(:,i) + d(1:n);
   mu(:,i+1) = mu(:,i) + d(n+1:n+p);
   fz(:,i+1) = f(z(:,i+1));
   i=i+1;

end



end


%Alles was benoetigt wird
function y=f(x)

y=2*(x(1)^2+x(2)^2-1)-x(1);

end


function y=h(x)

y=x(1)^2+x(2)^2-1;

end

function y=grad_f(x)

y=[4*x(1)-1;4*x(2)];

end


function y=grad_h(x)

y=[2*x(1);2*x(2)];

end


function A=Hess_f(x)

A=[4 0;0 4];

end


function A=Hess_h(x)

A=[2 0;0 2];

end


function y=grad_L(x,mu)

y = grad_f(x) + mu'*grad_h(x);

end


function A = Hess_L(x,mu)

A=Hess_f(x) + mu'*Hess_h(x);

end