function [lsg,fx0,gx0]=startsqp(fktname,x0,gradname,glzahl,nlines);
% ======================================================================
% Startroutine fuer das SQP-Verfahren. 
%
% Erstellt fuer die Veranstaltung
% Weiterbildenes Studium fuer Ingenieure WS 98/99
%
% Helge Baumann, 17.8.1998, geaendert S.F.
% ----------------------------------------------------------------------
% Aurufparameter:
%  
% fktname : Name der Funktion, die einen Vektor [F,G] zurueckliefert,
%   wobei F (R^n->R) unter den m Nebenbedingungen 
%   G <=0 (G==0) (G:R^n ->R^m) minimiert werden soll.
%   An <fktname> wird ein Spaltenvektor(!) uebergeben. 
% x0 : Spaltenvektor der Dimension n, Startvektor
% gradname : Name der Gradientenfunktion zu F aus fktname.
% glzahl : Zahl der Gleichungsnebenbedingungen.
% nlines : Anzahl der zu zeichnenden Hoehenlinien. Standard: 40

%% Initialisierung 
% xmin haelt fuer jede Dimension die kleinste waehrend der Iteration
%      auftretende Zahl fest.
% xmax haelt fuer jede Dimension die groeste waehrend der Iteration
%      auftretende Zahl fest.
% xmin und xmax werden benoetigt um den richtigen Bildschirmausschnitt
% bei der Darstellung des Iterationsverlaufes anzuzeigen.
global xmin
global xmax

% Anzahl der Hoehenlinien
if nargin < 5, nlines = 40; end;

% Anzahl der Hoehenlinien
if nargin < 4, glzahl=0; end;
if isempty(glzahl), glzahl=0;end;

% Startvektor
if nargin < 2, x0 = [-2; 2]; end
[n,m] = size(x0); % n = Dimension des Startvektors/ der Startvektoren

% Name der zu minimierenden Funktion
if nargin < 1, fktname = 'rosenbrock'; end

% Name der Gradientenfunktion zu fktname
if nargin < 3, 
  if length(fktname)==length('rosenbrock')&fktname=='rosenbrock'; 
	gradname = 'grosenbrock'; 
  else 
        gradname=[];
  end;
end;

% Name der Ausgaberoutine
global outputname
outputname = '';

% Optionen vorgeben
options = [];
% Genauigkeit einstellen
%options(2) = 1e-6;
%options(3) = 1e-6;
% maximale Anzahl der Funktionsauswertungen:
%options(14) = 10000;
% Standardmaessig werden die m Nebenbedingungen in G als Ungleichungs-
% bedingungen aufgefasst (G<=0). 
% Sind unter den m Nebenbedingung aber k Gleichungsnebenbedingungen (G==0),
% so muessen diese in den ersten k Komponenten von G definiert werden und
% mittels der Einstellung option(13)=k als Gleichungsnebenbedingungen 
% deklariert werden.
options(13) = glzahl;

% Konstanten
xpoints = 41;   % Anzahl der Gitterpunkte in x Richtung fuer contour-plot
ypoints = 41;   % Anzahl der Gitterpunkte in y Richtung fuer contour-plot

% Im Fall n==2 kann eine graphische Ausgabe erfolgen. Dafuer sind einige
% Vorarbeiten noetig, die fuer die eigentliche Minimierung keine Rolle
% spielen. Der gesamte 'if'-Block dient ausschliesslich diesen Vorbereitungen.
if (n==2)  
   % Fuer xmin (xmax) werden die groestmoeglichen (kleinstmoeglichen) Zahlen
   % vorgegeben. Dient nur dazu einen Vektor in entsprechender
   % Groesse bereitzustellen.
   xmin = ones(n,1)*realmax;
   xmax = -ones(n,1)*realmax;
   
   outputname = 'minmax';
   % Der erste Durchlauf ist nur zur Bestimmung von xmin und xmax gedacht.
   [lsg,options] = constr(fktname,x0,options,[],[],gradname);

   % Beobachtungsfenster festlegen, um 10% groesser als ermittelte Werte
   xlaenge = (xmax-xmin)/2; xmitte = (xmax+xmin)/2;
   xmin = xmitte -(xlaenge* 1.1);
   xmax = xmitte + (xlaenge*1.1);

   % Einheitliches Koordinatensystem und Graphikausgabe in ein Fenster
   clf;
   axis([xmin(1), xmax(1), xmin(2), xmax(2)]);
   hold on
   %% Hoehenlinien einzeichnen
   % Gitter erzeugen; zu den Knoten wird der Funktionswert ermittelt.
   xrange = xmin(1):(xmax(1)-xmin(1))/(xpoints-1):xmax(1);
   yrange = xmin(2):(xmax(2)-xmin(2))/(ypoints-1):xmax(2);
  [xgitter,ygitter] = meshgrid(xrange,yrange);
  for i=1:xpoints
     [zgitter(i,1:ypoints),ghilfe] = ...
        feval(fktname,[xrange(i)*ones(1,ypoints);yrange(1:ypoints)]);     
	ghilfedim=size(ghilfe);
	if ghilfedim(1)==1
           ggitter(i,1:ypoints)=ghilfe;
        else
           ggitter(i,1:ypoints)=max(ghilfe);
        end;
   end 

   % Erzeugzung der Hoehenlinien
   fx0=feval(fktname,lsg);
   zmax=max(max(zgitter));
   zmin=min(min(min(zgitter)),fx0);
   zlmax=zmax;
   zlmin=zmin;
   ZZ=linspace(zlmin,zlmax,nlines);
   contour(xgitter,ygitter,zgitter',ZZ);
   
   if options(10) < options(14)
      plot(lsg(1),lsg(2),'r*','Erasemode','none','Markersize',10);
drawnow;
	% Minimum einzeichnen
%     legend('berechnetes Minimum');
   end


   % Restriktionen einzeichnen
   contour(xgitter,ygitter,ggitter',[0,0.1,0.25,0.5,1,2,4,8,16,32,64],'k');
%   contour(xgitter,ygitter,ggitter',10,'k');

   % Erlaeuterungen in die Zeichnung einfuegen
   title('SQP - Verfahren');
   % nur wenn ein Minimum gefunden wurde, wird es eingezeichnet
   
   outputname = 'zeichnung'; % Name der Ausgaberoutine
end

% Der Durchlauf, bei dem in jedem Iterationsschritt eine Ausgabe
% erfolgen kann.
[lsg,options] = constr(fktname,x0,options,[],[],gradname);

if n==2,

   if options(10) < options(14)
      plot(lsg(1),lsg(2),'r*','Markersize',10);
%drawnow;		% Minimum einzeichnen
      legend('berechnetes Minimum');
   end
    hold off;
%colorbar
end;

[fx0,gx0]=feval(fktname,lsg);
% Ergebnisse anzeigen
format long
disp('Loesungsvektor')
disp(lsg)
disp(sprintf('%s\n%s %d %s', 'Fuer einen Durchlauf des SQP-Verfahrens',' wurden',options(10),'Funktionsauswertungen benoetigt.'))
disp(sprintf('%s\n%s %d %s', 'Fuer einen Durchlauf des SQP-Verfahrens',' wurden',options(11),'Gradientenauswertungen benoetigt.'))












