Why does this code give error?

When I run the main, it gives me the following error:

Unable to perform assignment because the indices on the left side are not compatible with

the size of the right side.

Error in main (line 55)

two(run_idx,:) = temp(run_idx,ix1);

>>

Answer 1

[Torsten]

Force the array ix1 to be of the same length as ix by explicitly allocating it with the size of ix:

    % Run the GBO algorithm multiple times
    for run_idx = 1:num_runs
        % Run the GBO Algorithm
        [Best_fitness, BestPositions, time1] = GBO(nP, MaxIt, lb, ub, dim, @(b) myfun(b, u, P_far, P_near));
        
        % Store output of GBO in desired variables
        one(run_idx) = Best_fitness; % Fitness value
        temp(run_idx, :) = BestPositions; % Best estimated vector
        time(run_idx) = time1; % Elapsed time
        
        % Perform swapping
        [~, ix] = sort(u); % u is my desired vector
        ix1 = zeros(size(ix));
        [~, ix1(ix)] = sort(1:length(temp(run_idx,:)));  
        two(run_idx,:) = temp(run_idx,ix1)     
    end

[Sadiq Akbar]

@Torsten

Thanks a lot for your kind response. I ran it and it works too but it doesn't work as want. I want that the arrangement of elements in vectors two and u must be the same but it changes the arrangement of elements in vector two from the that of u.

[Torsten]

Maybe like this - I don't know exactly what you want.

clear; clc;

% Define the set of desired vectors and their corresponding P_far and P_near
u_sets = {
    [30 60 80 1.2],   
    [90 30 80 110 0.3 1.2],   
    [20 31 61 70 1.5],   
    [20 31 61 70 35 45 4.1 1.5] 
};

% Corresponding P_far and P_near values
P_values = [
    2, 1;  % For u = [30 60 80 1.2]
    2, 2;  % For u = [90 30 80 110 0.3 1.2]
    3, 1;  % For u = [20 31 61 70 1.5]
    4, 2   % For u = [20 31 61 70 35 45 4.1 1.5]
];

% Define GBO parameters
nP = 20; % Number of population members
MaxIt = 100; % Maximum number of iterations
num_runs = 5; % Number of runs per vector

% Loop over each desired vector
for idx = 1:length(u_sets)
    % Get the current vector and corresponding P_far and P_near
    u = u_sets{idx};
    [~, ix] = sort(u); 
    P_far = P_values(idx, 1);
    P_near = P_values(idx, 2);

    % Define bounds based on P_far and P_near
    dim = P_far + 2 * P_near; % Dimension of the problem
    lb = zeros(1, dim);
    ub = [180 * ones(1, P_far), 180 * ones(1, P_near), 5 * ones(1, P_near)];
    
    % Initialize result storage for the current vector
    one = zeros(num_runs, 1); % Fitness values
    temp = zeros(num_runs, dim); % Best estimated vectors
    time = zeros(num_runs, 1); % Elapsed times
    two = zeros(num_runs, dim); % Swapped best estimated vectors

    % Run the GBO algorithm multiple times
    for run_idx = 1:num_runs
        % Run the GBO Algorithm
        [Best_fitness, BestPositions, time1] = GBO(nP, MaxIt, lb, ub, dim, @(b) myfun(b, u, P_far, P_near));
        
        % Store output of GBO in desired variables
        one(run_idx) = Best_fitness; % Fitness value
        temp(run_idx, :) = BestPositions; % Best estimated vector
        time(run_idx) = time1; % Elapsed time
         
        two(run_idx,:) = temp(run_idx,ix);
        
    end
    
    two
    
    % Find the best fitness and corresponding estimated vector
    [best_fitness, best_idx] = min(one);
    best_positions = temp(best_idx, :);
    
    % Determine the file name based on P_far and P_near
    fileName = sprintf('%dF_and_%dJ.mat', P_far, P_near);
    
    % Save the workspace data for the current vector
    save(fileName, 'best_fitness', 'best_positions', 'time', 'two', 'u', 'P_far', 'P_near');
end

two = 5x4

5.0000 86.2511 32.3919 0.0000 5.0000 86.2510 32.3919 0 1.2026 60.0007 30.0016 80.0015 2.0142 60.0354 30.0625 80.2433 1.2194 30.0003 59.9999 80.0079

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x6

0.2932 1.2839 29.9713 79.8428 90.0052 110.1599 0.3262 2.2912 89.8344 80.6461 30.0182 111.0829 2.5512 3.9190 89.6414 82.7801 30.2667 111.3979 0.7265 2.7079 100.9539 89.0844 30.1995 82.9226 0.3547 4.9883 29.9967 81.1165 101.0871 89.3145

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x5

1.6803 30.7550 20.1940 60.8768 70.0652 1.3824 30.9597 20.1198 60.9735 69.7866 1.9491 20.0025 30.8888 60.9394 70.4321 0.8689 179.9927 64.6170 34.3261 73.9395 4.8041 30.7451 20.0053 60.8541 71.1950

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x8

1.5659 4.9999 71.7245 145.0666 3.1133 25.1807 30.5631 17.3745 1.6442 3.2021 66.0182 76.7475 4.3100 16.1569 22.3659 29.3521 2.8297 1.7876 29.6651 65.4358 26.3171 1.0665 21.0578 77.3154 3.3692 1.6603 64.8210 74.8839 34.5288 0 21.6582 29.5972 1.4346 2.4476 62.1805 30.5778 20.0715 48.3679 20.1977 70.2608

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

function e = myfun(b, u, P_far, P_near)
    M = length(b);

    
    u_far = u(1:P_far); 
    u_near = u(P_far+1:P_far+P_near); 
    u_near_ranges = u(P_far+P_near+1:end); 

    b_far = b(1:P_far);
    b_near = b(P_far+1:P_far+P_near); 
    b_near_ranges = b(P_far+P_near+1:end);

    
    k = (0:M-1)'; % Column vector

    % calculations1
    exp_far_u = exp(-1i * k * pi .* cosd(u_far)); 
    exp_far_b = exp(-1i * k * pi .* cosd(b_far)); 
    xo_far = sum(exp_far_u, 2); 
    xe_far = sum(exp_far_b, 2); 

    
    k_sq = k.^2; % Square of indices
    exp_near_u = zeros(M, P_near); 
    exp_near_b = zeros(M, P_near); 

    
    for i = 1:P_near
    
        exp_near_u(:, i) = exp(1i * (k * (-pi/2) .* cosd(u_near(i)) + ...
            (k_sq .* pi ./ (16 * u_near_ranges(i))) .* cosd(u_near(i)).^2)); 
    
        exp_near_b(:, i) = exp(1i * (k * (-pi/2) .* cosd(b_near(i)) + ...
            (k_sq .* pi ./ (16 * b_near_ranges(i))) .* cosd(b_near(i)).^2)); 
    end

    
    xo_near = sum(exp_near_u, 2); 
    xe_near = sum(exp_near_b, 2); 

   
    xo = xo_far + xo_near;
    xe = xe_far + xe_near;

    % Calculate Mean Squared Error (MSE)
    e = mean(abs(xo - xe).^2);
end

%___________________________________________________________________%                                                           %
%  Gradien-Based Optimizer source code (Developed in MATLAB R2017a) %
%                                                                   %
%  programming: Iman Ahmadianfar                                    %
%                                                                   %
%         e-Mail: [e-mail]                             %
%                 [e-mail]                         %
%                                                                   %
% Source codes demo version 1.0
% ------------------------------------------------------------------------------------------------------------
% Main paper (Please refer to the main paper):
% Gradient-Based Optimizer: A New Metaheuristic Optimization Algorithm
% Iman Ahmadianfar, Omid Bozorg-Haddad, Xuefeng Chu
% Information Sciences,2020
% DOI: https://doi.org/10.1016/j.ins.2020.06.037
% https://www.sciencedirect.com/science/article/pii/S0020025520306241
% ------------------------------------------------------------------------------------------------------------
% Website of GBO: http://imanahmadianfar.com/
% You can find and run the GBO code online at http://imanahmadianfar.com/

% You can find the GBO paper at https://doi.org/10.1016/j.ins.2020.06.037
% Please follow the paper for related updates in researchgate: https://www.researchgate.net/profile/Iman_Ahmadianfar
% ------------------------------------------------------------------------------------------------------------
%  Co-author:
%             Omid Bozorg-Haddad([e-mail])
%             Xuefeng Chu([e-mail])           
% _____________________________________________________
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% This function initialize the first population 
function X=initialization(nP,dim,ub,lb)
Boundary_no= size(ub,2); % numnber of boundaries

% If the boundaries of all variables are equal and user enter a signle
% number for both ub and lb

if Boundary_no==1
    X=rand(nP,dim).*(ub-lb)+lb;
end

% If each variable has a different lb and ub
if Boundary_no>1
    for i=1:dim
        X(:,i)=rand(nP,1).*(ub(i)-lb(i))+lb(i);
    end
end
end

%---------------------------------------------------------------------------------------------------------------------------
%  Author, inventor and programmer: Iman Ahmadianfar,
%  Assistant Professor, Department of Civil Engineering, Behbahan Khatam Alanbia University of Technology, Behbahan, Iran

%  Researchgate: https://www.researchgate.net/profile/Iman_Ahmadianfar

%  e-Mail: [e-mail], [e-mail],
%---------------------------------------------------------------------------------------------------------------------------
%  Co-author:
%             Omid Bozorg-Haddad([e-mail])
%             Xuefeng Chu([e-mail]) 
%---------------------------------------------------------------------------------------------------------------------------
% Please refer to the main paper:
% Gradient-Based Optimizer: A New Metaheuristic Optimization Algorithm
% SIman Ahmadianfar, Omid Bozorg-Haddad, Xuefeng Chu
% Information Sciences,2020
% DOI: https://doi.org/10.1016/j.ins.2020.06.037
% https://www.sciencedirect.com/science/article/pii/S0020025520306241
% ------------------------------------------------------------------------------------------------------------
% Website of GBO: http://imanahmadianfar.com/
% You can find and run the GBO code online at http://imanahmadianfar.com/

% You can find the GBO paper at https://doi.org/10.1016/j.ins.2020.06.037
% Please follow the paper for related updates in researchgate: https://www.researchgate.net/profile/Iman_Ahmadianfar
%--------------------------------------------------------------------------------------------------------------------------- 
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%function [Best_Cost,Best_X,Convergence_curve]=GBO(nP,MaxIt,lb,ub,dim,fobj)
function [Best_Cost,Best_X,time1]=GBO(nP,MaxIt,lb,ub,dim,fobj)
tic; % By Me
%% Initialization
nV = dim;                        % Number f Variables
pr = 0.51;%0.5;                  % Probability Parameter
%Note: Results on pr=0.5 was worse but best on pr=0.51 and above

lb = ones(1,dim).*lb;            % lower boundary 
ub = ones(1,dim).*ub;            % upper boundary
Cost = zeros(nP,1);
X = initialization(nP,nV,ub,lb); %Initialize the set of random solutions
Convergence_curve = zeros(1,MaxIt);

for i=1:nP
    Cost(i) = fobj(X(i,:));      % Calculate the Value of Objective Function
end

[~,Ind] = sort(Cost);     
Best_Cost = Cost(Ind(1));        % Determine the vale of Best Fitness
Best_X = X(Ind(1),:);
Worst_Cost=Cost(Ind(end));       % Determine the vale of Worst Fitness
Worst_X=X(Ind(end),:);
%% Main Loop
for it=1:MaxIt
    
    beta = 0.2+(1.2-0.2)*(1-(it/MaxIt)^3)^2;                        % Eq.(14.2)
    alpha = abs(beta.*sin((3*pi/2+sin(3*pi/2*beta))));              % Eq.(14.1)

    for i=1:nP                
        A1 = fix(rand(1,nP)*nP)+1;                                  % Four positions randomly selected from population        
        r1 = A1(1);r2 = A1(2);   
        r3 = A1(3);r4 = A1(4);        
        
        Xm = (X(r1,:)+X(r2,:)+X(r3,:)+X(r4,:))/4;                   % Average of Four positions randomly selected from population        
        ro = alpha.*(2*rand-1);ro1 = alpha.*(2*rand-1);        
        eps = 5e-3*rand;                                            % Randomization Epsilon
        
        DM = rand.*ro.*(Best_X-X(r1,:));Flag = 1;                   % Direction of Movement Eq.(18)
        GSR=GradientSearchRule(ro1,Best_X,Worst_X,X(i,:),X(r1,:),DM,eps,Xm,Flag);      
        DM = rand.*ro.*(Best_X-X(r1,:));
        X1 = X(i,:) - GSR + DM;                                     % Eq.(25)
        
        DM = rand.*ro.*(X(r1,:)-X(r2,:));Flag = 2;
        GSR=GradientSearchRule(ro1,Best_X,Worst_X,X(i,:),X(r1,:),DM,eps,Xm,Flag); 
        DM = rand.*ro.*(X(r1,:)-X(r2,:));
        X2 = Best_X - GSR + DM;                                     % Eq.(26)            
        
        Xnew=zeros(1,nV);
        for j=1:nV                                                  
            ro=alpha.*(2*rand-1);                       
            X3=X(i,j)-ro.*(X2(j)-X1(j));           
            ra=rand;rb=rand;
            Xnew(j) = ra.*(rb.*X1(j)+(1-rb).*X2(j))+(1-ra).*X3;     % Eq.(27)          
        end
        
        % Local escaping operator(LEO)                              % Eq.(28)
        if rand<pr           
            k = fix(rand*nP)+1;
            f1 = -1+(1-(-1)).*rand();f2 = -1+(1-(-1)).*rand();         
            ro = alpha.*(2*rand-1);
            Xk = unifrnd(lb,ub,1,nV);%lb+(ub-lb).*rand(1,nV);       % Eq.(28.8)

            L1=rand<0.5;u1 = L1.*2*rand+(1-L1).*1;u2 = L1.*rand+(1-L1).*1;
            u3 = L1.*rand+(1-L1).*1;                                    
            L2=rand<0.5;            
            Xp = (1-L2).*X(k,:)+(L2).*Xk;                           % Eq.(28.7)
                                                 
            if u1<0.5
                Xnew = Xnew + f1.*(u1.*Best_X-u2.*Xp)+f2.*ro.*(u3.*(X2-X1)+u2.*(X(r1,:)-X(r2,:)))/2;     
            else
                Xnew = Best_X + f1.*(u1.*Best_X-u2.*Xp)+f2.*ro.*(u3.*(X2-X1)+u2.*(X(r1,:)-X(r2,:)))/2;   
            end
        end
        
        % Check if solutions go outside the search space and bring them back
        Flag4ub=Xnew>ub;
        Flag4lb=Xnew<lb;
        Xnew=(Xnew.*(~(Flag4ub+Flag4lb)))+ub.*Flag4ub+lb.*Flag4lb;                
        Xnew_Cost = fobj(Xnew);
        % Update the Best Position        
        if Xnew_Cost<Cost(i)
            X(i,:)=Xnew;
            Cost(i)=Xnew_Cost;
            if Cost(i)<Best_Cost
                Best_X=X(i,:);
                Best_Cost=Cost(i);
            end            
        end
       % Update the Worst Position 
        if Cost(i)>Worst_Cost
            Worst_X= X(i,:);
            Worst_Cost= Cost(i);
        end               
    end
    
    Convergence_curve(it) = Best_Cost;
    % Show Iteration Information
    %disp(['Iteration ' num2str(it) ': Best Fitness = ' num2str(Convergence_curve(it))]);
end
time1=toc; % By Me
end
% _________________________________________________
% Gradient Search Rule
function GSR=GradientSearchRule(ro1,Best_X,Worst_X,X,Xr1,DM,eps,Xm,Flag)
    nV = size(X,2);
    Delta = 2.*rand.*abs(Xm-X);                            % Eq.(16.2)
    Step = ((Best_X-Xr1)+Delta)/2;                         % Eq.(16.1)
    DelX = rand(1,nV).*(abs(Step));                        % Eq.(16)
    
    GSR = randn.*ro1.*(2*DelX.*X)./(Best_X-Worst_X+eps);   % Gradient search rule  Eq.(15)
    if Flag == 1
      Xs = X - GSR + DM;                                   % Eq.(21)
    else
      Xs = Best_X - GSR + DM;
    end    
    yp = rand.*(0.5*(Xs+X)+rand.*DelX);                    % Eq.(22.6)
    yq = rand.*(0.5*(Xs+X)-rand.*DelX);                    % Eq.(22.7)
    GSR = randn.*ro1.*(2*DelX.*X)./(yp-yq+eps);            % Eq.(23)   
end

[Sadiq Akbar]

@Torsten

Thanks a lot for your kind response. Actually I want like this:

Look at the u_set at top. The 1st vector in u_set [30 60 80 1.2]. When the code runs, it becomes:

u = [30 60 80 1.2];

Now the algorithm will give me an estimated vector "BestPositions" having approximately the same values as in above u. But usually the positions of these elements in the estimated "BestPositions" are not the same as is in above u. Say for example the algorithm gives me the "BestPositions" as:

BestPositions = [60.01 30.05 79.99 1.20];

Now if you compare this "BestPositions" with above u, the arrangement of the values must be like that of u i.e., the "BestPositions" must look like this:

BestPositions = [30.05 60.01 79.99 1.20];

Each time the algorithm runs, it sends a random vector in which the positions of the elements are not the same as is in u. So I want to do swapping of those positions in the elements of "BestPositions" so that it becomes like u.

[Torsten]

But does it make sense that the algorithm returns the BestPositions vector in a different order than the u-vector ? If yes: what do you suggest as rule for ordering the BestPositions vector ?

[Sadiq Akbar]

@Torsten

Thanks a lot for your kind response. Yes, I have to take the data in the same order which I have to process further. If it's not in the same order, then the huge data is difficult to process further.

what do you suggest as rule for ordering the BestPositions vector ? That is what I have tried to implement in my code but it was giving error. you have removed the error but that doesn't meet my requirement. If you look at my swapping piece of code, it works in simple situation i.e., if you have a desired vector u and you have to arrange any other temp vector in the same order as is u, it works but I don't know why it doesn't work here. For example:

%%%%%%%%%%%%%%%%%%%%
% Swapping vector b
%%%%%%%%%%%%%%%%%%%%
u=[1 2 3 4 5 6];% Desired vector
b=[3 5 1 4 2 6];
[~, ix] = sort(u); 
[~, ix1(ix)] = sort(b);  
b = b(ix1)

If you run this piece of code, you will see that b is always arranged like vector u whatever arrangement of elements is taken in b, but it will always be arranged like u.

But I don't know what this logic doesn't work here.

[Torsten]

This code uses the ordering of "u" to order the rows of the array "two".

If there are bigger discrepancies between "u" and "two", the reason is in the optimization, not in the ordering of the elements of "two" with respect to "u".

clear; clc;
% Define the set of desired vectors and their corresponding P_far and P_near
u_sets = {
    [30 60 80 1.2],   
    [90 30 80 110 0.3 1.2],   
    [20 31 61 70 1.5],   
    [20 31 61 70 35 45 4.1 1.5] 
};
% Corresponding P_far and P_near values
P_values = [
    2, 1;  % For u = [30 60 80 1.2]
    2, 2;  % For u = [90 30 80 110 0.3 1.2]
    3, 1;  % For u = [20 31 61 70 1.5]
    4, 2   % For u = [20 31 61 70 35 45 4.1 1.5]
];
% Define GBO parameters
nP = 20; % Number of population members
MaxIt = 100; % Maximum number of iterations
num_runs = 5; % Number of runs per vector
% Loop over each desired vector
for idx = 1:length(u_sets)
    % Get the current vector and corresponding P_far and P_near
    u = u_sets{idx}
    [~, iu] = sort(u); 
    P_far = P_values(idx, 1);
    P_near = P_values(idx, 2);
    % Define bounds based on P_far and P_near
    dim = P_far + 2 * P_near; % Dimension of the problem
    lb = zeros(1, dim);
    ub = [180 * ones(1, P_far), 180 * ones(1, P_near), 5 * ones(1, P_near)];
    
    % Initialize result storage for the current vector
    one = zeros(num_runs, 1); % Fitness values
    temp = zeros(num_runs, dim); % Best estimated vectors
    time = zeros(num_runs, 1); % Elapsed times
    two = zeros(num_runs, dim); % Swapped best estimated vectors
    % Run the GBO algorithm multiple times
    itemp = zeros(size(iu));
    for run_idx = 1:num_runs
        % Run the GBO Algorithm
        [Best_fitness, BestPositions, time1] = GBO(nP, MaxIt, lb, ub, dim, @(b) myfun(b, u, P_far, P_near));
        
        % Store output of GBO in desired variables
        one(run_idx) = Best_fitness; % Fitness value
        temp(run_idx, :) = BestPositions; % Best estimated vector
        time(run_idx) = time1; % Elapsed time
        
        [~,itemp(iu)] = sort(temp(run_idx,:));
        two(run_idx,:) = temp(run_idx,itemp);        
    end
    
    two
    
    % Find the best fitness and corresponding estimated vector
    [best_fitness, best_idx] = min(one);
    best_positions = temp(best_idx, :);
    
    % Determine the file name based on P_far and P_near
    fileName = sprintf('%dF_and_%dJ.mat', P_far, P_near);
    
    % Save the workspace data for the current vector
    save(fileName, 'best_fitness', 'best_positions', 'time', 'two', 'u', 'P_far', 'P_near');
end

u = 1x4

30.0000 60.0000 80.0000 1.2000

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x4

51.7038 67.8083 179.9983 0.3997 5.0000 32.3919 86.2510 0.0000 30.0150 60.0099 80.2332 2.0858 30.0246 59.9941 80.1404 1.6926 5.0000 32.3919 86.2510 0.0008

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

u = 1x6

90.0000 30.0000 80.0000 110.0000 0.3000 1.2000

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x6

89.8147 30.2715 82.6205 110.8155 1.8021 2.1298 89.1124 30.1177 82.0499 101.0197 0.5196 4.4192 89.2303 30.0569 81.5234 101.0437 0.4138 4.3086 89.0464 30.2145 81.8429 100.9386 0.5079 1.8858 89.6163 30.3250 86.4456 110.3150 1.6628 3.2659

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

u = 1x5

20.0000 31.0000 61.0000 70.0000 1.5000

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x5

20.0247 30.7237 60.8448 71.2082 4.9954 5.0000 16.4220 32.8678 82.0236 0.0001 5.0000 16.4229 32.8673 82.0237 0.0004 19.9936 31.0032 61.0127 70.2146 1.6133 2.9148 33.6758 83.1080 179.9998 0.0488

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

u = 1x8

20.0000 31.0000 61.0000 70.0000 35.0000 45.0000 4.1000 1.5000

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

two = 5x8

19.5796 30.5865 62.2965 76.7366 33.5695 38.9464 4.8948 3.4709 3.3904 20.8595 64.4505 74.5641 29.9665 38.0050 1.7394 1.4672 9.6896 21.2386 65.2747 77.1490 24.5271 29.7454 2.8984 1.7447 20.4511 30.9837 60.9132 69.9424 36.4743 44.3011 4.5772 1.4903 20.4882 20.8783 64.0927 74.0697 30.4107 41.1541 4.6266 1.8831

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

function e = myfun(b, u, P_far, P_near)
    M = length(b);
    
    u_far = u(1:P_far); 
    u_near = u(P_far+1:P_far+P_near); 
    u_near_ranges = u(P_far+P_near+1:end); 
    b_far = b(1:P_far);
    b_near = b(P_far+1:P_far+P_near); 
    b_near_ranges = b(P_far+P_near+1:end);
    
    k = (0:M-1)'; % Column vector
    % calculations1
    exp_far_u = exp(-1i * k * pi .* cosd(u_far)); 
    exp_far_b = exp(-1i * k * pi .* cosd(b_far)); 
    xo_far = sum(exp_far_u, 2); 
    xe_far = sum(exp_far_b, 2); 
    
    k_sq = k.^2; % Square of indices
    exp_near_u = zeros(M, P_near); 
    exp_near_b = zeros(M, P_near); 
    
    for i = 1:P_near
    
        exp_near_u(:, i) = exp(1i * (k * (-pi/2) .* cosd(u_near(i)) + ...
            (k_sq .* pi ./ (16 * u_near_ranges(i))) .* cosd(u_near(i)).^2)); 
    
        exp_near_b(:, i) = exp(1i * (k * (-pi/2) .* cosd(b_near(i)) + ...
            (k_sq .* pi ./ (16 * b_near_ranges(i))) .* cosd(b_near(i)).^2)); 
    end
    
    xo_near = sum(exp_near_u, 2); 
    xe_near = sum(exp_near_b, 2); 
   
    xo = xo_far + xo_near;
    xe = xe_far + xe_near;
    % Calculate Mean Squared Error (MSE)
    e = mean(abs(xo - xe).^2);
end
%___________________________________________________________________%                                                           %
%  Gradien-Based Optimizer source code (Developed in MATLAB R2017a) %
%                                                                   %
%  programming: Iman Ahmadianfar                                    %
%                                                                   %
%         e-Mail: [e-mail]                             %
%                 [e-mail]                         %
%                                                                   %
% Source codes demo version 1.0
% ------------------------------------------------------------------------------------------------------------
% Main paper (Please refer to the main paper):
% Gradient-Based Optimizer: A New Metaheuristic Optimization Algorithm
% Iman Ahmadianfar, Omid Bozorg-Haddad, Xuefeng Chu
% Information Sciences,2020
% DOI: https://doi.org/10.1016/j.ins.2020.06.037
% https://www.sciencedirect.com/science/article/pii/S0020025520306241
% ------------------------------------------------------------------------------------------------------------
% Website of GBO: http://imanahmadianfar.com/
% You can find and run the GBO code online at http://imanahmadianfar.com/
% You can find the GBO paper at https://doi.org/10.1016/j.ins.2020.06.037
% Please follow the paper for related updates in researchgate: https://www.researchgate.net/profile/Iman_Ahmadianfar
% ------------------------------------------------------------------------------------------------------------
%  Co-author:
%             Omid Bozorg-Haddad([e-mail])
%             Xuefeng Chu([e-mail])           
% _____________________________________________________
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This function initialize the first population 
function X=initialization(nP,dim,ub,lb)
Boundary_no= size(ub,2); % numnber of boundaries
% If the boundaries of all variables are equal and user enter a signle
% number for both ub and lb
if Boundary_no==1
    X=rand(nP,dim).*(ub-lb)+lb;
end
% If each variable has a different lb and ub
if Boundary_no>1
    for i=1:dim
        X(:,i)=rand(nP,1).*(ub(i)-lb(i))+lb(i);
    end
end
end
%---------------------------------------------------------------------------------------------------------------------------
%  Author, inventor and programmer: Iman Ahmadianfar,
%  Assistant Professor, Department of Civil Engineering, Behbahan Khatam Alanbia University of Technology, Behbahan, Iran
%  Researchgate: https://www.researchgate.net/profile/Iman_Ahmadianfar
%  e-Mail: [e-mail], [e-mail],
%---------------------------------------------------------------------------------------------------------------------------
%  Co-author:
%             Omid Bozorg-Haddad([e-mail])
%             Xuefeng Chu([e-mail]) 
%---------------------------------------------------------------------------------------------------------------------------
% Please refer to the main paper:
% Gradient-Based Optimizer: A New Metaheuristic Optimization Algorithm
% SIman Ahmadianfar, Omid Bozorg-Haddad, Xuefeng Chu
% Information Sciences,2020
% DOI: https://doi.org/10.1016/j.ins.2020.06.037
% https://www.sciencedirect.com/science/article/pii/S0020025520306241
% ------------------------------------------------------------------------------------------------------------
% Website of GBO: http://imanahmadianfar.com/
% You can find and run the GBO code online at http://imanahmadianfar.com/
% You can find the GBO paper at https://doi.org/10.1016/j.ins.2020.06.037
% Please follow the paper for related updates in researchgate: https://www.researchgate.net/profile/Iman_Ahmadianfar
%--------------------------------------------------------------------------------------------------------------------------- 
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%function [Best_Cost,Best_X,Convergence_curve]=GBO(nP,MaxIt,lb,ub,dim,fobj)
function [Best_Cost,Best_X,time1]=GBO(nP,MaxIt,lb,ub,dim,fobj)
tic; % By Me
%% Initialization
nV = dim;                        % Number f Variables
pr = 0.51;%0.5;                  % Probability Parameter
%Note: Results on pr=0.5 was worse but best on pr=0.51 and above
lb = ones(1,dim).*lb;            % lower boundary 
ub = ones(1,dim).*ub;            % upper boundary
Cost = zeros(nP,1);
X = initialization(nP,nV,ub,lb); %Initialize the set of random solutions
Convergence_curve = zeros(1,MaxIt);
for i=1:nP
    Cost(i) = fobj(X(i,:));      % Calculate the Value of Objective Function
end
[~,Ind] = sort(Cost);     
Best_Cost = Cost(Ind(1));        % Determine the vale of Best Fitness
Best_X = X(Ind(1),:);
Worst_Cost=Cost(Ind(end));       % Determine the vale of Worst Fitness
Worst_X=X(Ind(end),:);
%% Main Loop
for it=1:MaxIt
    
    beta = 0.2+(1.2-0.2)*(1-(it/MaxIt)^3)^2;                        % Eq.(14.2)
    alpha = abs(beta.*sin((3*pi/2+sin(3*pi/2*beta))));              % Eq.(14.1)
    for i=1:nP                
        A1 = fix(rand(1,nP)*nP)+1;                                  % Four positions randomly selected from population        
        r1 = A1(1);r2 = A1(2);   
        r3 = A1(3);r4 = A1(4);        
        
        Xm = (X(r1,:)+X(r2,:)+X(r3,:)+X(r4,:))/4;                   % Average of Four positions randomly selected from population        
        ro = alpha.*(2*rand-1);ro1 = alpha.*(2*rand-1);        
        eps = 5e-3*rand;                                            % Randomization Epsilon
        
        DM = rand.*ro.*(Best_X-X(r1,:));Flag = 1;                   % Direction of Movement Eq.(18)
        GSR=GradientSearchRule(ro1,Best_X,Worst_X,X(i,:),X(r1,:),DM,eps,Xm,Flag);      
        DM = rand.*ro.*(Best_X-X(r1,:));
        X1 = X(i,:) - GSR + DM;                                     % Eq.(25)
        
        DM = rand.*ro.*(X(r1,:)-X(r2,:));Flag = 2;
        GSR=GradientSearchRule(ro1,Best_X,Worst_X,X(i,:),X(r1,:),DM,eps,Xm,Flag); 
        DM = rand.*ro.*(X(r1,:)-X(r2,:));
        X2 = Best_X - GSR + DM;                                     % Eq.(26)            
        
        Xnew=zeros(1,nV);
        for j=1:nV                                                  
            ro=alpha.*(2*rand-1);                       
            X3=X(i,j)-ro.*(X2(j)-X1(j));           
            ra=rand;rb=rand;
            Xnew(j) = ra.*(rb.*X1(j)+(1-rb).*X2(j))+(1-ra).*X3;     % Eq.(27)          
        end
        
        % Local escaping operator(LEO)                              % Eq.(28)
        if rand<pr           
            k = fix(rand*nP)+1;
            f1 = -1+(1-(-1)).*rand();f2 = -1+(1-(-1)).*rand();         
            ro = alpha.*(2*rand-1);
            Xk = unifrnd(lb,ub,1,nV);%lb+(ub-lb).*rand(1,nV);       % Eq.(28.8)
            L1=rand<0.5;u1 = L1.*2*rand+(1-L1).*1;u2 = L1.*rand+(1-L1).*1;
            u3 = L1.*rand+(1-L1).*1;                                    
            L2=rand<0.5;            
            Xp = (1-L2).*X(k,:)+(L2).*Xk;                           % Eq.(28.7)
                                                 
            if u1<0.5
                Xnew = Xnew + f1.*(u1.*Best_X-u2.*Xp)+f2.*ro.*(u3.*(X2-X1)+u2.*(X(r1,:)-X(r2,:)))/2;     
            else
                Xnew = Best_X + f1.*(u1.*Best_X-u2.*Xp)+f2.*ro.*(u3.*(X2-X1)+u2.*(X(r1,:)-X(r2,:)))/2;   
            end
        end
        
        % Check if solutions go outside the search space and bring them back
        Flag4ub=Xnew>ub;
        Flag4lb=Xnew<lb;
        Xnew=(Xnew.*(~(Flag4ub+Flag4lb)))+ub.*Flag4ub+lb.*Flag4lb;                
        Xnew_Cost = fobj(Xnew);
        % Update the Best Position        
        if Xnew_Cost<Cost(i)
            X(i,:)=Xnew;
            Cost(i)=Xnew_Cost;
            if Cost(i)<Best_Cost
                Best_X=X(i,:);
                Best_Cost=Cost(i);
            end            
        end
       % Update the Worst Position 
        if Cost(i)>Worst_Cost
            Worst_X= X(i,:);
            Worst_Cost= Cost(i);
        end               
    end
    
    Convergence_curve(it) = Best_Cost;
    % Show Iteration Information
    %disp(['Iteration ' num2str(it) ': Best Fitness = ' num2str(Convergence_curve(it))]);
end
time1=toc; % By Me
end
% _________________________________________________
% Gradient Search Rule
function GSR=GradientSearchRule(ro1,Best_X,Worst_X,X,Xr1,DM,eps,Xm,Flag)
    nV = size(X,2);
    Delta = 2.*rand.*abs(Xm-X);                            % Eq.(16.2)
    Step = ((Best_X-Xr1)+Delta)/2;                         % Eq.(16.1)
    DelX = rand(1,nV).*(abs(Step));                        % Eq.(16)
    
    GSR = randn.*ro1.*(2*DelX.*X)./(Best_X-Worst_X+eps);   % Gradient search rule  Eq.(15)
    if Flag == 1
      Xs = X - GSR + DM;                                   % Eq.(21)
    else
      Xs = Best_X - GSR + DM;
    end    
    yp = rand.*(0.5*(Xs+X)+rand.*DelX);                    % Eq.(22.6)
    yq = rand.*(0.5*(Xs+X)-rand.*DelX);                    % Eq.(22.7)
    GSR = randn.*ro1.*(2*DelX.*X)./(yp-yq+eps);            % Eq.(23)   
end

[Sadiq Akbar]

@Torsten

Thanks a lot for your consistent help. But this code also doesn't meet my requirements. I mean again, the arrangement of two and u is not the same.

If there are bigger discrepancies between "u" and "two", the reason is in the optimization:

As far as I know every meteheuristic algorithm generates random population and that's why the arrangement is different from u and my this logic helps in every case but not here. I don't know why? So, this is the best platefrom to ask for. You have solved my other issues several times in the past. May be you or one of your colleagues will do it.

Regards,

[Torsten]

The arrangement of "u" and "two" now is the same - meaning that if u has its i-th biggest element at position j, "two" will also have its i-th biggest element at position j. So the ordering of the two arrays is consistent.

I don't know why you get a different arrangement of BestPositions and u - in my opinion, this shouldn't be the case so that no post-ordering should be necessary (e.g. if you use an objective function like sum(BestPositions-u).^2).

[Sadiq Akbar]

@Torsten

Thanks a lot for your kind help. Yes, you are right. Actually I was comparing "best_positions" with "u" which was not same as u. Then I replaced the 56th line which was as:

best_positions = temp(best_idx, :);

by the following

best_positions = two(best_idx, :);

and re-ran it. Then the results were the same. Thanks a lot.