function [xmin,ymin,nS,vetxmin,vetymin] = ag_r(fun, nvi, valmin, valmax, prec, tolx, ...
   nind, sind, nger, reptx, muttx, rep, nfilhos, maxd, felitismo, nelite, ...
   freset, nindot, greset)

% Função AG_R - Acha o mínimo de uma função utilizando algoritmo genético
% para funções contínuas - codificação real
%
%		xmin = ag_r(fun,nvi,valmin,valmax,prec)
%
%	onde:   fun - função objetivo (string)
%		nvi - número de variáveis independentes (inteiro, >0)
%		valmin - vetor ou escalar com os valores mínimos para variáveis 
%					independentes
%		valmax - vetor ou escalar com os valores máximos para variáveis 
%					independentes
%		prec - passo para geração das variáveis independentes (vetor ou escalar)
%				 para gerar a primeira geração (pais)	 
% 
%  sendo: xmin - o valor da(s) variável(is) independente(s) onde o
% 					  valor mínimo foi encontrado
%
% Na formulação anterior uma série de parâmetros de ajuste foi suprimida.
% Subsequentemente mostra-se a chamada completa da função:
% [xmin,ymin,vetxmin,vetymin] = ag_r(fun, nvi, valmin, valmax, prec, tolx, nind, sind, nger, ...
%					reptx, muttx, rep, nfilhos, maxd, felitismo, nelite, ...
%					freset, nindot, greset)
%		
% Onde:
%		tolx - tolerância relativa na variável de decisão (default 1e-4)
%		nind - Número de indivíduos (default: 200)	
%		sind - Algoritmo de escolha da geração inicial (pais) 
%				 1 - roleta (default)
%				 2 - torneio
%		nger - Número máximo de gerações (default: 20)
%		reptx - Taxa de cruzamento (default: 70 %)
%		muttx - Taxa de mutação: (default: 5 %)
%		rep - Tipo de cruzamento: 
%			1 - Igual aos pais
%			2 - Média aritimética entre os pais
%			3 - Média ponderada com pesos aleatórios entre pai e mãe
%		nfilhos - Número de filhos por cruzamento: (default: 2)
%		maxd -  maxima porcentagem de deslocamento para mutacao 		
%		felitismo - flag de ativação do elitismo 
%				1 - ligado (default) 
%				0 - desligado 
%		nelite - Número de indivíduos da elite: (default: 5)
%		freset - Flag relativo ao Drift(reset periódico da população menos apta)
%				1 - ligado  
%				0 - desligado (default)
%		nindot - Número de individuos carregados para a geracao seguinte (default: 10)
%		greset - Número de geração a cada reset da população (default: 5)

% Marcelo Farenzena - 12/2004
% GIMSCOP - DEQUI - UFRGS

% ==========================================================================
% Verificação da bobagem feita pelo usuário
% ==========================================================================

if nargin == 0
   help ag_r
    error('Incluir os argumentos')
end

if nargin < 5
   error('A chamada mínima da função é: ag_r(fun,nvi,valmin,valmax,prec)')
end

if isempty(fun)
   error('A função deve ser informada!')
end
   
if isempty(nvi)
   error('O número de variáveis independentes deve ser informado!') 
end

if isempty(valmin)
   error('O valor mínimo deve ser informado!')  
end

if isempty(valmax)
   error('O valor máximo deve ser informado!')  
end

if isempty(prec)
   error('A precisão deve ser informada') 
end

switch nargin
case 5
   tolx = 1e-4;
   nind = 200;		% número inicial de individuos
   sind = 1;		% algoritmo de escolha da geração inicial (pais)
   nger = 20;		% número máximo de gerações
   reptx = 70;		% taxa de reprodução (%)
   muttx = 5;		% taxa de mutação    (%)
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população

case 6    
   nind = 200;		% número inicial de individuos
   sind = 1;		% algoritmo de escolha da geração inicial (pais)
   nger = 20;		% número máximo de gerações
   reptx = 70;		% taxa de reprodução (%)
   muttx = 5;		% taxa de mutação    (%)
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 7
   sind = 1;		% algoritmo de escolha da geração inicial (pais)
   nger = 20;		% número máximo de gerações
   reptx = 70;		% taxa de reprodução (%)
   muttx = 5;		% taxa de mutação    (%)
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 8
   nger = 20;		% número máximo de gerações
   reptx = 70;		% taxa de reprodução (%)
   muttx = 5;		% taxa de mutação    (%)
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 9
   reptx = 70;		% taxa de reprodução (%)
   muttx = 5;		% taxa de mutação    (%)
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 10
   muttx = 5;		% taxa de mutação    (%)
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 11
   rep = 3;			% tipo de reprodução (cruzamento)
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 12
   nfilhos = 2;   % número de filhos por cruzamento
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 13
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 14
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 15
   nelite = 5;		% Número de indivíduos da elite 
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população

case 16
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 17
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
   greset = 5;		% Número de geração a cada reset da população
   
case 18
   greset = 5;		% Número de geração a cada reset da população
end

if isempty(nind)
nind = 200;		% número inicial de individuos   
end

if isempty(sind)
   sind = 1;		% algoritmo de escolha da geração inicial (pais)
end

if isempty(nger)
   nger = 20;		% número máximo de gerações
end

if isempty(reptx)
   reptx = 70;		% taxa de reprodução (%)
end

if isempty(muttx)
   muttx = 5;		% taxa de mutação    (%)
end

if isempty(rep)
   rep = 3;			% tipo de reprodução (cruzamento)
end

if isempty(nfilhos)
   nfilhos = 2;   % número de filhos por cruzamento
end

if isempty(maxd)
   maxd = 10;		% maxima porcentagem de deslocamento para mutacao
end

if isempty(felitismo)
   felitismo = 1; % Flag do elitismo (0 - ativado, 1 - ligado)
end

if isempty(nelite)
   nelite = 5;		% Número de indivíduos da elite 
end

if isempty(freset)
   freset = 0;		% Flag do reset periódico da população (drift) (0 - ativado, 1 - ligado)
end

if isempty(nindot)
   nindot = 10;	% Número de individuos carregados para a geracao seguinte
end

if isempty(greset)
   greset = 5;		% Número de geração a cada reset da população
end

% ==========================================================================
% Não fez bobagem, agora vai
% ==========================================================================

% =================================================
% selecao da primeira geracao (pais)
% =================================================

ntind = abs(valmax - valmin)./prec;			% calculo do numero total de individuos

% verifica se ha reset 
if freset
   nciclos = ceil(nger/greset);			% quantas vezes ocorrera reset
   nger = greset;
else
   nciclos = 1;
end

nS = 0;
for nc = 1:nciclos   
   for i = 1 : nind
      if sind == 1                     % roleta
         npind = 5; 							% numero de individuos da roleta
         pind = rand(npind,nvi);			% escolhe os indivíduos 
         
         pind = adeqscal(pind,valmax-valmin,valmin);            
         
         for j = 1:npind						% avalia a aptidao dos provaveis indivíduos
            aptpi(j) = feval(fun,pind(j,:));                  
         end
         nS = nS + npind;
         
         if min(aptpi > 0)
            aptpi = 1./(aptpi);				% defini aptidão como inverso do val da func objetivo
         else
            aptpi = 1./(aptpi - 2*min(aptpi));
         end
         
         prob = aptpi./sum(aptpi);			% normaliza as aptidoes, para obter-se a faixa de probabilidades
         
         for j = 1:npind
            proba(j) = sum(prob(1:j));		% calcula as probabilidades acumuladas
         end
         
         indesc = rand(1);					% individuo escolhido
         
         indmen = find(proba < indesc);
         if isempty(indmen)
            ind(i,:) = pind(1,:);
            apt(i) = aptpi(1);
         else
            ind(i,:) = pind(max(indmen) + 1,:);
            apt(i) = aptpi(max(indmen) + 1);
         end      
         
      elseif sind == 2
         npind = 5; 							% numero de individuos da roleta
         pind = rand(npind,nvi);			% escolhe os indivíduos 
         
         pind = adeqscal(pind,valmax-valmin,valmin);            
         for j = 1:npind						% avalia a aptidao dos provaveis indivíduos
            aptpi(j) = feval(fun,pind(j,:));                  
         end
         nS = nS + npind;
         
         if min(aptpi > 0)
            aptpi = 1./(aptpi);				% defini aptidão como inverso do val da func objetivo
         else
            aptpi = 1./(aptpi - 2*min(aptpi));
         end
         
         % acha o individuo mais apto
         
         ind(i,:) = pind(find(aptpi == min(aptpi)),:);
         apt(i) = min(aptpi);      
      else
         error('Defina corretamente a seleção da primeira geração')
         break
      end         
   end
   
   if nc > 1
      ind(end-nindot+1:end,:) = indold;   
   	apt(end-nindot+1:end) = aptold;        
   end      
   
   % =================================================
   % geraçao das próximas gerações
   % =================================================
   % calculo das novas aptidoes
   for j=1:size(ind,1)
      apt(j) = feval(fun,ind(j,:));
   end
   nS = nS + size(ind,1);
   
   if min(apt > 0)
      apt = 1./(apt);				% defini aptidão como inverso do val da func objetivo
   else
      apt = 1./(apt - 2*min(apt));
   end      
   
   ger = 1;
   aptnew = 1000;
   tolc = 10;
   
   while ger <= nger && nind > 4 %& tolc > tolx
    
      % determinacao dos cruzamentos
      apto = sort(apt);   						% ordena as aptidões
      nindrep = ceil(reptx/100*nind);		% determina o numero de ind que irao ficar
      pcorte = apto(end - nindrep);
      
      [ninapt,inapt] = find(apt > pcorte);  
      
      if length(inapt)<16
         [ninapt,inapt] = find(apt >= pcorte);   
        inapt = inapt(1:nindrep);
      end           
            
      newind = [];
      for j = 1:length(inapt)/2
         switch rep
         case 1			% igual ao pai (coisinha do papai)
            switch nfilhos 
            case 1
               newind(j,:) = ind(inapt(1 + (j-1)*2),:);
            case 2
               newind(2*j-1,:) = ind(inapt(1 + (j-1)*2),:);
               newind(2*j,:) = ind(inapt(2 + (j-1)*2),:);               
            otherwise
               newind(nfilhos*(j-1)+1,:) = ind(inapt(1 + (j-1)*2),:);
               newind(nfilhos*(j-1)+2,:) = ind(inapt(2 + (j-1)*2),:);               
               
               for k = 1 : nfilhos - 2
                  p = rand(1);
                  newind(nfilhos*(j-1) + k + 2,:) = p*ind(inapt(1 + (j-1)*2),:) + (1-p)*ind(inapt(1 + (j-1)*2),:);
               end
            end            
            
         case 2			% média aritimética entre o pai e mãe
            switch nfilhos 
            case 1            
               newind(j,:) = 0.5*(ind(inapt(1 + (j-1)*2),:) + ind(inapt(2 + (j-1)*2),:));                  
            otherwise
               newind(nfilhos*(j-1)+1,:) = 0.5*(ind(inapt(1 + (j-1)*2),:) + ind(inapt(2 + (j-1)*2),:));                              
               for k = 1 : nfilhos - 1
                  p = rand(1);
                  newind(nfilhos*(j-1) + k + 1,:) = p*ind(inapt(1 + (j-1)*2),:) + (1-p)*ind(inapt(2 + (j-1)*2),:);
               end  
            end   
         case 3			% média ponderada (mais parecido com a mãe)
            for k = 1 : nfilhos
               p = rand(1);
               newind(nfilhos*(j-1) + k ,:) = p*ind(inapt(1 + (j-1)*2),:) + (1-p)*ind(inapt(2 + (j-1)*2),:);
            end  
         end      
      end
                    
      
      % ============================================================================
      % determinacao das mutações 
      % altera-se uma porcentagem do indivíduo, para mais ou menos, aleatoriamente)
      % ============================================================================
      
      nindm = floor(muttx/100*nind); % numero de individuos mutados
      
      if nindm > 0
         for k = 1 : nindm
            pnewind = ind(ceil(nind*rand(1)),:);		%escolha do provavel indivíduo
            pnewind = pnewind + 2*(rand(1)-0.5)*pnewind*maxd/100;
            while any(pnewind > valmax) || any(pnewind < valmin)
               pnewind = pnewind + 2*(rand(1)-0.5)*pnewind*maxd/100;      % verifica se não é um individuo expurio   
            end      
            newind(nfilhos*j+k,:) = pnewind;
         end
      end
      
      % ============================================================================   
      % verifica se será utilizado elitismo
      % ou seja, clona-se os individuos mais aptos
      % ============================================================================
      
      if felitismo
         [ninapt,inapt] = find(apt > apto(end - nelite));       
         if length(inapt) == nelite
            newind(nfilhos*j+k+1:nfilhos*j+k+nelite,:) = ind(inapt,:);
         else
            [ninapt,inapt] = find(apt >= apto(end - nelite));       % Caso haja aptidoes repetidas
            newind(nfilhos*j+k+1:nfilhos*j+k+nelite,:) = ind(inapt(1:nelite),:);            
         end         
      end
      
      clear ind  apt
      ind = newind;
      
      % calculo das novas aptidoes
      for j=1:size(ind,1)
         apt(j) = feval(fun,ind(j,:));
      end
      nS = nS + size(ind,1);
      
      if min(apt > 0)
         apt = 1./(apt);				% defini aptidão como inverso do val da func objetivo
      else
         apt = 1./(apt - 2*min(apt));
      end      
      
      nind = size(newind,1);      
      clear newind
      
%      disp(['Geração numero ',int2str(ger)])
      ger = ger + 1;  
      
      aptold = aptnew;
      aptnew = max(apt);
      
      % calcula a diferença relativa entre as aptidões
      tolc = abs(aptnew/aptold - 1);
      
   end
   
   if freset
      apto = sort(apt);
      aptold = apto(end-nindot+1:end);
      for k=1:length(aptold)
         aptf(k) = find(apt == aptold(k));
      end      
      indold = ind(aptf,:);			% guarda os individuos q passarão pelo drift
   end         
end

[val,pos] = max(apt);
xmin = ind(pos,:);
ymin = feval(fun,xmin); 
vetxmin = ind;

for i=1:size(ind,1)
	vetymin(i) = feval(fun,vetxmin(i,:));      
end

nS = nS + size(ind,1) + 1;

% ===============================================
% Funções auxiliares
% ===============================================

function matres = adeqscal(mat,coefa,coefl) 

for i = 1:length(coefa)
   matres(:,i) = mat(:,i) * coefa(i) + coefl(i);
end




