function G=createlattice(n,k,p)
%CREATELATTICE.m--reoganizes a regular lattice from Strogatz-"Small World" paper
%
%G=createlattice(n,k,p)
%
%Creates a lattice by defining an n-by-n matrix G.  G is upper 
%triangular and sparse G(i,j)=1 if an edge exists between vertices
% i and j. The integer k specifies the number of nearest neighbors
%to each vertex.  p specifies the probability of changing an edge during each
%reworking pass.
%
%From Watts and Strogatz (1998) "Collective dynamics of 'small-world'
%networks", Nature, 393:440:442
%


%1) Create initial G
%
%The n vertices lie on a circle.  k=1 implies each is connected
%only to the vertices immediately to the right and left.  k=2
%adds connections to the vertices 2 to the right and left.  Thus,
%each vertex will be connected to 2*k others.  However, because
%the graph is undirected, we will only represent the upper-triangular
%portion.
%
%
I=(1:n)'*ones(1,k); %the rows of G
P=ones(n,1)*(1:k); 
J=I+P;		%the columns of G are the rows plus 1:k
I=I(:);
J=J(:);
K=find(J>n);  %Some entries in J  are > n. For these edges
				  %we'll represent them going to the left
tmp=mod(J(K),n);
J(K)=I(K);
I(K)=tmp;

G=sparse(I,J,ones(n*k,1), n,n);


P=rand(n*k,1);
K=find(P<p);	%Edges with P<p must be changed
if(~isempty(K))
   m=length(K);
	N=rand(m,1);
   N=N*n+0.5;
   N=round(N); %new edges
	
   I=I(K);
   origJ=J(K);
   J(K)=N;
   for j=1:m
   	p=I(j);
      q=J(j);
      if(p>q)
      	q=p;
         p=J(j);
      end
      if(p~=q)
   		if(~G(p,q)) %new edge doesn't exist
         	G(I(j),origJ(j))=0; %delete old
      		G(p,q)=1;%add new
         end
      end
   end
end