function ks=kshell(net)
%net records the undirected links of the network. Each link is only
%recorded once. The node index is from 1 to n.
n=max(max(net));
deg=zeros(n,1);
L=length(net);
for i=1:L
    deg(net(i,1))=deg(net(i,1))+1;
    deg(net(i,2))=deg(net(i,2))+1;
end

nlsize=sum(deg);
nl=zeros(nlsize,1);
part=ones(n+1,1);
for i=2:n+1
    part(i)=part(i-1)+deg(i-1);
end
pos=part;
for i=1:L
    node1=net(i,1);
    node2=net(i,2);
    nl(pos(node1))=node2;
    pos(node1)=pos(node1)+1;
    nl(pos(node2))=node1;
    pos(node2)=pos(node2)+1;
end


md=max(deg);
bin=zeros(md+1,1);
for v=1:n
    bin(deg(v)+1)=bin(deg(v)+1)+1;
end
start=1;
for d=0:md
    num=bin(d+1);
    bin(d+1)=start;
    start=start+num;
end
pos=zeros(n,1);
vert=zeros(n,1);
for v=1:n
    pos(v)=bin(deg(v)+1);
    vert(pos(v))=v;
    bin(deg(v)+1)=bin(deg(v)+1)+1;
end
for d=md:-1:1
    bin(d+1)=bin(d);
end
bin(1)=1;
for i=1:n
    v=vert(i);
    neighbor=nl(part(v):part(v+1)-1);
    for j=1:part(v+1)-part(v)
        u=neighbor(j);
        if pos(u)>i
            if deg(u)>deg(v)
                du=deg(u);
                pu=pos(u);
                pw=bin(du+1);
                w=vert(pw);
                if u~=w
                    pos(u)=pw;
                    vert(pu)=w;
                    pos(w)=pu;
                    vert(pw)=u;
                end
                bin(du+1)=bin(du+1)+1;
                deg(u)=deg(u)-1;
            end
        end
    end
end
ks=deg;

