一、理论基础

二、核心程序

三、仿真测试结果

作者ID ：fpga和matlab
CSDN主页：https://blog.csdn.net/ccsss22?type=blog
擅长技术：
1.无线基带,无线图传,编解码
2.机器视觉,图像处理,三维重建
3.人工智能,深度学习
4.智能控制,智能优化
5.其他

一、理论基础

随着计算机和通信系统的发展，大规模信息系统的协调问题成为新兴的研究热点。基于代理的协调算法，可以应用于竞争系统和合作系统。在竞争的环境中，每个自私的代理的目标是使自己的利益最大化。其中，关键问题是如何在保证个体利益的同时，使系统的利益最大化。

随着微型多功能传感器网络的快速发展，无线传感器网络的规模得到了极大的发展，在一个无线传感器网络中，每个网络节点都可以通过其无线传感器以及通信信道进行通信和信号的传输。传统的算法已经无法满足需求，目前应用较多的算法——和积算法(Sum-Product Algorithm，SPA)就是信息传递算法中比较高效的一种。

SPA是Judea Pearl提出的基于多环可信传播模型的立体匹配模型，解决图形模型中的边际分析，可以广泛应用到组合优化、信息编码、人工智能和计算机视觉等多个领域。有关该算法的理论基础研究，文献给出了有环因子图的收敛的充分条件。Rogers等人用SPA协调无线传感器网络中代理的工作。他的工作显示了SPA在大规模系统协调问题中的应用潜力。然而，对于大规模协调系统的应用，算法的收敛速度、资源消耗和鲁棒性是必须考虑的3个重要问题。本文在该算法中引入混沌机制对节点信息进行优化，并且优化了算法的信息传递模型，不仅显著提高了收敛速度，而且大大降低了算法的资源消耗，改善了算法的鲁棒性。

使用MAX-SUM算法在无线传感器网络中使系统的应用获得最佳的性能，一个无线传感器网络的主要结构如下所示：

图1无线传感器网络的基本结构

从图1可以看到，无线传感器网络的基本组成可以分为传感器，接口，控制器，系统驱动等几个组成。 MAX-SUM算法是一个基于因子图的算法，整个算法的节点表示其中的变量和约束。每个节点表示各个系统的变量或者是函数节点。每个变量节点连接着所有的函数节点，同理，每个函数节点连接着所有的变量节点。在MAX-SUM算法中，不同的节点表示的不同的角色。

图2 DCOP变换的因子图

这里假设节点连接不同的变量节点和函数节点，每个节点均包含对应的节点约束表示其对应的作用。其中，变量节点和函数节点在MAX-SUM算法中均表示代理节点。通过这些节点可以发送消息，读取消息和性能计算。图2.1所示的是三个节点的情况，通过DCOP变化的因子图的变化效果。

通过总结，整个MAX-SUM算法的基本流程如下所示：

 其中第5和第6步骤就是执行从变量节点到函数节点的消息传递过程，第7和第8步骤就是执行函数节点到变量节点的消息传递过程。

二、核心程序

clc;
clear;
% close all;
warning off;

%每次产生的随机数固定
s = RandStream.create('mt19937ar','seed',1);
RandStream.setDefaultStream(s);
%首先顶一个无线传感器网络环境
noOfNodes = 30;%网络节点个数,仿真的时候，网络节点不易过多，会out of memory
L = 1000;
R = 300; %相邻规则，小于R，表示相邻
MAX_INTER = 50;%最大迭代次数
coff = 0.9;%收敛系数，表示信号的接收能量
rand('state', 3); %节点状态初始化

figure;
clf;grid on;
hold on;
netXloc = rand(1,noOfNodes)*L;
netYloc = rand(1,noOfNodes)*L;
for i = 1:noOfNodes
plot(netXloc(i), netYloc(i), 'r.');
text(netXloc(i)+8, netYloc(i)+8, num2str(i));
for j = 1:noOfNodes
distance = sqrt((netXloc(i) - netXloc(j))^2 + (netYloc(i) - netYloc(j))^2);
%标志相邻节点，使用1表示
if distance <= R
matrix(i, j) = 1; % there is a link;
line([netXloc(i) netXloc(j)], [netYloc(i) netYloc(j)], 'LineStyle', ':');
else
matrix(i, j) = 0;
end
end
end
%matrix为节点连接图，如果matrix(i,j)=1，说明节点i和节点j相邻，否则不相邻
%初始化节点信息
for i=1:noOfNodes
for j=1:noOfNodes
if matrix(i, j) == 1;
q(i,j,1:2) = [rand(1,1),rand(1,1)];
r(i,j,1:2) = [rand(1,1),rand(1,1)];
else
q(i,j,1:2) = 0;
r(i,j,1:2) = 0;
end
end
end
%初始化结点的效用值
for i=1:noOfNodes
F(i,1:2) = [rand(1,1),rand(1,1)];
end
tmp_r = r;
tmp_q = q;
%利用消息传递规则进行节点信息的更新，知道节点值收敛
z = zeros(noOfNodes,2,MAX_INTER);%该变量保存每个节点的边界值
States = zeros(noOfNodes,1);%该变量保存每个节点的状态值
tmp1 = 0;
tmp2 = 0;
tmp3 = 0;
tmp4 = 0;
U0 = zeros(noOfNodes,noOfNodes);
U1 = zeros(noOfNodes,noOfNodes);
U2 = zeros(noOfNodes,noOfNodes);
U3 = zeros(noOfNodes,noOfNodes);
for ii = 1:MAX_INTER
if ii>1
load data.mat
end
%计算每个节点的边界值
for nodes = 1:noOfNodes
z(nodes,1,ii) = sum(r(:,nodes,1));
z(nodes,2,ii) = sum(r(:,nodes,2));
if z(nodes,1,ii) >= z(nodes,2,ii)
States(nodes) = 1;
else
States(nodes) = 0;
end
end
%对各个节点的边界值的更新，更新规则按照MAX-SUM算法的消息传递规则进行
%A:变化->函数更新，qi-j
for i = 1:noOfNodes
for j = 1:noOfNodes
if matrix(i, j) == 1
if i==j
tmp1 = sum(r(:,i,1)) - r(j,i,1);
tmp2 = sum(r(:,i,2)) - r(j,i,2);
q(i,j,1) = tmp1 - (tmp1+tmp2)/2;%使和满足0
q(i,j,2) = tmp2 - (tmp1+tmp2)/2;%使和满足0
else
q(i,j,1) = coff*sqrt(2)*(q(i,i,1)/(sum(matrix(i,:))^2) + tmp_q(i,j,1));%使和满足0
q(i,j,2) = coff*sqrt(2)*(q(i,i,2)/(sum(matrix(i,:))^2) + tmp_q(i,j,2));%使和满足0
end
end
end
end
%B:函数->变量更新，rj-i
%1）计算对应的效用值Uij
for nodes = 1:noOfNodes
for j = 1:noOfNodes
if (matrix(nodes, j) == 1) & (j ~= nodes)
U0(nodes,j) = F(nodes,1) - 1 + xor(0,0);%四个状态
U1(nodes,j) = F(nodes,2) - 1 + xor(0,1);%四个状态
U2(nodes,j) = F(nodes,2) - 1 + xor(1,0);%四个状态
U3(nodes,j) = F(nodes,1) - 1 + xor(1,1);%四个状态
end
end
end
%2）函数->变量更新
for j = 1:noOfNodes
for i = 1:noOfNodes
if matrix(j, i) == 1
if i==j
tmp1 = 0.5*log((2*pi*exp(1))^9) * (U0(j,i) + sum(q(j,:,1)) - q(j,i,1));%统计输入的r值
tmp2 = 0.5*log((2*pi*exp(1))^9) * (U1(j,i) + sum(q(j,:,1)) - q(j,i,1));%统计输入的r值
tmp3 = 0.5*log((2*pi*exp(1))^9) * (U2(j,i) + sum(q(j,:,1)) - q(j,i,1));%统计输入的r值
tmp4 = 0.5*log((2*pi*exp(1))^9) * (U3(j,i) + sum(q(j,:,1)) - q(j,i,1));%统计输入的r值
r(j,i,1) = max([tmp1 tmp2 tmp3 tmp4]);
tmp1 = 0.5*log((2*pi*exp(1))^9) * (U0(j,i) + sum(q(j,:,2)) - q(j,i,2));%统计输入的r值
tmp2 = 0.5*log((2*pi*exp(1))^9) * (U1(j,i) + sum(q(j,:,2)) - q(j,i,2));%统计输入的r值
tmp3 = 0.5*log((2*pi*exp(1))^9) * (U2(j,i) + sum(q(j,:,2)) - q(j,i,2));%统计输入的r值
tmp4 = 0.5*log((2*pi*exp(1))^9) * (U3(j,i) + sum(q(j,:,2)) - q(j,i,2));%统计输入的r值
r(j,i,2) = max([tmp1 tmp2 tmp3 tmp4]);
else
r(j,i,1) = coff*sqrt(2)*(r(j,j,1)/(sum(matrix(j,:))^2) + tmp_r(i,j,1));
r(j,i,2) = coff*sqrt(2)*(r(j,j,2)/(sum(matrix(j,:))^2) + tmp_r(i,j,2));
end
end
end
end
save data.mat r q
end


%画出节点收敛曲线图，节点更新后的节点矩阵
%画出每个节点的值收敛情况，取其中三个节点作为测试节点
figure
for hh = 1:MAX_INTER
plot_z11(hh) = z(1,1,hh);
plot_z21(hh) = z(1,2,hh);

plot_z12(hh) = z(15,1,hh);
plot_z22(hh) = z(15,2,hh);

plot_z13(hh) = z(30,1,hh);
plot_z23(hh) = z(30,2,hh);
end
plot(plot_z11,'b','LineWidth',2);hold on
plot(plot_z21,'b--','LineWidth',2);hold on
plot(plot_z12,'k','LineWidth',2);hold on
plot(plot_z22,'k--','LineWidth',2);hold on
plot(plot_z13,'r','LineWidth',2);hold on
plot(plot_z23,'r--','LineWidth',2);hold on
legend('节点1','','节点2','','节点3','');
grid on
ylabel('节点边界值收敛');
xlabel('迭代次数');

三、仿真测试结果

  A16-10