Bug #337
generating pairs causing trouble when other nonbon_params are available
Description
If there are any nonbond_params are given in the topology file and setting gen-pairs to yes. Gromacs is messing around with these values.
Given the following topolgy file:
[ defaults ]
; nbfunc comb-rule gen-pairs fudgeLJ fudgeQQ
1 3 yes 0.5 0.5 [ atomtypes ]
;name mass charge ptype c6 c12
ow 16.00000 -0.84760 A 0.31656 0.65019
hw 1.00000 0.42380 A 0.00000 0.00000
c3_1 15.02400 0.25000 A 0.37500 0.81482
o_2 16.00000 -0.50000 A 0.28000 0.45730
c2_3 14.01600 0.25000 A 0.39500 0.38247[ nonbond_params ]
; i j func c6 c12
hw ow 1 0.000000E+00 0.000000E+00
hw c3_1 1 0.000000E+00 0.000000E+00
hw o_2 1 0.000000E+00 0.000000E+00
hw c2_3 1 0.000000E+00 0.000000E+00
ow c3_1 1 0.497636E-02 0.850573E-05
ow c2_3 1 0.404546E-02 0.820460E-05
c3_1 o_2 1 0.301272E-02 0.371729E-05
o_2 c2_3 1 0.247229E-02 0.365378E-05
c3_1 c2_3 1 0.727198E-02 0.236819E-04
ow o_2 1 0.153612E-02 0.108186E-05[ moleculetype ]
; molname nrexcl
SOL 1[ atoms ]
; nr type resnr residue atom cgnr charge mass
1 ow 1 SOL ow 1 -0.8476 15.99940
2 hw 1 SOL hw 1 0.4238 1.00800
3 hw 1 SOL hw 1 0.4238 1.00800[ settles ]
; ow funct doh dhh
1 1 0.1 0.16333[ exclusions ]
1 2 3
2 1 3
3 1 2[ moleculetype ]
; Name nrexcl
poe3 3[ atoms ]
; nr type resnr residue atom cgnr charge mass
1 c3_1 1 poe3 c3_1 1 0.2500 15.0240
2 o_2 1 poe3 o_2 2 -0.5000 16.0000
3 c2_3 1 poe3 c2_3 3 0.2500 14.0160
4 c2_3 1 poe3 c2_3 4 0.2500 14.0160
5 o_2 1 poe3 o_2 5 -0.5000 16.0000
6 c2_3 1 poe3 c2_3 6 0.2500 14.0160
7 c2_3 1 poe3 c2_3 7 0.2500 14.0160
8 o_2 1 poe3 o_2 8 -0.5000 16.0000
9 c3_1 1 poe3 c3_1 9 0.2500 15.0240[ constraints ]
; ai aj
2 1 1 0.141000E+00
3 2 1 0.141000E+00
4 3 1 0.154000E+00
5 4 1 0.141000E+00
6 5 1 0.141000E+00
7 6 1 0.154000E+00
8 7 1 0.141000E+00
9 8 1 0.141000E+00[ angles ]
; ai aj ak
1 2 3 1 0.112000E+03 0.519600E+03
2 3 4 1 0.112000E+03 0.418200E+03
3 4 5 1 0.112000E+03 0.418200E+03
4 5 6 1 0.112000E+03 0.519600E+03
5 6 7 1 0.112000E+03 0.418200E+03
6 7 8 1 0.112000E+03 0.418200E+03
7 8 9 1 0.112000E+03 0.519600E+03[ pairs ]
; ai aj
4 1 1 ; 0.363599E-02 0.118410E-04 0.500000E+00
5 2 1 ; 0.440735E-03 0.212386E-06 0.500000E+00
6 3 1 ; 0.290541E-02 0.110354E-04 0.500000E+00
7 4 1 ; 0.290541E-02 0.110354E-04 0.500000E+00
8 5 1 ; 0.440735E-03 0.212386E-06 0.500000E+00
9 6 1 ; 0.363599E-02 0.118410E-04 0.500000E+00[ dihedrals ]
; ai aj ak al
4 3 2 1 1 0.0000 -5.1600 1
4 3 2 1 1 0.0000 -0.6971 2
4 3 2 1 1 0.0000 5.3501 3
4 3 2 1 1 0.0000 0.8031 4
4 3 2 1 1 0.0000 0.2831 5
4 3 2 1 1 0.0000 0.0953 6
4 3 2 1 1 0.0000 -0.0580 7
5 4 3 2 1 0.0000 7.5853 1
5 4 3 2 1 0.0000 6.7052 2
5 4 3 2 1 0.0000 8.4007 3
5 4 3 2 1 0.0000 0.6322 4
5 4 3 2 1 0.0000 0.1106 5
5 4 3 2 1 0.0000 0.3596 6
5 4 3 2 1 0.0000 0.0168 7
6 5 4 3 1 0.0000 -5.1600 1
6 5 4 3 1 0.0000 -0.6971 2
6 5 4 3 1 0.0000 5.3501 3
6 5 4 3 1 0.0000 0.8031 4
6 5 4 3 1 0.0000 0.2831 5
6 5 4 3 1 0.0000 0.0953 6
6 5 4 3 1 0.0000 -0.0580 7
7 6 5 4 1 0.0000 -5.1600 1
7 6 5 4 1 0.0000 -0.6971 2
7 6 5 4 1 0.0000 5.3501 3
7 6 5 4 1 0.0000 0.8031 4
7 6 5 4 1 0.0000 0.2831 5
7 6 5 4 1 0.0000 0.0953 6
7 6 5 4 1 0.0000 -0.0580 7
8 7 6 5 1 0.0000 7.5853 1
8 7 6 5 1 0.0000 6.7052 2
8 7 6 5 1 0.0000 8.4007 3
8 7 6 5 1 0.0000 0.6322 4
8 7 6 5 1 0.0000 0.1106 5
8 7 6 5 1 0.0000 0.3596 6
8 7 6 5 1 0.0000 0.0168 7
9 8 7 6 1 0.0000 -5.1600 1
9 8 7 6 1 0.0000 -0.6971 2
9 8 7 6 1 0.0000 5.3501 3
9 8 7 6 1 0.0000 0.8031 4
9 8 7 6 1 0.0000 0.2831 5
9 8 7 6 1 0.0000 0.0953 6
9 8 7 6 1 0.0000 -0.0580 7[ exclusions ]
; ai aj
2 1
3 2
4 3
5 4
6 5
7 6
8 7
9 8
1 3
2 4
3 5
4 6
5 7
6 8
7 9
4 1
5 2
6 3
7 4
8 5
9 6[ system ]
; Name
TEST[ molecules ]
; Compound #mols
SOL 924
poe3 76gmxdump gives following results
[snip]
ffparams:
atnr=5
ntypes=47
functype[0]=LJ_SR, c6= 2.61719758e-03, c12= 2.63373886e-06
functype[1]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[2]=LJ_SR, c6= 5.16704677e-19, c12= 7.84717270e-33
functype[3]=LJ_SR, c6= 5.68567269e-23, c12= 7.47021001e-40
functype[4]=LJ_SR, c6= 1.43854954e-19, c12= 6.30568414e-34
functype[5]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[6]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[7]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[8]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[9]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[10]=LJ_SR, c6= 5.16704677e-19, c12= 7.84717270e-33
functype[11]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[12]=LJ_SR, c6= 9.06377845e-03, c12= 2.52055906e-05
functype[13]=LJ_SR, c6= 1.11183180e-20, c12= 8.31365008e-36
functype[14]=LJ_SR, c6= 1.40085303e-17, c12= 2.07161281e-30
functype[15]=LJ_SR, c6= 5.68567269e-23, c12= 7.47021001e-40
functype[16]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[17]=LJ_SR, c6= 1.11183180e-20, c12= 8.31365008e-36
functype[18]=LJ_SR, c6= 8.81473767e-04, c12= 4.24773674e-07
functype[19]=LJ_SR, c6= 3.33732824e-21, c12= 7.62071036e-37
functype[20]=LJ_SR, c6= 1.43854954e-19, c12= 6.30568414e-34
functype[21]=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype[22]=LJ_SR, c6= 1.40085303e-17, c12= 2.07161281e-30
functype[23]=LJ_SR, c6= 3.33732824e-21, c12= 7.62071036e-37
functype[24]=LJ_SR, c6= 5.81085496e-03, c12= 2.20710335e-05
[snip]whereas setting the gen-pair option to know results in something completley different :
[snip]
ffparams:
atnr=5
ntypes=47
functype0=LJ_SR, c6= 2.61719758e-03, c12= 2.63373886e-06
functype1=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype2=LJ_SR, c6= 4.87049343e-03, c12= 8.14769919e-06
functype3=LJ_SR, c6= 1.51887815e-03, c12= 1.05770641e-06
functype4=LJ_SR, c6= 3.89976287e-03, c12= 7.62425771e-06
functype5=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype6=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype7=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype8=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype9=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype10=LJ_SR, c6= 4.87049343e-03, c12= 8.14769919e-06
functype11=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype12=LJ_SR, c6= 9.06377845e-03, c12= 2.52055906e-05
functype13=LJ_SR, c6= 2.82656797e-03, c12= 3.27210637e-06
functype14=LJ_SR, c6= 7.25729251e-03, c12= 2.35862972e-05
functype15=LJ_SR, c6= 1.51887815e-03, c12= 1.05770641e-06
functype16=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype17=LJ_SR, c6= 2.82656797e-03, c12= 3.27210637e-06
functype18=LJ_SR, c6= 8.81473767e-04, c12= 4.24773674e-07
functype19=LJ_SR, c6= 2.26320908e-03, c12= 3.06189349e-06
functype20=LJ_SR, c6= 3.89976287e-03, c12= 7.62425771e-06
functype21=LJ_SR, c6= 0.00000000e+00, c12= 0.00000000e+00
functype22=LJ_SR, c6= 7.25729251e-03, c12= 2.35862972e-05
functype23=LJ_SR, c6= 2.26320908e-03, c12= 3.06189349e-06
functype24=LJ_SR, c6= 5.81085496e-03, c12= 2.20710335e-05
functype25=SETTLE, doh= 1.00000001e-01, dhh= 1.63330004e-01
[snip]
Diff between the two files:
184,186c184,186
< functype2=LJ_SR, c6= 4.87049343e-03, c12= 8.14769919e-06
< functype3=LJ_SR, c6= 1.51887815e-03, c12= 1.05770641e-06
< functype4=LJ_SR, c6= 3.89976287e-03, c12= 7.62425771e-06
---
functype2=LJ_SR, c6= 5.16704677e-19, c12= 7.84717270e-33
functype3=LJ_SR, c6= 5.68567269e-23, c12= 7.47021001e-40
functype4=LJ_SR, c6= 1.43854954e-19, c12= 6.30568414e-34
192c192
< functype10=LJ_SR, c6= 4.87049343e-03, c12= 8.14769919e-06
---
functype10=LJ_SR, c6= 5.16704677e-19, c12= 7.84717270e-33
195,197c195,197
< functype13=LJ_SR, c6= 2.82656797e-03, c12= 3.27210637e-06
< functype14=LJ_SR, c6= 7.25729251e-03, c12= 2.35862972e-05
< functype15=LJ_SR, c6= 1.51887815e-03, c12= 1.05770641e-06
---
functype13=LJ_SR, c6= 1.11183180e-20, c12= 8.31365008e-36
functype14=LJ_SR, c6= 1.40085303e-17, c12= 2.07161281e-30
functype15=LJ_SR, c6= 5.68567269e-23, c12= 7.47021001e-40
199c199
< functype17=LJ_SR, c6= 2.82656797e-03, c12= 3.27210637e-06
---
functype17=LJ_SR, c6= 1.11183180e-20, c12= 8.31365008e-36
201,202c201,202
< functype19=LJ_SR, c6= 2.26320908e-03, c12= 3.06189349e-06
< functype20=LJ_SR, c6= 3.89976287e-03, c12= 7.62425771e-06
---
functype19=LJ_SR, c6= 3.33732824e-21, c12= 7.62071036e-37
functype20=LJ_SR, c6= 1.43854954e-19, c12= 6.30568414e-34
204,205c204,205
< functype22=LJ_SR, c6= 7.25729251e-03, c12= 2.35862972e-05
< functype23=LJ_SR, c6= 2.26320908e-03, c12= 3.06189349e-06
---
functype22=LJ_SR, c6= 1.40085303e-17, c12= 2.07161281e-30
functype23=LJ_SR, c6= 3.33732824e-21, c12= 7.62071036e-37
224c224
< functype42=LJ14, c6A= 3.62864626e-03, c12A= 1.17931486e-05, c6B= 3.62864626e-03, c12B= 1.17931486e-05
---
functype42=LJ14, c6A= 7.00426516e-18, c12A= 1.03580641e-30, c6B= 7.00426516e-18, c12B= 1.03580641e-30
****************************
taking out the nonbond sections and letting do gromacs all the work the .tpr file is identical to the one with gen-pairs set to no and given C6 and C12 Values
History
#1 Updated by Berk Hess about 10 years ago
Could you attach all the input files required to reproduce this?
BTW I don't know what the three commented out parameters
are in your pairs section, but with comb. rule you should have
sigma and epsilon.
Berk
#2 Updated by Sascha Hempel about 10 years ago
Created an attachment (id=381)
.mdp file
#3 Updated by Sascha Hempel about 10 years ago
Created an attachment (id=382)
input .gro file
#4 Updated by Sascha Hempel about 10 years ago
Created an attachment (id=383)
input .top file
#5 Updated by Berk Hess about 10 years ago
I can't reproduce this.
I get identical tpr files with and without gen-pairs.
Are you really using Gromacs 4.0?
Berk
#6 Updated by Sascha Hempel about 10 years ago
First of all thanks for noticing my mistake in the pairs sections. I was messing around a little to much with my files :)
Shouldn't there be an error created by grompp if the wrong number of parameters is given?
The values in the pairs section are sigma and epsilon and the last value is just the fudge factor.
If done several tests of possible combinations of parameters. Lets say the first gmxdump i posted gives bad results (with the e-33) and the second one gives good results.
X means commented out (or no in case of gen-pairs)
O means the opposite
gen-pairs | nonbond-section | pairs-section (with top file from the attachment)
1. X O X bad results
2. X O O bad results
3. O O X bad results
4. O O O bad results
5. X X X good results
6. O X X good results
7. O X O good results - Additional Terms
8. X X O good results - Additional Terms
To further complicate the situation:
diff between #5 & #6 shows that they are identical, as well as #7 & #8. between #5 & #7 there some differences.
When pairs is active some additional terms are generated
**********************
diff gen_no-nonbond_no-pairs_yes gen_no-nonbond_no-pairs_no
181c181
< ntypes=47
---
ntypes=44
224,228c224,225
< functype42=LJ14, c6A= 1.09442095e-19, c12A= 2.52883503e-34, c6B= 1.09442095e-19, c12B= 2.52883503e-34
< functype43=LJ14, c6A= 6.22660161e-27, c12A= 0.00000000e+00, c6B= 6.22660161e-27, c12B= 0.00000000e+00
< functype44=LJ14, c6A= 2.65517206e-20, c12A= 1.59711886e-35, c6B= 2.65517206e-20, c12B= 1.59711886e-35
< functype45=CONSTR, dA= 1.41000003e-01, dB= 1.41000003e-01
< functype46=CONSTR, dA= 1.53999999e-01, dB= 1.53999999e-01
---
functype42=CONSTR, dA= 1.41000003e-01, dB= 1.41000003e-01
functype43=CONSTR, dA= 1.53999999e-01, dB= 1.53999999e-01
379,387d375
< LJ-14:
< nr: 18
< iatoms:
< 0 type=42 (LJ14) 3 0
< 1 type=43 (LJ14) 4 1
< 2 type=44 (LJ14) 5 2
< 3 type=44 (LJ14) 6 3
< 4 type=43 (LJ14) 7 4
< 5 type=42 (LJ14) 8 5
391,398c379,386
< 0 type=45 (CONSTR) 1 0
< 1 type=45 (CONSTR) 2 1
< 2 type=46 (CONSTR) 3 2
< 3 type=45 (CONSTR) 4 3
< 4 type=45 (CONSTR) 5 4
< 5 type=46 (CONSTR) 6 5
< 6 type=45 (CONSTR) 7 6
< 7 type=45 (CONSTR) 8 7
---
0 type=42 (CONSTR) 1 0
1 type=42 (CONSTR) 2 1
2 type=43 (CONSTR) 3 2
3 type=42 (CONSTR) 4 3
4 type=42 (CONSTR) 5 4
5 type=43 (CONSTR) 6 5
6 type=42 (CONSTR) 7 6
7 type=42 (CONSTR) 8 7
At this moment i am getting lost totally in a bunch of numbers... I hope you understand what i am trying to explain here.
So long,
Sascha
#7 Updated by Berk Hess about 10 years ago
Now we seem to be talking about something completely different.
Looking at your nonbond_params section,
you also seem to have entered C6 and C12 iso sigma and epsilon.
Berk
#8 Updated by Berk Hess about 10 years ago
I assume your problem has been resolved.
Can we close this bug?
Berk
#9 Updated by Sascha Hempel about 10 years ago
Of course,
I am sorry for the inconveniences i caused. I promise to look up my data better next time.
Thanks for the help anyway
Sascha
#10 Updated by Berk Hess about 10 years ago
This was not a bug, but the result of an incorrect interpretation
of the topology file format.
Berk