antinface.m

蚁群算法采用matlab开发的仿真平台：算法实现

%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end



function edit_optNO_Callback(hObject, eventdata, handles)
% hObject    handle to edit_optNO (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user state (see GUIDATA)

% Hints: get(hObject,'String') returns contents of edit_optNO as text
%        str2double(get(hObject,'String')) returns contents of edit_optNO as a double


% --- Executes during object creation, after setting all properties.
function edit_optNO_CreateFcn(hObject, eventdata, handles)
% hObject    handle to edit_optNO (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end


% --- Executes on button press in calbutton.
function calbutton_Callback(calbutton, eventdata, handles)
% hObject    handle to calbutton (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user state (see GUIDATA)
%%%%%%%%蚁群算法用于路径规划求解%%%%%%%%%%
%%%%%%%%%Initialization phase%%%%%%%%%%
m = str2double(get(handles.edit_antsum,'String'));      % set ant number by using Matlab GUI
initao = str2double(get(handles.edit_initao,'String')); %信息量
alpha = str2double(get(handles.edit_alpha,'String'));   %信息启发因子，大趋向于选择其他蚂蚁经过的路径
beta = str2double(get(handles.edit_beta,'String'));     %期望启发因子，大趋向于选择其他蚂蚁经过的路径
Q = str2double(get(handles.edit_q,'String'));           %信息强度
rou = str2double(get(handles.edit_rou,'String'));       %信息挥发系数
NCmax = str2double(get(handles.edit_ncmax,'String'));   %最大循环次数
n = 0;  %初始化城市个数
%%%%%%%%从读取数据处传递过来的%%%%%%%%%%%%%%%%%%%%%
Ai0 = get(handles.obtext1,'UserData');      %城市障碍物地图
Bi0 = get(handles.openfile,'UserData');     %城市地图
B=reshape(Bi0,4,[]);
%disp(B);
%disp(Ai0);
ng=max(B(1,:));        %第一行最大数即为路径条数,为计算方便每条路径上的城市个数相等
n = str2double(get(handles.edit_citysum,'String'));           %每条路径上的长度
%disp(n);   
    for k=1:m
            for i=1:n
               g(k,i)= i; %按最大循环次数和蚂蚁数量产生城市索引值,即路径
            end
    end;
disp(g);
%初始化路径显示
plot(B(3,1:n),B(4,1:n),'-r');
hold on;
plot(B(3,n+1:2*n),B(4,n+1:2*n),'-r');
hold on;
plot(B(3,2*n+1:3*n),B(4,2*n+1:3*n),'-r');
hold on;
plot(B(3,3*n+1:4*n),B(4,3*n+1:4*n),'-r');
hold on;
plot(B(3,4*n+1:5*n),B(4,4*n+1:5*n),'-r');
hold on;
plot(B(3,5*n+1:6*n),B(4,5*n+1:6*n),'-r');
hold on;
plot(B(3,6*n+1:7*n),B(4,6*n+1:7*n),'-r');
hold on;
plot(B(3,7*n+1:8*n),B(4,7*n+1:8*n),'-r');
hold on;
plot(B(3,8*n+1:9*n),B(4,8*n+1:9*n),'-r');
hold on;
plot(B(3,9*n+1:10*n),B(4,9*n+1:10*n),'-r');
hold on;
plot(B(3,10*n+1:11*n),B(4,10*n+1:11*n),'-r');
hold on;
plot(B(3,11*n+1:12*n),B(4,11*n+1:12*n),'-r');
hold on;
disp('Variable initiate ok!');

initao = initao.*ones(1,ng);
for NC =1:NCmax
  tao(:,NC)=initao;           %初始化信息量
end;
detatao = zeros(ng);          %增加的信息量 
%disp(tao(:,2));
bestsolution = 12000;         %最后选最优路径时用到(20*600)
%%%%%%%%%初始化完毕%%%%%%%%
%%%%%%%%%可能路径%%%%%%%%%%%%%
for NC = 1:NCmax 
  maxp(1,NC)=0;   
  tabu = zeros(m,ng,NC);        %各条路径上的禁忌表  ????
  for j = 1:ng
            tao(j,NC)=(1-rou)*tao(j,NC)+rou*initao(1,j); %Local updating for the first iteration 公式2.4.5信息素初始化
  end;
  for l=1:ng                           %产生每次循环可能的路径
      xc(1)=B(3,1);
      yc(1)=B(4,1);
      xc(n)=B(3,31);
      yc(n)=B(4,31);
      for i=2:n-1  
        xc(i)=B(3,i+n*(l-1))+rand(1)*20;
        yc(i)=B(4,i+n*(l-1))+rand(1)*20;        
    end;
    C(:,:,l,NC)=cat(1,xc,yc);         %保存每次的各条路径的座标
    disp(C);
    %plot(xc,yc,'-r');
    %disp(C(1,:,l,NC));   
    %disp(C(2,:,l,NC));
    for i=1:n-1                        %分别计算每条路径上各城市坐标点之间的距离
          distance(i,l,NC)=sqrt((C(1,i,l,NC)-C(1,i+1,l,NC))^2+(C(2,i,l,NC)-C(2,i+1,l,NC))^2);  
    end;
    disp(distance);
    distancel(l)=sum(distance(:,l,NC));  %计算出每次每条路径的长度
  end;
  distl(:,NC)=distancel(:);              %保存每次的各条路径长度
  %disp(distl);
  %%%%%%%%%%%%%修改各路径上的信息量%%%%%%%%%%%%%%%%
  for l=1:ng  
   %disp(yita);
   for k=1:m                    %第k只蚂蚁      
        %disp(tabu); 
        for j=1:ng
          yita(j,NC)=1/distl(j,NC);     %分别计算每条路径的启发函数,表示第NC次选择路径l的期望程度 
          psum_medium0(j,NC)=(tao(j,NC)^alpha).*(yita(j,NC)^beta); 
        end;
        disp(psum_medium0);  
        %disp(psum(l,NC));  
        for j=1:ng
            psum(j,NC)=sum(psum_medium0(:,NC)); %???没有明白  
            p(j,NC)=(tao(j,NC)^alpha).*(yita(j,NC)^beta)/psum(j,NC); %2.4.3公式                      
        end;
        %disp(p); 
        maxp(1,NC)=max(p(:,NC));   %找出最大概率的路径                
        %disp(maxp);
        [Linear_index(l,NC)]=find(maxp(1,NC)==p(:,NC));  %在p的第一行中查找和maxp(1)相等元素的索引值        
        plot(C(1,:,Linear_index(l,NC),NC),C(2,:,Linear_index(l,NC),NC),'-r');
        %disp(C);
        %%%%%%%%%%%%%%%%%%%NC次选择的最优路径%%%%%%%%%%%%%%%%%%%%%%%%
        solution_NC(l,NC)=Linear_index(l,NC);
        %disp(solution_NC);               
        tabu(k,solution_NC(l,NC),NC)=1;  %每次循环的修改禁忌索引表
        %disp(tabu);             
    end;   
    %%%%%%%%%%%%%% In this phase global updating occurs%%%%%%%%%%%%%%%
    %更新第条路径上的信息量    
    for j=1:ng  
        if j==Linear_index(l,NC)
           detatao(l,NC)=Q/distl(l,NC); 
        else 
           detatao(j,NC)=0;
        end;
        tao(j,NC)=(1-rou).*tao(j,NC)+rou.*detatao(j,NC);  
    end;  
        disp(tao); %m只蚂蚁走完一遍                    
  end;             %所有的蚂蚁都走完
end;
disp(distl);
lmin=min(distl(:,:));
disp(lmin);
[solution_lindex]=find(min(lmin)==lmin);
disp(solution_lindex);
[solution_NCindex]=find(min(distl(:,solution_lindex))==distl(:,solution_lindex));
disp(solution_NCindex);
%plot(C(1,:,solution_lindex,solution_NCindex),C(2,:,solution_lindex,solution_NCindex),'-y','LineWidth',4);
%%%%%%%%%%%%%%%%%%%%%%效率和能耗输出%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
v=5;      %v=5m/s
pow=0.7;    %power=700mAh
for NC=1:NCmax    
   for l=1:ng
    rtime(l,NC)=distl(l,NC)/v;
    otime(NC)=distl(solution_NCindex,NC)/v;
    rpowsum(l,NC)=pow*rtime(l,NC);
    opowsum(NC)=pow*otime(NC);   
   end
end

    h1=figure(1);
    %subplot(2,1,1);   
    x1=1:NC;
    y1=rtime(:,:);    
    plot(x1,y1,'-g');
    hold on;
    x2=1:NC;
    y2=otime(:);    
    plot(x2,y2,'-*c');     
    xlabel('Cycle No.');
    ylabel('Time(s)');
    
    h2=figure(2);
    %subplot(2,1,2);    
    x3=1:NC;
    y3=rpowsum(:,:);
    plot(x3,y3,'-r');
    hold on;
    x4=1:NC;
    y4=opowsum(:);
    plot(x4,y4,'-*b'); 
    xlabel('Cycle No.');
    ylabel('Power consume(mAh)');





% --- Executes on button press in stopbutton.
function stopbutton_Callback(hObject, eventdata, handles)
% hObject    handle to stopbutton (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user state (see GUIDATA)
clf;


% --- Executes on button press in loadmap.
function loadmap_Callback(hObject, eventdata, handles)
% hObject    handle to loadmap (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

%load e-map
%mapshow boston.tif
%axis image off

% The orthophoto is in survey feet, the roads are in meters;
% convert the road units to feet before overlaying them.
%S = shaperead('boston_roads.shp');
%5surveyFeetPerMeter = unitsratio('sf','meter');
%x = surveyFeetPerMeter * [S.X];
%y = surveyFeetPerMeter * [S.Y]; 
%mapshow(x,y)
%%%%%%%%%%%%%%% creat map %%%%%%%%%%%%%%%%
%边界初始点与目标点
%初始化障碍物
obstacle=zeros(600,600);
xi=[10];
yi=[10];
plot(xi,yi,'ob');
hold on;
xo=[565];
yo=[565];
plot(xo,yo,'ob');
hold on;

x1=[50];
y1=[500];
%obtxt1(1,:)=x1;
%obtxt1(2,:)=y1;
for i=1:50
    for j=1:30    
     obstacle(x1+i,y1+j)=1;
     %plot(x1+i,y1+j,'.m');
    end
end

%测试物体
x2=[220];
y2=[420];
%obtxt1(1,:)=x1;
%obtxt1(2,:)=y1;
for i=1:80
    for j=1:80
     obstacle(x2+i,y2+j)=1;     
     %plot(x2+i,y2+j,'.m');
    end
end

x3=[300];
y3=[200];
%obtxt3(1,:)=x1;
%obtxt3(2,:)=y1;
for i=1:100
    for j=1:80
     obstacle(x3+i,y3+j)=1;
     %plot(x3+i,y3+j,'.m');
    end
end

x4=[350];
y4=[280];
%obtxt3(1,:)=x1;
%obtxt3(2,:)=y1;
for i=1:50
    for j=1:80
     obstacle(x4+i,y4+j)=1;
     %plot(x4+i,y4+j,'.m');
    end
end

%边界
x5=[100];
y5=[120];
%obtxt4(1,:)=x1;
%obtxt4(2,:)=y1;
for i=1:60
    for j=1:50
     obstacle(x5+i,y5+j)=1;     
     %plot(x5+i,y5+j,'.m');
    end
end

x6=[500];
y6=[100];
%obtxt5(1,:)=x1;
%obtxt5(2,:)=y1;
for i=1:50
    for j=1:70
     obstacle(x6+i,y6+j)=1;     
     %plot(x6+i,y6+j,'.m');
    end
end
hold on;
grid on;
for i=1:600
    for j=1:600
        if obstacle(i,j)==1
           plot(i,j,'.m');
        end
    end
end

%disp(obstacle);
%%%%地图信息数据%%%%%
set(handles.obtext1,'UserData',[obstacle]); %把产生的g放入state文本中
disp('map ok!');

%%%%%%%%%%%%%%%%%    动态部署节点算法仿真V1.0    %%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%　　2008.3.27　　%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 500＊500m无障碍物区中进行移动节点部署
% 初始化区域为50＊50m，然后，经计算展开到300＊300m场
% 节点传感半径25m；通讯半径50m
% 节点能感知GPS信息
lon=10;          %GPS经度基点(需为整数）
lan=10;          %GPS纬度基点(需为整数）
n=500;              % the square of mobile  
x=rand(1,n);        % the X position of mobile nodes 
y=rand(1,n);        % the Y position of mobile nodes 
rs=25;              % the sensor radius
rc=50;              % the communication radius
rd=2;               % the detection obstacle radius
robs=1;                % the voiding obstacle R
rspan=50;              % the span between two nodes
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% initialization position of node
for i=1:0.1*n       %0.1n=30个移动节点
    x0(i)=lon+0.1*n*x(i);  %初始化在最小方格50＊50中
    y0(i)=lan+0.1*n*y(i);  
end
%　显示初始化节点位置，随机分布
for i=1:0.1*n
    plot(x0(i),y0(i),'s','LineWidth',2,'MarkerEdgeColor','k','MarkerFaceColor','g','MarkerSize',10);       
    hold on;
end
grid on;
%  初始化时的通讯链路
for i=1:0.1*n
    for j=2:0.1*n
     if i~=j
      dij=sqrt((x0(i)-x0(j))^2+(y0(i)-y0(j))^2);       %i节点和j节点之间的距离
      if (dij<=rc)                                     %小于通讯距离则可以通讯，存在通讯链路
         x=[x0(i),x0(j)];
         y=[y0(i),y0(j)];       
         line(x,y);                                    %短线表示通讯链路
         hold on;
      end
     end
     %text(x0(i)+0.5,y0(i),num2str(i),'color','r');     %循环标注节点号,调试时关闭
    end
end
xlabel('X(m)');
ylabel('Y(m)');