clear;clear;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%PARAMETERS%%%%%%%%%%%%%%%%%%%%%%%%%%%%w=4;%FieldDimensions-xandymaximum(inmeters)100*100的区域xm=100;ym=100;%xandyCoordinatesoftheSinksink节点的坐标sink.x=100;sink.y=200;%NumberofNodesinthefield总共200个节点n=50*w;%OptimalElectionProbabilityofanode%tobecomeclusterhead簇节点选择概率p=0.05;packetLength=4000;%数据包长度ctrPacketLength=100;%控制包长度%EnergyModel(allvaluesinJoules)%InitialEnergy初始能量Eo=0.5;%Eelec=Etx=Erxenode0=0;ETX=50*0.000000001;%50nj/bitERX=50*0.000000001;%50nj/bit%TransmitAmplifiertypesEfs=10*0.000000000001;%10pj/bit/m2Emp=0.0013*0.000000000001;%0.0013pj/bit/m4%DataAggregationEnergyEDA=5*0.000000001;%5nj/bit/signal%ValuesforHetereogeneity%Percentageofnodesthanareadvanced%m=1;%\alpha%a=1;INFINITY=999999999999999;%maximumnumberofrounds最大轮询次数rmax=3000;node_1=0;node_2=0;node_3=0;node_1_x=0;node_2_x=0;node_3_x=0;node_1_z=0;node_2_z=0;node_3_z=0;%%%%%%%%%%%%%%%%%%%%%%%%%ENDOFPARAMETERS%%%%%%%%%%%%%%%%%%%%%%%%%Computationofdo一个距离的计算根据能量公式参考wendi的phd论文do=sqrt(Efs/Emp);%取平方根%CreationoftherandomSensorNetwork%figure(1);fori=1:1:n%n=200200个节点数设置节点位置及初始能量和状态S(i).xd=rand(1,1)*xm;%坐标x,rand(1,1)产生一个0,1之间的随机数s(i)表示的是数组XR(i)=S(i).xd;S(i).yd=rand(1,1)*ym;YR(i)=S(i).yd;%将随机产生的点S(i).G=0;%初始化的时候,没有簇头产生的个数集合%initiallytherearenoclusterheadsonlynodesS(i).type='N';%普通节点%temp_rnd0=i;%RandomElectionofNormalNodes%if(temp_rnd0=m*n+1)S(i).E=Eo;%每个节点的能量保持为0.05该数组表示的是能量的存贮%S(i).ENERGY=0;%初始消耗的能量???未知ing%end%RandomElectionofAdvancedNodespriority(i)=randint(1,1,[1,3]);figure(1);plot(S(i).xd,S(i).yd,'+');holdon;endnode_1_x=node_1;node_2_x=node_2;node_3_x=node_3;S(n+1).xd=sink.x;%基站位置x坐标S(n+1).yd=sink.y;%基站位置y坐标%FirstIteration%figure(1);%counterforCHstotal簇头个数countCHs=0;%counterforCHsperround每轮的簇头个数rcountCHs=0;cluster=1;rcountCHs=rcountCHs+countCHs;flag_first_dead_old=0;flag_first_dead=0;node_p1=0;node_p2=0;node_p3=0;node_p1_x=0;node_p2_x=0;node_p3_x=0;forr=0:1:rmax%主循环,每次1轮r;ENOLD(r+1)=0;enode0=0;fori=1:1:n;ENOLD(r+1)=S(i).E+ENOLD(r+1);%第r轮所有节点的总能量endEAOLD(r+1)=ENOLD(r+1)/n;%这个是第r轮的平均能量%Operationforepochif(mod(r,round(1/p))==0)%20的整数倍轮数时fori=1:1:nS(i).G=0;%还未当选簇头的节点集合,因为是要判断前1/p轮中是否有节点时簇节点S(i).cl=0;%plot(S(i).xd,S(i).yd,'*');%这是绘制节点位置endendholdoff;%Numberofdeadnodes每一轮开始的时候将该轮的统计参数清零dead=0;dead_old=0;%NumberofdeadAdvancedNodes%dead_a=0;%NumberofdeadNormalNodesdead_n=0;dead_n_old=0;node_p1=0;node_p2=0;node_p3=0;%counterforbittransmittedtoBasesStationandtoClusterHeads%counterforbittransmittedtoBasesStationandtoClusterHeads%perround%figure(1);fori=1:1:n%checkingifthereisadeadnodeifS(i).E0S(i).type='N';%如果节点能量大于0,则类型为正常endend%countCHs=0;cluster=1;fori=1:1:nif(S(i).E0)temp_rand=rand;if((S(i).G)=0)%如果该节点在候选集合中保证钱1/p轮没该节点不是簇节点%ElectionofClusterHeadsif(temp_rand=(p/(1-p*mod(r,round(1/p)))))%判断一个节点的概率是否小于阈值countCHs=countCHs+1;%簇节点数加一S(i).type='C';%该节点的类型变为簇节点S(i).G=round(1/p)-1;%簇头节点的位置C(cluster).xd=S(i).xd;C(cluster).yd=S(i).yd;%plot(S(i).xd,S(i).yd,'k*');distance=sqrt((S(i).xd-(S(n+1).xd))^2+(S(i).yd-(S(n+1).yd))^2);%到sink的距离,因为是簇节点所以直接到sink节点C(cluster).distance=distance;C(cluster).id=i;%表示第i个节点时簇节点X(cluster)=S(i).xd;Y(cluster)=S(i).yd;cluster=cluster+1;%广播自成为簇头distanceBroad=50;%sqrt((S(i).xd)^2+(S(i).yd)^2);???????????????广播包的距离是50%%%需要判断所剩的能量是否足以发送if(S(i).E=(ETX*ctrPacketLength+Emp*ctrPacketLength*(distanceBroad*distanceBroad*distanceBroad*distanceBroad))||S(i).E=(ETX*ctrPacketLength+Emp*ctrPacketLength*(distanceBroad*distanceBroad)))if(distanceBroaddo)%大于极限距离就是4次方,小于则是2次方发送控制比特也就是广播包,广播给其他节点,不是给sinkS(i).E=S(i).E-(ETX*ctrPacketLength+Emp*ctrPacketLength*(distanceBroad*distanceBroad*distanceBroad*distanceBroad));%广播自成为簇头elseS(i).E=S(i).E-(ETX*ctrPacketLength+Efs*ctrPacketLength*(distanceBroad*distanceBroad));endelseS(i).E=0;endif(S(i).E=((ETX+EDA)*(packetLength)+Emp*packetLength*(distance*distance*distance*distance))||S(i).E=((ETX+EDA)*(packetLength)+Emp*packetLength*(distance*distance)))%CalculationofEnergydissipated簇头自己发送数据包能量消耗distance;%这是到sink节点的距离if(distancedo)%发送数据到sink节点所剩的能量S(i).E=S(i).E-((ETX+EDA)*(packetLength)+Emp*packetLength*(distance*distance*distance*distance));elseS(i).E=S(i).E-((ETX+EDA)*(packetLength)+Efs*packetLength*(distance*distance));endelseS(i).E=0;endendendendendSTATISTICSOLD(r+1).CLUSTERHEADS=cluster-1;%统计第r轮簇头数目,r是从0开始的,所以加1;cluster最后要-1,是应为上面的循环多加了1CLUSTERHSOLD(r+1)=cluster-1;%ElectionofAssociatedClusterHeadforNormalNodesfori=1:1:nif(S(i).type=='N'&&S(i).E0)%普通节点min_dis=sqrt((S(i).xd-S(n+1).xd)^2+(S(i).yd-S(n+1).yd)^2);%默认距离是到sink的距离%min_dis=INFINITY;if(cluster-1=1)%如果有簇头存在min_dis_cluster=1;%加入最近的簇头forc=1:1:cluster-1%簇头数量一共是cluster-1%temp=min(min_dis,sqrt((S(i).xd-C(c).xd)^2+(S(i).yd-C(c).yd)^2));temp=sqrt((S(i).xd-C(c).xd)^2+(S(i).yd-C(c).yd)^2);%判断到已有簇节点的距离if(tempmin_dis)min_dis=temp;min_dis_cluster=c;%标记簇节点标号(节点所属的簇节点)通过循环找到到簇节点距离最小的簇节点end%接收簇头发来的广播的消耗S(i).E=S(i).E-ETX*ctrPacketLength;endenergy_1=ETX*(ctrPacketLength)+Emp*ctrPacketLength*(min_dis*min_dis*min_dis*min_dis)+ETX*(packetLength)+Emp*packetLength*(min_dis*min_dis*min_dis*min_dis);energy_2=ETX*(ctrPacketLength)+Emp*