- Posts: 210
Simulation blows up when the protein crosses periodic boundary
- h_and
- Topic Author
- Visitor
6 years 5 months ago #7412
by h_and
Simulation blows up when the protein crosses periodic boundary was created by h_and
Hello,
I simulate modelled proteins in bilayer using elnedyn 2.2p. I energy minimize my system and do backbone-restrained NVT and NPT equilibrations with gradually increasing time steps. When I start MD production there is no problem with my system but blows up after several nanoseconds (usually around 200-400ns).
I realized that the LINCS warnings preceeding the explosion happen when my protein starts to cross the periodic boundary. I also observed that at least one of the affected beads are S-beads (histidine or triptophan).
I tried to change several options to circumvent this problem: longer energy minimization, equilibration, decreasing time steps. When I martinize my protein I receive the warning that elnedyn topologies are stable only if the mass of S-beads are set to 72. I did that but did not help. The only way I could avoid blowing up is using Berendsen pressure coupling instead of Parrinello-Rahman. I experienced this on several kind of membrane-embed proteins, no matter how it is generated (Charmm-GUI or martinize+insane).
What could be wrong with my system? Is there anything that I am missing?
MDP options are below. I can share my itp and tpr files if necessary.
I always run my multi-threaded simulations with the parameter -rdd 1.4.
Gromacs version: 5.1.4
martinize.py version: 2.6
particle definition: martini_v2.2refP.itp
I simulate modelled proteins in bilayer using elnedyn 2.2p. I energy minimize my system and do backbone-restrained NVT and NPT equilibrations with gradually increasing time steps. When I start MD production there is no problem with my system but blows up after several nanoseconds (usually around 200-400ns).
I realized that the LINCS warnings preceeding the explosion happen when my protein starts to cross the periodic boundary. I also observed that at least one of the affected beads are S-beads (histidine or triptophan).
I tried to change several options to circumvent this problem: longer energy minimization, equilibration, decreasing time steps. When I martinize my protein I receive the warning that elnedyn topologies are stable only if the mass of S-beads are set to 72. I did that but did not help. The only way I could avoid blowing up is using Berendsen pressure coupling instead of Parrinello-Rahman. I experienced this on several kind of membrane-embed proteins, no matter how it is generated (Charmm-GUI or martinize+insane).
What could be wrong with my system? Is there anything that I am missing?
MDP options are below. I can share my itp and tpr files if necessary.
I always run my multi-threaded simulations with the parameter -rdd 1.4.
Gromacs version: 5.1.4
martinize.py version: 2.6
particle definition: martini_v2.2refP.itp
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
; PREPROCESSING
; RUN CONTROL
integrator = md
dt = 0.02
nsteps = 10000000 ;200 ns
comm-mode = Linear
comm-grps = SYSTEM
; OUTPUT CONTROL
nstlog = 25000
nstenergy = 25000
nstxout = 25000
nstxout_compressed = 25000
nstvout = 0
nstfout = 0
energygrps = PROTEIN MEMBRANE SOLVENT
; NEIGHBOUR SEARCHING
cutoff-scheme = Verlet
nstlist = 20
; ELECTROSTATICS
coulombtype = Cut-off
coulomb-modifier = Potential-shift
rcoulomb = 1.1
epsilon-r = 2.5
; VAN DER VAALS
vdwtype = Cut-off
vdw-modifier = Potential-shift
rvdw = 1.1
DispCorr = No
; TEMPERATURE COUPLING
tcoupl = v-rescale
tc-grps = PROTEIN_MEMBRANE SOLVENT
tau-t = 1.0 1.0
ref-t = 310 310
; PRESSURE COUPLING
pcoupl = Parrinello-Rahman
pcoupltype = Semiisotropic
compressibility = 3e-4 3e-4
ref-p = 1 1
tau-p = 12.0
refcoord-scaling = com
; BONDS
constraints = none
lincs-order = 6
lincs-iter = 8
lincs-warnangle = 90
; RUN CONTROL
integrator = md
dt = 0.02
nsteps = 10000000 ;200 ns
comm-mode = Linear
comm-grps = SYSTEM
; OUTPUT CONTROL
nstlog = 25000
nstenergy = 25000
nstxout = 25000
nstxout_compressed = 25000
nstvout = 0
nstfout = 0
energygrps = PROTEIN MEMBRANE SOLVENT
; NEIGHBOUR SEARCHING
cutoff-scheme = Verlet
nstlist = 20
; ELECTROSTATICS
coulombtype = Cut-off
coulomb-modifier = Potential-shift
rcoulomb = 1.1
epsilon-r = 2.5
; VAN DER VAALS
vdwtype = Cut-off
vdw-modifier = Potential-shift
rvdw = 1.1
DispCorr = No
; TEMPERATURE COUPLING
tcoupl = v-rescale
tc-grps = PROTEIN_MEMBRANE SOLVENT
tau-t = 1.0 1.0
ref-t = 310 310
; PRESSURE COUPLING
pcoupl = Parrinello-Rahman
pcoupltype = Semiisotropic
compressibility = 3e-4 3e-4
ref-p = 1 1
tau-p = 12.0
refcoord-scaling = com
; BONDS
constraints = none
lincs-order = 6
lincs-iter = 8
lincs-warnangle = 90
Please Log in or Create an account to join the conversation.
- peterkroon
- Offline
- Gold Boarder
Less
More
6 years 5 months ago #7415
by peterkroon
Replied by peterkroon on topic Simulation blows up when the protein crosses periodic boundary
Heyhey,
Nothing looks really wrong at first glance, so that's good. Since you're using 2.2P you should probably consider using PME (or maybe RF) for the electrostatics and polarizable water.
Being close to the periodic boundary should not make any difference (otherwise it's a Gromacs bug). Are you sure the pressure is converged after your equilibration? You could try increasing tau-p to, e.g. 24. And personally I'd set -rdd 1.6, but I don't think that'll matter here.
Peter
Nothing looks really wrong at first glance, so that's good. Since you're using 2.2P you should probably consider using PME (or maybe RF) for the electrostatics and polarizable water.
Being close to the periodic boundary should not make any difference (otherwise it's a Gromacs bug). Are you sure the pressure is converged after your equilibration? You could try increasing tau-p to, e.g. 24. And personally I'd set -rdd 1.6, but I don't think that'll matter here.
Peter
Please Log in or Create an account to join the conversation.
- h_and
- Topic Author
- Visitor
6 years 4 months ago #7445
by h_and
Replied by h_and on topic Simulation blows up when the protein crosses periodic boundary
Setting tau-p to 24 helped, thank you very much.
To answer your question: Yes, the pressure converged, the crash actually happened long-long time after equilibration (even microseconds sometimes). The only common pattern in my crashed simulations was the protein crossing the periodic boundary. It seems that higher tau-p solves this problem.
To answer your question: Yes, the pressure converged, the crash actually happened long-long time after equilibration (even microseconds sometimes). The only common pattern in my crashed simulations was the protein crossing the periodic boundary. It seems that higher tau-p solves this problem.
Please Log in or Create an account to join the conversation.
Time to create page: 0.092 seconds