Yes, the above section is included in the topol.itp files.
So, is it better to remove define=-DPOSRES from the .mdp file?
Yes, the above section is included in the topol.itp files.
So, is it better to remove define=-DPOSRES from the .mdp file?
If that section (#ifdef POSRES
…) is present in the protein itp files then you should not keep -DPOSRES
in the mdp file, unless you want to restrain the protein. However, you might want to check that all restraints you want are really present in the posre_dna.itp
file.
You might also want to look through the topol_DNA_chain*.itp
files and consider if you should change any #ifdef POSRES
sections (if any) to #ifdef POSRES_DNA
. Exactly which atoms you want to restrain is up to you to decide.
@MagnusL Thank you for your suggestions. Here I only want to restrain DNA and for that, I generated an index file for DNA and using genrestr generated the restrain file for DNA . And included in the topol.top file.
Is the method which I have done is correct or any suggestions?
As it says in the manual (https://manual.gromacs.org/documentation/current/onlinehelp/gmx-genrestr.html#gmx-genrestr), position restraints are defined per molecule. So, you should generate one for each of your DNA molecule types (.itp files) and include the restraint file from the .itp file. You must make sure that the atom numbering starts at 1 (where 1 is the first atom in the molecule type, i.e., you do not have to restrain atom number 1).
@ MagnusL I am facing a problem, the restrain on DNA is not working. I am using the mdp file
title = KALP15-DPPC Production MD
define = -DPOSRES_DNA
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 50000000 ; 2 * 50000000 = 100000 ps (100 ns)
dt = 0.002 ; 2 fs
; Output control
nstxout = 1000000 ; save coordinates every 2 ps
nstvout = 1000000 ; save velocities every 2 ps
nstxtcout = 1000000 ; xtc compressed trajectory output every 2 ps
nstenergy = 1000000 ; save energies every 2 ps
nstlog = 1000000 ; update log file every 2 ps
; Bond parameters
continuation = yes ; Restarting after NPT
constraint_algorithm = lincs ; holonomic constraints
constraints = h-bonds ; all bonds (even heavy atom-H bonds) constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
ns_type = grid ; search neighboring grid cels
nstlist = 5 ; 10 fs
rlist = 1.2 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm)
rvdw = 1.2 ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
; Temperature coupling is on
tcoupl = Nose-Hoover ; More accurate thermostat
tc-grps = DNA_DFPE Water_and_ions ; two coupling groups - more accurate
tau_t = 0.5 0.5 ; time constant, in ps
ref_t = 300 300 ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl = C-rescale ; Pressure coupling on in NPT
pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 1.0 ; reference pressure, x-y, z (in bar)
compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = no ; Velocity generation is off
; COM motion removal
; These options remove motion of the protein/bilayer relative to the solvent/ions
nstcomm = 100
comm-mode = Linear
comm-grps = DNA_DFPE Water_and_ions
; Scale COM of reference coordinates
refcoord_scaling = all
and the topol.top
File ‘topol.top’ was generated
; By user: nabanita (1001)
; On host: kashi
; At date: Wed Apr 17 10:32:08 2024
;
; This is a standalone topology file
;
; Created by:
; :-) GROMACS - gmx pdb2gmx, 2023.3 (-:
;
; Executable: /usr/local/gromacs/bin/gmx
; Data prefix: /usr/local/gromacs
; Working dir: /data/nabanita/13th_demo/hj_test/17th_apr
; Command line:
; gmx pdb2gmx -f new.pdb -o dna_mem.gro -ter
; Force field was read from current directory or a relative path - path added.
;
; Include forcefield parameters
#include “./charmm36-jul2022.ff/forcefield.itp”
; Include chain topologies
#include “topol_DNA_chain_D.itp”
#include “topol_DNA_chain_D2.itp”
#include “topol_DNA_chain_D3.itp”
#include “topol_DNA_chain_D4.itp”
#include “topol_Other_chain_P.itp”
#include “topol_Other_chain_P2.itp”
#include “topol_Other_chain_P3.itp”
#include “topol_Other_chain_P4.itp”
#include “topol_Other_chain_P5.itp”
#include “topol_Other_chain_P6.itp”
; Include position restrain file
#ifdef POSRES_DNA
#include “posre_dna.itp”
#endif
; Include water topology
#include “./charmm36-jul2022.ff/tip3p.itp”
#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
; i funct fcx fcy fcz
1 1 1000 1000 1000
#endif
; Include topology for ions
#include “./charmm36-jul2022.ff/ions.itp”
[ system ]
; Name
Protein in water
[ molecules ]
; Compound #mols
DNA_chain_D 1
DNA_chain_D2 1
DNA_chain_D3 1
DNA_chain_D4 1
Other_chain_P 1
Other_chain_P2 1
Other_chain_P3 1
Other_chain_P4 1
Other_chain_P5 1
Other_chain_P6 1
SOL 87531
CL 47
MG 47
CL 47
NA 155
CL 63
Can you tell me where exactly I am going wrong
As I said in my post above, the
#ifdef POSRES_DNA
#include “posre_dna.itp”
#endif
section must be moved to each separate topol_DNA_chain_*.itp
and the posre_dna.itp
must be specific to each molecule type (i.e., you will need four such files).
There is more information in the manual link I posted.
@MagnusL the topol_DNA_chain_*.itp contains the following
[ dihedrals ]
; ai aj ak al funct c0 c1 c2 c3
21 16 12 22 2
22 21 23 24 2
45 46 49 47 2
45 50 48 49 2
52 50 54 53 2
79 76 81 80 2
83 81 85 84 2
85 83 77 86 2
118 113 109 119 2
119 118 120 121 2
142 143 146 144 2
142 147 145 146 2
149 147 151 150 2
176 173 178 177 2
180 178 182 181 2
182 180 174 183 2
215 210 206 216 2
216 215 217 218 2
240 237 242 241 2
244 242 246 245 2
246 244 238 247 2
274 269 276 275 2
277 276 272 278 2
278 277 279 280 2
304 299 306 305 2
307 306 302 308 2
308 307 309 310 2
339 334 330 340 2
340 339 341 342 2
364 361 366 365 2
368 366 370 369 2
370 368 362 371 2
418 421 420 419 2
; Include Position restraint file
#ifdef POSRES
#include “posre_DNA_chain_D.itp”
#endif
The modified topol.top is
File ‘topol.top’ was generated
; By user: nabanita (1001)
; On host: kashi
; At date: Wed Apr 17 10:32:08 2024
;
; This is a standalone topology file
;
; Created by:
; :-) GROMACS - gmx pdb2gmx, 2023.3 (-:
;
; Executable: /usr/local/gromacs/bin/gmx
; Data prefix: /usr/local/gromacs
; Working dir: /data/nabanita/13th_demo/hj_test/17th_apr
; Command line:
; gmx pdb2gmx -f new.pdb -o dna_mem.gro -ter
; Force field was read from current directory or a relative path - path added.
;
; Include forcefield parameters
#include “./charmm36-jul2022.ff/forcefield.itp”
; Include chain topologies
#include “topol_DNA_chain_D.itp”
#include “topol_DNA_chain_D2.itp”
#include “topol_DNA_chain_D3.itp”
#include “topol_DNA_chain_D4.itp”
#include “topol_Other_chain_P.itp”
#include “topol_Other_chain_P2.itp”
#include “topol_Other_chain_P3.itp”
#include “topol_Other_chain_P4.itp”
#include “topol_Other_chain_P5.itp”
#include “topol_Other_chain_P6.itp”
; Include position restrain file
#ifdef POSRES_DNAD
#include “posre_DNA_chain_D.itp”
#endif
; Include position restrain file
#ifdef POSRES_DNAD2
#include “posre_DNA_chain_D2.itp”
#endif
; Include position restrain file
#ifdef POSRES_DNAD3
#include “posre_DNA_chain_D3.itp”
#endif
; Include position restrain file
#ifdef POSRES_DNAD4
#include “posre_DNA_chain_D4.itp”
#endif
; Include water topology
#include “./charmm36-jul2022.ff/tip3p.itp”
#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
; i funct fcx fcy fcz
1 1 1000 1000 1000
#endif
; Include topology for ions
#include “./charmm36-jul2022.ff/ions.itp”
[ system ]
; Name
Protein in water
[ molecules ]
; Compound #mols
DNA_chain_D 1
DNA_chain_D2 1
DNA_chain_D3 1
DNA_chain_D4 1
Other_chain_P 1
Other_chain_P2 1
Other_chain_P3 1
Other_chain_P4 1
Other_chain_P5 1
Other_chain_P6 1
SOL 87531
CL 47
MG 47
CL 47
NA 155
CL 63
and the modified .mdp file
title = KALP15-DPPC Production MD
define = -DPOSRES_DNAD -DPOSRES_DNAD2 -DPOSRES_DNAD3 -DPOSRES_DNAD4
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 50000000 ; 2 * 50000000 = 100000 ps (100 ns)
dt = 0.002 ; 2 fs
; Output control
nstxout = 1000000 ; save coordinates every 2 ps
nstvout = 1000000 ; save velocities every 2 ps
nstxtcout = 1000000 ; xtc compressed trajectory output every 2 ps
nstenergy = 1000000 ; save energies every 2 ps
nstlog = 1000000 ; update log file every 2 ps
; Bond parameters
continuation = yes ; Restarting after NPT
constraint_algorithm = lincs ; holonomic constraints
constraints = h-bonds ; all bonds (even heavy atom-H bonds) constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
ns_type = grid ; search neighboring grid cels
nstlist = 5 ; 10 fs
rlist = 1.2 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm)
rvdw = 1.2 ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
; Temperature coupling is on
tcoupl = Nose-Hoover ; More accurate thermostat
tc-grps = DNA_DFPE Water_and_ions ; two coupling groups - more accurate
tau_t = 0.5 0.5 ; time constant, in ps
ref_t = 300 300 ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl = C-rescale ; Pressure coupling on in NPT
pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 1.0 ; reference pressure, x-y, z (in bar)
compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = no ; Velocity generation is off
; COM motion removal
; These options remove motion of the protein/bilayer relative to the solvent/ions
nstcomm = 100
comm-mode = Linear
comm-grps = DNA_DFPE Water_and_ions
; Scale COM of reference coordinates
refcoord_scaling = all
@MagnusL Please let me know whether it’s okay now
No. That will still not work. The restraints must be defined in the molecule type (or before the next molecule type block Is started, not after.
Easiest for you is to just change
; Include Position restraint file
#ifdef POSRES
#include “posre_DNA_chain_D.itp”
#endif
to
; Include Position restraint file
#ifdef POSRES_DNA
#include “posre_DNA_chain_D.itp”
#endif
in the topol_DNA_chain_D.itp
file and do the corresponding in all other DNA files.
or you can do:
; Include chain topologies
#include “topol_DNA_chain_D.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D.itp”
#endif
#include “topol_DNA_chain_D2.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D2.itp”
#endif
#include “topol_DNA_chain_D3.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D3.itp”
#endif
#include “topol_DNA_chain_D4.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D4.itp”
#endif
#include “topol_Other_chain_P.itp”
#include “topol_Other_chain_P2.itp”
#include “topol_Other_chain_P3.itp”
#include “topol_Other_chain_P4.itp”
#include “topol_Other_chain_P5.itp”
#include “topol_Other_chain_P6.itp”
; Include water topology
#include “./charmm36-jul2022.ff/tip3p.itp”
...
in the topol.top
file.
Then you still only need define = -DPOSRES_DNA
in the mdp file.
I am getting an error
Invalid delimiter for filename in #include statement - File topol.top, line 24
Last line read:
‘#include “topol_DNA_chain_D.itp”’
The problem is almost certainly the quotation marks. The forum uses special more elegant left and right curved quotation marks (“”), when they are not in a block of preformatted text, whereas GROMACS can only handle typewriter quotation marks (""
). You’ll just have to replace those yourself.
@MagnusL I am using the following md.mdp file for running the simulation. But in that case, I see no significant change in the membrane. Is the membrane getting restrained somehow.
Here is the md.mdp file .
title = KALP15-DPPC Production MD
define = -DPOSRES_DNA
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 600000000 ; 1.5 * 600000000 = 900000 ps (90 ns)
dt = 0.0015 ; 1.5 fs
; Output control
nstxout = 1000000 ; save coordinates every 2 ps
nstvout = 1000000 ; save velocities every 2 ps
nstxtcout = 1000000 ; xtc compressed trajectory output every 2 ps
nstenergy = 1000000 ; save energies every 2 ps
nstlog = 1000000 ; update log file every 2 ps
; Bond parameters
continuation = yes ; Restarting after NPT
constraint_algorithm = lincs ; holonomic constraints
constraints = h-bonds ; all bonds (even heavy atom-H bonds) constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
ns_type = grid ; search neighboring grid cels
nstlist = 5 ; 10 fs
rlist = 1.2 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm)
rvdw = 1.2 ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
; Temperature coupling is on
tcoupl = Nose-Hoover ; More accurate thermostat
tc-grps = DNA_DFPE Water_and_ions ; two coupling groups - more accurate
tau_t = 0.5 0.5 ; time constant, in ps
ref_t = 300 300 ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl = C-rescale ; Pressure coupling on in NPT
pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 1.0 ; reference pressure, x-y, z (in bar)
compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = no ; Velocity generation is off
; COM motion removal
; These options remove motion of the protein/bilayer relative to the solvent/ions
nstcomm = 100
comm-mode = Linear
comm-grps = DNA_DFPE Water_and_ions
; Scale COM of reference coordinates
refcoord_scaling = all
The only thing that controls restraints are the settings in your topology files (.top and the .itp files included from you .top file). The restraints are usually set by using define statements, such as -DPOSRES
in the .mdp file. Unless the -DPOSRES_DNA
statement also affects lipids, by mistake, there is no reason to believe they are restrained. What do you mean with significant change? Do the atoms in the lipids not move? DPPC will be in a gel-state at 300K.
@MagnusL I am using the topol.top
File ‘topol.top’ was generated
; By user: nabanita (1001)
; On host: kashi
; At date: Wed Apr 17 10:32:08 2024
;
; This is a standalone topology file
;
; Created by:
; :-) GROMACS - gmx pdb2gmx, 2023.3 (-:
;
; Executable: /usr/local/gromacs/bin/gmx
; Data prefix: /usr/local/gromacs
; Working dir: /data/nabanita/13th_demo/hj_test/17th_apr
; Command line:
; gmx pdb2gmx -f new.pdb -o dna_mem.gro -ter
; Force field was read from current directory or a relative path - path added.
;
; Include forcefield parameters
#include “./charmm36-jul2022.ff/forcefield.itp”
; Include chain topologies
#include “topol_DNA_chain_D.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D.itp”
#endif
#include “topol_DNA_chain_D2.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D2.itp”
#endif
#include “topol_DNA_chain_D3.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D3.itp”
#endif
#include “topol_DNA_chain_D4.itp”
#ifdef POSRES_DNA
#include “posre_DNA_chain_D4.itp”
#endif
#include “topol_Other_chain_P.itp”
#include “topol_Other_chain_P2.itp”
#include “topol_Other_chain_P3.itp”
#include “topol_Other_chain_P4.itp”
#include “topol_Other_chain_P5.itp”
#include “topol_Other_chain_P6.itp”
; Include water topology
#include “./charmm36-jul2022.ff/tip3p.itp”
#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
; i funct fcx fcy fcz
1 1 1000 1000 1000
#endif
; Include topology for ions
#include “./charmm36-jul2022.ff/ions.itp”
[ system ]
; Name
Protein in water
[ molecules ]
; Compound #mols
DNA_chain_D 1
DNA_chain_D2 1
DNA_chain_D3 1
DNA_chain_D4 1
Other_chain_P 1
Other_chain_P2 1
Other_chain_P3 1
Other_chain_P4 1
Other_chain_P5 1
Other_chain_P6 1
SOL 87531
CL 47
MG 47
CL 47
NA 155
CL 63
And in the posre.itp for lipids it is:
; Include Position restraint file
#ifdef POSRES
#include “posre_Other_chain_P.itp”
#endif
With the name of the molecules it is not easy to guess which one is DPPC. I guess it is one of the “Other_chain_P” molecules, or? If there is no position restraints file included without an “#ifdef” statement or with an “#ifdef POSRES_DNA” statement it should be fine.
I think it’s a matter of temperature. For testing purposes, try increasing it to 325 K (then DPPC should have melted). I think you’ll see more movement then.
@MagnusL it seems my simulation is restarting everytime from initial configuration after running for few ns.
Is there any problem with my md.mdp file
title = KALP15-DPPC Production MD
define = -DPOSRES_DNA
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 600000000 ; 1.5 * 600000000 = 900000 ps (900 ns)
dt = 0.0015 ; 1.5 fs
; Output control
nstxout = 1000000 ; save coordinates every 2 ps
nstvout = 1000000 ; save velocities every 2 ps
nstxtcout = 1000000 ; xtc compressed trajectory output every 2 ps
nstenergy = 1000000 ; save energies every 2 ps
nstlog = 1000000 ; update log file every 2 ps
; Bond parameters
continuation = yes ; Restarting after NPT
constraint_algorithm = lincs ; holonomic constraints
constraints = h-bonds ; all bonds (even heavy atom-H bonds) constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
ns_type = grid ; search neighboring grid cels
nstlist = 5 ; 10 fs
rlist = 1.2 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm)
rvdw = 1.2 ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
; Temperature coupling is on
tcoupl = Nose-Hoover ; More accurate thermostat
tc-grps = DNA_DFPE Water_and_ions ; two coupling groups - more accurate
tau_t = 0.5 0.5 ; time constant, in ps
ref_t = 300 300 ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl = C-rescale ; Pressure coupling on in NPT
pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 1.0 ; reference pressure, x-y, z (in bar)
compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = no ; Velocity generation is off
; COM motion removal
; These options remove motion of the protein/bilayer relative to the solvent/ions
nstcomm = 100
comm-mode = Linear
comm-grps = DNA_DFPE Water_and_ions
; Scale COM of reference coordinates
refcoord_scaling = all