Hi all,
I want to test the g_dipoles tool and made a 20 ns simulation with SPC/E
water. Sadly the result for the dielectric constant (e_sim=109.713)
looks very different from experimental (e_exp=78) and reported values
for SPCE (e_paper=70). I put the command line and the output below.
Gromacs 4.5.3, SPC/E water, 20ns, NPT, 300K, 1bar
g_dipoles -corr total -c dipoles -f md_20ns.xtc -s md_20ns.tpr -b 0 -e 20000
#==========================================
Using 5 as mu_max and -1 as the dipole moment.
WARNING: EpsilonRF = 0.0, this really means EpsilonRF = infinity
Selected 0: 'System'
There are 1024 molecules in the selection
Using Volume from topology: 30.2049 nm3
Average volume over run is 30.2103
t0 0, t 20000, teller 10001
Dipole moment (Debye)
---------------------
Average = 2.4883 Std. Dev. = 0.1342 Error = 0.0000
The following averages for the complete trajectory have been calculated:
Total< M_x> = 1.45986 Debye
Total< M_y> = 5.71202 Debye
Total< M_z> = 3.39196 Debye
Total< M_x2 > = 10964.8 Debye2
Total< M_y2 > = 10794.4 Debye2
Total< M_z2 > = 10762.4 Debye2
Total< |M|^2> = 32521.6 Debye2
Total |< M>|^2 = 46.2638 Debye2
< |M|^2> - |< M>|^2 = 32475.4 Debye2
Finite system Kirkwood g factor G_k = 5.12208
Infinite system Kirkwood g factor g_k = 3.43028
Epsilon = 109.713
#==========================================
I would be really grateful if anybody could point me to my error. I
attached my .mdp file in case the error lies there.
Best regards
Tom
;================== run control ===================================
integrator = md
tinit = 0
dt = 0.002 ; time step [ps]
nsteps = 10000000 ; number of steps
comm_mode = Linear ; Remove center of mass translation
nstcomm = 10 ; reset c.o.m. motion
;================== Output control ===================================
nstxout = 1000 ; write coords
nstvout = 1000 ; write velocities
nstfout = 0 ; do NOT write forces
nstlog = 1000 ; print to logfile
nstenergy = 1000
nstxtcout = 1000 ; write coords to xtc-trajectory file
xtc-precision = 1000
;================== Neighborsearching and short-range nonbonded interactions
nstlist = 10 ; update pairlist
ns_type = grid ; pairlist method
pbc = xyz ; periodic boundary conditions
rlist = 1.3 ; cut-off for neighborsearching - rlist
>= rvdw+ (0.1 till 0.3)
;================== Electrostatics
coulombtype = PME
rcoulomb = 1.3 ; cut-off for coulomb - rcoulomb = rlist
fourierspacing = 0.112 ; nm
pme_order = 4
ewald_rtol = 1e-05
optimize_fft = no
epsilon_surface = 0
;================== dispersion correction
dispcorr = EnerPres
;================== Van der Waals
vdw-type = shift ; switch
rvdw_switch = 0.9
rvdw = 1.0 ; cut-off for vdw
;================== OPTIONS FOR WEAK COUPLING ALGORITHMS ===============
Tcoupl = v-rescale ; temperature coupling
(v-rescale)
tc-grps = System
ref_t = 298.15
tau_t = 1.0
Pcoupl = berendsen ; pressure bath (berendsen)
Pcoupltype = isotropic ; pressure geometry
tau_p = 1.0 ; p-coupoling time
compressibility = 4.5e-5
ref_p = 1.01325
;================== GENERATE VELOCITIES FOR STARTUP RUN =================
gen_vel = yes ; generate initial velocities
gen_temp = 298.15 ; initial temperature
gen_seed = -1 ; random seed
continuation = no ; since we came from nvt simulation
;================== OPTIONS FOR BONDS ===================================
constraints = h-bonds
constraint_algorithm = lincs
shake_tol = 0.00001
lincs_order = 4
lincs_warnangle = 30
morse = no
lincs_iter = 2
--
gmx-users mailing list gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
