I am running a system with protein and carbon nanotube using GROMACS 5.1.4. I used position restraints for carbon nanotube however, during NPT run (both protein and carbon nanotube fixed),it appears that the carbon nanotube drifts from its original position (not parallel to z axis; image attached).
My npt.mdp is as below:
title = NPT equilibration
define = -DPROSRES -DPOSRES_CNT ; position restrain the protein
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 75000 ; 2 * 50000 = 100 ps
dt = 0.002 ; 2 fs
; Output control
nstxout = 5000 ; save coordinates every 1.0 ps
nstvout = 5000 ; save velocities every 1.0 ps
nstenergy = 5000 ; save energies every 1.0 ps
nstlog = 5000 ; update log file every 1.0 ps
energygrps = Protein CNT Water NA ;Group(s) to write to energy file
; Bond parameters
continuation = yes ; Restarting after NVT
constraint_algorithm = lincs ; holonomic constraints
constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
cutoff-scheme = Verlet
ns_type = grid ; search neighboring grid cells
nstlist = 20 ; 20 fs, largely irrelevant with Verlet scheme
rcoulomb = 1.4 ; short-range electrostatic cutoff (in nm)
rvdw = 1.4 ; 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 = V-rescale ; modified Berendsen thermostat
tc-grps = Protein CNT Water_and_ions ; two coupling groups - more accurate
tau_t = 0.1 0.1 0.1 ; time constant, in ps
ref_t = 310 310 310 ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT
pcoupltype = isotropic ; uniform scaling of box vectors
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 ; reference pressure, in bar
compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1
refcoord_scaling = com
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
periodic_molecules = yes
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = no ; Velocity generation is off
The position restraints in cnt.itp file are defined as
#ifdef POSRES_CNT
; Position restraint for each carbon
[ position_restraints ]
; i funct fcx fcy fcz
1 1 1000 1000 1000
#endif
Should I use freeze groups or use position restraints? Feedback from people doing CNT simulations in GROMACS is appreciated.
Hi Neha,
I guess you can use the following:
; Non-equilibrium MD stuff
freezegrps = Protein
freezedim = Y Y Y
You can define the groups which you wanna freeze in an index file.
http://manual.gromacs.org/online/mdp_opt.html#neq
You might want to have a look at the above link.
Hi Vishal,
I went through GROMACS mailing list. They don't recommend using freeze groups or position restraints with pressure coupling. Most cases and papers on CNT use position restraints. Should I use freezegroups to get the density right and then switch to NVT with position restraints?
Summary from CNT tutorial
Don't apply pressure coupling along the axis of the nanotube. In fact, for debugging purposes, it might be better to turn off pressure coupling altogether until you figure out if anything is going wrong, and if so, what.
Post from gromacs mailing list:
http://thread.gmane.org/gmane.science.biology.gromacs.user/63316
Hey Neha,
I am not sure whether you can do that or not. I guess you might wanna have some words with Zak. He might be able to help you.
- Badges
- Science topic
- Similar topics
- Thermodynamics
- Pressure