[gmx-users] Is anyone also using lammps?s

2009-10-27 Thread Peng Yi


I am trying to simulate alkane melt and found out that gromacs and lammps 
gave different results, particularly the bonded interaction energy.

I wonder if anyone has such experience.  Thanks,
___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php


Re: [gmx-users] Is anyone also using lammps?s

2009-10-27 Thread Mark Abraham

Peng Yi wrote:


I am trying to simulate alkane melt and found out that gromacs and 
lammps gave different results, particularly the bonded interaction energy.

I wonder if anyone has such experience.  Thanks,


Even two installations of the same version of GROMACS can give different 
results. The question is whether when using comparable model physics you 
observe the same ensemble averages.


Mark
___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php


Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread Peng Yi


On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:


I am trying to simulate alkane melt and found out that gromacs and lammps 
gave different results, particularly the bonded interaction energy.

I wonder if anyone has such experience.  Thanks,


Even two installations of the same version of GROMACS can give different 
results. The question is whether when using comparable model physics you 
observe the same ensemble averages.


Mark


Hi, Mark,

Thanks for reply!  The difference is statistically significant.  And I am
wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison please
see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  There are
bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+corr:-9238-9350
Volume(nm^3): 66.7 66.5

where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
Elj+corr is the total LJ energy including tail correction.

However, if I use integration time 2fs, Lammps results do not change 
much, but Gromacs results changed a lot:


Lammps   Gromacs
Ebond(kJ/mol):2133 2700 
Eangle:   1799 1640

Etors:2552 2200
Elj+corr:-9292-9886 
Volume:   66.7 64.0


The results given by Lammps is more reasonable because the Ebond should
be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
240 n-octanes have total 1680 bonds and 1440 angles.

The bond and angle interactions are both harmonic functions.  Bond
interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
ossilation period 16 fs.

Is there something related to the integrator?

Here I attached my grompp.mdp and topol.top files.

##
grompp.mdp
##

; VARIOUS PREPROCESSING OPTIONS
title= Yo
cpp  = /usr/bin/cpp
include  = 
define   =


; RUN CONTROL PARAMETERS
integrator   = md
tinit= 0
dt   = 0.001
nsteps   = 200
init_step= 0
comm-mode= Linear
nstcomm  = 1
comm-grps=

; OUTPUT CONTROL OPTIONS
nstxout  = 5000
nstvout  = 5000
nstfout  = 5000
nstcheckpoint= 1
nstlog   = 1000
nstenergy= 1000
nstxtcout= 5000
xtc-precision= 1000
xtc-grps = 
energygrps   =


; NEIGHBORSEARCHING PARAMETERS
nstlist  = 10
ns_type  = grid
pbc  = xyz
rlist= 1.0025
domain-decomposition = no

; OPTIONS FOR ELECTROSTATICS AND VDW
coulombtype  = Cut-off
rcoulomb-switch  = 0
rcoulomb = 1.0025
epsilon-r= 1
vdw-type = Cut-off
rvdw-switch  = 0		; default 
rvdw = 1.0025	; default 1 nm

DispCorr = EnerPres
;table-extension  = 1.5
fourierspacing   = 0.12
fourier_nx   = 0
fourier_ny   = 0
fourier_nz   = 0
pme_order= 4
ewald_rtol   = 1e-05
ewald_geometry   = 3d
epsilon_surface  = 0
optimize_fft = no


; OPTIONS FOR WEAK COUPLING ALGORITHMS
Tcoupl   = nose-hoover
tc-grps  = System
tau_t= 0.02
ref_t= 300.0
Pcoupl   = Parrinello-Rahman
Pcoupltype   = isotropic
tau_p= 2.0
compressibility  = 4.5e-5
ref_p= 1.0
andersen_seed= 815131

; GENERATE VELOCITIES FOR STARTUP RUN
gen_vel  = yes
gen_temp = 300
gen_seed = 2009

; OPTIONS FOR BONDS 
constraints  = none

constraint-algorithm = Lincs
unconstrained-start  = no
Shake-SOR= no
shake-tol= 1e-04
lincs-order  = 4
lincs-iter   = 1
lincs-warnangle  = 30
morse= no

; ENERGY GROUP EXCLUSIONS
; Pairs of energy groups for which all non-

Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread aherz
Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs can be
rather bad for single precision.

Alex


Peng Yi schrieb:
>
> On Wed, 28 Oct 2009, Mark Abraham wrote:
>
>> Peng Yi wrote:
>>>
>>> I am trying to simulate alkane melt and found out that gromacs and
>>> lammps gave different results, particularly the bonded interaction
>>> energy.
>>> I wonder if anyone has such experience.  Thanks,
>>
>> Even two installations of the same version of GROMACS can give
>> different results. The question is whether when using comparable
>> model physics you observe the same ensemble averages.
>>
>> Mark
>
> Hi, Mark,
>
> Thanks for reply!  The difference is statistically significant.  And I am
> wondering if it is caused by the integrator: Leap-frog for Gromacs and
> Velocity-verlet for Lammps.  Detail description of the comparison please
> see below:
>
> It is an NPT simulation of a melt of 240 n-octane molecules using
> united-atom model, i.e., CHx group is considered as one atom.  There are
> bond, angle, torsion and LJ interactions.  T=300K and P=1atm.
>
> Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
> nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
> for
> thermostat and barostat are 0.02ps and 2.0ps, respectively.
>
> If I use integration time 1fs, Lammps and Gromacs gave consistent
> results:
> Lammps   Gromacs
> Ebond(kJ/mol):2092 2146
> Eangle:   1757 1760
> Etors:2510 2500
> Elj+corr:-9238-9350
> Volume(nm^3): 66.7 66.5
>
> where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
> Elj+corr is the total LJ energy including tail correction.
>
> However, if I use integration time 2fs, Lammps results do not change
> much, but Gromacs results changed a lot:
>
> Lammps   Gromacs
> Ebond(kJ/mol):2133 2700 Eangle:  
> 1799 1640
> Etors:2552 2200
> Elj+corr:-9292-9886 Volume:  
> 66.7 64.0
>
> The results given by Lammps is more reasonable because the Ebond should
> be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
> to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
> 240 n-octanes have total 1680 bonds and 1440 angles.
>
> The bond and angle interactions are both harmonic functions.  Bond
> interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
> ossilation period 16 fs.
>
> Is there something related to the integrator?
>
> Here I attached my grompp.mdp and topol.top files.
>
> ##
> grompp.mdp
> ##
>
> ; VARIOUS PREPROCESSING OPTIONS
> title= Yo
> cpp  = /usr/bin/cpp
> include  = define   =
>
> ; RUN CONTROL PARAMETERS
> integrator   = md
> tinit= 0
> dt   = 0.001
> nsteps   = 200
> init_step= 0
> comm-mode= Linear
> nstcomm  = 1
> comm-grps=
>
> ; OUTPUT CONTROL OPTIONS
> nstxout  = 5000
> nstvout  = 5000
> nstfout  = 5000
> nstcheckpoint= 1
> nstlog   = 1000
> nstenergy= 1000
> nstxtcout= 5000
> xtc-precision= 1000
> xtc-grps = energygrps   =
>
> ; NEIGHBORSEARCHING PARAMETERS
> nstlist  = 10
> ns_type  = grid
> pbc  = xyz
> rlist= 1.0025
> domain-decomposition = no
>
> ; OPTIONS FOR ELECTROSTATICS AND VDW
> coulombtype  = Cut-off
> rcoulomb-switch  = 0
> rcoulomb = 1.0025
> epsilon-r= 1
> vdw-type = Cut-off
> rvdw-switch  = 0; default rvdw
> = 1.0025; default 1 nm
> DispCorr = EnerPres
> ;table-extension  = 1.5
> fourierspacing   = 0.12
> fourier_nx   = 0
> fourier_ny   = 0
> fourier_nz   = 0
> pme_order= 4
> ewald_rtol   = 1e-05
> ewald_geometry   = 3d
> epsilon_surface  = 0
> optimize_fft = no
>
>
> ; OPTIONS FOR WEAK COUPLING ALGORITHMS
> Tcoupl   = nose-hoover
> tc-grps  = System
> tau_t= 0.02
> ref_t= 300.0
> Pcoupl   = Parrinello-Rahman
> Pcoupltype   = isotropic
> tau_p= 2.0
> compressibility  = 4.5e-5
> ref_p= 1.0
> andersen_seed= 815131
>
> ; GENERATE VELOCITIES FOR STARTUP RUN
> gen_vel   

Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread David van der Spoel

aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs can be
rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you do 
not have, you may need to run NVT for that.


Note that your integration time step is quite large, and the temperature 
coupling constant is very small.


You could try a shifted LJ + dispersion correction, it is not clear to 
me how LAMMPS treats cutoffs, couldn't find it in the manual.




Alex


Peng Yi schrieb:

On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:

I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,

Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark

Hi, Mark,

Thanks for reply!  The difference is statistically significant.  And I am
wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison please
see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  There are
bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent
results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+corr:-9238-9350
Volume(nm^3): 66.7 66.5

where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
Elj+corr is the total LJ energy including tail correction.

However, if I use integration time 2fs, Lammps results do not change
much, but Gromacs results changed a lot:

Lammps   Gromacs
Ebond(kJ/mol):2133 2700 Eangle:  
1799 1640

Etors:2552 2200
Elj+corr:-9292-9886 Volume:  
66.7 64.0


The results given by Lammps is more reasonable because the Ebond should
be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
240 n-octanes have total 1680 bonds and 1440 angles.

The bond and angle interactions are both harmonic functions.  Bond
interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
ossilation period 16 fs.

Is there something related to the integrator?

Here I attached my grompp.mdp and topol.top files.

##
grompp.mdp
##

; VARIOUS PREPROCESSING OPTIONS
title= Yo
cpp  = /usr/bin/cpp
include  = define   =

; RUN CONTROL PARAMETERS
integrator   = md
tinit= 0
dt   = 0.001
nsteps   = 200
init_step= 0
comm-mode= Linear
nstcomm  = 1
comm-grps=

; OUTPUT CONTROL OPTIONS
nstxout  = 5000
nstvout  = 5000
nstfout  = 5000
nstcheckpoint= 1
nstlog   = 1000
nstenergy= 1000
nstxtcout= 5000
xtc-precision= 1000
xtc-grps = energygrps   =

; NEIGHBORSEARCHING PARAMETERS
nstlist  = 10
ns_type  = grid
pbc  = xyz
rlist= 1.0025
domain-decomposition = no

; OPTIONS FOR ELECTROSTATICS AND VDW
coulombtype  = Cut-off
rcoulomb-switch  = 0
rcoulomb = 1.0025
epsilon-r= 1
vdw-type = Cut-off
rvdw-switch  = 0; default rvdw
= 1.0025; default 1 nm

DispCorr = EnerPres
;table-extension  = 1.5
fourierspacing   = 0.12
fourier_nx   = 0
fourier_ny   = 0
fourier_nz   = 0
pme_order= 4
ewald_rtol   = 1e-05
ewald_geometry   = 3d
epsilon_surface  = 0
optimize_fft = no


; OPTIONS FOR WEAK COUPLING ALGORITHMS
Tcoupl   = nose-hoover
tc-grps  = System
tau_t= 0.02
ref_t= 300.0
Pcoupl   = Parrinello-Rahman
Pcoupltype   = isotropic
t

Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread Ran Friedman
Hi Peng,

Note that you're not using any bond constraints in Gromacs and a
timestep of 2fs may be too long.
Also, tau_t=0.02 seems too short for me.

With 1fs timescale the agreement seem good enough, but you didn't
include estimated errors so it's hard to tell. Also, I assume you run
GMX in single and LAMMPS in double precision. Did you check for convergence?

Ran

Peng Yi wrote:
>
> On Wed, 28 Oct 2009, Mark Abraham wrote:
>
>> Peng Yi wrote:
>>>
>>> I am trying to simulate alkane melt and found out that gromacs and
>>> lammps gave different results, particularly the bonded interaction
>>> energy.
>>> I wonder if anyone has such experience.  Thanks,
>>
>> Even two installations of the same version of GROMACS can give
>> different results. The question is whether when using comparable
>> model physics you observe the same ensemble averages.
>>
>> Mark
>
> Hi, Mark,
>
> Thanks for reply!  The difference is statistically significant.  And I am
> wondering if it is caused by the integrator: Leap-frog for Gromacs and
> Velocity-verlet for Lammps.  Detail description of the comparison please
> see below:
>
> It is an NPT simulation of a melt of 240 n-octane molecules using
> united-atom model, i.e., CHx group is considered as one atom.  There are
> bond, angle, torsion and LJ interactions.  T=300K and P=1atm.
>
> Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
> nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
> for
> thermostat and barostat are 0.02ps and 2.0ps, respectively.
>
> If I use integration time 1fs, Lammps and Gromacs gave consistent
> results:
> Lammps   Gromacs
> Ebond(kJ/mol):2092 2146
> Eangle:   1757 1760
> Etors:2510 2500
> Elj+corr:-9238-9350
> Volume(nm^3): 66.7 66.5
>
> where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
> Elj+corr is the total LJ energy including tail correction.
>
> However, if I use integration time 2fs, Lammps results do not change
> much, but Gromacs results changed a lot:
>
> Lammps   Gromacs
> Ebond(kJ/mol):2133 2700 Eangle:  
> 1799 1640
> Etors:2552 2200
> Elj+corr:-9292-9886 Volume:  
> 66.7 64.0
>
> The results given by Lammps is more reasonable because the Ebond should
> be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
> to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
> 240 n-octanes have total 1680 bonds and 1440 angles.
>
> The bond and angle interactions are both harmonic functions.  Bond
> interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
> ossilation period 16 fs.
>
> Is there something related to the integrator?
>
> Here I attached my grompp.mdp and topol.top files.
>
> ##
> grompp.mdp
> ##
>
> ; VARIOUS PREPROCESSING OPTIONS
> title= Yo
> cpp  = /usr/bin/cpp
> include  = define   =
>
> ; RUN CONTROL PARAMETERS
> integrator   = md
> tinit= 0
> dt   = 0.001
> nsteps   = 200
> init_step= 0
> comm-mode= Linear
> nstcomm  = 1
> comm-grps=
>
> ; OUTPUT CONTROL OPTIONS
> nstxout  = 5000
> nstvout  = 5000
> nstfout  = 5000
> nstcheckpoint= 1
> nstlog   = 1000
> nstenergy= 1000
> nstxtcout= 5000
> xtc-precision= 1000
> xtc-grps = energygrps   =
>
> ; NEIGHBORSEARCHING PARAMETERS
> nstlist  = 10
> ns_type  = grid
> pbc  = xyz
> rlist= 1.0025
> domain-decomposition = no
>
> ; OPTIONS FOR ELECTROSTATICS AND VDW
> coulombtype  = Cut-off
> rcoulomb-switch  = 0
> rcoulomb = 1.0025
> epsilon-r= 1
> vdw-type = Cut-off
> rvdw-switch  = 0; default rvdw
> = 1.0025; default 1 nm
> DispCorr = EnerPres
> ;table-extension  = 1.5
> fourierspacing   = 0.12
> fourier_nx   = 0
> fourier_ny   = 0
> fourier_nz   = 0
> pme_order= 4
> ewald_rtol   = 1e-05
> ewald_geometry   = 3d
> epsilon_surface  = 0
> optimize_fft = no
>
>
> ; OPTIONS FOR WEAK COUPLING ALGORITHMS
> Tcoupl   = nose-hoover
> tc-grps  = System
> tau_t= 0.02
> ref_t= 300.0
> Pcoupl   = Parrinello-Rahman
> Pcoupltype   = isotropic
> tau

Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread Peng Yi


Hi, Ran,

I didn't use bond restraints.  I checked that the bond length had a
Gaussian-like distributes, and the length range looked normal.

I estimated the fastest timescale in the system, which is the bond
ossilation period, around 16fs.  Would that require as integration
timestep much smaller than 1fs?

With the parameters I have, could you recommend a set of tau_t and tau_p?
I did mention the fluctuation, 100 kJ/mol for energy and 1 nm^3 for
volume.  And I ran GMX in single.  Not sure about Lammps, should be
double.  All measured physical quantities converged well.  Would you
expect differece if I compile GMX in double?  Would that be much slower?
-Peng

On Thu, 29 Oct 2009, Ran Friedman wrote:


Hi Peng,

Note that you're not using any bond constraints in Gromacs and a
timestep of 2fs may be too long.
Also, tau_t=0.02 seems too short for me.

With 1fs timescale the agreement seem good enough, but you didn't
include estimated errors so it's hard to tell. Also, I assume you run
GMX in single and LAMMPS in double precision. Did you check for convergence?

Ran

Peng Yi wrote:


On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:


I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,


Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark


Hi, Mark,

Thanks for reply!  The difference is statistically significant.  And I am
wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison please
see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  There are
bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent
results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+corr:-9238-9350
Volume(nm^3): 66.7 66.5

where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
Elj+corr is the total LJ energy including tail correction.

However, if I use integration time 2fs, Lammps results do not change
much, but Gromacs results changed a lot:

Lammps   Gromacs
Ebond(kJ/mol):2133 2700 Eangle:
1799 1640
Etors:2552 2200
Elj+corr:-9292-9886 Volume:
66.7 64.0

The results given by Lammps is more reasonable because the Ebond should
be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
240 n-octanes have total 1680 bonds and 1440 angles.

The bond and angle interactions are both harmonic functions.  Bond
interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
ossilation period 16 fs.

Is there something related to the integrator?

Here I attached my grompp.mdp and topol.top files.

##
grompp.mdp
##

; VARIOUS PREPROCESSING OPTIONS
title= Yo
cpp  = /usr/bin/cpp
include  = define   =

; RUN CONTROL PARAMETERS
integrator   = md
tinit= 0
dt   = 0.001
nsteps   = 200
init_step= 0
comm-mode= Linear
nstcomm  = 1
comm-grps=

; OUTPUT CONTROL OPTIONS
nstxout  = 5000
nstvout  = 5000
nstfout  = 5000
nstcheckpoint= 1
nstlog   = 1000
nstenergy= 1000
nstxtcout= 5000
xtc-precision= 1000
xtc-grps = energygrps   =

; NEIGHBORSEARCHING PARAMETERS
nstlist  = 10
ns_type  = grid
pbc  = xyz
rlist= 1.0025
domain-decomposition = no

; OPTIONS FOR ELECTROSTATICS AND VDW
coulombtype  = Cut-off
rcoulomb-switch  = 0
rcoulomb = 1.0025
epsilon-r= 1
vdw-type = Cut-off
rvdw-switch  = 0; default rvdw
= 1.0025; default 1 nm
DispCorr = EnerPres
;table-extension  = 1.5
fourierspacing   = 0.12
fourier_nx   = 0
fourier_ny   = 0
fourier_nz   = 0
pme

Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread Peng Yi




On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs can be
rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you do not 
have, you may need to run NVT for that.


Note that your integration time step is quite large, and the temperature 
coupling constant is very small.


You could try a shifted LJ + dispersion correction, it is not clear to me how 
LAMMPS treats cutoffs, couldn't find it in the manual.




I tried NPT simulation with only LJ interaction.  With tau_t=0.2ps for 
both Gromacs and Lammps, all the other parameters same as I mentioned

before.  At integration time 2fs, Gromacs and Lammps produce same results
and agree with published data. (N=4000, T=2.0, P=0.5 -> rou=0.3).  That's
why I focused my attention to the difference on bonded energy.

I think I used a simple cutoff plus tail correction for energy and
pressure for both Gromacs and Lammps.



Alex


Peng Yi schrieb:

On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:

I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,

Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark

Hi, Mark,

Thanks for reply!  The difference is statistically significant.  And I am
wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison please
see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  There are
bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent
results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+corr:-9238-9350
Volume(nm^3): 66.7 66.5

where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
Elj+corr is the total LJ energy including tail correction.

However, if I use integration time 2fs, Lammps results do not change
much, but Gromacs results changed a lot:

Lammps   Gromacs
Ebond(kJ/mol):2133 2700 Eangle:  1799 
1640

Etors:2552 2200
Elj+corr:-9292-9886 Volume:  66.7 
64.0


The results given by Lammps is more reasonable because the Ebond should
be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
240 n-octanes have total 1680 bonds and 1440 angles.

The bond and angle interactions are both harmonic functions.  Bond
interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
ossilation period 16 fs.

Is there something related to the integrator?

Here I attached my grompp.mdp and topol.top files.

##
grompp.mdp
##

; VARIOUS PREPROCESSING OPTIONS
title= Yo
cpp  = /usr/bin/cpp
include  = define   =

; RUN CONTROL PARAMETERS
integrator   = md
tinit= 0
dt   = 0.001
nsteps   = 200
init_step= 0
comm-mode= Linear
nstcomm  = 1
comm-grps=

; OUTPUT CONTROL OPTIONS
nstxout  = 5000
nstvout  = 5000
nstfout  = 5000
nstcheckpoint= 1
nstlog   = 1000
nstenergy= 1000
nstxtcout= 5000
xtc-precision= 1000
xtc-grps = energygrps   =

; NEIGHBORSEARCHING PARAMETERS
nstlist  = 10
ns_type  = grid
pbc  = xyz
rlist= 1.0025
domain-decomposition = no

; OPTIONS FOR ELECTROSTATICS AND VDW
coulombtype  = Cut-off
rcoulomb-switch  = 0
rcoulomb = 1.0025
epsilon-r= 1
vdw-type = Cut-off
rvdw-switch  = 0; default rvdw= 
1.0025; default 1 nm

DispCorr = EnerPres
;table-extension  = 1.5
fourierspacing   = 0.12
fourier_

Re: [gmx-users] Is anyone also using lammps?s

2009-10-29 Thread Ran Friedman
Hi Peng,

The time scale should be much shorter than the fastest vibration. A rule
of thumb from the reference below is a factor of ten, but it would
depend on the precision. Running with double precision is shorter but I
didn't make benchmarks (perhaps other users have).

Appropriate values of tau_t and tau_p have also been discussed in this
list (search for references by Berk). I tend to use something like
tau_t=0.2 and tau_p=1.0.

I advise you to follow David's suggestion as well and run with NVT.

Ran.

Reference:
@book{Becker2001,
   Author = {Becker, O. M. and MacKerell, A. D. Jr. and Roux, B.  and
Watanabe, M.},
   Title = {Computational biochemistry and biophysics},
   Publisher = {Dekker, M.},
   Address = {New York},
 Year = {2001} }

Peng Yi wrote:
>
> Hi, Ran,
>
> I didn't use bond restraints.  I checked that the bond length had a
> Gaussian-like distributes, and the length range looked normal.
>
> I estimated the fastest timescale in the system, which is the bond
> ossilation period, around 16fs.  Would that require as integration
> timestep much smaller than 1fs?
>
> With the parameters I have, could you recommend a set of tau_t and tau_p?
> I did mention the fluctuation, 100 kJ/mol for energy and 1 nm^3 for
> volume.  And I ran GMX in single.  Not sure about Lammps, should be
> double.  All measured physical quantities converged well.  Would you
> expect differece if I compile GMX in double?  Would that be much slower?
> -Peng
>
> On Thu, 29 Oct 2009, Ran Friedman wrote:
>
>> Hi Peng,
>>
>> Note that you're not using any bond constraints in Gromacs and a
>> timestep of 2fs may be too long.
>> Also, tau_t=0.02 seems too short for me.
>>
>> With 1fs timescale the agreement seem good enough, but you didn't
>> include estimated errors so it's hard to tell. Also, I assume you run
>> GMX in single and LAMMPS in double precision. Did you check for
>> convergence?
>>
>> Ran
>>
>> Peng Yi wrote:
>>>
>>> On Wed, 28 Oct 2009, Mark Abraham wrote:
>>>
 Peng Yi wrote:
>
> I am trying to simulate alkane melt and found out that gromacs and
> lammps gave different results, particularly the bonded interaction
> energy.
> I wonder if anyone has such experience.  Thanks,

 Even two installations of the same version of GROMACS can give
 different results. The question is whether when using comparable
 model physics you observe the same ensemble averages.

 Mark
>>>
>>> Hi, Mark,
>>>
>>> Thanks for reply!  The difference is statistically significant.  And
>>> I am
>>> wondering if it is caused by the integrator: Leap-frog for Gromacs and
>>> Velocity-verlet for Lammps.  Detail description of the comparison
>>> please
>>> see below:
>>>
>>> It is an NPT simulation of a melt of 240 n-octane molecules using
>>> united-atom model, i.e., CHx group is considered as one atom.  There
>>> are
>>> bond, angle, torsion and LJ interactions.  T=300K and P=1atm.
>>>
>>> Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
>>> nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
>>> for
>>> thermostat and barostat are 0.02ps and 2.0ps, respectively.
>>>
>>> If I use integration time 1fs, Lammps and Gromacs gave consistent
>>> results:
>>> Lammps   Gromacs
>>> Ebond(kJ/mol):2092 2146
>>> Eangle:   1757 1760
>>> Etors:2510 2500
>>> Elj+corr:-9238-9350
>>> Volume(nm^3): 66.7 66.5
>>>
>>> where energy fluctuation is 100 kJ/mol and volume fluctuation is 1
>>> nm^3,
>>> Elj+corr is the total LJ energy including tail correction.
>>>
>>> However, if I use integration time 2fs, Lammps results do not change
>>> much, but Gromacs results changed a lot:
>>>
>>> Lammps   Gromacs
>>> Ebond(kJ/mol):2133 2700 Eangle:
>>> 1799 1640
>>> Etors:2552 2200
>>> Elj+corr:-9292-9886 Volume:
>>> 66.7 64.0
>>>
>>> The results given by Lammps is more reasonable because the Ebond should
>>> be equal to the total # of bonds times 1/2k_BT and Eangle should be
>>> equal
>>> to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25
>>> kJ/mol.
>>> 240 n-octanes have total 1680 bonds and 1440 angles.
>>>
>>> The bond and angle interactions are both harmonic functions.  Bond
>>> interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
>>> ossilation period 16 fs.
>>>
>>> Is there something related to the integrator?
>>>
>>> Here I attached my grompp.mdp and topol.top files.
>>>
>>> ##
>>> grompp.mdp
>>> ##
>>>
>>> ; VARIOUS PREPROCESSING OPTIONS
>>> title= Yo
>>> cpp  = /usr/bin/cpp
>>> include  = define   =
>>>
>>> ; RUN CONTROL PARAMETERS
>>> integrator   = md
>>> tinit= 0
>>> dt 

Re: [gmx-users] Is anyone also using lammps?s

2009-11-01 Thread Peng Yi


Thank for your reply!  I have done some NVT runs per your suggestion, and
the results are similar to NPT runs, i.e., Gromacs results is more
affected by changing integration timestep than Lammps.  Details below:

A melt of 240 octane chains by united-atom model. T=300K, V=55.46 nm^3.
Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.
All other parameters in .top and .mdp files are the same as previously
attached..

If I use integration time 1fs, Lammps and Gromacs produce
consistent results:
 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21332160 100
Eangle:   17571780  80
Etors:25312510  80
Elj+corr:   -10711  -10767  90
P(atm):   35003250 500

if I use integration time 2fs, Lammps results remain unchanged, but
Gromacs results change significantly, particularly bonded energy:

 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21752710 100
Eangle:   17991640  70
Etors:25732230  80
Elj+corr:   -10711  -11007 100 
P(atm):   32002730 700


Would that be a result of using different integrator between Lammps and 
Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses Leap-frog.

Thanks,
-Peng

On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs can be
rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you do not 
have, you may need to run NVT for that.


Note that your integration time step is quite large, and the temperature 
coupling constant is very small.


You could try a shifted LJ + dispersion correction, it is not clear to me how 
LAMMPS treats cutoffs, couldn't find it in the manual.




Alex


Peng Yi schrieb:

On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:

I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,

Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark

Hi, Mark,

Thanks for reply!  The difference is statistically significant.  And I am
wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison please
see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  There are
bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent
results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+corr:-9238-9350
Volume(nm^3): 66.7 66.5

where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 nm^3,
Elj+corr is the total LJ energy including tail correction.

However, if I use integration time 2fs, Lammps results do not change
much, but Gromacs results changed a lot:

Lammps   Gromacs
Ebond(kJ/mol):2133 2700 Eangle:  1799 
1640

Etors:2552 2200
Elj+corr:-9292-9886 Volume:  66.7 
64.0


The results given by Lammps is more reasonable because the Ebond should
be equal to the total # of bonds times 1/2k_BT and Eangle should be equal
to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 kJ/mol.
240 n-octanes have total 1680 bonds and 1440 angles.

The bond and angle interactions are both harmonic functions.  Bond
interaction constant kl=292880 kJ/mol/nm^2, corresponding to a bond
ossilation period 16 fs.

Is there something related to the integrator?

Here I attached my grompp.mdp and topol.top files.

##
grompp.mdp
##

; VARIOUS PREPROCESSING OPTIONS
title= Yo
cpp  = /usr/bin/cpp
include  = define   =

; RUN CONTROL PARAMETERS
integrator   = md
tinit= 0
dt   = 0.001
nsteps   = 200
init_step= 0
comm-mo

Re: [gmx-users] Is anyone also using lammps?s

2009-11-01 Thread David van der Spoel

Peng Yi wrote:


Thank for your reply!  I have done some NVT runs per your suggestion, and
the results are similar to NPT runs, i.e., Gromacs results is more
affected by changing integration timestep than Lammps.  Details below:

A melt of 240 octane chains by united-atom model. T=300K, V=55.46 nm^3.
Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.


As some people have note the tau_t is short for Nose Hoover. Are you 
sure this means the same in Lammps and in Gromacs? For one thing, there 
is no tau_t in the NH algorithm as far as I know, and Gromacs converts 
it to an appropriate weight or whatever that is called. What does Lammps 
do with this tau_t. To be a the safe side you could run both with 
Berendsen as well. Is the std err identical in both packages? And in the 
2 fs run, are both simulation equlibrated with this time step as well?



All other parameters in .top and .mdp files are the same as previously
attached..

If I use integration time 1fs, Lammps and Gromacs produce
consistent results:
 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21332160 100
Eangle:   17571780  80
Etors:25312510  80
Elj+corr:   -10711  -10767  90
P(atm):   35003250 500

if I use integration time 2fs, Lammps results remain unchanged, but
Gromacs results change significantly, particularly bonded energy:

 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21752710 100
Eangle:   17991640  70
Etors:25732230  80
Elj+corr:   -10711  -11007 100 
P(atm):   32002730 700


Would that be a result of using different integrator between Lammps 
and Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses Leap-frog.

Thanks,
-Peng

On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs 
can be

rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you 
do not have, you may need to run NVT for that.


Note that your integration time step is quite large, and the 
temperature coupling constant is very small.


You could try a shifted LJ + dispersion correction, it is not clear 
to me how LAMMPS treats cutoffs, couldn't find it in the manual.




Alex


Peng Yi schrieb:

On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:

I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,

Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark

Hi, Mark,

Thanks for reply!  The difference is statistically significant.  
And I am

wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison 
please

see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  
There are

bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent
results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+corr:-9238-9350
Volume(nm^3): 66.7 66.5

where energy fluctuation is 100 kJ/mol and volume fluctuation is 1 
nm^3,

Elj+corr is the total LJ energy including tail correction.

However, if I use integration time 2fs, Lammps results do not change
much, but Gromacs results changed a lot:

Lammps   Gromacs
Ebond(kJ/mol):2133 2700 Eangle:  
1799 1640

Etors:2552 2200
Elj+corr:-9292-9886 Volume:  
66.7 64.0


The results given by Lammps is more reasonable because the Ebond 
should
be equal to the total # of bonds times 1/2k_BT and Eangle should be 
equal
to the total # of angles times 1/2k_BT.  At T=300K, 1/2k_BT=1.25 
kJ/mol.

240 n-octanes have total 1680 bonds and 1440 angles.

The bond and angle interactions are both harmonic functions.  Bond
interaction constant kl=292880 kJ/mo

Re: [gmx-users] Is anyone also using lammps?s

2009-11-03 Thread Peng Yi


Hi, David,

I used Berendsen thermostat and a bigger tau_t=1ps to redo the simulations.
The general conclusion is the same.  The std err is the same in both
packages.  And during each simulation, the integration time does not
change.  Details below:

Integration time step 1fs:
Lammps  GromacsStd.Err. (for both)
Ebond(kJ/mol):   21122170   100
Eangle:  17991770   100
Etors:   25522490   100
Elj+corr:  -10711  -10777   100
P(atm):  32503216   500

Integration time step 2fs:
Lammps  GromacsStd.Err. (for both)
Ebond:   21542654   120
Eangle:  18401645   120
Etors:   25732236   120
Elj+corr:  -10711  -11019   100 
P(atm):  32502590   600


-Peng

On Sun, 1 Nov 2009, David van der Spoel wrote:


Peng Yi wrote:


Thank for your reply!  I have done some NVT runs per your suggestion, and
the results are similar to NPT runs, i.e., Gromacs results is more
affected by changing integration timestep than Lammps.  Details below:

A melt of 240 octane chains by united-atom model. T=300K, V=55.46 nm^3.
Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.


As some people have note the tau_t is short for Nose Hoover. Are you sure 
this means the same in Lammps and in Gromacs? For one thing, there is no 
tau_t in the NH algorithm as far as I know, and Gromacs converts it to an 
appropriate weight or whatever that is called. What does Lammps do with this 
tau_t. To be a the safe side you could run both with Berendsen as well. Is 
the std err identical in both packages? And in the 2 fs run, are both 
simulation equlibrated with this time step as well?



All other parameters in .top and .mdp files are the same as previously
attached..

If I use integration time 1fs, Lammps and Gromacs produce
consistent results:
 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21332160 100
Eangle:   17571780  80
Etors:25312510  80
Elj+corr:   -10711  -10767  90
P(atm):   35003250 500

if I use integration time 2fs, Lammps results remain unchanged, but
Gromacs results change significantly, particularly bonded energy:

 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21752710 100
Eangle:   17991640  70
Etors:25732230  80
Elj+corr:   -10711  -11007 100 P(atm): 
32002730 700


Would that be a result of using different integrator between Lammps and 
Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses Leap-frog.

Thanks,
-Peng

On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs can be
rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you do 
not have, you may need to run NVT for that.


Note that your integration time step is quite large, and the temperature 
coupling constant is very small.


You could try a shifted LJ + dispersion correction, it is not clear to me 
how LAMMPS treats cutoffs, couldn't find it in the manual.




Alex


Peng Yi schrieb:

On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:

I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,

Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark

Hi, Mark,

Thanks for reply!  The difference is statistically significant.  And I 
am

wondering if it is caused by the integrator: Leap-frog for Gromacs and
Velocity-verlet for Lammps.  Detail description of the comparison please
see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  There are
bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time constants
for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Gromacs gave consistent
results:
Lammps   Gromacs
Ebond(kJ/mol):2092 2146
Eangle:   1757 1760
Etors:2510 2500
Elj+cor

Re: [gmx-users] Is anyone also using lammps?s

2009-11-03 Thread David van der Spoel

Peng Yi wrote:


Hi, David,

I used Berendsen thermostat and a bigger tau_t=1ps to redo the simulations.
The general conclusion is the same.  The std err is the same in both
packages.  And during each simulation, the integration time does not
change.  Details below:

Integration time step 1fs:
Lammps  GromacsStd.Err. (for both)
Ebond(kJ/mol):   21122170   100
Eangle:  17991770   100
Etors:   25522490   100
Elj+corr:  -10711  -10777   100
P(atm):  32503216   500

Integration time step 2fs:
Lammps  GromacsStd.Err. (for both)
Ebond:   21542654   120
Eangle:  18401645   120
Etors:   25732236   120
Elj+corr:  -10711  -11019   100 P(atm):  
32502590   600




How about the temperature in both systems? Was Lammps also run with 
Berendsen? It could also still be a topology error. Maybe you can turn 
off the torsion potential to test this.



-Peng

On Sun, 1 Nov 2009, David van der Spoel wrote:


Peng Yi wrote:


Thank for your reply!  I have done some NVT runs per your suggestion, 
and

the results are similar to NPT runs, i.e., Gromacs results is more
affected by changing integration timestep than Lammps.  Details below:

A melt of 240 octane chains by united-atom model. T=300K, V=55.46 nm^3.
Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.


As some people have note the tau_t is short for Nose Hoover. Are you 
sure this means the same in Lammps and in Gromacs? For one thing, 
there is no tau_t in the NH algorithm as far as I know, and Gromacs 
converts it to an appropriate weight or whatever that is called. What 
does Lammps do with this tau_t. To be a the safe side you could run 
both with Berendsen as well. Is the std err identical in both 
packages? And in the 2 fs run, are both simulation equlibrated with 
this time step as well?



All other parameters in .top and .mdp files are the same as previously
attached..

If I use integration time 1fs, Lammps and Gromacs produce
consistent results:
 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21332160 100
Eangle:   17571780  80
Etors:25312510  80
Elj+corr:   -10711  -10767  90
P(atm):   35003250 500

if I use integration time 2fs, Lammps results remain unchanged, but
Gromacs results change significantly, particularly bonded energy:

 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21752710 100
Eangle:   17991640  70
Etors:25732230  80
Elj+corr:   -10711  -11007 100 P(atm): 
32002730 700


Would that be a result of using different integrator between Lammps 
and Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses Leap-frog.

Thanks,
-Peng

On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs 
can be

rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you 
do not have, you may need to run NVT for that.


Note that your integration time step is quite large, and the 
temperature coupling constant is very small.


You could try a shifted LJ + dispersion correction, it is not clear 
to me how LAMMPS treats cutoffs, couldn't find it in the manual.




Alex


Peng Yi schrieb:

On Wed, 28 Oct 2009, Mark Abraham wrote:


Peng Yi wrote:

I am trying to simulate alkane melt and found out that gromacs and
lammps gave different results, particularly the bonded interaction
energy.
I wonder if anyone has such experience.  Thanks,

Even two installations of the same version of GROMACS can give
different results. The question is whether when using comparable
model physics you observe the same ensemble averages.

Mark

Hi, Mark,

Thanks for reply!  The difference is statistically significant.  
And I am
wondering if it is caused by the integrator: Leap-frog for Gromacs 
and
Velocity-verlet for Lammps.  Detail description of the comparison 
please

see below:

It is an NPT simulation of a melt of 240 n-octane molecules using
united-atom model, i.e., CHx group is considered as one atom.  
There are

bond, angle, torsion and LJ interactions.  T=300K and P=1atm.

Lammps uses nose-hoover thermostat and barostat, and Gromacs uses
nose-hoover thermostat and Parranello-Rahman barostat.  Time 
constants

for
thermostat and barostat are 0.02ps and 2.0ps, respectively.

If I use integration time 1fs, Lammps and Grom

Re: [gmx-users] Is anyone also using lammps?s

2009-11-03 Thread Peng Yi


I turned off the torsion interaction.  The difference between Lammps and 
Gromacs at integration time step 2fs was reduced.  Details below:


A melt of 240 n-octane (united-atom model), NVT, T=300K, V=55.46nm^3.
Both Lammps and Gromacs use berendsen thermostat with tau_t=1ps.

Integration time step 1fs:
 LammpsGromacsStd. Err. (for both)
Ebond(kJ/mol):2092  2109   100
Eangle:   1778  175480
Elj+corr:   -10501-10553   100
T(K):  300   299 5
P(atm):   3188  3016   700

integration time step 2fs:
 LammpsGromacsStd.Err. (for both)
Ebond:2133  2232   100
Eangle:   1803  173780
Elj+corr:   -10501-10623   100
T: 300   298 5
P:3133  2955   600

Lammps results remain almost unchanged when dt increases from 1fs to 2fs,
and Ebond : Eangle = 7 : 6, which is the ratio between # of bonds and
# of angles.

Gromacs results change more significantly when dt goes from 1fs to 2fs.
and the trend of Ebond and Eangle are opposite.  It is more significant
when torsion interaction is present.


On Tue, 3 Nov 2009, David van der Spoel wrote:


Peng Yi wrote:


Hi, David,

I used Berendsen thermostat and a bigger tau_t=1ps to redo the simulations.
The general conclusion is the same.  The std err is the same in both
packages.  And during each simulation, the integration time does not
change.  Details below:

Integration time step 1fs:
Lammps  GromacsStd.Err. (for both)
Ebond(kJ/mol):   21122170   100
Eangle:  17991770   100
Etors:   25522490   100
Elj+corr:  -10711  -10777   100
P(atm):  32503216   500

Integration time step 2fs:
Lammps  GromacsStd.Err. (for both)
Ebond:   21542654   120
Eangle:  18401645   120
Etors:   25732236   120
Elj+corr:  -10711  -11019   100 P(atm):  3250 
2590   600




How about the temperature in both systems? Was Lammps also run with 
Berendsen? It could also still be a topology error. Maybe you can turn off 
the torsion potential to test this.



-Peng

On Sun, 1 Nov 2009, David van der Spoel wrote:


Peng Yi wrote:


Thank for your reply!  I have done some NVT runs per your suggestion, and
the results are similar to NPT runs, i.e., Gromacs results is more
affected by changing integration timestep than Lammps.  Details below:

A melt of 240 octane chains by united-atom model. T=300K, V=55.46 nm^3.
Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.


As some people have note the tau_t is short for Nose Hoover. Are you sure 
this means the same in Lammps and in Gromacs? For one thing, there is no 
tau_t in the NH algorithm as far as I know, and Gromacs converts it to an 
appropriate weight or whatever that is called. What does Lammps do with 
this tau_t. To be a the safe side you could run both with Berendsen as 
well. Is the std err identical in both packages? And in the 2 fs run, are 
both simulation equlibrated with this time step as well?



All other parameters in .top and .mdp files are the same as previously
attached..

If I use integration time 1fs, Lammps and Gromacs produce
consistent results:
 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21332160 100
Eangle:   17571780  80
Etors:25312510  80
Elj+corr:   -10711  -10767  90
P(atm):   35003250 500

if I use integration time 2fs, Lammps results remain unchanged, but
Gromacs results change significantly, particularly bonded energy:

 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21752710 100
Eangle:   17991640  70
Etors:25732230  80
Elj+corr:   -10711  -11007 100 P(atm): 3200 
2730 700


Would that be a result of using different integrator between Lammps and 
Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses Leap-frog.

Thanks,
-Peng

On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in gromacs can 
be

rather bad for single precision.


Please refer to the gromacs 4.0 paper for a discussion of the drift.
If you want to compare energies you need the same density, which you do 
not have, you may need to run NVT for that.


Note that your integration time step is quite large, and the temperature 
coupling constant is very small.


You could try a shifted LJ + dispersion correction,

Re: [gmx-users] Is anyone also using lammps?s

2009-11-04 Thread Ran Friedman
Dear Peng,

Did you also try to run GMX in double precision at some point?

Ran

Peng Yi wrote:
>
> I turned off the torsion interaction.  The difference between Lammps
> and Gromacs at integration time step 2fs was reduced.  Details below:
>
> A melt of 240 n-octane (united-atom model), NVT, T=300K, V=55.46nm^3.
> Both Lammps and Gromacs use berendsen thermostat with tau_t=1ps.
>
> Integration time step 1fs:
>  LammpsGromacsStd. Err. (for both)
> Ebond(kJ/mol):2092  2109   100
> Eangle:   1778  175480
> Elj+corr:   -10501-10553   100
> T(K):  300   299 5
> P(atm):   3188  3016   700
>
> integration time step 2fs:
>  LammpsGromacsStd.Err. (for both)
> Ebond:2133  2232   100
> Eangle:   1803  173780
> Elj+corr:   -10501-10623   100
> T: 300   298 5
> P:3133  2955   600
>
> Lammps results remain almost unchanged when dt increases from 1fs to 2fs,
> and Ebond : Eangle = 7 : 6, which is the ratio between # of bonds and
> # of angles.
>
> Gromacs results change more significantly when dt goes from 1fs to 2fs.
> and the trend of Ebond and Eangle are opposite.  It is more significant
> when torsion interaction is present.
>
>
> On Tue, 3 Nov 2009, David van der Spoel wrote:
>
>> Peng Yi wrote:
>>>
>>> Hi, David,
>>>
>>> I used Berendsen thermostat and a bigger tau_t=1ps to redo the
>>> simulations.
>>> The general conclusion is the same.  The std err is the same in both
>>> packages.  And during each simulation, the integration time does not
>>> change.  Details below:
>>>
>>> Integration time step 1fs:
>>> Lammps  GromacsStd.Err. (for both)
>>> Ebond(kJ/mol):   21122170   100
>>> Eangle:  17991770   100
>>> Etors:   25522490   100
>>> Elj+corr:  -10711  -10777   100
>>> P(atm):  32503216   500
>>>
>>> Integration time step 2fs:
>>> Lammps  GromacsStd.Err. (for both)
>>> Ebond:   21542654   120
>>> Eangle:  18401645   120
>>> Etors:   25732236   120
>>> Elj+corr:  -10711  -11019   100 P(atm): 
>>> 3250 2590   600
>>>
>>
>> How about the temperature in both systems? Was Lammps also run with
>> Berendsen? It could also still be a topology error. Maybe you can
>> turn off the torsion potential to test this.
>>
>>> -Peng
>>>
>>> On Sun, 1 Nov 2009, David van der Spoel wrote:
>>>
 Peng Yi wrote:
>
> Thank for your reply!  I have done some NVT runs per your
> suggestion, and
> the results are similar to NPT runs, i.e., Gromacs results is more
> affected by changing integration timestep than Lammps.  Details
> below:
>
> A melt of 240 octane chains by united-atom model. T=300K, V=55.46
> nm^3.
> Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.

 As some people have note the tau_t is short for Nose Hoover. Are
 you sure this means the same in Lammps and in Gromacs? For one
 thing, there is no tau_t in the NH algorithm as far as I know, and
 Gromacs converts it to an appropriate weight or whatever that is
 called. What does Lammps do with this tau_t. To be a the safe side
 you could run both with Berendsen as well. Is the std err identical
 in both packages? And in the 2 fs run, are both simulation
 equlibrated with this time step as well?

> All other parameters in .top and .mdp files are the same as
> previously
> attached..
>
> If I use integration time 1fs, Lammps and Gromacs produce
> consistent results:
>  Lammps  Gromacs  std. err.
> Ebond(kJ/mol):21332160 100
> Eangle:   17571780  80
> Etors:25312510  80
> Elj+corr:   -10711  -10767  90
> P(atm):   35003250 500
>
> if I use integration time 2fs, Lammps results remain unchanged, but
> Gromacs results change significantly, particularly bonded energy:
>
>  Lammps  Gromacs  std. err.
> Ebond(kJ/mol):21752710 100
> Eangle:   17991640  70
> Etors:25732230  80
> Elj+corr:   -10711  -11007 100 P(atm):
> 3200 2730 700
>
> Would that be a result of using different integrator between
> Lammps and Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses
> Leap-frog.
> Thanks,
> -Peng
>
> On Thu, 29 Oct 2009, David van der Spoel 

Re: [gmx-users] Is anyone also using lammps?s

2009-11-04 Thread Peng Yi


Hi, Ran,

No, I haven't.  I still have to find out how to install in double
precision.  Would double precision be slower than single?  If so,
how much?  Or just double the memory used?  Thanks,

-Peng


On Wed, 4 Nov 2009, Ran Friedman wrote:


Dear Peng,

Did you also try to run GMX in double precision at some point?

Ran

Peng Yi wrote:


I turned off the torsion interaction.  The difference between Lammps
and Gromacs at integration time step 2fs was reduced.  Details below:

A melt of 240 n-octane (united-atom model), NVT, T=300K, V=55.46nm^3.
Both Lammps and Gromacs use berendsen thermostat with tau_t=1ps.

Integration time step 1fs:
 LammpsGromacsStd. Err. (for both)
Ebond(kJ/mol):2092  2109   100
Eangle:   1778  175480
Elj+corr:   -10501-10553   100
T(K):  300   299 5
P(atm):   3188  3016   700

integration time step 2fs:
 LammpsGromacsStd.Err. (for both)
Ebond:2133  2232   100
Eangle:   1803  173780
Elj+corr:   -10501-10623   100
T: 300   298 5
P:3133  2955   600

Lammps results remain almost unchanged when dt increases from 1fs to 2fs,
and Ebond : Eangle = 7 : 6, which is the ratio between # of bonds and
# of angles.

Gromacs results change more significantly when dt goes from 1fs to 2fs.
and the trend of Ebond and Eangle are opposite.  It is more significant
when torsion interaction is present.


On Tue, 3 Nov 2009, David van der Spoel wrote:


Peng Yi wrote:


Hi, David,

I used Berendsen thermostat and a bigger tau_t=1ps to redo the
simulations.
The general conclusion is the same.  The std err is the same in both
packages.  And during each simulation, the integration time does not
change.  Details below:

Integration time step 1fs:
Lammps  GromacsStd.Err. (for both)
Ebond(kJ/mol):   21122170   100
Eangle:  17991770   100
Etors:   25522490   100
Elj+corr:  -10711  -10777   100
P(atm):  32503216   500

Integration time step 2fs:
Lammps  GromacsStd.Err. (for both)
Ebond:   21542654   120
Eangle:  18401645   120
Etors:   25732236   120
Elj+corr:  -10711  -11019   100 P(atm):
3250 2590   600



How about the temperature in both systems? Was Lammps also run with
Berendsen? It could also still be a topology error. Maybe you can
turn off the torsion potential to test this.


-Peng

On Sun, 1 Nov 2009, David van der Spoel wrote:


Peng Yi wrote:


Thank for your reply!  I have done some NVT runs per your
suggestion, and
the results are similar to NPT runs, i.e., Gromacs results is more
affected by changing integration timestep than Lammps.  Details
below:

A melt of 240 octane chains by united-atom model. T=300K, V=55.46
nm^3.
Both Gromacs and Lammps use Nose-Hoover thermostat with tau_t=0.2 ps.


As some people have note the tau_t is short for Nose Hoover. Are
you sure this means the same in Lammps and in Gromacs? For one
thing, there is no tau_t in the NH algorithm as far as I know, and
Gromacs converts it to an appropriate weight or whatever that is
called. What does Lammps do with this tau_t. To be a the safe side
you could run both with Berendsen as well. Is the std err identical
in both packages? And in the 2 fs run, are both simulation
equlibrated with this time step as well?


All other parameters in .top and .mdp files are the same as
previously
attached..

If I use integration time 1fs, Lammps and Gromacs produce
consistent results:
 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21332160 100
Eangle:   17571780  80
Etors:25312510  80
Elj+corr:   -10711  -10767  90
P(atm):   35003250 500

if I use integration time 2fs, Lammps results remain unchanged, but
Gromacs results change significantly, particularly bonded energy:

 Lammps  Gromacs  std. err.
Ebond(kJ/mol):21752710 100
Eangle:   17991640  70
Etors:25732230  80
Elj+corr:   -10711  -11007 100 P(atm):
3200 2730 700

Would that be a result of using different integrator between
Lammps and Gromacs?  Lammps uses Velocity-Verlet, and Gromacs uses
Leap-frog.
Thanks,
-Peng

On Thu, 29 Oct 2009, David van der Spoel wrote:


aherz wrote:

Hey,

are you running single or double precision gromacs?
Afaik, depending on the circumstances the energy drift in
gromacs can be
rather bad for single precision.


Please refer to the gromacs 4.0 paper 

Re: [gmx-users] Is anyone also using lammps?s

2009-11-04 Thread Ran Friedman
Hi Peng,

It would be slower - never made any benchmarks on how much slower. But
you don't run a very long simulation, do you?
Installing it isn't a problem. If you use PME you also need fftw in
double precision though.

Ran

Peng Yi wrote:
>
> Hi, Ran,
>
> No, I haven't.  I still have to find out how to install in double
> precision.  Would double precision be slower than single?  If so,
> how much?  Or just double the memory used?  Thanks,
>
> -Peng

___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php


Re: [gmx-users] Is anyone also using lammps?s

2009-11-04 Thread Peng Yi


Hi, Ran,

Do you mean the simulation I used here for testing purpose?  That is
not long, maybe a few hours.  But my research will require simulations
for days.  If I will not use PME, do I still need a fftw?

-Peng

On Wed, 4 Nov 2009, Ran Friedman wrote:


Hi Peng,

It would be slower - never made any benchmarks on how much slower. But
you don't run a very long simulation, do you?
Installing it isn't a problem. If you use PME you also need fftw in
double precision though.

Ran

Peng Yi wrote:


Hi, Ran,

No, I haven't.  I still have to find out how to install in double
precision.  Would double precision be slower than single?  If so,
how much?  Or just double the memory used?  Thanks,

-Peng


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php


Re: [gmx-users] Is anyone also using lammps?s

2009-11-04 Thread Ran Friedman
Hi Peng,

AFAIK GMX only uses fftw for PME.

If you have no reason to prefer flexible bonds you can use LINCS or
SHAKE to constrain the bond lengths and run with a timestep of 2fs at
room temperature.

You may want to try LAMMPS and GMX with SHAKE and a timestep of 2fs and
see if the results are similar before recompiling GMX. If I understood
you correctly you didn't constrain the bonds.

Ran

Peng Yi wrote:
>
> Hi, Ran,
>
> Do you mean the simulation I used here for testing purpose?  That is
> not long, maybe a few hours.  But my research will require simulations
> for days.  If I will not use PME, do I still need a fftw?
>
> -Peng
>
> On Wed, 4 Nov 2009, Ran Friedman wrote:
>
>> Hi Peng,
>>
>> It would be slower - never made any benchmarks on how much slower. But
>> you don't run a very long simulation, do you?
>> Installing it isn't a problem. If you use PME you also need fftw in
>> double precision though.
>>
>> Ran
>>
>> Peng Yi wrote:
>>>
>>> Hi, Ran,
>>>
>>> No, I haven't.  I still have to find out how to install in double
>>> precision.  Would double precision be slower than single?  If so,
>>> how much?  Or just double the memory used?  Thanks,
>>>
>>> -Peng
>>
>> ___
>> gmx-users mailing listgmx-users@gromacs.org
>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>> Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
>>
> ___
> gmx-users mailing listgmx-users@gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php

___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php