% ======================================================== % 
% Algorithm of the Innovative Gunner - AIG
% ======================================================== %   

% Basic version
% AIG optimization algorithm for nonlinear constrained optimization     %
% Programmed by Pawe³ Pijarski and Piotr Kacejko, Lublin University of Technology 2018, Poland   %
% p.pijarski@pollub.pl, p.kacejko@pollub.pl


% fobj - objective function
% n_iter - number of iteration
% ns - number of shots
% lb - lower bounds of variables
% ub - upper bounds of variables
% nvars - number of variables
% a - coefficiont of decrease angles 
% angle - max angle alpha and beta

% initial solution
 for i=1:ns
     shots(i,:)=lb+(ub-lb).*rand(size(lb));
 end
fitness=10^10*ones(ns,1);
[fmin,bestshot,shots,fitness]=best_shot(fobj,shots,shots,fitness);

%% Starting
angle=[alpha beta];
k=1;
while k<=n_iter        
    angle_c=decrease(angle,a,k);
    new_shot=get_shots(shots,lb,ub,angle_c);   
    [fnew,best,shots,fitness]=best_shot(fobj,shots,new_shot,fitness); 
        if fnew<fmin
            fmin=fnew;
            bestshot=best;
        end
k=k+1;
end

best_f=fmin;
best_x=bestshot;

%% --------------- functions ------------------
function shots=get_shots(shots,lb,ub,angle)

    [ns,nvars]=size(shots);
    w(1:ns,1:nvars)=random('Normal',0,angle(1,1)/3,ns,nvars);
    u(1:ns,1:nvars)=random('Normal',0,angle(1,2)/3,ns,nvars);
%     w(1:ns,1:nvars)=random('Uniform',-angle(1,1),angle(1,1),ns,nvars);
%     u(1:ns,1:nvars)=random('Uniform',-angle(1,2),angle(1,2),ns,nvars);

    delta1=w;
    delta2=u;
    d1u=delta1<0;
    d1d=delta1>=0;
    d2u=delta2<0;
    d2d=delta2>=0;
    wsp(1:ns,1:nvars)=1;
    wsp(d1u)=wsp(d1u).*cos(delta1(d1u)*pi/180);
    wsp(d1d)=wsp(d1d)./cos(delta1(d1d)*pi/180);
    wsp(d2u)=wsp(d2u).*cos(delta2(d2u)*pi/180);
    wsp(d2d)=wsp(d2d)./cos(delta2(d2d)*pi/180);
    shots=wsp.*shots;
    
i=1;
while i<=ns
    shots(i,:)=bounds(shots(i,:),lb,ub);
i=i+1;
end
end

function [fmin,best,shots,fitness]=best_shot(fobj,shots,newshot,fitness)

for j=1:size(shots,1)
    fnew=feval(fobj,newshot(j,:));
    if fnew<=fitness(j)
       fitness(j)=fnew;
       shots(j,:)=newshot(j,:);
    end
end
[fmin,t]=min(fitness);
best=shots(t,:);
end

function s=bounds(s,lb,ub)
  s1=s;
  i=s1<lb;
  s1(i)=lb(i);
  j=s1>ub;
  s1(j)=ub(j); 
  s=s1;
end
 function angle_c = decrease(angle,a,t)
 angle_c(1) = angle(1)*a(1)^t;
 angle_c(2) = angle(2)*a(2)^t;
 end
 
 % ----------------- end ------------------------------