function fig2

% These are graphic display commands
% The computation engine is in the functions opDyn and boundConf

subplot(1,2,1);
t = 20;
n= 100;
e = 17;
X = [0:n]'; X(:,2)=1;
Y = opDyn(X,e+10^-10,t); % e+10^-10 must be chosen due to numerical vulnerablity 
                         % of abs(X-Xnew)<=e when abs(X-Xnew)==e
c=colormap(hsv(n));
box on
set(gca,'Layer','top','nextplot','replacechildren','ColorOrder',c,'Position',[0.1 0.145 0.7/5 0.77])
axis([0 t 0 n])
plot(0:t,Y');
xlabel('t')
ylabel('agents')
title('101 agents, \epsilon = 17')

subplot(1,2,2);
t = 100;
n= 400;
e = 68;
X = [0:n]'; X(:,2)=1;
Y = opDyn(X,e+10^-10,t);
c=colormap(hsv(n));
box on
set(gca,'Layer','top','nextplot','replacechildren','ColorOrder',c,'Position',[0.2+0.7/5 0.145 4*0.7/5 0.77])
axis([0 t 0 n])
plot(0:t,Y');
xlabel('t')
ylabel('agents')
title('401 agents, \epsilon = 68')
rectangle('Position',[1 100 10 120]);
text(14,165,'\rightarrow')

axes('Position',[0.45 0.155 0.2 0.35])
box on
set(gca,'Layer','top','nextplot','replacechildren','ColorOrder',c,'XTick',[],'YTick',[])
plot(0:t,Y')
axis([1 10 100 220])


set(gcf,'Position',[105   251 600 280],'PaperPositionMode','auto')
print -dtiff -r96 fig1




function Y = opDyn(X,e,t)

% Y = opDyn(X,e)
% X    : Matrix n x 2, 
%        first column : n opinions, 
%        second column : number of agents with this opinion
% e    : bound of confidence
% t    : time step until we compute the dynamic
%
% Y(:,k) is the opinion profile of time step k

Y = X(:,1);

for i=1:t
    Xnew = boundConf(X,e);
    X = reduceOp(Xnew); % to spare computation time
	Y = [Y allOp(X)];
end


function X = boundConf(X,e)

% X = boundConf(X,e)
% X    : Matrix n x 2, 
%        first column : n opinions, 
%        second column : number of agents with this opinion
% e    : bound of confidence
%
% boundConf calculates one step in the HK process of OD
% 
% 02/2003 Jan Lorenz

n = size(X,1); % number of opinions

A=[]; % will get the confidence matrix

% produce rows of the confidence matrix
for i = 1:n
   d = abs(X(:,1)-X(i,1));
   a = [d <= e ]'; % a is confidence set of opinion i
   a = a.*X(:,2)'; % a weighted with number of agents of opinions
   a = a/sum(a); 
   A = [A; a]; 
end
X(:,1)=A*X(:,1);


function X = reduceOp(X)

% X = reduceOp(X)
% X    : Matrix n x 2, 
%        first column : n opinions, 
%        second column : number of agents with this opinion
%
% reduceOp collects all opinions with the same value in one row
% and adds the numbers of agents with this opinion
%
% 09/2002 Jan Lorenz

k=1; % row counter

while size(X,1) > k
   I = ( find(X(:,1) == X(k,1)) ); % row indices equal to X(k,1)
   X(k,2) = sum(X(I,2));
   X = X( sort([k setdiff([1:size(X,1)],I)]) ,:); % delete rows
   k=k+1; 
end




function X = allOp(X)

% X = allOp(X)
%
% X (in): Matrix n x 2, 
%        first column : n opinions, 
%        second column : number of agents with this opinion
% X (out): column Vector
%
% allOp blows up an opinion profile with collected equal opinions
% to a vector with single opinions, for plotting reasons

M = [];
for i = 1:size(X,1)
    m = zeros(X(i,2),size(X,1));
    m(:,i) = ones(X(i,2),1);
    M = [M;m];
end
X = M*X;
X = X(:,1);
    