Link Information not getting reflected in itp file

[majumderS]

Hi, So I was defining link information for backbone bonds, angles and dihedrals under the link directive as follows but this information do not show up in the molecule.itp file. Is this a bug or am I doing something wrong with the definitions:

;;;; Links

;; Links for COIL. We apply them first as coil is the default.
[ link ]
resname $na_resnames
[ bonds ]
BB1   BB2   1      0.363             20000 {"group": "Backbone bonds"}
BB2   BB3   1      0.202             40000 {"group": "Backbone bonds"}
BB3   +BB1  1      0.354             10000 {"group": "Backbone bonds"}
-BB3   BB1  1      0.354             10000 {"group": "Backbone bonds"}

[ link ]
resname $na_resnames
[ angles ]
#meta {"group": "BBB angles"}
BB1   BB2   BB3 2         117   175 {"group": "BBB angles", "comment": "BB1-BB2-BB3"}
BB2   BB3   +BB1    2         95       105 {"group": "BBB angles", "comment": "BB2-BB3-BB1"}
-BB2  -BB3   BB1  2        95    105 {"group": "BBB angles", "comment": "BB2-BB3-BB1"} 
BB3   +BB1   +BB2   2         93       75  {"group": "BBB angles", "comment": "BB3-BB1-BB2"}
-BB3   BB1   BB2    2         93       75  {"group": "BBB angles", "comment": "BB3-BB1-BB2"}

[ link ]
resname $na_resnames
[dihedrals]
#meta  {"group": "BB-BB-BB-BB"}
 BB1   BB2  BB3  +BB1   2         0       3.5   1   {"group": "BB-BB-BB-BB"}
 BB2   BB3  +BB1  +BB2  1         0.0      1     4  {"group": "BB-BB-BB-BB"}
 BB3   +BB1  +BB2  +BB3 9         -10.0    1.5   2  {"group": "BB-BB-BB-BB"}
 BB3   +BB1  +BB2  +BB3 9         10.0     1.5   2  {"group": "BB-BB-BB-BB"}
-BB3   BB1   BB2   BB3  9         10.0     1.5   2  {"group": "BB-BB-BB-BB"}
-BB3   BB1   BB2   BB3  9         -10.0    1.5   2  {"group": "BB-BB-BB-BB"}
-BB2  -BB3   BB1   BB2  1         0.0      1     4  {"group": "BB-BB-BB-BB"}
-BB1  -BB2   -BB3  BB1  2         0       3.5   1   {"group": "BB-BB-BB-BB"}

[pckroon]

Hello!

In this snippet you define 3 links. You probably want/need to split them up into smaller pieces. Also take a look at #401. Let's take your first link as example:

[ link ]
resname $na_resnames
[ bonds ]
BB1 BB2 1 0.363 20000 {"group": "Backbone bonds"}
BB2 BB3 1 0.202 40000 {"group": "Backbone bonds"}
BB3 +BB1 1 0.354 10000 {"group": "Backbone bonds"}
-BB3 BB1 1 0.354 10000 {"group": "Backbone bonds"}

This defines 1 Link, which matches a sequence of 5 BB beads, named BB3-BB1-BB2-BB3-BB1, where the middle three are of the same residue, and the outer two are part of the previous and next residue (defined by residue ID), respectively. I would advice putting the "internal" bonds (BB1-BB2-BB3) in the respective Block definition, and using links only for the interactions bridging across residue boundaries. Then you will only need to make a Link between BB3 and +BB1, since that will be the same as between -BB3 and BB1. Similarly for your angle link, BB2 BB3 +BB1 and -BB2 -BB3 BB1 would apply to the same fragment.

Finally, since you are working with triplets it's good to note that we use induced subgraph isomorphisms to check where links apply: the path/line A-B-C does not match the triangle A-B-C(-back to A).

HTH

[majumderS]

Hi @pckroon Thanks a ton for the suggestions. I did shift the intra residue definitions to the block section while kept the inter-residue links under 'Links' . I also removed the redundant bonds/angles/dihedrals. I also defined the links as separate fragments instead of putting them all together. Heres the snippet.

[ link ]
resname $na_resnames
[ angles ]
#meta {"group": "BBB angles"}
BB3   +BB1   +BB2   2 93 75 

[ link ]
resname $na_resnames
[dihedrals]
#meta  {"group": "BB-BB-BB-BB"}
 BB1   BB2  BB3  +BB1   2 0 3.5 1   

[ link ]
resname $na_resnames
[dihedrals]
#meta  {"group": "BB-BB-BB-BB"}
 BB2   BB3  +BB1  +BB2  1 0.0 1 4

[ link ]
resname $na_resnames
[dihedrals]
#meta  {"group": "BB-BB-BB-BB"}
 BB3   +BB1  +BB2  +BB3 9 -10.0 1.5 2 

But no luck so far. The itp files didnt reflect the links section at all.

Also while going through your suggestion, "Finally, since you are working with triplets it's good to note that we use induced subgraph isomorphisms to check where links apply: the path/line A-B-C does not match the triangle A-B-C(-back to A)." --------I wanted to tell you that the inter-residue connections such as BB2 BB3 +BB1 goes in a straight line from the current residue to the next, so I am not sure if the links are getting rejected for this reason.

Also I was going through the codebase to try to understand where exactly the algorithm is defined that rejects the links section. The ffinput.py , molecule.py all seems to read the links section, but when the itp.py file writes the itp file, it ignores the links section completely [I printed out the molecule.interactions that gives a list of all the block elements and none from links :( ].

If you could direct me to the files / workflow in vermouth that helps to process the links, or if there was some error statements coming up whenever links were getting rejected such as "Bad Link Definition" , it would be great!

[pckroon]

Links are applied by https://github.com/marrink-lab/vermouth-martinize/blob/master/vermouth/processors/do_links.py. Not the easiest part of the codebase, debugging links is not simple.

Did you check that ffinput.py reads your links (just sprinkle in some print statements)? Otherwise your file(s) may be in the wrong place

Also while going through your suggestion, "Finally, since you are working with triplets it's good to note that we use induced subgraph isomorphisms to check where links apply: the path/line A-B-C does not match the triangle A-B-C(-back to A)." --------I wanted to tell you that the inter-residue connections such as BB2 BB3 +BB1 goes in a straight line from the current residue to the next, so I am not sure if the links are getting rejected for this reason.

Maybe to be more explicit: The link BB2 BB3 +BB1 needs edges BB2 - BB3 and BB3 - +BB1, and there cannot be an edge BB2 +BB1 in the molecule you want to apply this link to.

[majumderS]

Yes I did go through the codebase using print statements and yup ffinput.py did read link data and when I tried to print the interactions being written using itp.py file, except the [links] portion everything else got printed out. I will check the edges which are not supposed to be connected . Thanks again.

[majumderS]

Thanks @pckroon I was able to resolve all the issues based on your suggestions and the link information is getting reflected in the itp files finally .

[pckroon]

For posterity, what was the issue in the end?

[majumderS]

So yes, the issue was mainly connectivity of beads as defined in the mapping section, so after doing sanity check with the beads and edges , I was able to resolve the issue. Also corrected few issues for example a few inter-residue connections were placed in [Blocks] instead of [links]. I had also defined the links separately , removing all redundant links . Also there is a specific syntax/order in which the links need to be defined. For example, I had defined the #meta elsewhere as a result one of the links was not getting read:

So a correct syntax that worked for me was: [ link ] resname $resnames [ angles ]

meta {"group": "BBB angles"}

BB2 BB3 +BB1 2 95 105