simu_for_pdmin.asv

无线ADHOC网络中

%**************************************************************************
%**  File name : interference_con.m
%**  Author  : Liufengwei (English name:andy)
%**  Creation date : 2006,6,13
%**  Modified date: 4,22,2007,Sunday
%**  Parameter explanation : See the note in the program
%**************************************************************************
% ** This simulation is done to investigate the probability of minimun network degree greater than n0
%      analyzed in paper<On the minimum node degree and connectivity of a wireless multihop network>
% ** Model: Distribute n=500 nodes in 1000*1000m^2 accoring to uniform distribution;
%      there is an edge between nodes if their distance is less than r0;We calculate the degree of every nodes
% ** Goal: simulate the formula in paper<On the minimum node degree and connectivity of a wireless multihop network>
% ** Explanation :    
%** Program process: initialize simulation parameter, then 
%************************************************************************

%function 
%***********first initialize simulaiton parameters***************
tic 
run_times = 500 ; % simulation times
Nodenum = 500 ; % number of nodes
square_edge = 1000 ; % the square edge length
N0 = [1 2 3] ; % minimum degree 
Len_N0 = length(N0) ; % length of N0
Radio_down = 10 ; % the minimum of radio range
Radio_step = 10 ; % increment step
Radio_up = 150 ; % the maximum of radio range
Radio = [Radio_down: Radio_step: Radio_up]; % radio range
Len = length(Radio) ; % number of radio range

%*****************first initialize simulaiton matrix ***********************

NodesXY = zeros(Nodenum, 2) ;% store the location of nodes 
Nodes_dis = zeros(1, (Nodenum-1)); % store the node distance
Neighbor = zeros((Nodenum-1), Len) ; % store the node distance
Ndeg = zeros(Nodenum, Len) ; % record degree of every node
G_than_N0 = zeros(Len_N0, Len); % record the number of degree of every node greater than N0
simu_Pdmin = zeros(Len_N0, Len); % the probability of minimun network degree greater than N0
ana_Pdmin = zeros(Len_N0, Len);

for run =1 : run_times  
    run
    
   % rand('state', sum(100*clock))   
    NodesXY = square_edge*(rand(Nodenum, 2)) ; % generate coordinates of node
    k = ones(Len,1) ; % index for links for every radio range
    radio_index =1 ;
    
   % ****************** calculate the degree of every node *****************
   %for radio = Radio_down:Radio_step:Radio_up
    for n =1:Nodenum
         %n
         if n==1
            branch_0 = NodesXY(2:Nodenum,:)-ones((Nodenum-1), 1)*NodesXY(1,:);
            Nodes_dis = sqrt(sum(branch_0'.^2)) ;
         else
             branch_1 = NodesXY(1:(n-1),:)-ones((n-1), 1)*NodesXY(n,:) ;% branch_1:(n-1)*2
             branch_2 = NodesXY((n+1):Nodenum,:)-ones((Nodenum-n), 1)*NodesXY(n, :); %branch_2 : (Nodenum-n)*2
             Nodes_dis = sqrt(sum([branch_1; branch_2]'.^2)) ; %  Nodes_dis: 1*(Nodenum-1)       
         end %  end if n==1   
         %**alternate method to calculate distance*********
         Nodes_dis = NodesXY - ones(Nodenum, 1)*Nodes
         Neighbor = Nodes_dis'*ones(1, Len) <= ones((Nodenum-1), 1)*Radio ; % store the nodes connecting node i :(Nodenum-1)*len  
         Ndeg(n , :) = sum(Neighbor) ; % the degree of node  i--Ndeg: Nodenum*len

    end % end for  i =1:(Nodenum-1)
        
   %*******Calculate the prob of degree of every node greater than N0*********
   for i = 1: Len_N0
       
        G_than_N0(i, :) = (min(Ndeg)>=N0(i))+ G_than_N0(i, :) ; % G_than_N0: Len*Len_N0,  min(Ndeg):1*len
   end
    %     radio_index = radio_index+1 ;
   %end % for radio   
end % for run=1 : run_times

simu_Pdmin = G_than_N0/run_times ; %calculate the prob of degree of every node greater than N0
% Pdmin_N=prod(G_than_N0);

%**********Calculate the analytical prob********************

A = square_edge^2 ; % the area of system plane
for i = 1: Len_N0
   j = 1;
   for radio = Radio_down: Radio_step: Radio_up
         for N=0 : N0(i)-1 
            ana_Pdmin(i,j)=((Nodenum/A)*pi*radio^2)^N/factorial(N)+ana_Pdmin(i,j);
         end
            ana_Pdmin(i,j)=(1-ana_Pdmin(i,j)*exp(-(Nodenum/A)*pi*radio^2))^Nodenum;
         j = j+1;
   end
end
%*******************draw figure************************

figure(1)
plot(Radio,simu_Pdmin(1,:),'ro',Radio,ana_Pdmin(1,:),'k-',Radio,simu_Pdmin(2,:),'r*',Radio,ana_Pdmin(2,:),'k--',Radio,simu_Pdmin(3,:),'rsquare',Radio,ana_Pdmin(3,:),'k-.');
grid on;
legend('\it simP(d_min>=1)','\ita analytical P(d_min>=1)','\it simP(d_min>=2)','\ita analytical P(d_min>=2)','\it simP(d_min>=3)','\ita analytical P(d_min>=3)'); 
hold on
axis([Radio_down Radio_up 0 1])
xlabel('range(in m)','color','b')
ylabel('prob for d_min','color','b')
title(' P(dmin>=N0 with Euclidian distance metric)')

toc