lagrange_newton.m

lagrange_newton.m — Objective-C source code, 2 kB (2.315 bytes)

Dateiinhalt

%--------------------------------------------------------------------------------------
%
% Lagrange-Newton-Verfahren zur Loesung des nichtlinearen
% Optimierungsproblems
% 
%     Minimiere f(x) u.d.N. g(x) = 0
% 
% Eingabe: 
% --------
% x0          : Startvektor fuer x
% lambda0 : Startvektor fuer Lagrange-Multiplikator
% tol          : Optimalitaetstoleranz (Verfahren endet, falls 
%                 Verletzung der KKT Bedingungen <= tol)
% maxiter  : maximale Iterationszahl
%
% Ausgabe:
% --------
% x            : Loesungsvektor x
% lambda   : Lagrange-Multiplikator
% fval         : Zielfunktionswert f(x)
%
%---------------------------------------------------------------------------------------
function [x,lambda,fval] = lagrange_newton(x0,lambda0,tol,maxiter);
% Initialisierung
  x      = x0;
  lambda = lambda0;
  n      = size(x,1);
  m      = size(lambda,1);
  %  
  % Schleife
  %
  fprintf('Iter          f(x)                 ||g(x)||                KKT \n');
  fprintf('======================================================================\n');
  for i=0:1:maxiter
     [fval,gval,grad,gjac,hess] = f(x,lambda);
     kkt   = norm(grad);
     gnorm = norm(gval);
     %
     % Ausgabe
     %
     fprintf('%4d  %20.16f  %20.16f  %20.16f\n',i,fval,gnorm,kkt)
     %
     % Toleranz erfuellt? 
     %
     if (max(kkt,gnorm) <= tol)        
       disp('Konvergenz!');
       return;
     end     
     %
     % Berechne Schrittweite
     %
     A = [ hess , gjac' ; 
           gjac , zeros(m,m) ];
     b = [ -grad ; 
           -gval ];
     d = A\b;
     % Update
     x      = x + d(1:n);
     lambda = lambda + d(n+1:n+m);
  end;
  fprintf('Maximale Anzahl von Iterationen!')



%--------------------------------------------------------
%
% Berechne Zielfunktionswert, 
% Nebenbedingungen,
% Gradient der Lagrangefunktion und 
% Hessematrix der Lagrangefunktion
% 
%--------------------------------------------------------
   function [fval,gval,grad,gjac,hess] = f(x,lambda);
   fval = 2.0*x(1)^4 + x(2)^2 + 4.0*x(1)^2 - x(1)*x(2) + 6.0*x(2)^2;
   gval = [ 2.0*x(1) - x(2) + 4.0 ];
   grad = [ 8.0*x(1)^3 + 8.0*x(1) - x(2) + lambda(1)*2.0;
	     4.0*x(2)^3 - x(1) + 12.0*x(2) - lambda(1) ];
   gjac = [ 2.0 , -1.0 ];
   hess = [24.0*x(1)^2 + 8.0, -1.0 ; 
	     -1.0 , 12.0*x(2)^2 + 12.0];