Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
On 2013-05-13 21:57, Reza wrote: So I think I figured out what was causing the discrepancy of Charmm27 energies between gromacs and NAMD. It appears that it's related to the gromacs version: energies from 4.5.7 and NAMD match very well while 4.6.1 gives energies that are different from both 4.5.7 and NAMD. Here are the results form my tests on POPC membrane: === BONDANGLEDIH IMP COULLJ POT gromacs-4.5.7/gromacs.log 1539.1203111.9021250.540 16.284 -1705.306 -1219.345 2993.188 gromacs-4.6.1/gromacs.log 1539.1203111.9021250.540 16.284 -5997.992 350.084 269.945 namd-2.9/namd.out 1539.1163111.9201250.542 16.284 -1705.307 -1219.343 2993.211 === I think this could be due to either 1) the simulation setup should be different between 4.5.7 and 4.6.1 when running single point calculations, or 2) my 4.6.1 compilation is not correct (although all regression tests passed), or 3) there is a bug somewhere. I'll attach a script that automates the comparison between GROMACS and NAMD in a separate post (it will make this post 50kB). If you run it and got similar/different results as mine, I'll be happy to know. Also any comment on the simulations setup is highly appreciated. Reza Thanks for reporting. I have now created a redmine issue for this http://redmine.gromacs.org/issues/1249 Could you please upload the datafiles needed to reproduce this there? Since gmx 4.5.7 is correct we can use that as a reference. On May 3, 2013, at 11:57 AM, Reza resa...@gmail.com wrote: Actually pdb which is similarly doesn't have velocity information. In this case however I'm mainly interested to see how the potential terms compare between the two packages with Charmm27 ff. Reza On 5/3/13 11:33 AM, Reza wrote: Thanks Mark! No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. The rerun option is interesting - in my case however the potential energy terms stayed identical but the kinetic term became zero. If you are using an .xtc file for the rerun, this makes sense since the .xtc does not store velocities and hence kinetic energy cannot be calculated. -Justin Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. I totally agree :) So far I found out that for no cut-off simulation in Gromacs, rather that specifying a large cut-off, it needs rlist=rvdw=rcoulomb=0 and pbc=no along with ns_type=simple and nstlist=0 (according to the manual). I am running various tests and will update if I find out what is causing the discrepancy. Reza On May 1, 2013, at 5:46 PM, Mark Abraham mark.j.abra...@gmail.com wrote: On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 The intra-molecule terms look fine. Since this is a lipid, there are non-bonded interactions that are intra-molecule, so the non-bondeds also seem fine. For POPC Memb (kcal/mol) Gromacs NAMD Diff
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
Sure - I added the input files to the Redmine issue. Reza On May 14, 2013, at 9:21 AM, David van der Spoel sp...@xray.bmc.uu.se wrote: On 2013-05-13 21:57, Reza wrote: So I think I figured out what was causing the discrepancy of Charmm27 energies between gromacs and NAMD. It appears that it's related to the gromacs version: energies from 4.5.7 and NAMD match very well while 4.6.1 gives energies that are different from both 4.5.7 and NAMD. Here are the results form my tests on POPC membrane: === BONDANGLEDIH IMP COULLJ POT gromacs-4.5.7/gromacs.log 1539.1203111.9021250.540 16.284 -1705.306 -1219.345 2993.188 gromacs-4.6.1/gromacs.log 1539.1203111.9021250.540 16.284 -5997.992 350.084 269.945 namd-2.9/namd.out 1539.1163111.9201250.542 16.284 -1705.307 -1219.343 2993.211 === I think this could be due to either 1) the simulation setup should be different between 4.5.7 and 4.6.1 when running single point calculations, or 2) my 4.6.1 compilation is not correct (although all regression tests passed), or 3) there is a bug somewhere. I'll attach a script that automates the comparison between GROMACS and NAMD in a separate post (it will make this post 50kB). If you run it and got similar/different results as mine, I'll be happy to know. Also any comment on the simulations setup is highly appreciated. Reza Thanks for reporting. I have now created a redmine issue for this http://redmine.gromacs.org/issues/1249 Could you please upload the datafiles needed to reproduce this there? Since gmx 4.5.7 is correct we can use that as a reference. On May 3, 2013, at 11:57 AM, Reza resa...@gmail.com wrote: Actually pdb which is similarly doesn't have velocity information. In this case however I'm mainly interested to see how the potential terms compare between the two packages with Charmm27 ff. Reza On 5/3/13 11:33 AM, Reza wrote: Thanks Mark! No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. The rerun option is interesting - in my case however the potential energy terms stayed identical but the kinetic term became zero. If you are using an .xtc file for the rerun, this makes sense since the .xtc does not store velocities and hence kinetic energy cannot be calculated. -Justin Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. I totally agree :) So far I found out that for no cut-off simulation in Gromacs, rather that specifying a large cut-off, it needs rlist=rvdw=rcoulomb=0 and pbc=no along with ns_type=simple and nstlist=0 (according to the manual). I am running various tests and will update if I find out what is causing the discrepancy. Reza On May 1, 2013, at 5:46 PM, Mark Abraham mark.j.abra...@gmail.com wrote: On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
So I think I figured out what was causing the discrepancy of Charmm27 energies between gromacs and NAMD. It appears that it's related to the gromacs version: energies from 4.5.7 and NAMD match very well while 4.6.1 gives energies that are different from both 4.5.7 and NAMD. Here are the results form my tests on POPC membrane: === BONDANGLEDIH IMP COULLJ POT gromacs-4.5.7/gromacs.log 1539.1203111.9021250.540 16.284 -1705.306 -1219.345 2993.188 gromacs-4.6.1/gromacs.log 1539.1203111.9021250.540 16.284 -5997.992 350.084 269.945 namd-2.9/namd.out 1539.1163111.9201250.542 16.284 -1705.307 -1219.343 2993.211 === I think this could be due to either 1) the simulation setup should be different between 4.5.7 and 4.6.1 when running single point calculations, or 2) my 4.6.1 compilation is not correct (although all regression tests passed), or 3) there is a bug somewhere. I'll attach a script that automates the comparison between GROMACS and NAMD in a separate post (it will make this post 50kB). If you run it and got similar/different results as mine, I'll be happy to know. Also any comment on the simulations setup is highly appreciated. Reza On May 3, 2013, at 11:57 AM, Reza resa...@gmail.com wrote: Actually pdb which is similarly doesn't have velocity information. In this case however I'm mainly interested to see how the potential terms compare between the two packages with Charmm27 ff. Reza On 5/3/13 11:33 AM, Reza wrote: Thanks Mark! No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. The rerun option is interesting - in my case however the potential energy terms stayed identical but the kinetic term became zero. If you are using an .xtc file for the rerun, this makes sense since the .xtc does not store velocities and hence kinetic energy cannot be calculated. -Justin Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. I totally agree :) So far I found out that for no cut-off simulation in Gromacs, rather that specifying a large cut-off, it needs rlist=rvdw=rcoulomb=0 and pbc=no along with ns_type=simple and nstlist=0 (according to the manual). I am running various tests and will update if I find out what is causing the discrepancy. Reza On May 1, 2013, at 5:46 PM, Mark Abraham mark.j.abra...@gmail.com wrote: On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 The intra-molecule terms look fine. Since this is a lipid, there are non-bonded interactions that are intra-molecule, so the non-bondeds also seem fine. For POPC Memb (kcal/mol) Gromacs NAMD Diff Bond 1539.1181 1539.1162 -0.0019 Angle 3111.9264 3111.9197 -0.0067 Dih 1250.5425 1250.5421 -0.0004 Imp 16.2837 16.2837 0. Coul -1837.8585
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
Here is the script:
#!/bin/bash
# note: you must have par_all27_prot_lipid.prm in the starting directory
set -e
PROGRAMS=/g1/home/resal/Programs
GMX_BIN_457=$PROGRAMS/gromacs/4.5.7/thread/bin
GMX_BIN_461=$PROGRAMS/gromacs/4.6.1/thread/bin
NMD_BIN=$PROGRAMS/namd/NAMD_2.9_Linux-x86_64-multicore
VMD_BIN=$PROGRAMS/vmd/1.9.1/bin
if [ -d prep ]; then
rm -rf prep
fi
if [ -d gromacs-4.5.7 ]; then
rm -rf gromacs-4.5.7
fi
if [ -d gromacs-4.6.1 ]; then
rm -rf gromacs-4.6.1
fi
if [ -d namd-2.9 ]; then
rm -rf namd-2.9
fi
mkdir prep
cd prep
cat psfgen EOF
# create the membrane
package require membrane
membrane -l POPC -x 40 -y 40 -o membrane
mol delete all
# delete water molecules
mol load psf membrane.psf pdb membrane.pdb
package require psfgen
readpsf membrane.psf
coordpdb membrane.pdb
set wat_sel [atomselect top water]
set wat_segs[lsort -unique [\$wat_sel get segid]]
foreach wseg \$wat_segs {
set sel [atomselect top segid \$wseg and water]
set wat_res [lsort -unique [\$sel get resid]]
foreach res \$wat_res {
delatom \$wseg \$res
}
}
writepsf memb_nowat.psf
writepdb memb_nowat.pdb
quit
EOF
$VMD_BIN/vmd -dispdev none -e psfgen
cd ..
mkdir gromacs-4.5.7
cd gromacs-4.5.7
cat pdb2gmx.in EOF
8
6
EOF
cat mdpfile.mdp EOF
integrator= md
nsteps= 0
nstlog= 1
nstlist = 0
ns_type = simple
rlist = 0
coulombtype = cut-off
rcoulomb = 0
rvdw = 0
pbc = no
EOF
GMXBIN=$GMX_BIN_457
source $GMXBIN/GMXRC
$GMXBIN/pdb2gmx -f ../prep/memb_nowat.pdb # pdb2gmx.in # choose CHARMM27 for
ff, None for water
$GMXBIN/grompp -f mdpfile.mdp -p topol.top -c ../prep/memb_nowat.pdb -o
topol.tpr
$GMXBIN/mdrun -nt 1 -s topol.tpr -g gromacs.log
cd ..
mkdir gromacs-4.6.1
cd gromacs-4.6.1
GMXBIN=$GMX_BIN_461
source $GMXBIN/GMXRC
$GMXBIN/pdb2gmx -f ../prep/memb_nowat.pdb # ../gromacs-4.5.7/pdb2gmx.in
$GMXBIN/grompp -f ../gromacs-4.5.7/mdpfile.mdp -p topol.top -c
../prep/memb_nowat.pdb -o topol.tpr
$GMXBIN/mdrun -nt 1 -s topol.tpr -g gromacs.log
cd ..
mkdir namd-2.9
cd namd-2.9
cat conf EOF
structure../prep/memb_nowat.psf
coordinates ../prep/memb_nowat.pdb
paratypecharmm on
parameters ../par_all27_prot_lipid.prm
exclude scaled1-4
1-4scaling 1.0
switching off
cutoff 1000
pairlistdist1000
timestep1.0
outputenergies 1
outputtiming1
binaryoutputno
outputname namd
dcdfreq 1
temperature 300
run 0
EOF
$NMD_BIN/namd2 conf namd.out
cd ..
cat analyze_energies.py EOF
import sys
lines = open(sys.argv[1]).readlines()
if sys.argv[1].endswith('log'):
match = lambda line, fields: all([f in line for f in fields])
for i, line in enumerate(lines):
if match(line, ('Bond', 'U-B', 'LJ-14')):
B, A, D, I, L1 = map(float, lines[i+1].split())
elif match(line, ('Coulomb-14','LJ (SR)')):
C1,L2, C2,P, K = map(float, lines[i+1].split())
elif match(line, ('Total Energy', 'Temperature', 'Pressure
(bar)')):
break
print sys.argv[1]
print '%10.3f %10.3f %10.3f %10.3f %10.3f %10.3f %10.3f' % (
B/4.184, A/4.184,D/4.184,I/4.184,
(C1+C2)/4.184, (L1+L2)/4.184, P/4.184)
elif sys.argv[1].endswith('out'):
for line in lines:
if line.startswith('ENERGY:'):
f = line.split()
if f[1] == '0':
E = map(float, f[2:])
B, A, D, I, C, L = E[0:6]
P = E[11]
break
print sys.argv[1]
print '%10.3f %10.3f %10.3f %10.3f %10.3f %10.3f %10.3f' % (
B, A,D,I, C, L, P)
EOF
echo
echo ===
echo BONDANGLEDIH IMP COULLJ
POT
python analyze_energies.py gromacs-4.5.7/gromacs.log
python analyze_energies.py gromacs-4.6.1/gromacs.log
python analyze_energies.py namd-2.9/namd.out
echo ===
echo
On May 13, 2013, at 3:57 PM, Reza resa...@gmail.com wrote:
So I think I figured out what was causing the discrepancy of Charmm27
energies between gromacs and NAMD. It appears that it's related to the
gromacs version: energies from 4.5.7 and NAMD match very well while 4.6.1
gives energies that are different from both 4.5.7 and NAMD.
Here are the results form my tests on POPC membrane:
===
BONDANGLEDIH IMP COULLJ POT
gromacs-4.5.7/gromacs.log
1539.1203111.9021250.540 16.284 -1705.306 -1219.345
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
Thanks Mark! No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. The rerun option is interesting - in my case however the potential energy terms stayed identical but the kinetic term became zero. Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. I totally agree :) So far I found out that for no cut-off simulation in Gromacs, rather that specifying a large cut-off, it needs rlist=rvdw=rcoulomb=0 and pbc=no along with ns_type=simple and nstlist=0 (according to the manual). I am running various tests and will update if I find out what is causing the discrepancy. Reza On May 1, 2013, at 5:46 PM, Mark Abraham mark.j.abra...@gmail.com wrote: On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 The intra-molecule terms look fine. Since this is a lipid, there are non-bonded interactions that are intra-molecule, so the non-bondeds also seem fine. For POPC Memb (kcal/mol) Gromacs NAMD Diff Bond 1539.1181 1539.1162 -0.0019 Angle 3111.9264 3111.9197 -0.0067 Dih 1250.5425 1250.5421 -0.0004 Imp 16.2837 16.2837 0. Coul -1837.8585 -1705.3075 132.5510 LJ -995.0311 -1219.3432 -224.3121 Pot 3084.9904 2993.2110 -91.7794 Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. Mark On Tue, Apr 30, 2013 at 8:33 PM, Justin Lemkul jalem...@vt.edu wrote: On 4/30/13 4:19 PM, Reza Salari wrote: Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) It appears you hit send too early. Please provide the entirety of the details you intended. Complete .mdp files and actual quantification of the differences you are observing are also very important. -Justin -- ==**== Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.**vt.edu/Pages/Personal/justin http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin ==**== -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-users http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/**
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
On 5/3/13 11:33 AM, Reza wrote: Thanks Mark! No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. The rerun option is interesting - in my case however the potential energy terms stayed identical but the kinetic term became zero. If you are using an .xtc file for the rerun, this makes sense since the .xtc does not store velocities and hence kinetic energy cannot be calculated. -Justin Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. I totally agree :) So far I found out that for no cut-off simulation in Gromacs, rather that specifying a large cut-off, it needs rlist=rvdw=rcoulomb=0 and pbc=no along with ns_type=simple and nstlist=0 (according to the manual). I am running various tests and will update if I find out what is causing the discrepancy. Reza On May 1, 2013, at 5:46 PM, Mark Abraham mark.j.abra...@gmail.com wrote: On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 The intra-molecule terms look fine. Since this is a lipid, there are non-bonded interactions that are intra-molecule, so the non-bondeds also seem fine. For POPC Memb (kcal/mol) Gromacs NAMD Diff Bond 1539.1181 1539.1162 -0.0019 Angle 3111.9264 3111.9197 -0.0067 Dih 1250.5425 1250.5421 -0.0004 Imp 16.2837 16.2837 0. Coul -1837.8585 -1705.3075 132.5510 LJ -995.0311 -1219.3432 -224.3121 Pot 3084.9904 2993.2110 -91.7794 Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. Mark On Tue, Apr 30, 2013 at 8:33 PM, Justin Lemkul jalem...@vt.edu wrote: On 4/30/13 4:19 PM, Reza Salari wrote: Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) It appears you hit send too early. Please provide the entirety of the details you intended. Complete .mdp files and actual quantification of the differences you are observing are also very important. -Justin -- ==**== Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.**vt.edu/Pages/Personal/justin http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin ==**== -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-users http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
Actually pdb which is similarly doesn't have velocity information. In this case however I'm mainly interested to see how the potential terms compare between the two packages with Charmm27 ff. Reza On 5/3/13 11:33 AM, Reza wrote: Thanks Mark! No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. The rerun option is interesting - in my case however the potential energy terms stayed identical but the kinetic term became zero. If you are using an .xtc file for the rerun, this makes sense since the .xtc does not store velocities and hence kinetic energy cannot be calculated. -Justin Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. I totally agree :) So far I found out that for no cut-off simulation in Gromacs, rather that specifying a large cut-off, it needs rlist=rvdw=rcoulomb=0 and pbc=no along with ns_type=simple and nstlist=0 (according to the manual). I am running various tests and will update if I find out what is causing the discrepancy. Reza On May 1, 2013, at 5:46 PM, Mark Abraham mark.j.abra...@gmail.com wrote: On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 The intra-molecule terms look fine. Since this is a lipid, there are non-bonded interactions that are intra-molecule, so the non-bondeds also seem fine. For POPC Memb (kcal/mol) Gromacs NAMD Diff Bond 1539.1181 1539.1162 -0.0019 Angle 3111.9264 3111.9197 -0.0067 Dih 1250.5425 1250.5421 -0.0004 Imp 16.2837 16.2837 0. Coul -1837.8585 -1705.3075 132.5510 LJ -995.0311 -1219.3432 -224.3121 Pot 3084.9904 2993.2110 -91.7794 Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. Mark On Tue, Apr 30, 2013 at 8:33 PM, Justin Lemkul jalem...@vt.edu wrote: On 4/30/13 4:19 PM, Reza Salari wrote: Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) It appears you hit send too early. Please provide the entirety of the details you intended. Complete .mdp files and actual quantification of the differences you are observing are also very important. -Justin -- ==**== Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 For POPC Memb (kcal/mol) Gromacs NAMD Diff Bond 1539.1181 1539.1162 -0.0019 Angle 3111.9264 3111.9197 -0.0067 Dih 1250.5425 1250.5421 -0.0004 Imp 16.2837 16.2837 0. Coul -1837.8585 -1705.3075 132.5510 LJ -995.0311 -1219.3432 -224.3121 Pot 3084.9904 2993.2110 -91.7794 On Tue, Apr 30, 2013 at 8:33 PM, Justin Lemkul jalem...@vt.edu wrote: On 4/30/13 4:19 PM, Reza Salari wrote: Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) It appears you hit send too early. Please provide the entirety of the details you intended. Complete .mdp files and actual quantification of the differences you are observing are also very important. -Justin -- ==**== Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.**vt.edu/Pages/Personal/justinhttp://www.bevanlab.biochem.vt.edu/Pages/Personal/justin ==**== -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-usershttp://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/** Support/Mailing_Lists/Searchhttp://www.gromacs.org/Support/Mailing_Lists/Searchbefore 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_Listshttp://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-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
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
On Wed, May 1, 2013 at 2:32 PM, Reza Salari resa...@gmail.com wrote: Hi Justin, I actually did :) but it ended up being bigger than 50 kb so it needed moderator approval to show up. I was hoping it would've been released by now. I'll attach a the details below. Any help/hint is highly appreciated. Reza Details: *1)* Both systems were prepared using VMD membrane package and then waters were removed. I used Gromacs 4.6.1 and NAMD 2.9. *2)* Simulations were run in vacuum as a single-point energy calculations (0 step). No, they weren't. See http://www.gromacs.org/Documentation/How-tos/Single-Point_Energy You cannot hope to reproduce accurately the energy of a configuration if you let the coordinates be manipulated. PME was not used. *3) *For Gromacs runs, pdb2gmx was used to prepare the system and the output was saved as the pdb format. The mdp file: integrator= md nsteps= 0 nstlog = 1 nstlist= 1 ns_type= grid rlist= 100.0 coulombtype= cut-off rcoulomb= 100.0 rvdw= 100.0 pbc = no *4) *NAMD input file: structure ../0_prep/memb_nowat.psf paratypecharmm on parameterspar_all27_prot_lipid.prm exclude scaled1-4 1-4scaling 1.0 switching off switchdist 8 cutoff 1000 pairlistdist 1000 margin 0 timestep 1.0 outputenergies 1 outputtiming 1 binaryoutput no coordinates ../0_prep/memb_nowat.pdb outputname out dcdfreq 10 temperature 300 run 0 *5)* Energies: For Single POPC (kcal/mol) Gromacs NAMD Diff Bond 43.0022 43.0015 -0.0007 Angle 80.6568 80.6571 0.0003 Dih 29.8083 29.8083 0. Imp 0.8452 0.8452 0. Coul -17.2983 -17.2983 0. LJ -7.0798 -7.0798 0. Pot 129.9343 129.9340 -0.0003 The intra-molecule terms look fine. Since this is a lipid, there are non-bonded interactions that are intra-molecule, so the non-bondeds also seem fine. For POPC Memb (kcal/mol) Gromacs NAMD Diff Bond 1539.1181 1539.1162 -0.0019 Angle 3111.9264 3111.9197 -0.0067 Dih 1250.5425 1250.5421 -0.0004 Imp 16.2837 16.2837 0. Coul -1837.8585 -1705.3075 132.5510 LJ -995.0311 -1219.3432 -224.3121 Pot 3084.9904 2993.2110 -91.7794 Something you think is equivalent is not :-) Move to testing a system with two lipids. Inspect all the logfile outputs very carefully for clues. Mark On Tue, Apr 30, 2013 at 8:33 PM, Justin Lemkul jalem...@vt.edu wrote: On 4/30/13 4:19 PM, Reza Salari wrote: Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) It appears you hit send too early. Please provide the entirety of the details you intended. Complete .mdp files and actual quantification of the differences you are observing are also very important. -Justin -- ==**== Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.**vt.edu/Pages/Personal/justin http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin ==**== -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-users http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/** Support/Mailing_Lists/Search http://www.gromacs.org/Support/Mailing_Lists/Searchbefore 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 http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-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 -- gmx-users mailing listgmx-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
[gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) -- gmx-users mailing listgmx-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
Re: [gmx-users] Comparing Charmm 27 energies from GROMACS and NAMD
On 4/30/13 4:19 PM, Reza Salari wrote: Hi I have set up two small systems, one with a single POPC lipid, and another system with 23 POPC's arranged as a lipid bilayer. I am hoping to do a similar comparison as in Table 1 in the Par Bjelkmar et al paper (porting charmm ff to gromacs) for my systems. My main question is that for the single POPC, all the potential energy terms match very well, but for the membrane system the non-bonding terms differ significantly. I am providing the full details below and greatly appreciate any hint for better comparison of the energies. Thanks, Reza Salari Details: 1) Both systems were prepared using VMD membrane package: 2) It appears you hit send too early. Please provide the entirety of the details you intended. Complete .mdp files and actual quantification of the differences you are observing are also very important. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-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
