- Posts: 11
Incorrect density
- edroaldo
-
Topic Author
- Offline
- Fresh Boarder
Less
More
12 years 4 months ago - 12 years 4 months ago #859
by edroaldo
Incorrect density was created by edroaldo
Dear,
I did a simulation in a system with a lipid bilayer (1200 lipids), ions, water and one nanoparticle. The box dimensions are 20x20x18 and I add about 49000 CG water beads in this system. After a 50ns simulation, I use the gromacs tool g_density to calc the density of some molecules in my system. The density for water was above of 1000 kg/m³ and the papers that I am following show results below to this value. One of this papers says that about 60000 CG water beads was added in a system that contains 1152 lipids. Someone please give me some direction? Do I need remove water from my system?
Thank you very much.
Best regards,
Edroaldo.
I did a simulation in a system with a lipid bilayer (1200 lipids), ions, water and one nanoparticle. The box dimensions are 20x20x18 and I add about 49000 CG water beads in this system. After a 50ns simulation, I use the gromacs tool g_density to calc the density of some molecules in my system. The density for water was above of 1000 kg/m³ and the papers that I am following show results below to this value. One of this papers says that about 60000 CG water beads was added in a system that contains 1152 lipids. Someone please give me some direction? Do I need remove water from my system?
Thank you very much.
Best regards,
Edroaldo.
Please Log in or Create an account to join the conversation.
- djurre
-
- Offline
- Admin
Less
More
- Posts: 272
12 years 4 months ago #862
by djurre
Replied by djurre on topic Incorrect density
Dear Edroaldo,
I'm guessing your taking the maximum peak hight from the density profile to be the water density. Are you averaging over the the whole 50ns or just using the last snapshot?
In the first case, the system might not be relaxed at first, you should leave some nano seconds of. In the second case you might be looking at a fluctuation.
Is the pressure in the system correct?
Groetnis,
Djurre
I'm guessing your taking the maximum peak hight from the density profile to be the water density. Are you averaging over the the whole 50ns or just using the last snapshot?
In the first case, the system might not be relaxed at first, you should leave some nano seconds of. In the second case you might be looking at a fluctuation.
Is the pressure in the system correct?
Groetnis,
Djurre
Please Log in or Create an account to join the conversation.
- edroaldo
-
Topic Author
- Offline
- Fresh Boarder
Less
More
- Posts: 11
12 years 4 months ago #864
by edroaldo
Replied by edroaldo on topic Incorrect density
Dear Djurre,
Thank you for your replay!
I am using g_density as: g_density -f dyn.xtc -s dyn.tpr -n index.ndx -o denW.xvg. So, I averaging over the whole 50nm, right? The curve that I get for density profile looks like a "U" format. I have some regions above 1000kg/m³ and other regions below this value. I set up the pressure coupling with the parameters
Pcoupl = Berendsen
pcoupltype = anisotropic
tau_p = 5.0 5.0 5.0 0 0 0
compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0
ref_p = 1.0 1.0 1.0 0 0 0
I look at the pressure with g_energy and looks like the pressure is wrong because the curve varies around zero. Any suggestion?
Thank you very much.
Thank you for your replay!
I am using g_density as: g_density -f dyn.xtc -s dyn.tpr -n index.ndx -o denW.xvg. So, I averaging over the whole 50nm, right? The curve that I get for density profile looks like a "U" format. I have some regions above 1000kg/m³ and other regions below this value. I set up the pressure coupling with the parameters
Pcoupl = Berendsen
pcoupltype = anisotropic
tau_p = 5.0 5.0 5.0 0 0 0
compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0
ref_p = 1.0 1.0 1.0 0 0 0
I look at the pressure with g_energy and looks like the pressure is wrong because the curve varies around zero. Any suggestion?
Thank you very much.
Please Log in or Create an account to join the conversation.
- djurre
-
- Offline
- Admin
Less
More
- Posts: 272
12 years 4 months ago #865
by djurre
Replied by djurre on topic Incorrect density
Dear Edroaldo,
The pressure might fluctuate, a lot actually, that normal, as long as the average is 1 bar.
I don't now how much the pressure is above 1000 kg/m3. If it is within a few percent, it is could be due to lack of sampling.
The pressure might fluctuate, a lot actually, that normal, as long as the average is 1 bar.
I don't now how much the pressure is above 1000 kg/m3. If it is within a few percent, it is could be due to lack of sampling.
Please Log in or Create an account to join the conversation.
- edroaldo
-
Topic Author
- Offline
- Fresh Boarder
Less
More
- Posts: 11
12 years 4 months ago #866
by edroaldo
Replied by edroaldo on topic Incorrect density
Hi, the average value of the pressure is 1.02bar, so the pressue is ok. The density is around 1150 kg/m³ in my system.
In relation to sampling, can you please tell me where I change this? I am using a configuration file similar to the MARTINI web site for DPPC membrane.
Thank you very much!
In relation to sampling, can you please tell me where I change this? I am using a configuration file similar to the MARTINI web site for DPPC membrane.
Thank you very much!
Please Log in or Create an account to join the conversation.
- djurre
-
- Offline
- Admin
Less
More
- Posts: 272
12 years 4 months ago #869
by djurre
Replied by djurre on topic Incorrect density
With sampling problem I mean you have to simulate (sample) longer. Simulate 1000 ns instead of 50.
Please Log in or Create an account to join the conversation.
- skavyani
-
- Offline
- Fresh Boarder
Less
More
- Posts: 3
11 years 3 months ago #1330
by skavyani
Replied by skavyani on topic Incorrect density
Dears,
I have this poroblem too.
Recently i simulated a system that contain a dendrimer molecule in a box of water. when i was trying to calculate density with g_density, the calculated density has a "U" shape but it's value is more than expected !.
After that, i simulate a single box with pure water but i gave the same result !!!.
the pressure average is 0.9999 bar and everything is fine except density !. i need correct density value because it can define my dendric system behavior.The spoilers are my .mdp and .top files for pure water system.
(NOTE:
-em tow times with two different methods.
-mdrun with 8 nodes.
)
.top
nvt :
npt :
md
I turned to CG martini force.f for less calculation time. unfortunately i can't do 1000ns simulation. this is so much. because best t.step that gave me a fine result with no system belowing up is 0.008ps for at least 50000 beads, with 8 nodes (more nodes make mdrun "DD" crash..i test many em setting and run methods ( as i can ) to use more nodes but i can't).
I would appreciate if you could help me to handle the problem.
Thank you very much.
Best regards,
S. Kavyani
I have this poroblem too.
Recently i simulated a system that contain a dendrimer molecule in a box of water. when i was trying to calculate density with g_density, the calculated density has a "U" shape but it's value is more than expected !.
After that, i simulate a single box with pure water but i gave the same result !!!.
the pressure average is 0.9999 bar and everything is fine except density !. i need correct density value because it can define my dendric system behavior.The spoilers are my .mdp and .top files for pure water system.
(NOTE:
-em tow times with two different methods.
-mdrun with 8 nodes.
)
.top
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
#include "martini_v2.2.itp"
[ system ]
; name
water
[ molecules ]
; name number
W 3175
[ system ]
; name
water
[ molecules ]
; name number
W 3175
nvt :
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
;
; 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.
;
define =-DPOSRES
;
;
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.0080
nsteps = 80000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f):
nstxout = 500
nstvout = 500
nstfout = 0
; Output frequency for energies to log (.log) file and energy (.edr) file:
nstlog = 100
nstenergy = 100
; 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:
;
;
;
;;continuation = no ; first dynamics run
;;constraint_algorithm = lincs ; holonomic constraints
;constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
;
;
;
;
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
;
;
;
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 = switch
; 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 = system
; Time constant (ps) and reference temperature (K):
tau_t = 1.0
ref_t = 310
; Pressure coupling:
Pcoupl = no
;Pcoupltype = isotropic
; Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):
;tau_p = 3.0
;compressibility = 3e-5
;ref_p = 1.0
; tau_p = 3.0 3.0
; compressibility = 3e-5 3e-5
; ref_p = 1.0 1.0
; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 320
gen_seed = 473529
; 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
; 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.
;
define =-DPOSRES
;
;
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.0080
nsteps = 80000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f):
nstxout = 500
nstvout = 500
nstfout = 0
; Output frequency for energies to log (.log) file and energy (.edr) file:
nstlog = 100
nstenergy = 100
; 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:
;
;
;
;;continuation = no ; first dynamics run
;;constraint_algorithm = lincs ; holonomic constraints
;constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
;
;
;
;
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
;
;
;
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 = switch
; 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 = system
; Time constant (ps) and reference temperature (K):
tau_t = 1.0
ref_t = 310
; Pressure coupling:
Pcoupl = no
;Pcoupltype = isotropic
; Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):
;tau_p = 3.0
;compressibility = 3e-5
;ref_p = 1.0
; tau_p = 3.0 3.0
; compressibility = 3e-5 3e-5
; ref_p = 1.0 1.0
; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 320
gen_seed = 473529
; 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
npt :
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
;
; 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.
;
define =-DPOSRES
;
;
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.0080
nsteps = 300000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f):
nstxout = 500
nstvout = 500
nstfout = 0
; Output frequency for energies to log (.log) file and energy (.edr) file:
nstlog = 100
nstenergy = 100
; 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:
;
;
;
;continuation = no ; Restarting after NVT
;constraint_algorithm = lincs ; holonomic constraints
;constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
;
;
;
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
;
;
;
;
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 = system
; Time constant (ps) and reference temperature (K):
tau_t = 1
ref_t = 310
; Pressure coupling:
Pcoupl = Parrinello-Rahman
Pcoupltype = isotropic
;Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):
tau_p = 3.0
compressibility = 3e-5
ref_p = 1.0
; tau_p = 3.0 3.0
; compressibility = 3e-5 3e-5
; ref_p = 1.0 1.0
; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 320
gen_seed = 473529
; 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
; 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.
;
define =-DPOSRES
;
;
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.0080
nsteps = 300000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f):
nstxout = 500
nstvout = 500
nstfout = 0
; Output frequency for energies to log (.log) file and energy (.edr) file:
nstlog = 100
nstenergy = 100
; 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:
;
;
;
;continuation = no ; Restarting after NVT
;constraint_algorithm = lincs ; holonomic constraints
;constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
;
;
;
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
;
;
;
;
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 = system
; Time constant (ps) and reference temperature (K):
tau_t = 1
ref_t = 310
; Pressure coupling:
Pcoupl = Parrinello-Rahman
Pcoupltype = isotropic
;Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):
tau_p = 3.0
compressibility = 3e-5
ref_p = 1.0
; tau_p = 3.0 3.0
; compressibility = 3e-5 3e-5
; ref_p = 1.0 1.0
; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 320
gen_seed = 473529
; 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
md
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
;
; 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.
;
define =-DPOSRES
;
;
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.0050
nsteps = 2000000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
; 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:
;
;
;
;continuation = yes ; Restarting after NVT
;constraint_algorithm = lincs ; holonomic constraints
;constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
;
;
;
;
;
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
;
;
;
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 = switch
; 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 = system
; Time constant (ps) and reference temperature (K):
tau_t = 1.0
ref_t = 310
; Pressure coupling:
Pcoupl = Parrinello-Rahman
Pcoupltype = isotropic
;Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):
tau_p = 3.0
compressibility = 3e-5
ref_p = 1.0
; tau_p = 3.0 3.0
; compressibility = 3e-5 3e-5
; ref_p = 1.0 1.0
; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 320
gen_seed = 473529
; 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
; 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.
;
define =-DPOSRES
;
;
integrator = md
; Start time and timestep in ps:
tinit = 0.0
dt = 0.0050
nsteps = 2000000
; nsteps = 500000
; Number of steps for center of mass motion removal:
nstcomm = 1
; 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:
;
;
;
;continuation = yes ; Restarting after NVT
;constraint_algorithm = lincs ; holonomic constraints
;constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained
;
;
;
;
;
coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
;
;
;
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 = switch
; 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 = system
; Time constant (ps) and reference temperature (K):
tau_t = 1.0
ref_t = 310
; Pressure coupling:
Pcoupl = Parrinello-Rahman
Pcoupltype = isotropic
;Pcoupltype = semiisotropic
; Time constant (ps), compressibility (1/bar) and reference P (bar):
tau_p = 3.0
compressibility = 3e-5
ref_p = 1.0
; tau_p = 3.0 3.0
; compressibility = 3e-5 3e-5
; ref_p = 1.0 1.0
; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no
gen_temp = 320
gen_seed = 473529
; 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
I turned to CG martini force.f for less calculation time. unfortunately i can't do 1000ns simulation. this is so much. because best t.step that gave me a fine result with no system belowing up is 0.008ps for at least 50000 beads, with 8 nodes (more nodes make mdrun "DD" crash..i test many em setting and run methods ( as i can ) to use more nodes but i can't).
I would appreciate if you could help me to handle the problem.
Thank you very much.
Best regards,
S. Kavyani
Please Log in or Create an account to join the conversation.
- djurre
-
- Offline
- Admin
Less
More
- Posts: 272
11 years 3 months ago #1334
by djurre
Replied by djurre on topic Incorrect density
If you're simulating a box of water you should not need two EM steps and you do not need to do a NVT relaxation. Also if you can not simulated at a timestep longer then 0.008 ps, something is wrong.
When you generate the box of water, do you use genbox? Do you set -vdwd 0.21 ? If you open the initial box of water in, for example, VMD do you see any difference in density in a certain area? Genbox has the tendency to make two overlapping cubes inside large CG water systems.
Fix it by first making a much larger box using genconf and use that as input for genbox.
When you generate the box of water, do you use genbox? Do you set -vdwd 0.21 ? If you open the initial box of water in, for example, VMD do you see any difference in density in a certain area? Genbox has the tendency to make two overlapping cubes inside large CG water systems.
Fix it by first making a much larger box using genconf and use that as input for genbox.
Please Log in or Create an account to join the conversation.
- skavyani
-
- Offline
- Fresh Boarder
Less
More
- Posts: 3
11 years 3 months ago #1335
by skavyani
Replied by skavyani on topic Incorrect density
Dear Djurre,
Many thanks for your reply.
For two em and nvt steps you are right. i was scrupulous!
I'm sorry, i was saying "same result" didn't mean that the pure water system has "U" shape. i just mean that the calculated density has a value more than expected.for pw sys density profile has a flat shape with but with wrong value, a value around 1200.
For miss understanding > i'm sorry. that's my fault
For p.w system i can use longer t.s but for dendric sys, i can't.
something is wrong with "W" beads.
why for such a simple system like pure water, density has such a ERROR?
I use "water.gro" that has been uploaded in matrini w.site, with 400 "w" beads (is for water.gro in martini-w.site -vdwd set to 0.21 ?! ) and then use genbox to fill larger box with this command:NOTE:
in vmd the box is uniform. water_box and empty_box.gro are in size of 7 x 7 x 7 nm dim. the output has 3175 water beads in 7 7 7 box dim.
NOTE:
Best regards,
S. Kavyani
Many thanks for your reply.
For two em and nvt steps you are right. i was scrupulous!
I'm sorry, i was saying "same result" didn't mean that the pure water system has "U" shape. i just mean that the calculated density has a value more than expected.for pw sys density profile has a flat shape with but with wrong value, a value around 1200.
For miss understanding > i'm sorry. that's my fault
For p.w system i can use longer t.s but for dendric sys, i can't.
something is wrong with "W" beads.
why for such a simple system like pure water, density has such a ERROR?
I use "water.gro" that has been uploaded in matrini w.site, with 400 "w" beads (is for water.gro in martini-w.site -vdwd set to 0.21 ?! ) and then use genbox to fill larger box with this command:
genbox -cs water.gro -cp empty_box.gro -o water_box.gro -nmol 10000
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
the "water_box.gro" with or without -nmol option is same.
in vmd the box is uniform. water_box and empty_box.gro are in size of 7 x 7 x 7 nm dim. the output has 3175 water beads in 7 7 7 box dim.
NOTE:
Warning: Spoiler!
[ Click to expand ]
[ Click to hide ]
As you said, i try genconf to make water system for this command :and the output has 3200 water beads in 7.28856 7.28856 7.28856 box dim. fortunately for this system seems that the genbox output has no overlapping.
genconf -f water.gro -nbox 2 2 2 -o w_system.groBest regards,
S. Kavyani
Please Log in or Create an account to join the conversation.
- skavyani
-
- Offline
- Fresh Boarder
Less
More
- Posts: 3
11 years 3 months ago #1337
by skavyani
Replied by skavyani on topic Incorrect density
Dears,
I apologize to anyone who wasted his time reading my first post.
I ignored that the tau_p parameter can affect my density profile. the value that has been used was 3.0 for "parrinello-rahman " and it better be around 0.5.
S. Kavyani
I apologize to anyone who wasted his time reading my first post.
I ignored that the tau_p parameter can affect my density profile. the value that has been used was 3.0 for "parrinello-rahman " and it better be around 0.5.
S. Kavyani
Please Log in or Create an account to join the conversation.
Time to create page: 0.172 seconds