[gmx-users] Force Calculation - SBLOCKS
Justin M. Shorb
shorb at wisc.edu
Fri Jun 22 00:36:50 CEST 2007
Greetings!
I have been running gromacs to simulate a chromophore (NMAd) in water.
I have verified that the coulombic force is calculated correctly by
doing a rerun of the simulation trajectory and extracting the forces
from a forcefield of vanishing force constants (thus removing all vdW
and bonded terms). This works, although it is fairly brute-force.
Now, I am trying to scale up and look at a single peptide Amide-I band
in a transmembrane protein. Unfortunately, if I use the same mdp file,
I get different results than if I use the mdp file with which the
production run was generated. By gmxdump-ing the tpr files, I can
figure out the differences between the run files. I verified that
varying PME parameters cause a negligible effect, and most of the basic
I/O parameters of course shouldn't matter.
The biggest difference between the two files is that one reports zero
for the multinr and uses nr=0 and nra=0, whereas the other tpr file
reports 123 for the multinr, nr=123, nra=6326, and then specifies
blocks[SBLOCKS][...]
The coulombic forces generated are not simply smaller or larger. In
some cases the vector goes in a different direction entirely (different
signs on the components). Is there something in the way these groups
are used to calculate the coulombic force that would account for such
wild discrepancy? How is an energy group used in PME? Is the lack of
vdW parameters causing superfluous exclusion rules? (I know that the
-check14 flag should be used to get accurate coulombic forces in the
absence of vdW forces). Is there somewhere in the manual that discusses
this? What is a VCM?
I have checked the archives, and note that Yang Ye and Dongsheng Zhang
have discussed many aspects of the tpr file, but nobody had any
information about SBLOCKS. I have included a diff of the two tpr files.
I truncated out the list of groups and kept only two representative
Sblocks definitions. I also checked the edr file, and found no other
potential energy outside of coulombic (so there shouldn't be any hidden
forces from bonded/nb terms). Also: ignore userint/user1_grps - they
are used in my on-the-fly electric field extraction code.
Thanks in advance,
Justin
====
1c1
< PM/memb.tpr: Run with the T-coupling groups, etc.
---
> cd3z_mdr2.tpr: Run with nothing fancy, mdr file from the single
chromophore in water
15c15
< nsteps = 20000
---
> nsteps = 0
18c18
< nstlist = 10
---
> nstlist = 1
22c22
< nstcomm = 10
---
> nstcomm = 1
24,30c24,30
< nstcheckpoint = 10
< nstlog = 10
< nstxout = 10
< nstvout = 10
< nstfout = 10
< nstenergy = 10
< nstxtcout = 10
---
> nstcheckpoint = 1000
> nstlog = 10000
> nstxout = 0
> nstvout = 0
> nstfout = 0
> nstenergy = 10000
> nstxtcout = 50000
32,38c32,38
< delta_t = 0.0005
< xtcprec = 10
< nkx = 54
< nky = 54
< nkz = 54
< pme_order = 4
< ewald_rtol = 1e-05
---
> delta_t = 0.002
> xtcprec = 1000
> nkx = 32
> nky = 32
> nkz = 32
> pme_order = 6
> ewald_rtol = 1e-06
41c41
< optimize_fft = TRUE
---
> optimize_fft = FALSE
45,48c45,48
< etc = Nose-Hoover
< epc = Parrinello-Rahman
< epctype = Anisotropic
< tau_p = 3
---
> etc = No
> epc = No
> epctype = Isotropic
> tau_p = 1
50,52c50,52
< ref_p[ 0]={ 1.00000e+00, 1.00000e+00, 1.00000e+00}
< ref_p[ 1]={ 1.00000e+00, 1.00000e+00, 1.00000e+00}
< ref_p[ 2]={ 1.00000e+00, 1.00000e+00, 1.00000e+00}
---
> ref_p[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> ref_p[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> ref_p[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
54,56c54,56
< compress[ 0]={ 4.50000e-05, 0.00000e+00, 0.00000e+00}
< compress[ 1]={ 0.00000e+00, 4.50000e-05, 0.00000e+00}
< compress[ 2]={ 0.00000e+00, 0.00000e+00, 4.50000e-05}
---
> compress[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> compress[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> compress[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
58c58
< rlist = 1
---
> rlist = 1.4
61c61
< rcoulomb = 1
---
> rcoulomb = 1.4
64c64
< rvdw = 1
---
> rvdw = 1.4
110,112c110,112
< userint1 = 0
< userint2 = 0
< userint3 = 0
---
> userint1 = 1
> userint2 = 4
> userint3 = 40
119,123c119,123
< nrdf: 1689 21915 11064
< ref_t: 310 310 310
< tau_t: 0.3 0.3 0.3
< anneal: No No No
< ann_npoints: 0 0 0
---
> nrdf: 40893
> ref_t: 0
> tau_t: 0
> anneal: No
> ann_npoints: 0
126,128c126
< energygrp_flags[ 0]: 0 0 0
< energygrp_flags[ 1]: 0 0 0
< energygrp_flags[ 2]: 0 0 0
---
> energygrp_flags[ 0]: 0
159c157
< nosehoover_xi: 0 0 0
---
> nosehoover_xi: 0
17445,34727c17443,34725
---
34733,34742c34731,34740
---
34747,52017c34745,52015
< grp[0] nr=3, name=[ Protein Solvent DMPC]
< grp[1] nr=3, name=[ Protein Solvent DMPC]
---
> grp[0] nr=1, name=[ rest]
> grp[1] nr=1, name=[ rest]
52020c52018
< grp[4] nr=1, name=[ rest]
---
> grp[4] nr=41, name=[ LEU31_1 LEU34_1 ISO38_1 LEU39_1 ISO41_1 GLY43_1
VAL44_1 ISO45_1 LEU46_1 LEU49_1 LEU31_2 LEU34_2 ISO38_2 LEU39_2 ISO41_2
GLY43_2 VAL44_2 ISO45_2 LEU46_2 LEU49_2 LEU31_3 LEU34_3 ISO38_3 LEU39_3
ISO41_3 GLY43_3 VAL44_3 ISO45_3 LEU46_3 LEU49_3 LEU31_4 LEU34_4 ISO38_4
LEU39_4 ISO41_4 GLY43_4 VAL44_4 ISO45_4 LEU46_4 LEU49_4 rest]
52022c52020
< grp[6] nr=3, name=[ Protein Solvent DMPC]
---
> grp[6] nr=1, name=[ rest]
90461c90459
< grpname (24):
---
> grpname (64):
90485c90483,90523
< grpname[23]={name="rest"}
---
> grpname[23]={name="LEU31_1"}
...
> grpname[62]={name="LEU49_4"}
> grpname[63]={name="rest"}
124641,125427c124679,124681
< multinr[division over processors]: 123 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
< nr=123
< nra=6326
< blocks[SBLOCKS][0][0..212]={0, 1, 2, 3, 4, 5, 6, 7, 8, 9,
10,
< 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24,
< 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38,
< 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52,
< 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66,
< 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
< 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94,
< 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106,
107,
< 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118,
119,
< 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130,
131,
< 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142,
143,
< 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154,
155,
< 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166,
167,
< 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178,
179,
< 180, 181, 182, 183, 184, 185, 186, 187, 188, 189, 190,
191,
< 192, 193, 194, 195, 196, 197, 198, 199, 200, 201, 202,
203,
< 204, 205, 206, 207, 208, 209, 210, 211, 212}
...
< blocks[SBLOCKS][4][852..897]={852, 853, 854, 855, 856,
857,
< 858, 859, 860, 861, 862, 863, 864, 865, 866, 867, 868,
869,
< 870, 871, 872, 873, 874, 875, 876, 877, 878, 879, 880,
881,
< 882, 883, 884, 885, 886, 887, 888, 889, 890, 891, 892,
893,
< 894, 895, 896, 897}
< 6325}
...
---
> multinr[division over processors]: 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
> nr=0
> nra=0
125431c124685
< ntypes=9
---
> ntypes=22
125433,125440c124687,124707
< functype[1]=ANGLES, thA= 0.00000e+00, ctA= 0.00000e+00,
thB= 0.00000e+00, ctB= 0.00000e+00
< functype[2]=PDIHS, phiA= 0.00000000e+00, cpA=
0.00000000e+00, phiB= 0.00000000e+00, cpB= 0.00000000e+00, mult=2
< functype[3]=PDIHS, phiA= 0.00000000e+00, cpA=
0.00000000e+00, phiB= 0.00000000e+00, cpB= 0.00000000e+00, mult=6
< functype[4]=PDIHS, phiA= 0.00000000e+00, cpA=
0.00000000e+00, phiB= 0.00000000e+00, cpB= 0.00000000e+00, mult=3
< functype[5]=IDIHS, xiA= 0.00000e+00, cxA= 0.00000e+00,
xiB= 0.00000e+00, cxB= 0.00000e+00
< functype[6]=LJ14, c6A= 0.00000000e+00, c12A=
0.00000000e+00, c6B= 0.00000000e+00, c12B= 0.00000000e+00
< functype[7]=CONSTR, dA= 0.00000000e+00, dB= 0.00000000e+00
< functype[8]=SETTLE, doh= 1.00000001e-01, dhh=
1.63299993e-01
---
> functype[1]=BONDS, b0A= 0.00000e+00, cbA= 0.00000e+00,
b0B= 0.00000e+00, cbB= 0.00000e+00
> functype[2]=ANGLES, thA= 1.21000e+02, ctA= 0.00000e+00,
thB= 1.21000e+02, ctB= 0.00000e+00
> functype[3]=ANGLES, thA= 1.15000e+02, ctA= 0.00000e+00,
thB= 1.15000e+02, ctB= 0.00000e+00
> functype[4]=ANGLES, thA= 1.24000e+02, ctA= 0.00000e+00,
thB= 1.24000e+02, ctB= 0.00000e+00
> functype[5]=ANGLES, thA= 1.23000e+02, ctA= 0.00000e+00,
thB= 1.23000e+02, ctB= 0.00000e+00
> functype[6]=ANGLES, thA= 1.22000e+02, ctA= 0.00000e+00,
thB= 1.22000e+02, ctB= 0.00000e+00
> functype[7]=ANGLES, thA= 1.09500e+02, ctA= 0.00000e+00,
thB= 1.09500e+02, ctB= 0.00000e+00
> functype[8]=ANGLES, thA= 1.11000e+02, ctA= 0.00000e+00,
thB= 1.11000e+02, ctB= 0.00000e+00
> functype[9]=ANGLES, thA= 1.20000e+02, ctA= 0.00000e+00,
thB= 1.20000e+02, ctB= 0.00000e+00
> functype[10]=ANGLES, thA= 1.17000e+02, ctA= 0.00000e+00,
thB= 1.17000e+02, ctB= 0.00000e+00
> functype[11]=ANGLES, thA= 1.26000e+02, ctA= 0.00000e+00,
thB= 1.26000e+02, ctB= 0.00000e+00
> functype[12]=ANGLES, thA= 0.00000e+00, ctA= 0.00000e+00,
thB= 0.00000e+00, ctB= 0.00000e+00
> functype[13]=PDIHS, phiA= 1.80000000e+02, cpA=
0.00000000e+00, phiB= 1.80000000e+02, cpB= 0.00000000e+00, mult=2
> functype[14]=PDIHS, phiA= 1.80000000e+02, cpA=
0.00000000e+00, phiB= 1.80000000e+02, cpB= 0.00000000e+00, mult=6
> functype[15]=PDIHS, phiA= 0.00000000e+00, cpA=
0.00000000e+00, phiB= 0.00000000e+00, cpB= 0.00000000e+00, mult=3
> functype[16]=PDIHS, phiA= 0.00000000e+00, cpA=
0.00000000e+00, phiB= 0.00000000e+00, cpB= 0.00000000e+00, mult=6
> functype[17]=PDIHS, phiA= 0.00000000e+00, cpA=
0.00000000e+00, phiB= 0.00000000e+00, cpB= 0.00000000e+00, mult=2
> functype[18]=IDIHS, xiA= 0.00000e+00, cxA= 0.00000e+00,
xiB= 0.00000e+00, cxB= 0.00000e+00
> functype[19]=IDIHS, xiA= 3.52640e+01, cxA= 0.00000e+00,
xiB= 3.52640e+01, cxB= 0.00000e+00
> functype[20]=LJ14, c6A= 0.00000000e+00, c12A=
0.00000000e+00, c6B= 0.00000000e+00, c12B= 0.00000000e+00
> functype[21]=SETTLE, doh= 1.00000001e-01, dhh=
1.63299993e-01
125442c124709,130930
< nr: 0
---
> nr: 18657
> multinr[division over processors]: 18657
> iatoms:
125459,133037c130947,138525
---
133052,135351c138540,140839
---
135360,136204c140848,141692
---
136211,141245c141699,146733
---
141291,147512c146779
< nr: 18657
< multinr[division over processors]: 18657
< iatoms:
---
> nr: 0
147519,151169c146786,150436
---
151199,151200c150466,150467
< T-Coupling : 852 10957 5474 (total 17283 atoms)
< Energy Mon. : 852 10957 5474 (total 17283 atoms)
---
> T-Coupling : 17283 (total 17283 atoms)
> Energy Mon. : 17283 (total 17283 atoms)
151203c150470
< User1 : 17283 (total 17283 atoms)
---
> User1 : 4 4 4 4 4 4 4
4 4 4 4 4 4 4 4 4 4 4
4 4 4 4 4 4 4 4 4 4
4 4 4 4 4 4 4 4 4 4
4 4 17123 (total 17283 atoms)
151205c150472
< VCM : 852 10957 5474 (total 17283 atoms)
---
> VCM : 17283 (total 17283 atoms)
_________________________________________________________________
Justin M. Shorb Phone: (608) 262-0483
Skinner Group shorb at wisc.edu
University of Wisconsin-Madison, Department of Chemistry
1101 University Ave., Madison, WI, 53706
More information about the gromacs.org_gmx-users
mailing list