Constraining terminal basepairs in DNA

I just wanted to give a final update on this.

I never did get the “pull” function to work. As @genie mentioned, it may have been due to my Gromacs install, but as I am running on someone else’s cluster, I was not able to modify that. I tried multiple different definitions for my groups and nothing seemed to work.

However, I recently found out about the “merge” function in pdb2gmx. I used this function to combine my two DNA strands into one topology. (Previously, they were two separate topologies that were included as .itp files in the top-level topology file.)

By merging the two DNA strands into one topology, I was able to create a bonded (restraint) interaction between the first and last basepairs (using the bond “10” function). Despite the fact that this should be identical (at least as far as I understand the documentation), this approach worked perfectly where the “pull” function had failed. The DNA did not melt and I did not receive any errors in the simulation.

Just in case anyone is trying to do the same thing, this is the code I added to the bottom of the [bonds] section of my topology (topol.top):

;	atom_i atom_j	func	r0	r1	r2 k(kJ mol^-1 nm^-2)
   22	618		10	0.0	0.5	0.6	10000
   16	615		10	0.0	0.5	0.6	10000
  300	331		10	0.0	0.5	0.6	10000
  303 	337		10	0.0	0.5	0.6	10000

The atoms indicated are the heavy atoms (nitrogen and oxygens) that form the hydrogen bonds in my terminal A-T bases.

Hope this helps somebody!

1 Like