Add elements to tuples or arrays

I didn’t look at the details, but from the discussion above you should use for (_, cyc) in pairs(cycs) here.

Thank you, have you an alternative for using ase with pycall? I’d like to use a julia package to extract the neighborlists from xyz files

I had never heard of ase before, thus no idea.

Not really, the code does not seem to be appending to cycs while traversing it, so I am not sure anymore if this is the reason of the infinite loop.

I did not inspect all the code but:

    for (i_v, j_v) in zip(nl[1], nl[2])
        push!(new[i_v + 1], j_v)
    end

Is probably wrong, because you are changing the identifier of the nodes but not the references to them. This is, you are saying that now node 0 is node 1, node 1 is node 2 and so on, but the j_v is referring to the nodes by their old numbers, starting at 0. I see some +1 and -1 in some variables that may be trying to fix this inconsistency, however, the ideal solution is to either change all indexes to base 1 at the start and then reduce 1 at after the end, or do not change the indexes at all (you are using an OrderedDict not a Vector, you do not need to start at base 1 if this is were you will keep all the references).

You should probably make the .xyz files available to us if you want us to run the code trying to find the problem. Even if most people (I included) will probably not install Python and some Python library to run the code. Instead, a code with the object ase_atoms_no_h hardcoded would be ideal for us to check for problems.

1 Like

The xyz file

168
Properties=species:S:1:pos:R:3 a=T pbc="F F F"
C      -12.39796606       5.75557728      -0.00041200
C      -11.19877446       5.02970593      -0.00032500
C      -12.39394703       7.14159987      -0.00041400
C      -11.19681055       7.84012337      -0.00032400
C       -9.96744844       7.16661422      -0.00023300
C       -9.96819366       5.74383274      -0.00024000
C      -11.18330163       3.57177356      -0.00031300
C       -9.94599750       2.87057938      -0.00023200
C       -8.72064338       3.60828242      -0.00016200
C       -8.72310290       5.02638964      -0.00016200
C       -7.49496079       5.73540164      -0.00008200
C       -7.47119819       7.16548340      -0.00006000
C       -8.69832579       7.88433709      -0.00012800
C       -8.64244256       9.28197414      -0.00008900
C       -7.43982971       9.95371555       0.00002900
C       -6.21247989       9.27715558       0.00010500
C       -6.22187177       7.86472647       0.00003800
C       -4.97636979       7.12169735       0.00007500
C       -5.03553510       5.75021698       0.00006300
C       -6.26161667       5.01942371      -0.00000900
C       -7.48271125       2.90026830      -0.00008500
C       -6.28203360       3.61980199      -0.00001300
C       -7.50012903       1.47302251      -0.00008200
C       -8.65702774       0.73407711      -0.00014500
C       -9.92446254       1.43904343      -0.00021700
C      -11.14174854       0.72262286      -0.00026600
C      -12.34257130       1.44521981      -0.00034900
C      -12.36447512       2.82254806      -0.00037400
C      -11.18344950      -7.85916442      -0.00005400
C       -9.95523438      -7.18356547      -0.00004800
C      -12.38177302      -7.16267824      -0.00009000
C      -12.38815194      -5.77666563      -0.00012300
C      -11.19019759      -5.04875581      -0.00012000
C       -9.95840276      -5.76078798      -0.00008100
C       -8.68488954      -7.89912934      -0.00001100
C       -7.45898928      -7.17818687      -0.00000900
C       -7.48518665      -5.74814774      -0.00003800
C       -8.71453531      -5.04122633      -0.00007300
C       -8.71448956      -3.62311770      -0.00009900
C       -9.94109821      -2.88750101      -0.00014100
C      -11.17720663      -3.59079847      -0.00015600
C      -12.35965442      -2.84358337      -0.00020600
C      -12.34009569      -1.46621958      -0.00024100
C      -11.14050425      -0.74158212      -0.00022600
C       -9.92200007      -1.45593064      -0.00016900
C       -8.65576609      -0.74880892      -0.00014100
C       -7.49761137      -1.48578514      -0.00010700
C       -7.47776502      -2.91299915      -0.00008800
C       -6.25306491      -5.03007198      -0.00003100
C       -6.27586411      -3.63048832      -0.00005500
C       -5.02574116      -5.75877619       0.00000000
C       -4.96423994      -7.13015173       0.00002500
C       -6.20847340      -7.87530143       0.00002300
C       -6.19667499      -9.28771224       0.00005400
C       -7.42286943      -9.96636153       0.00005600
C       -8.62662507      -9.29666966       0.00002400
C        1.21450823     -13.61472667       0.00022800
C        1.24353338     -12.21325957       0.00019700
C        0.01216833     -14.30425732       0.00023100
C       -1.19134122     -13.61676991       0.00020300
C       -1.22274459     -12.21535646       0.00017100
C        0.00979104     -11.50461068       0.00016900
C        2.49840294     -11.47089481       0.00019400
C        2.48700183     -10.04876049       0.00015900
C        1.23545432      -9.35642485       0.00013300
C        0.00856839     -10.06760901       0.00013700
C       -1.21952599      -9.35851261       0.00011100
C       -2.46989304     -10.05297615       0.00011200
C       -2.47887384     -11.47512729       0.00014000
C       -3.71720311     -12.12555146       0.00013700
C       -4.90025658     -11.41993218       0.00010900
C       -4.92801419     -10.01873685       0.00008300
C       -3.70011978      -9.32065305       0.00008300
C       -3.67939065      -7.87050088       0.00005500
C       -2.46207174      -7.23599761       0.00005600
C       -1.21614479      -7.93241823       0.00008200
C        1.22964557      -7.93033819       0.00010200
C        0.00617072      -7.25028844       0.00007900
C        2.47438565      -7.23180047       0.00009300
C        3.69278161      -7.86423016       0.00011300
C        3.71598009      -9.31434484       0.00015200
C        4.94506338     -10.01033595       0.00018300
C        4.91968930     -11.41157691       0.00021900
C        3.73783892     -12.11920941       0.00022400
C       12.39796607      -5.75557729       0.00006900
C       11.19877447      -5.02970594       0.00003000
C       12.39394714      -7.14159888       0.00013300
C       11.19681155      -7.84012350       0.00015900
C        9.96744844      -7.16661423       0.00012300
C        9.96819366      -5.74383275       0.00006200
C       11.18330164      -3.57177357      -0.00004200
C        9.94599650      -2.87057927      -0.00006400
C        8.72064339      -3.60828242      -0.00002400
C        8.72310291      -5.02638964       0.00003200
C        7.49496079      -5.73540164       0.00005800
C        7.47119819      -7.16548341       0.00011100
C        8.69832579      -7.88433709       0.00014500
C        8.64244256      -9.28197414       0.00019700
C        7.43982971      -9.95371556       0.00021100
C        6.21247977      -9.27715658       0.00017500
C        6.22187277      -7.86472660       0.00012800
C        4.97636979      -7.12169735       0.00009600
C        5.03553511      -5.75021699       0.00004900
C        6.26161657      -5.01942471       0.00002800
C        7.48271126      -2.90026831      -0.00004800
C        6.28203249      -3.61980289      -0.00002200
C        7.50012804      -1.47302240      -0.00009600
C        8.65702674      -0.73407700      -0.00012700
C        9.92446254      -1.43904344      -0.00012900
C       11.14174854      -0.72262287      -0.00019800
C       12.34257131      -1.44521982      -0.00017100
C       12.36447512      -2.82254807      -0.00009200
C       11.18345049       7.85916430      -0.00054900
C        9.95523438       7.18356546      -0.00041400
C       12.38177390       7.16267712      -0.00065100
C       12.38815183       5.77666463      -0.00061500
C       11.19019759       5.04875580      -0.00048200
C        9.95840277       5.76078797      -0.00038900
C        8.68488966       7.89913033      -0.00029400
C        7.45898939       7.17818786      -0.00018300
C        7.48518666       5.74814773      -0.00017700
C        8.71453543       5.04122732      -0.00026800
C        8.71448858       3.62311780      -0.00023600
C        9.94109721       2.88750112      -0.00031000
C       11.17720663       3.59079846      -0.00043100
C       12.35965442       2.84358336      -0.00049100
C       12.34009569       1.46621958      -0.00042100
C       11.14050325       0.74158222      -0.00029400
C        9.92199909       1.45593074      -0.00025500
C        8.65576510       0.74880902      -0.00015400
C        7.49761050       1.48578624      -0.00008500
C        7.47776402       2.91299926      -0.00011800
C        6.25306492       5.03007198      -0.00006300
C        6.27586311       3.63048842      -0.00003700
C        5.02574017       5.75877631       0.00002100
C        4.96424006       7.13015272       0.00000800
C        6.20847340       7.87530143      -0.00006800
C        6.19667511       9.28771322      -0.00002000
C        7.42286943       9.96636153      -0.00014100
C        8.62662507       9.29666965      -0.00027900
C       -1.21450933      13.61472577       0.00072000
C       -1.24353337      12.21325957       0.00055900
C       -0.01216931      14.30425744       0.00079100
C        1.19134024      13.61677001       0.00069700
C        1.22274459      12.21535645       0.00053600
C       -0.00979104      11.50461067       0.00047600
C       -2.49840294      11.47089480       0.00047000
C       -2.48700182      10.04876048       0.00032700
C       -1.23545432       9.35642485       0.00026700
C       -0.00856839      10.06760901       0.00033000
C        1.21952599       9.35851261       0.00024600
C        2.46989315      10.05297714       0.00028200
C        2.47887385      11.47512729       0.00042100
C        3.71720312      12.12555145       0.00043800
C        4.90025670      11.41993317       0.00030300
C        4.92801419      10.01873684       0.00015300
C        3.70011990       9.32065404       0.00017000
C        3.67939076       7.87050186       0.00006700
C        2.46207284       7.23599848       0.00002800
C        1.21614479       7.93241823       0.00010800
C       -1.22964458       7.93033807       0.00012700
C       -0.00616972       7.25028833       0.00005400
C       -2.47438477       7.23179936       0.00006600
C       -3.69278061       7.86423003       0.00012300
C       -3.71598008       9.31434484       0.00023600
C       -4.94506239      10.01033583       0.00025100
C       -4.91968830      11.41157679       0.00040300
C       -3.73783892      12.11920941       0.00051400

@DNF, @Henrique_Becker, The variables

> ase_atoms_no_h 
PyObject Atoms(symbols='C168', pbc=False)
> neighborlist = pyimport("ase.neighborlist")
> neighbor_list = neighborlist.neighbor_list("ij", ase_atoms_no_h, 2.0)
> print(neighbor_list)
([0, 0, 1, 1, 1, 2, 2, 3, 3, 4, 4, 4, 5, 5, 5, 6, 6, 6, 7, 7, 7, 8, 8, 8, 9, 9, 9, 10, 10, 10, 11, 11, 11, 12, 12, 12, 13, 13, 14, 14, 15, 15, 15, 16, 16, 16, 17, 17, 17, 18, 18, 19, 19, 19, 20, 20, 20, 21, 21, 22, 22, 23, 23, 23, 24, 24, 24, 25, 25, 25, 26, 26, 27, 27, 28, 28, 29, 29, 29, 30, 30, 31, 31, 32, 32, 32, 33, 33, 33, 34, 34, 34, 35, 35, 35, 36, 36, 36, 37, 37, 37, 38, 38, 38, 39, 39, 39, 40, 40, 40, 41, 41, 42, 42, 43, 43, 43, 44, 44, 44, 45, 45, 45, 46, 46, 47, 47, 47, 48, 48, 48, 49, 49, 50, 50, 51, 51, 51, 52, 52, 52, 53, 53, 53, 54, 54, 55, 55, 56, 56, 57, 57, 57, 58, 58, 59, 59, 60, 60, 60, 61, 61, 61, 62, 62, 62, 63, 63, 63, 64, 64, 64, 65, 65, 65, 66, 66, 66, 67, 67, 67, 68, 68, 68, 69, 69, 70, 70, 71, 71, 71, 72, 72, 72, 73, 73, 73, 74, 74, 75, 75, 75, 76, 76, 76, 77, 77, 78, 78, 79, 79, 79, 80, 80, 80, 81, 81, 81, 82, 82, 83, 83, 84, 84, 85, 85, 85, 86, 86, 87, 87, 88, 88, 88, 89, 89, 89, 90, 90, 90, 91, 91, 91, 92, 92, 92, 93, 93, 93, 94, 94, 94, 95, 95, 95, 96, 96, 96, 97, 97, 98, 98, 99, 99, 99, 100, 100, 100, 101, 101, 101, 102, 102, 103, 103, 103, 104, 104, 104, 105, 105, 106, 106, 107, 107, 107, 108, 108, 108, 109, 109, 109, 110, 110, 111, 111, 112, 112, 113, 113, 113, 114, 114, 115, 115, 116, 116, 116, 117, 117, 117, 118, 118, 118, 119, 119, 119, 120, 120, 120, 121, 121, 121, 122, 122, 122, 123, 123, 123, 124, 124, 124, 125, 125, 126, 126, 127, 127, 127, 128, 128, 128, 129, 129, 129, 130, 130, 131, 131, 131, 132, 132, 132, 133, 133, 134, 134, 135, 135, 135, 136, 136, 136, 137, 137, 137, 138, 138, 139, 139, 140, 140, 141, 141, 141, 142, 142, 143, 143, 144, 144, 144, 145, 145, 145, 146, 146, 146, 147, 147, 147, 148, 148, 148, 149, 149, 149, 150, 150, 150, 151, 151, 151, 152, 152, 152, 153, 153, 154, 154, 155, 155, 155, 156, 156, 156, 157, 157, 157, 158, 158, 159, 159, 159, 160, 160, 160, 161, 161, 162, 162, 163, 163, 163, 164, 164, 164, 165, 165, 165, 166, 166, 167, 167], [2, 1, 5, 6, 0, 3, 0, 4, 2, 5, 3, 12, 1, 4, 9, 7, 1, 27, 8, 6, 24, 9, 7, 20, 5, 10, 8, 9, 11, 19, 12, 16, 10, 4, 11, 13, 12, 14, 13, 15, 14, 16, 165, 11, 15, 17, 16, 18, 163, 19, 17, 21, 10, 18, 8, 22, 21, 19, 20, 20, 23, 24, 22, 45, 7, 23, 25, 24, 26, 43, 25, 27, 6, 26, 29, 30, 28, 33, 34, 28, 31, 30, 32, 40, 31, 33, 29, 32, 37, 29, 35, 55, 34, 36, 52, 35, 37, 48, 36, 38, 33, 37, 39, 47, 40, 38, 44, 39, 32, 41, 42, 40, 41, 43, 25, 44, 42, 39, 45, 43, 23, 46, 44, 47, 45, 38, 46, 49, 36, 49, 50, 47, 48, 48, 51, 50, 73, 52, 53, 51, 35, 71, 54, 52, 55, 53, 54, 34, 57, 58, 61, 56, 62, 56, 59, 58, 60, 59, 68, 61, 60, 57, 65, 57, 63, 83, 80, 62, 64, 63, 76, 65, 66, 61, 64, 65, 67, 75, 66, 68, 72, 60, 67, 69, 68, 70, 69, 71, 53, 70, 72, 67, 71, 73, 74, 72, 51, 73, 75, 66, 74, 77, 77, 78, 64, 76, 75, 76, 79, 78, 80, 101, 79, 81, 63, 80, 82, 99, 83, 81, 62, 82, 85, 86, 84, 90, 89, 87, 84, 86, 88, 87, 89, 96, 88, 93, 85, 91, 111, 85, 90, 92, 108, 91, 93, 104, 89, 92, 94, 103, 93, 95, 94, 96, 100, 88, 95, 97, 96, 98, 97, 99, 81, 98, 100, 95, 99, 101, 79, 100, 102, 103, 101, 94, 102, 105, 106, 105, 92, 104, 103, 107, 104, 129, 106, 108, 91, 107, 109, 108, 110, 127, 109, 111, 90, 110, 113, 114, 112, 117, 118, 112, 115, 114, 116, 124, 115, 117, 113, 116, 121, 113, 119, 139, 118, 120, 136, 119, 121, 132, 117, 120, 122, 121, 123, 131, 124, 122, 128, 123, 116, 125, 126, 124, 127, 125, 126, 128, 109, 123, 129, 127, 107, 128, 130, 129, 131, 122, 130, 133, 120, 133, 134, 131, 132, 132, 135, 134, 136, 157, 137, 135, 119, 136, 138, 155, 137, 139, 118, 138, 141, 142, 140, 146, 145, 140, 143, 142, 144, 143, 145, 152, 149, 141, 144, 141, 147, 167, 164, 148, 146, 160, 149, 147, 150, 145, 148, 151, 149, 159, 150, 152, 156, 144, 153, 151, 152, 154, 153, 155, 137, 154, 156, 151, 155, 157, 135, 156, 158, 157, 159, 161, 158, 150, 148, 162, 161, 159, 160, 163, 160, 164, 162, 17, 163, 165, 147, 15, 166, 164, 167, 165, 166, 146])
> a = convert_neighbor_list(neighbor_list)
> print(a)
OrderedDict{Any, Any}(1 => [2, 1], 2 => [5, 6, 0], 3 => [3, 0], 4 => [4, 2], 5 => [5, 3, 12], 6 => [1, 4, 9], 7 => [7, 1, 27], 8 => [8, 6, 24], 9 => [9, 7, 20], 10 => [5, 10, 8], 11 => [9, 11, 19], 12 => [12, 16, 10], 13 => [4, 11, 13], 14 => [12, 14], 15 => [13, 15], 16 => [14, 16, 165], 17 => [11, 15, 17], 18 => [16, 18, 163], 19 => [19, 17], 20 => [21, 10, 18], 21 => [8, 22, 21], 22 => [19, 20], 23 => [20, 23], 24 => [24, 22, 45], 25 => [7, 23, 25], 26 => [24, 26, 43], 27 => [25, 27], 28 => [6, 26], 29 => [29, 30], 30 => [28, 33, 34], 31 => [28, 31], 32 => [30, 32], 33 => [40, 31, 33], 34 => [29, 32, 37], 35 => [29, 35, 55], 36 => [34, 36, 52], 37 => [35, 37, 48], 38 => [36, 38, 33], 39 => [37, 39, 47], 40 => [40, 38, 44], 41 => [39, 32, 41], 42 => [42, 40], 43 => [41, 43], 44 => [25, 44, 42], 45 => [39, 45, 43], 46 => [23, 46, 44], 47 => [47, 45], 48 => [38, 46, 49], 49 => [36, 49, 50], 50 => [47, 48], 51 => [48, 51], 52 => [50, 73, 52], 53 => [53, 51, 35], 54 => [71, 54, 52], 55 => [55, 53], 56 => [54, 34], 57 => [57, 58], 58 => [61, 56, 62], 59 => [56, 59], 60 => [58, 60], 61 => [59, 68, 61], 62 => [60, 57, 65], 63 => [57, 63, 83], 64 => [80, 62, 64], 65 => [63, 76, 65], 66 => [66, 61, 64], 67 => [65, 67, 75], 68 => [66, 68, 72], 69 => [60, 67, 69], 70 => [68, 70], 71 => [69, 71], 72 => [53, 70, 72], 73 => [67, 71, 73], 74 => [74, 72, 51], 75 => [73, 75], 76 => [66, 74, 77], 77 => [77, 78, 64], 78 => [76, 75], 79 => [76, 79], 80 => [78, 80, 101], 81 => [79, 81, 63], 82 => [80, 82, 99], 83 => [83, 81], 84 => [62, 82], 85 => [85, 86], 86 => [84, 90, 89], 87 => [87, 84], 88 => [86, 88], 89 => [87, 89, 96], 90 => [88, 93, 85], 91 => [91, 111, 85], 92 => [90, 92, 108], 93 => [91, 93, 104], 94 => [89, 92, 94], 95 => [103, 93, 95], 96 => [94, 96, 100], 97 => [88, 95, 97], 98 => [96, 98], 99 => [97, 99], 100 => [81, 98, 100], 101 => [95, 99, 101], 102 => [79, 100, 102], 103 => [103, 101], 104 => [94, 102, 105], 105 => [106, 105, 92], 106 => [104, 103], 107 => [107, 104], 108 => [129, 106, 108], 109 => [91, 107, 109], 110 => [108, 110, 127], 111 => [109, 111], 112 => [90, 110], 113 => [113, 114], 114 => [112, 117, 118], 115 => [112, 115], 116 => [114, 116], 117 => [124, 115, 117], 118 => [113, 116, 121], 119 => [113, 119, 139], 120 => [118, 120, 136], 121 => [119, 121, 132], 122 => [117, 120, 122], 123 => [121, 123, 131], 124 => [124, 122, 128], 125 => [123, 116, 125], 126 => [126, 124], 127 => [127, 125], 128 => [126, 128, 109], 129 => [123, 129, 127], 130 => [107, 128, 130], 131 => [129, 131], 132 => [122, 130, 133], 133 => [120, 133, 134], 134 => [131, 132], 135 => [132, 135], 136 => [134, 136, 157], 137 => [137, 135, 119], 138 => [136, 138, 155], 139 => [137, 139], 140 => [118, 138], 141 => [141, 142], 142 => [140, 146, 145], 143 => [140, 143], 144 => [142, 144], 145 => [143, 145, 152], 146 => [149, 141, 144], 147 => [141, 147, 167], 148 => [164, 148, 146], 149 => [160, 149, 147], 150 => [150, 145, 148], 151 => [151, 149, 159], 152 => [150, 152, 156], 153 => [144, 153, 151], 154 => [152, 154], 155 => [153, 155], 156 => [137, 154, 156], 157 => [151, 155, 157], 158 => [135, 156, 158], 159 => [157, 159], 160 => [161, 158, 150], 161 => [148, 162, 161], 162 => [159, 160], 163 => [163, 160], 164 => [164, 162, 17], 165 => [163, 165, 147], 166 => [15, 166, 164], 167 => [167, 165], 168 => [166, 146])

That neighbor list is distance based (easy) or is it based on the chemical nature and valence of each atom?

Unfortunately, I found the original code to be too confusing, and in the end I rewrote the whole algorithm. One thing I was not sure during the implementation is if you are interested only on cycles of exactly length 6? Or any length smaller than 6 will do? My implementation only gets the cycles of size 6, but is easy to change this.

using LinearAlgebra
using Statistics
using DataStructures

const edges = ([0, 0, 1, 1, 1, 2, 2, 3, 3, 4, 4, 4, 5, 5, 5, 6, 6, 6, 7, 7, 7, 8, 8, 8, 9, 9, 9, 10, 10, 10, 11, 11, 11, 12, 12, 12, 13, 13, 14, 14, 15, 15, 15, 16, 16, 16, 17, 17, 17, 18, 18, 19, 19, 19, 20, 20, 20, 21, 21, 22, 22, 23, 23, 23, 24, 24, 24, 25, 25, 25, 26, 26, 27, 27, 28, 28, 29, 29, 29, 30, 30, 31, 31, 32, 32, 32, 33, 33, 33, 34, 34, 34, 35, 35, 35, 36, 36, 36, 37, 37, 37, 38, 38, 38, 39, 39, 39, 40, 40, 40, 41, 41, 42, 42, 43, 43, 43, 44, 44, 44, 45, 45, 45, 46, 46, 47, 47, 47, 48, 48, 48, 49, 49, 50, 50, 51, 51, 51, 52, 52, 52, 53, 53, 53, 54, 54, 55, 55, 56, 56, 57, 57, 57, 58, 58, 59, 59, 60, 60, 60, 61, 61, 61, 62, 62, 62, 63, 63, 63, 64, 64, 64, 65, 65, 65, 66, 66, 66, 67, 67, 67, 68, 68, 68, 69, 69, 70, 70, 71, 71, 71, 72, 72, 72, 73, 73, 73, 74, 74, 75, 75, 75, 76, 76, 76, 77, 77, 78, 78, 79, 79, 79, 80, 80, 80, 81, 81, 81, 82, 82, 83, 83, 84, 84, 85, 85, 85, 86, 86, 87, 87, 88, 88, 88, 89, 89, 89, 90, 90, 90, 91, 91, 91, 92, 92, 92, 93, 93, 93, 94, 94, 94, 95, 95, 95, 96, 96, 96, 97, 97, 98, 98, 99, 99, 99, 100, 100, 100, 101, 101, 101, 102, 102, 103, 103, 103, 104, 104, 104, 105, 105, 106, 106, 107, 107, 107, 108, 108, 108, 109, 109, 109, 110, 110, 111, 111, 112, 112, 113, 113, 113, 114, 114, 115, 115, 116, 116, 116, 117, 117, 117, 118, 118, 118, 119, 119, 119, 120, 120, 120, 121, 121, 121, 122, 122, 122, 123, 123, 123, 124, 124, 124, 125, 125, 126, 126, 127, 127, 127, 128, 128, 128, 129, 129, 129, 130, 130, 131, 131, 131, 132, 132, 132, 133, 133, 134, 134, 135, 135, 135, 136, 136, 136, 137, 137, 137, 138, 138, 139, 139, 140, 140, 141, 141, 141, 142, 142, 143, 143, 144, 144, 144, 145, 145, 145, 146, 146, 146, 147, 147, 147, 148, 148, 148, 149, 149, 149, 150, 150, 150, 151, 151, 151, 152, 152, 152, 153, 153, 154, 154, 155, 155, 155, 156, 156, 156, 157, 157, 157, 158, 158, 159, 159, 159, 160, 160, 160, 161, 161, 162, 162, 163, 163, 163, 164, 164, 164, 165, 165, 165, 166, 166, 167, 167], [2, 1, 5, 6, 0, 3, 0, 4, 2, 5, 3, 12, 1, 4, 9, 7, 1, 27, 8, 6, 24, 9, 7, 20, 5, 10, 8, 9, 11, 19, 12, 16, 10, 4, 11, 13, 12, 14, 13, 15, 14, 16, 165, 11, 15, 17, 16, 18, 163, 19, 17, 21, 10, 18, 8, 22, 21, 19, 20, 20, 23, 24, 22, 45, 7, 23, 25, 24, 26, 43, 25, 27, 6, 26, 29, 30, 28, 33, 34, 28, 31, 30, 32, 40, 31, 33, 29, 32, 37, 29, 35, 55, 34, 36, 52, 35, 37, 48, 36, 38, 33, 37, 39, 47, 40, 38, 44, 39, 32, 41, 42, 40, 41, 43, 25, 44, 42, 39, 45, 43, 23, 46, 44, 47, 45, 38, 46, 49, 36, 49, 50, 47, 48, 48, 51, 50, 73, 52, 53, 51, 35, 71, 54, 52, 55, 53, 54, 34, 57, 58, 61, 56, 62, 56, 59, 58, 60, 59, 68, 61, 60, 57, 65, 57, 63, 83, 80, 62, 64, 63, 76, 65, 66, 61, 64, 65, 67, 75, 66, 68, 72, 60, 67, 69, 68, 70, 69, 71, 53, 70, 72, 67, 71, 73, 74, 72, 51, 73, 75, 66, 74, 77, 77, 78, 64, 76, 75, 76, 79, 78, 80, 101, 79, 81, 63, 80, 82, 99, 83, 81, 62, 82, 85, 86, 84, 90, 89, 87, 84, 86, 88, 87, 89, 96, 88, 93, 85, 91, 111, 85, 90, 92, 108, 91, 93, 104, 89, 92, 94, 103, 93, 95, 94, 96, 100, 88, 95, 97, 96, 98, 97, 99, 81, 98, 100, 95, 99, 101, 79, 100, 102, 103, 101, 94, 102, 105, 106, 105, 92, 104, 103, 107, 104, 129, 106, 108, 91, 107, 109, 108, 110, 127, 109, 111, 90, 110, 113, 114, 112, 117, 118, 112, 115, 114, 116, 124, 115, 117, 113, 116, 121, 113, 119, 139, 118, 120, 136, 119, 121, 132, 117, 120, 122, 121, 123, 131, 124, 122, 128, 123, 116, 125, 126, 124, 127, 125, 126, 128, 109, 123, 129, 127, 107, 128, 130, 129, 131, 122, 130, 133, 120, 133, 134, 131, 132, 132, 135, 134, 136, 157, 137, 135, 119, 136, 138, 155, 137, 139, 118, 138, 141, 142, 140, 146, 145, 140, 143, 142, 144, 143, 145, 152, 149, 141, 144, 141, 147, 167, 164, 148, 146, 160, 149, 147, 150, 145, 148, 151, 149, 159, 150, 152, 156, 144, 153, 151, 152, 154, 153, 155, 137, 154, 156, 151, 155, 157, 135, 156, 158, 157, 159, 161, 158, 150, 148, 162, 161, 159, 160, 163, 160, 164, 162, 17, 163, 165, 147, 15, 166, 164, 167, 165, 166, 146])
function edges2neighbors(xs, ys)
    # The input is base-zero, this is the reason for the "+1".
    max_vertex = max(maximum(xs), maximum(ys)) + 1
    neighbors = Vector{Int}[Int[] for i in 1:max_vertex]
    for (i_v, j_v) in zip(xs, ys)
        # Again, we updated every node number from base zero to base one.
        push!(neighbors[i_v + 1], j_v + 1)
    end
    return neighbors
end

function inner_dfs_cycle_detection!(cycles, candidate, neighbors, cycle_length)
    @assert length(candidate) <= cycle_length
    for next_vertex in neighbors[last(candidate)]
        if next_vertex == first(candidate)
            # The cycle needs to have exactly length 6?
            if length(candidate) < cycle_length
                continue
            else
                push!(cycles, copy(candidate))
            end
        elseif next_vertex in candidate
            # Do not cycle back to the middle of candidate.
            continue
        elseif length(candidate) < cycle_length
            # If the candidate does not have a subcycle and it is smaller
            # than the desired cycle length, then recurse.
            push!(candidate, next_vertex)
            inner_dfs_cycle_detection!(cycles, candidate, neighbors, cycle_length)
        end
    end
    return
end

function dfs_cycle_detection(neighbors, cycle_length)
    cycles = Vector{eltype(neighbors)}()
    for initial_vertex in eachindex(neighbors)
        inner_dfs_cycle_detection!(cycles, [initial_vertex], neighbors, cycle_length)
    end
    return cycles
end

neighbors = edges2neighbors(edges[1], edges[2])
cycles_base_1 = dfs_cycle_detection(neighbors, 6)
cycles_base_0 = map(vertices -> vertices .- 1, cycles_base_1)
println(cycles_base_0)

and the output

[[0, 2, 3, 4, 5, 1], [0, 2, 3, 4, 5, 1], [1, 5, 4, 3, 2, 0], [2, 3, 4, 5, 1, 6], [3, 4, 5, 1, 6, 7], [4, 5, 1, 6, 7, 8], [5, 1, 6, 7, 8, 9], [6, 7, 8, 9, 5, 1], [6, 7, 8, 9, 5, 1], [7, 8, 9, 5, 1, 6], [7, 8, 9, 5, 1, 6], [8, 9, 5, 1, 6, 7], [8, 9, 5, 1, 6, 7], [9, 5, 1, 6, 7, 8], [9, 5, 1, 6, 7, 8], [10, 9, 5, 1, 6, 7], [11, 12, 4, 5, 1, 6], [12, 4, 5, 1, 6, 7], [13, 12, 4, 5, 1, 6], [14, 13, 12, 4, 5, 1], [15, 14, 13, 12, 4, 5], [16, 11, 12, 4, 5, 1], [17, 16, 11, 12, 4, 5], [18, 19, 21, 20, 8, 9], [20, 8, 9, 5, 1, 6], [22, 20, 8, 9, 5, 1], [23, 24, 7, 8, 9, 5], [24, 7, 8, 9, 5, 1], [25, 24, 7, 8, 9, 5], [26, 25, 24, 7, 8, 9], [27, 6, 7, 8, 9, 5], [30, 28, 29, 33, 32, 40], [31, 30, 28, 29, 33, 32], [31, 30, 28, 29, 33, 32], [32, 40, 39, 38, 37, 36], [33, 29, 28, 30, 31, 32], [33, 29, 28, 30, 31, 32], [34, 29, 28, 30, 31, 32], [35, 34, 29, 28, 30, 31], [36, 35, 34, 29, 28, 30], [37, 36, 35, 34, 29, 28], [38, 37, 36, 35, 34, 29], [43, 25, 24, 7, 8, 9], [44, 39, 40, 32, 31, 30], [45, 23, 24, 7, 8, 9], [46, 47, 38, 37, 36, 35], [47, 38, 37, 36, 35, 34], [48, 36, 35, 34, 29, 28], [49, 47, 38, 37, 36, 35], [50, 48, 36, 35, 34, 29], [51, 50, 48, 36, 35, 34], [52, 53, 71, 70, 69, 68], [56, 57, 61, 60, 59, 58], [56, 57, 61, 60, 59, 58], [57, 61, 60, 59, 58, 56], [57, 61, 60, 59, 58, 56], [58, 56, 57, 61, 60, 59], [58, 56, 57, 61, 60, 59], [59, 58, 56, 57, 61, 60], [59, 58, 56, 57, 61, 60], [60, 59, 58, 56, 57, 61], [60, 59, 58, 56, 57, 61], [61, 60, 59, 58, 56, 57], [61, 60, 59, 58, 56, 57], [62, 57, 61, 60, 59, 58], [63, 80, 79, 78, 76, 77], [64, 63, 80, 79, 78, 76], [64, 63, 80, 79, 78, 76], [67, 66, 65, 61, 60, 59], [68, 60, 59, 58, 56, 57], [69, 68, 60, 59, 58, 56], [70, 69, 68, 60, 59, 58], [72, 67, 66, 65, 61, 60], [75, 66, 65, 61, 60, 59], [78, 76, 77, 75, 66, 65], [79, 78, 76, 77, 75, 66], [80, 79, 78, 76, 77, 75], [81, 80, 79, 78, 76, 77], [82, 83, 62, 57, 61, 60], [83, 62, 57, 61, 60, 59], [85, 84, 86, 87, 88, 89], [88, 87, 86, 84, 85, 90], [89, 88, 87, 86, 84, 85], [89, 88, 87, 86, 84, 85], [91, 90, 111, 110, 109, 108], [92, 91, 90, 111, 110, 109], [93, 89, 88, 87, 86, 84], [95, 94, 103, 102, 101, 79], [96, 88, 87, 86, 84, 85], [97, 96, 88, 87, 86, 84], [98, 97, 96, 88, 87, 86], [99, 81, 80, 79, 78, 76], [100, 95, 94, 103, 102, 101], [100, 95, 94, 103, 102, 101], [101, 79, 78, 76, 77, 75], [102, 103, 94, 93, 89, 88], [104, 106, 107, 129, 128, 123], [105, 104, 106, 107, 129, 128], [106, 107, 129, 128, 123, 124], [108, 91, 90, 111, 110, 109], [108, 91, 90, 111, 110, 109], [109, 108, 91, 90, 111, 110], [109, 108, 91, 90, 111, 110], [110, 109, 108, 91, 90, 111], [110, 109, 108, 91, 90, 111], [111, 90, 91, 92, 93, 89], [114, 112, 113, 117, 116, 124], [115, 114, 112, 113, 117, 116], [115, 114, 112, 113, 117, 116], [116, 124, 123, 122, 121, 117], [116, 124, 123, 122, 121, 117], [117, 113, 112, 114, 115, 116], [117, 113, 112, 114, 115, 116], [118, 113, 112, 114, 115, 116], [119, 118, 113, 112, 114, 115], [120, 119, 118, 113, 112, 114], [121, 117, 113, 112, 114, 115], [122, 121, 117, 113, 112, 114], [125, 126, 127, 128, 123, 124], [125, 126, 127, 128, 123, 124], [128, 123, 124, 116, 115, 114], [130, 129, 107, 106, 104, 105], [131, 122, 121, 117, 113, 112], [132, 120, 119, 118, 113, 112], [133, 131, 122, 121, 117, 113], [134, 132, 120, 119, 118, 113], [135, 134, 132, 120, 119, 118], [138, 137, 136, 135, 134, 132], [139, 118, 113, 112, 114, 115], [141, 140, 142, 143, 144, 145], [142, 140, 141, 146, 147, 164], [143, 142, 140, 141, 146, 147], [144, 143, 142, 140, 141, 146], [145, 149, 150, 151, 152, 144], [145, 149, 150, 151, 152, 144], [146, 141, 140, 142, 143, 144], [147, 164, 163, 162, 160, 148], [147, 164, 163, 162, 160, 148], [149, 150, 151, 152, 144, 143], [152, 144, 143, 142, 140, 141], [153, 152, 144, 143, 142, 140], [154, 153, 152, 144, 143, 142], [155, 137, 136, 135, 134, 132], [156, 151, 150, 149, 145, 141], [157, 135, 134, 132, 120, 119], [158, 157, 135, 134, 132, 120], [159, 161, 160, 148, 149, 150], [162, 163, 164, 165, 15, 14], [165, 15, 14, 13, 12, 4]]

EDIT: I just perceived that I get all classes of symmetry for each cycle, but I think it is better to remove this in postprocessing. Transforming each cycle in a Set of edges and calling unique should reduce the cycle list to only one representation of each cycle.

1 Like

Hey, @Henrique_Becker thank you for rewriting the code:

using MolecularGraph
using LinearAlgebra
using Statistics
using DataStructures

function compound(SMILES::String) 
    mol = smilestomol(SMILES)
    mol_ed = [] # molecule_edges extracted from MolecularGraph
    i_vert = []
    j_vert = []
    for edge in mol.edges
        push!(mol_ed, edge, reverse(edge))
    end
    for ed in mol_ed
        push!(i_vert, ed[1])
        push!(j_vert, ed[2])
    end

    neighbour_list = tuple(i_vert, j_vert)
    return neighbour_list
end

function edges2neighbors(xs, ys)
    max_vertex = max(maximum(xs), maximum(ys))
    neighbours = Vector{Int}[Int[] for i in 1:max_vertex]
    for (i_v, j_v) in zip(xs, ys)
        push!(neighbours[i_v], j_v)
    end
    return neighbours
end

a = edges2neighbors(compound("C1=CC=CC=C1")[1], compound("C1=CC=CC=C1")[2])

print(a)

function inner_dfs_cycle_detection!(cycles, candidate, neighbours, cycle_length)
    @assert length(candidate) <= cycle_length
    for next_vertex in neighbours[last(candidate)]
        if next_vertex == first(candidate)
            if length(candidate) < cycle_length
                continue
            else
                push!(cycles, copy(candidate))
            end
        elseif next_vertex in candidate
            # Do not cycle back to the middle of candidate.
            continue
        elseif length(candidate) < cycle_length
            # If the candidate does not have a subcycle and it is smaller
            # than the desired cycle length, then recurse.
            push!(candidate, next_vertex)
            inner_dfs_cycle_detection!(cycles, candidate, neighbours, cycle_length)
        end
    end
    return
end

function dfs_cycle_detection(neighbours, cycle_length)
    cycles = Vector{eltype(neighbours)}()
    for initial_vertex in eachindex(neighbours)
        inner_dfs_cycle_detection!(cycles, [initial_vertex], neighbours, cycle_length)
    end
    return cycles
end

cycles_base_1 = dfs_cycle_detection(a, 6)
print(cycles_base_1)

# output: [[1, 2, 3, 4, 5, 6], [2, 1, 6, 5, 4, 3], [3, 2, 1, 6, 5, 4], [3, 2, 1, 6, 5, 4], [4, 3, 2, 1, 6, 5], [4, 3, 2, 1, 6, 5], [5, 4, 3, 2, 1, 6], [5, 4, 3, 2, 1, 6], [6, 5, 4, 3, 2, 1], [6, 5, 4, 3, 2, 1]]

I have added the compound function, that converts SMILES to a list of edges. But now, I have another issue: this code you wrote returns the same cycle with different orders, as you can see in the output:

[[1, 2, 3, 4, 5, 6], [2, 1, 6, 5, 4, 3], [3, 2, 1, 6, 5, 4], [3, 2, 1, 6, 5, 4], [4, 3, 2, 1, 6, 5], [4, 3, 2, 1, 6, 5], [5, 4, 3, 2, 1, 6], [5, 4, 3, 2, 1, 6], [6, 5, 4, 3, 2, 1], [6, 5, 4, 3, 2, 1]]

I would like to return an array with only different cycles, because I’m dealing with molecules, which are instances of undirected graphs.

And I would like to identify cycles with less than 4 vertices

I had already mentioned how to only have unique cycles.

The following code solve the uniqueness problem:

function unique_cycles(cycles)
    (x -> first.(x)).(unique!(map(c -> Set(c .=> circshift(c, 1)), cycles)))
end

unique_cycles_base_0 = unique_cycles(cycles_base_0)
println(unique_cycles_base_0)

Less than 6 vertices you mean? Or between 4 and 6 vertices?

The part of the code that needs to be edited is:

function inner_dfs_cycle_detection!(cycles, candidate, neighbors, cycle_length)
    @assert length(candidate) <= cycle_length
    for next_vertex in neighbors[last(candidate)]
        if next_vertex == first(candidate)
            # The cycle needs to have exactly length 6?
            if length(candidate) < cycle_length
                continue
            else
                push!(cycles, copy(candidate))
            end
[...]

The last if clause before continue indicates the condition for ignoring a cycle. If you want every cycle larger than a self-loop (i.e., single vertice loops, not sure if they are possible there) then just change from length(candidate) < cycle_length to length(candidate) < 1. The parameter cycle_length will keep defining the maximum cycle size.