library(igraph)

g = graph(c(2,6,8,2,2,4,8,5,1,2,9,4,9,2,9,5,6,7,1,10,3,8),directed = FALSE)
plot(g) #Fig. 1a

#finding maximum degree
degrees = 1:10
for (v in 1:10) {
degrees[v] = degree(g,v)
}
Delta = max(degrees)
 
#finding JDM
J = matrix(ncol = Delta, nrow = Delta,0)
E = get.edgelist(g)
for (i in 1:dim(E)[1]) {J[degree(g,E[i,1]),degree(g,E[i,2])] = J[degree(g,E[i,1]),degree(g,E[i,2])] +1}
J = J+t(J)
J
 
[,1] [,2] [,3] [,4] [,5]
[1,] 0 2 1 0 0
[2,] 2 0 3 0 3
[3,] 1 3 0 0 2
[4,] 0 0 0 0 0
[5,] 0 3 2 0 0
 
#walking in a fiber
steps = 3
g1 = g
for (i in 1:steps) {
A = get.adjacency(g1)
#move. finding a pair of edges e1 = (v11,v12) and e2 = (v21,v22), such that deg(v11)=deg(v21)
repeat {
quart = sample(1:10,4)
v1 = quart[1];v2=quart[2];v3 = quart[3];v4=quart[4]
if ((A[v1,v2]==1)&&(A[v3,v4]==1)&&(A[v1,v4]==0)&&(A[v3,v2]==0)&&degree(g1,v1)==degree(g1,v3)) {break}
}
Anew = A
Anew[v1,v2] = 0;Anew[v3,v4] = 0;Anew[v1,v4]=1;Anew[v3,v2]=1;Anew[v2,v1] = 0;Anew[v4,v3] = 0;Anew[v4,v1]=1;Anew[v2,v3]=1;
g1 = graph_from_adjacency_matrix(Anew,mode = 'undirected')
}
plot(g1)#Fig. 1b
 


#taking a difference
A = get.adjacency(g)
A1 = get.adjacency(g1)
Adiff= abs(A-A1)
gdiff = graph_from_adjacency_matrix(Adiff,mode = 'undirected')
plot(gdiff) #Fig. 1c