GROMACS version: 2021
GROMACS modification: No
Here post your question
I have two conserved water molecules, which should be considered as ligands interacting with the protein. The [ molecules ] section in the topology file is like this. The first SOL group has 2 molecules, which represents the 2 conserved water; the second SOL group is the solvent water. But these two SOL groups essentially share the same water topology (i.e. tip3p.itp).
[ molecules ]
; Compound #mols
Protein_chain_A 1
lig 1
SOL 2
SOL 7457
NA 14
CL 8
Justin’s protein-ligand tutorial tells us we can create a new index that includes protein and ligand, so as to separate it from the rest in the tc-grps.
But I think I also need to add a new [ moleculetype ] specifically for the two conserved water, and use the same POSRES position restraints as protein
#ifdef POSRES
; Position restraint for each water oxygen
[ position_restraints ]
; i funct fcx fcy fcz
1 1 1000 1000 1000
#endif
to constrain the 2 conserved water, so that define = -DPOSRES will constrain not only the protein but also the 2 conserved water. Can I ask what I need to modify to this new [ moleculetype ]? Is there a explanation for the default tip3p.itp (e.g. #ifndef FLEXIBLE, [ settles ]) ?
Hi,
If you want to position restrain the 2 conserved water molecule you can copy the water topology in your working dir, change the name of the itp file and change the molecular name in the file. Update the [molecules] directive in topol.top and include the new water topology in the correct order in the topol.top together with a new POSRES definition (with a different name from the water one). Then add the new POSRES name in mdp file.
But I wonder if you really need this, since (I guess) you will release the restrain in data production and the water molecules usually exchange with another water molecule.
Hi @alevilla, thank you very much. I think I am following your instructions:
I save the default tip3p.itp as tip3p1.itp and put into my working directory, replace all the SOL strings with FIX. But I am not sure how to deal with the [ settles ] and [ exclusions ] sections. Based on this thread, I replace the [ settles ] and [ exclusions ] sections with [ constraints ], is that correct? But that thread is for “flat-bottom restraint”, what if I want to use the normal “position restraint”?
[ moleculetype ]
; molname nrexcl
FIX 2
[ atoms ]
; id at type res nr res name at name cg nr charge mass
1 OW 1 FIX OW 1 -0.834 16.00000
2 HW 1 FIX HW1 1 0.417 1.00800
3 HW 1 FIX HW2 1 0.417 1.00800
#ifndef FLEXIBLE
;;[ settles ]
; OW funct doh dhh
;;1 1 0.09572 0.15139
;;
;;[ exclusions ]
;;1 2 3
;;2 1 3
;;3 1 2
[ constraints ]
1 2 1 0.09572
1 3 1 0.09572
2 3 1 0.15139
#else
[ bonds ]
; i j funct length force_constant
1 2 1 0.09572 502416.0 0.09572 502416.0
1 3 1 0.09572 502416.0 0.09572 502416.0
[ angles ]
; i j k funct angle force_constant
2 1 3 1 104.52 628.02 104.52 628.02
#endif
Yes, I have also updated the topology file to include the new water molecule. But I use the same POSRES for the conserved water as protein.
leave [ settles ] and [ exclusions ] how they were originally. You do not want to change the description of water, I understood that you want just position restrain the conserved water molecules.
Sorry I do not think so. If the [ settles ] and [ exclusions ] sections are kept in the conserved water, this error occurs:
Fatal error:
The [molecules] section of your topology specifies more than one block of
a [moleculetype] with a [settles] block. Only one such is allowed. If you
are trying to partition your solvent into different *groups* (e.g. for
freezing, T-coupling, etc.) then you are using the wrong approach. Index
files specify groups. Otherwise, you may wish to change the least-used
block of molecules with SETTLE constraints into 3 normal constraints.
Any time you wish to treat some subset of waters differently, yes, you do need to replace [settles] with [constraints]. The exclusions are the same either way.
Thank you, @jalemkul, do you mean I should only choose the words of [ settles ] to [ constraints ], or should I change the [ settles ]section to [ constraints ]section?
The below is only changing the words. If the section needs to change, what is the content for the [ constraints ] section?
[ moleculetype ]
; molname nrexcl
FIX 2
[ atoms ]
; id at type res nr res name at name cg nr charge mass
1 OW 1 FIX OW 1 -0.834 16.00000
2 HW 1 FIX HW1 1 0.417 1.00800
3 HW 1 FIX HW2 1 0.417 1.00800
#ifndef FLEXIBLE
[ constraints ]
; OW funct doh dhh
1 1 0.09572 0.15139
[ exclusions ]
1 2 3
2 1 3
3 1 2
#else
[ bonds ]
; i j funct length force_constant
1 2 1 0.09572 502416.0 0.09572 502416.0
1 3 1 0.09572 502416.0 0.09572 502416.0
[ angles ]
; i j k funct angle force_constant
2 1 3 1 104.52 628.02 104.52 628.02
#endif
Thank you Justin, would it be possible if you can show me an example for the “three constraints” for “position restraint”? Is it the same for “position restraint” and “flat-bottom restraint” (earlier thread)?
Replacing SETTLE by three constraints is not dependent upon the type of restraint being used (in fact it’s not technically related at all). What you have listed for constraints is correct.