function [best,fmin,N_iter]=bat_algorithm(para)
% Default parameters
if nargin<1, para=[10 0.25 0.5]; end
n=para(1); % Population size, typically 10 to 25
A=para(2); % Loudness (constant or decreasing)
r=para(3); % Pulse rate (constant or decreasing)
% This frequency range determines the scalings
Qmin=0; % Frequency minimum
Qmax=2; % Frequency maximum
% Iteration parameters
%% In order to obtain better/more accurate results, N_iter
%% should be increased to N_iter=2000 or more if necessary.
N_iter=1000; % Total number of function evaluations
% Dimension of the search variables
d=3;
% Initial arrays
Q=zeros(n,1); % Frequency
v=zeros(n,d); % Velocities
% Initialize the population/solutions
for i=1:n,
Sol(i,:)=randn(1,d);
Fitness(i)=Fun(Sol(i,:));
end
% Find the current best
[fmin,I]=min(Fitness);
best=Sol(I,:);% ====================================================== %
% Note: As this is a demo, here we did not implement the %
% reduction of loudness and increase of emission rates. %
% Interested readers can do some parametric studies %
% and also implementation various changes of A and r etc %
% ====================================================== %% Start the iterations -- Bat Algorithm
for i_ter=1:N_iter,
% Loop over all bats/solutions
for i=1:n,
Q(i)=Qmin+(Qmin-Qmax)*rand;
v(i,:)=v(i,:)+(Sol(i,:)-best)*Q(i);
S(i,:)=Sol(i,:)+v(i,:);
% Pulse rate
if rand>r
S(i,:)=best+0.01*randn(1,d);
end % Evaluate new solutions
Fnew=Fun(S(i,:));
% If the solution improves or not too loudness
if (Fnew<=Fitness(i)) & (rand<A) ,
Sol(i,:)=S(i,:);
Fitness(i)=Fnew;
end % Update the current best
if Fnew<=fmin,
best=S(i,:);
fmin=Fnew;
end
end