normal failed to run simulation of mix system

  • syahidah
  • syahidah's Avatar Topic Author
  • Offline
  • Senior Boarder
More
11 years 6 months ago #1184 by syahidah
failed to run simulation of mix system was created by syahidah
Dear all,

i am doing coarse-grained simulation of mix system whereby in my system there are 11 different molecules. and all of them are more than 1. the energy minimization for this system was done before the production run.

command for production run:
~grompp -f mdsemi.mdp -c em1.gro -p cg.top -o md.tpr -n ndx.ndx
~ mdrun -v -s md.tpr -deffnm md

the md.mdp file is as below;
;
; STANDARD MD INPUT OPTIONS FOR MARTINI 2.1
;
; for use with GROMACS 4
;

; RUN CONTROL PARAMETERS
; MARTINI - Most simulations are stable with dt=40 fs, some (especially rings)
; require 20-30 fs.
; The range of time steps used for parametrization is 20-40 fs, using smaller
; time steps is therefore not recommended.
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.030
nsteps = 1000000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
comm-grps = OC OCL OL OLN OM OO OP OS DPPC CMP W

; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f):
nstxout = 5000
nstvout = 5000
nstfout = 0
; Output frequency for energies to log (.log) file and energy (.edr) file:
nstlog = 1000
nstenergy = 1000
; Output frequency and precision for .xtc file:
nstxtcout = 1000
xtc_precision = 100

; NEIGHBORSEARCHING PARAMETERS
; MARTINI - No need for more frequent updates or larger neighborlist cut-off
; due to the use of shifted potential energy functions.
; Neighborlist update frequency:
nstlist = 10
; Neighbor searching algorithm (simple|grid):
ns_type = grid
; Periodic boundary conditions (xyz|none):
pbc = xyz
; Neighborlist cut-off:
rlist = 1.2

; OPTIONS FOR ELECTROSTATICS AND VDW
; MARTINI - VdW and electrostatic interactions are used in their shifted forms.
; Changing to other types of electrostatics will affect the general performance
; of the model.
; Method for doing electrostatics:
coulombtype = Shift
rcoulomb_switch = 0.0
rcoulomb = 1.2
; Dielectric constant (DC) for cut-off or DC of reaction field:
epsilon_r = 15
; Method for doing Van der Waals:
vdw_type = Shift
; Cut-off lengths:
rvdw_switch = 0.9
rvdw = 1.2
; Apply long range dispersion corrections for Energy and Pressure?
DispCorr = No

; OPTIONS FOR WEAK COUPLING ALGORITHMS
; MARTINI -Normal temperature and pressure coupling schemes can be used. It
; is recommended to couple individual groups in your system seperately.
; Temperature coupling:
tcoupl = berendsen
; Groups to couple separately:
tc-grps = OC OCL OL OLN OM OO OP OS DPPC CMP W
; Time constant (ps) and reference temperature (K):
tau_t = 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5
ref_t = 298 298 298 298 298 298 298 298 298 298 298
; Pressure coupling:
Pcoupl = berendsen
Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):

tau_p = 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0
compressibility = 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5
ref_p = 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 298
gen_seed = -1

; OPTIONS FOR BONDS
; MARTINI - For ring systems constraints are defined which are best handled
; using Lincs.
constraints = none
; Type of constraint algorithm:
constraint_algorithm = Lincs
; Do not constrain the start configuration:
unconstrained_start = no
; Highest order in the expansion of the constraint coupling matrix:
lincs_order = 4
; Lincs will write a warning to the stderr if in one step a bond rotates over
; more degrees than:
lincs_warnangle = 30


however, when it came to production run,this error was given,
step 109800, will finish Wed Oct 10 15:12:20 2012vol 0.77! imb F 68%
A list of missing interactions:
Proper Dih. of 4824 missing 1

Molecule type 'CMP'
the first 10 missing interactions, except for exclusions:
Proper Dih. atoms 79 80 81 22 global 2697 2698 2699 2640

Back Off! I just backed up dd_dump_err_0_n1.pdb to ./#dd_dump_err_0_n1.pdb.5#

Back Off! I just backed up dd_dump_err_0_n3.pdb to ./#dd_dump_err_0_n3.pdb.5#

Back Off! I just backed up dd_dump_err_0_n0.pdb to ./#dd_dump_err_0_n0.pdb.5#

Back Off! I just backed up dd_dump_err_0_n2.pdb to ./#dd_dump_err_0_n2.pdb.5#

Program mdrun, VERSION 4.5
Source code file: domdec_top.c, line: 356

Fatal error:
1 of the 11435 bonded interactions could not be calculated because some atoms involved moved further apart than the multi-body cut-off distance (2.84735 nm) or the two-body cut-off distance (2.84735 nm), see option -rdd, for pairs and tabulated bonds also see option -ddcheck
For more information and tips for troubleshooting, please check the GROMACS
website at www.gromacs.org/Documentation/Errors

i try to increase the -rdd value to 3.4 but still the same problem produced.

do you guys have any idea about this problem?

your reply is greatly appreciated,

Thank you so much..

cheers,
syahida

Please Log in or Create an account to join the conversation.

More
11 years 6 months ago #1186 by xavier
Replied by xavier on topic failed to run simulation of mix system
Well well well, It think the problem is that you have constructed your system in a very funny way. It is actually not clear what the exact source of the problem is but there are a few issues to solve in your mdp file first.

- for equilibration you may want to decrease your time step (dt) to 10 fs. You can increase it to 30 fs later on.
- the removal of the center of mass motion (comm-grps) should not be applied on each molecule of molecule type as you are doing it. If you have a bilayer with a aqueous phase it is strongly advised to apply the COM motion removal separately as they can slide one relative to another, but the molecules solvated in the aqueous phase of the lipid bilayer should by associated with the corresponding phase. And you could also start by removing the COM motion of the "system" and later on switch to the decomposed on.
- the same reasoning should apply to the temperature control. Although it should not be a problem in principle you do not want to couple the different molecules separately. as you do. You may want to group them by phase ... Note that on one hand the smaller the group of atoms you use the higher the fluctuations will be and thus the more the coupling will be used to control the temperature of the system; on the other hand the bigger the group the more chances you have to observe transfer of kinetic energy from low frequency motion to high frequency ones, which you want to avoid. So reduce the number of temperature coupling groups.
- your manner to define the pressure coupling is somehow disturbing or confusing. The semi-isotropic scheme need only two values! In case that helps it is a global variable, which cannot be controlled per molecule but should be applied to the complete system.

Once you have fixed all these issue things might go better ...

syahidah wrote: Dear all,

i am doing coarse-grained simulation of mix system whereby in my system there are 11 different molecules. and all of them are more than 1. the energy minimization for this system was done before the production run.

command for production run:
~grompp -f mdsemi.mdp -c em1.gro -p cg.top -o md.tpr -n ndx.ndx
~ mdrun -v -s md.tpr -deffnm md

the md.mdp file is as below;
;
; STANDARD MD INPUT OPTIONS FOR MARTINI 2.1
;
; for use with GROMACS 4
;

; RUN CONTROL PARAMETERS
; MARTINI - Most simulations are stable with dt=40 fs, some (especially rings)
; require 20-30 fs.
; The range of time steps used for parametrization is 20-40 fs, using smaller
; time steps is therefore not recommended.
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.030
nsteps = 1000000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
comm-grps = OC OCL OL OLN OM OO OP OS DPPC CMP W

; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f):
nstxout = 5000
nstvout = 5000
nstfout = 0
; Output frequency for energies to log (.log) file and energy (.edr) file:
nstlog = 1000
nstenergy = 1000
; Output frequency and precision for .xtc file:
nstxtcout = 1000
xtc_precision = 100

; NEIGHBORSEARCHING PARAMETERS
; MARTINI - No need for more frequent updates or larger neighborlist cut-off
; due to the use of shifted potential energy functions.
; Neighborlist update frequency:
nstlist = 10
; Neighbor searching algorithm (simple|grid):
ns_type = grid
; Periodic boundary conditions (xyz|none):
pbc = xyz
; Neighborlist cut-off:
rlist = 1.2

; OPTIONS FOR ELECTROSTATICS AND VDW
; MARTINI - VdW and electrostatic interactions are used in their shifted forms.
; Changing to other types of electrostatics will affect the general performance
; of the model.
; Method for doing electrostatics:
coulombtype = Shift
rcoulomb_switch = 0.0
rcoulomb = 1.2
; Dielectric constant (DC) for cut-off or DC of reaction field:
epsilon_r = 15
; Method for doing Van der Waals:
vdw_type = Shift
; Cut-off lengths:
rvdw_switch = 0.9
rvdw = 1.2
; Apply long range dispersion corrections for Energy and Pressure?
DispCorr = No

; OPTIONS FOR WEAK COUPLING ALGORITHMS
; MARTINI -Normal temperature and pressure coupling schemes can be used. It
; is recommended to couple individual groups in your system seperately.
; Temperature coupling:
tcoupl = berendsen
; Groups to couple separately:
tc-grps = OC OCL OL OLN OM OO OP OS DPPC CMP W
; Time constant (ps) and reference temperature (K):
tau_t = 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5 2.5
ref_t = 298 298 298 298 298 298 298 298 298 298 298
; Pressure coupling:
Pcoupl = berendsen
Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):

tau_p = 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0
compressibility = 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5 4.5e-5
ref_p = 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 298
gen_seed = -1

; OPTIONS FOR BONDS
; MARTINI - For ring systems constraints are defined which are best handled
; using Lincs.
constraints = none
; Type of constraint algorithm:
constraint_algorithm = Lincs
; Do not constrain the start configuration:
unconstrained_start = no
; Highest order in the expansion of the constraint coupling matrix:
lincs_order = 4
; Lincs will write a warning to the stderr if in one step a bond rotates over
; more degrees than:
lincs_warnangle = 30


however, when it came to production run,this error was given,
step 109800, will finish Wed Oct 10 15:12:20 2012vol 0.77! imb F 68%
A list of missing interactions:
Proper Dih. of 4824 missing 1

Molecule type 'CMP'
the first 10 missing interactions, except for exclusions:
Proper Dih. atoms 79 80 81 22 global 2697 2698 2699 2640

Back Off! I just backed up dd_dump_err_0_n1.pdb to ./#dd_dump_err_0_n1.pdb.5#

Back Off! I just backed up dd_dump_err_0_n3.pdb to ./#dd_dump_err_0_n3.pdb.5#

Back Off! I just backed up dd_dump_err_0_n0.pdb to ./#dd_dump_err_0_n0.pdb.5#

Back Off! I just backed up dd_dump_err_0_n2.pdb to ./#dd_dump_err_0_n2.pdb.5#


Program mdrun, VERSION 4.5
Source code file: domdec_top.c, line: 356

Fatal error:
1 of the 11435 bonded interactions could not be calculated because some atoms involved moved further apart than the multi-body cut-off distance (2.84735 nm) or the two-body cut-off distance (2.84735 nm), see option -rdd, for pairs and tabulated bonds also see option -ddcheck
For more information and tips for troubleshooting, please check the GROMACS
website at www.gromacs.org/Documentation/Errors

i try to increase the -rdd value to 3.4 but still the same problem produced.

do you guys have any idea about this problem?

your reply is greatly appreciated,

Thank you so much..

cheers,
syahida

Please Log in or Create an account to join the conversation.

Time to create page: 0.103 seconds