function [Fbest,BestChart,BestFitness,farmlayout,farmlayout_NA]=AGPSO(popsize,wf,iterations,t,NA_type,tn,wt,algorithmDir)

BestChart=zeros(iterations,1);
BestFitness=zeros(iterations,1);
farmlayout = zeros(iterations,wf.rows*wf.cols);
farmlayout_NA = zeros(iterations,wf.rows*wf.cols);
pop_power = zeros(popsize,wf.cols * wf.rows);

D = wf.turbine_num; down = 1;up = wf.rows *wf.cols;

% Population initialization
[popu,lu] = windfarm_init(popsize, wf.turbine_num,wf);


% Population evaluation
[popuFitness,power_order,lp_power_accum]= wf_fitness(wf,popu);
optimal = max(popuFitness);
for i = 1:popsize
    pop_power(i,power_order(i,:))  = pop_power(i,power_order(i,:)) + (lp_power_accum(i,:) / sum(lp_power_accum(i,:)));
end



%%
% Local parameters setting
vel = zeros(popsize, D);
pBest = popu;
pBestFit = popuFitness;
pBest_power = power_order;
[~, gBestId] = max(pBestFit);
gBest = pBest(gBestId,:);

omega = 0.7298;
c = 1.49618;
flag = zeros(popsize,1);
sg = 7;
pm = 0.01;
sum_count = 3;
count = 0;
for iter = 1:iterations
    %% Search machine

    popu = repalce_worst(popu, popsize, power_order,wf.NA_loc,pop_power);


    for i = 1:popsize
        %% Exemplar Update: Crossover
        offsPbest = zeros(1,D);
        for d = 1:D
            k = randperm(popsize,1);

            if pBestFit(i) > pBestFit(k)
                r = rand;
                offsPbest(d) = r * pBest(i,d) + (1 - r) * gBest(1,d);
            else
                offsPbest(d) = pBest(k,d);
            end
        end
        %% Exemplar Update: Mutation
        for d = 1:D
            if rand < pm
                offsPbest(d) = lu(1,d) + rand * (lu(2,d) - lu(1,d));
            end
        end

        %% Exemplar Update: Selection


        offsPbest = windfarm_constraint(offsPbest, wf.NA_loc, D,down,up);
        [offsPbestFitness,offs_power_order,offs_lp_power_accum]= wf_fitness(wf,offsPbest);
        optimal = max(popuFitness);

        if offsPbestFitness > optimal
            optimal = offsPbestFitness;

        end

        if offsPbestFitness > pBestFit(i)
            pBest(i,:) = offsPbest;
            pBestFit(i) = offsPbestFitness;
            pBest_power(i,:) = offs_power_order;
            pop_power(i,offs_power_order)  = pop_power(i,offs_power_order) +  (offs_lp_power_accum / sum(offs_lp_power_accum));
        end

        %% 20%M tournament
        if flag(i) == sg
            flag(i) = 0;
            competitor = randperm(popsize, 0.2 * popsize);
            [~, winId] = max(pBestFit(competitor));
            pBest(i,:) = pBest(competitor(winId),:);
            pBestFit(i) = pBestFit(competitor(winId));
            pBest_power(i,:) = pBest_power(competitor(winId),:);
        end

        %% Particle Update
        for d = 1:D
            vel(i,d) =   c * rand * (pBest(i,d) - popu(i,d)); %omega * vel(i,d) + omega * vel(i,d) +
            popu(i,d) = popu(i,d) + vel(i,d);
        end
    end
    %% Boundary detection
    popu = windfarm_constraint(popu, wf.NA_loc, D,down,up);
    [popuFitness,power_order,lp_power_accum] = wf_fitness(wf,popu);
    for i = 1:popsize
        pop_power(i,power_order(i,:))  = pop_power(i,power_order(i,:)) + (lp_power_accum(i,:) / sum(lp_power_accum(i,:)));
    end



    %% Evaluatition
    if max(popuFitness) > optimal
        optimal = max(popuFitness);
    else
        count = count +1;
    end

    %% Update pBest
    pos = popuFitness > pBestFit;
    flag(~pos) = flag(~pos) + 1;
    pBestFit(pos) = popuFitness(pos);
    pBest(pos,:) = popu(pos,:);
    [gBestFit, gBestId] = max(pBestFit);
    gBest = pBest(gBestId,:);
    %% Survival

    Fbest = gBestFit;
    BestChart(iter) = gBestFit / wf.power_total;
    BestFitness(iter) = gBestFit;
    [best_farmlayout,best_farmlayout_NA]  = gene_layout_by_indices_one(wf,gBest);
    farmlayout(iter,:) = best_farmlayout;
    farmlayout_NA(iter,:) = best_farmlayout_NA;
%     fprintf('%s NA %d Turbine Num:%d Wind %s run: %d iteration: %d  eta %f fitness %f\n',algorithmDir,NA_type,tn,wt,t,iter,(Fbest / wf.power_total),Fbest)
end

end