Hi! So glad to use this nice package!
I wonder how I can calculate all angles/dihedrals in an undirected graph?
For example, if I have a graph
i -- j -- k -- l
|
m
[(i, j, k), (j, k, l), (i, j, m)...];[(i, j, k, l), (m, j, k, l)][(j, i, m, k)]Are there any convenient methods to calculate those many-body? Currently I calculate it by for-loop each atom’s neighborlist.
It’s great that you gave some examples of what you want to compute, but this alone doesn’t make your question entirely clear.
Can you provide a mathematically precise explanation of what you are looking for?
What about this one? I need to think about where I can find the mathematic descriptions
https://www.charmm-gui.org/?doc=lecture&module=scientific&lesson=9
For example, what do you mean by “improper tension”?
Could you perhaps show the code that you currently use to calculate it? That might provide some understanding of what you are interested in.
Sorry, I left my office, can I provide some pseudo code?
for angle, I do like this:
angles = []
for j in graph:
for i, k in permutation(j.neighbors, 2):
angles.append((i, j, k))
for dihedral:
for j in graph:
for i, k in permutation(j.neighbors, 2):
for l in k.neighbors:
dihedrals.append((i, j, k, l))
for improper tension:
for i in graph:
for j, k, l in permuation(i, 3):
impropers.append((i, j, k, l))
As you can see, I want to enumerate all bond angle/dihedral and improper tension in the molecule graph, and know the total number. You also can know why we need to calculate those thing from aforementioned link :
So you want to count certain subgraphs? Are you looking for induced subgraphs only, or all subgraphs?
Sorry I dont know what is induced subgraph, and I think I want to find all subgraph with a special pattern. For example, the angles, not only I need to get all three-body, but also need the middle index should be the middle node.
Maybe this requirement is too special in graph theory…
This is standard terminology that you can quickly look up … as opposed to “improper tension”, which you did not explain … We need to have a common language to understand each other. Let’s try in a different way:
What is the output you expect for this graph?
Thanks for your patience!!
the angles I want to get are:
(2, 1, 3), (2, 1, 4), (3, 1, 4), # 1 is the vertex of an angle
(2, 3, 1), (2, 3, 4), (1, 3, 4), (2, 3, 5), (4, 3, 5), (1, 3, 5)
(1, 2, 3), (1, 4, 3),
dihedral:
(4, 1, 2, 3), # twist 1-2 by using 4 and 3, that's why it calls dihedral
(2, 1, 4, 3),
(1, 4, 3, 2),
(1, 2, 3, 4),
(2, 1, 3, 4), (2, 1, 3, 5), (4, 1, 3, 5),
# can not twist 3-5 since 5 doesn't have neighbor
improper dihedral angle(improper torsion(sorry for the typo before)) is how atom left the plane formed by other three atoms. Here is another tutorial https://wanglab.hosted.uark.edu/DLPOLY2/node50.html
(3, 2, 1, 5), (3, 2, 4, 5), (3, 1, 2, 4), (3, 1, 4, 5) # 3 left the plane formed by 2-1-5 etc.
(1, 2, 3, 4)
# 2, 4 , 5 doesn't have improper torsion because they don't have enough three neighbors
When we talk about “twist”, “plane”, those spital concept because we need to calculate it in 3D chemical space, but that chemical connection is topological and defined on a graph. So I need to infer those angle/dihedral/improper(which atoms form an angle, get their index and access position, then calculate degree of the angle) from the molecule graph, and then I can calculate intramolecular energy.
count_subisomorphisms_vf2() will count all mappings of a small graph onto a big one.
For example, let’s take the P_4 subgraph, which you called “dihedral”, and labels its vertices as follows:
Then we have the following mappings to the graph from my previous comment:
This is 20 in total. Notice that P_4 can me mapped to itself in two ways (1,2,3,4) → (1,2,3,4) and (1,2,3,4) → (4,3,2,1). Thus each occurrence of P_4 appears in two different ways in the above list. After eliminating these duplicates we are left with a count of 10.
Your list includes only 7. Perhaps you forgot some (easy to do this) or intended to skip some?
The complete list, after eliminating duplicates, would be:
{{1, 2, 3, 4}, {1, 2, 3, 5}, {1, 4, 3, 2}, {1, 4, 3, 5}, {2, 1, 3, 4},
{2, 1, 3, 5}, {2, 1, 4, 3}, {2, 3, 1, 4}, {3, 2, 1, 4}, {4, 1, 3, 5}}
With the 3-star S_3 subgraph (which you called “improper torsion”) we have 6 ways to map it onto itself: there’s a central vertex, and all 6 permutations of the other three are equivalent. count_subisomorphisms_vf2() counts 30 mappings, which is indeed 5 subgraphs (as you stated) after dividing by 6 to eliminate duplicates.
I hope this helps.
P.S. There’s also the motifs_randesu() function, but this effectively counts induced subgraphs, i.e. both connections and non-connections need to match. While this function is much faster, it’d take a bit of work to transform its output to the counts you are looking for.
Sooooo thanks for your kind help!!!
It seems your result is right, I do forget some cases. I will have a try with two functions. I knew that there must be some mathematic tools could solve my problem!
Thank you again for bothering you in your nice weekend!