  function [xnew,fnew,gnew] = SDaskel(xc,fc,gc,Lmax,auto)
% [xnew,fnew,gnew] = SDStep(xc,f,g,Lmax)
% Generates a steepest descent step.
%
%  Globaalit funktiot:
%
%  objf.m ja objg.m
%    z=objf(x)
%    Syöte: n-sarakevektori x
%    Tulos: kohdefunktion arvon x:ssä, skalaari  
%
%  v=objg(x)
%    Syöte: n-sarakevektori x
%    Tulos: kohdefunktion gradientin arvon x:ssä, n-sarakevektori
%
%  Argumentit:
%    input:
%     xc     a column-vector, an approximate minimizer
%     fc     f(xc)
%     gc     g(xc)
%   Lmax     max step length
%    output:
%   xnew     a column n-vector, an improved minimizer
%   fnew     f(xnew)
%   gnew     g(new)
%
% Esimerkkifunktiot: Tiedostoissa objf.m ja objg.m
%function z=objf(X)
%%  z=objf(x)
%%  Syöte: n-sarakevektori x
%%  Tulos: kohdefunktion arvon x:ssä, skalaari  
%x=X(1);y=X(2);
%z=(x-2)^4+(x-2*y)^2;

%function v=objg(X)
%%  z=objg(X)
%%  Syöte: n-sarakevektori X
%%  Tulos: kohdefunktion arvon X:ssä, n-sarakevektori
%x=X(1);y=X(2);
%v=[4*(x-2)^3+2*x-4*y;-4*x+8*y];


nplotvals = 20;

% Get the Steepest Descent Search Direction
sc = -gc;
% Line Search
% Manuaalinen, vuorovaikutteinen haku:
if auto==0
  % Sample across [0,L]
  L=Lmax;
  lambdavals = linspace(0,L,nplotvals);
  fvals = zeros(1,nplotvals);
  for k=1:nplotvals
    fvals(k) = objf(xc+lambdavals(k)*sc);
  end
  figure
  plot(lambdavals,fvals)

else
  % Automaattinen haku väliä puolittamalla
  % Try to get L<=Lmax so xc+L*sc is at least as good as xc.
  L = Lmax;
  Halvings = 0;
  fL = objf(xc+L*sc);
  while (fL>=fc) & Halvings<=10
    L = L/2;
    Halvings = Halvings+1;
    fL = objf(xc+L*sc);
   end  
% Sample across [0,L]
   lambdavals = linspace(0,L,nplotvals);
   fvals = zeros(1,nplotvals);
   for k=1:nplotvals
      fvals(k) = objf(xc+lambdavals(k)*sc);
   end
end

   % Line search.
[fnew,i] = min(fvals(1:nplotvals));
xnew = xc + lambdavals(i(1))*sc;
gnew = objg(xnew);

%% Alifunktiot - subfunctions, ei käytetä sittenkään %%%


%  Kohdefunktio --  f
%function y=f(X)
%x=X(1);y=X(2);
%y=(x-2)^4+(x-2*y)^2;

%  gradientti -- g
%function z=g(X)
%x=X(1);y=X(2);
%z=[4*(x-2)^3+2*x-4*y;-4*x+8*y];







