[gmx-users] Diagnosing + system blowing up

Mark Abraham Mark.Abraham at anu.edu.au
Mon Jul 30 15:35:41 CEST 2012


On 30/07/2012 8:49 PM, Shima Arasteh wrote:
>
>   
>> If you want further help, you will need to paste your .rtp entry(s), the N-terminal fragment of your pdb2gmx -f input, pdb2gmx command line, full pdb2gmx output, the N-terminal fragment of the grompp -c input, grompp >command line, and full grompp output. You may think your context is clear, but nobody else is paying your problem as much attention as you are. I won't help further if you only provide partial information.
> 1. rtp entry
>
> [ FVAL ]
>   [ atoms ]
>      CN    C     0.357    0
>      ON    O    -0.51    1
>      H1    HA    0.100    2
>      N    NH1    -0.423    3
>      HN    H    0.333    4
>      CA    CT1    0.034    5
>      HA    HB    0.09    6
>      CB    CT1    -0.093    7
>      HB    HA    0.09    8
>      CG1    CT3    -0.268    9
>      HG11    HA    0.09    10
>      HG12    HA    0.09    11
>      HG13    HA    0.09    12
>      CG2    CT3    -0.268    13
>      HG21    HA    0.09    14
>      HG22    HA    0.09    15
>      HG23    HA    0.09    16
>      C    C    0.528    17
>      O    O    -0.510    18
>   [ bonds ]
>      CN    H1
>      CN    ON
>      CN    N
>      N    HN
>      CA    N
>      CA    HA
>      CA    C
>      C    O
>      CA    CB
>      CB    HB
>      CB    CG1
>      CB    CG2
>      CG2    HG21
>      CG2    HG22
>      CG2    HG23
>      CG1    HG11
>      CG1    HG12
>      CG1    HG13grompp -f ions.mdp -c monomer_solv.gro -p topol.top -o ions.tpr
>       
>   [ impropers ]
>      CN     N    ON    H1
>
> 2. hdb entry

In which file?

> FVAL    6
> 1    1    H1    CN    N    ON

This should be generating H1 bonded to CN...

> 1    1    HN    N    C    CA
> 1    5    HA    CA    N    C    CB
> 1    5    HB    CB    CA    CG1    CG2
> 3    4    HG1    CG1    CB    CA
> 3    4    HG2    CG2    CB    CA
>
> 3. N-terminal fragment
> HETATM    1  CN  FVAL    1      -0.721   1.600   1.249
> HETATM    2  ON  FVAL    1      -0.839   2.806   1.453
> ATOM      3  N   FVAL    1      -1.227   0.728   2.125
> ATOM      4  CA  FVAL    1      -1.918   1.159   3.323
> ATOM      5  C   FVAL    1      -1.969   2.678   3.410
> ATOM      6  O   FVAL    1      -0.931   3.335   3.447
> ATOM      7  CB  FVAL    1      -1.219   0.644   4.576
> ATOM      8  CG1 FVAL    1       0.208   1.178   4.618
> ATOM      9  CG2 FVAL    1      -1.976   1.118   5.812
>
> 4. pdb2gmx command
> #pdb2gmx -f monomer.pdb -o monomer.gro -water tip3p -ter
>
> 5. pdb2gmx output
> Using the Charmm36-modified force field in directory ./charmm36-modified.ff
>
> Opening force field file ./charmm36-modified.ff/aminoacids.r2b
> Opening force field file ./charmm36-modified.ff/rna.r2b
> Reading monomer.pdb...
> Read 177 atoms
> Analyzing pdb file
> Splitting chemical chains based on TER records or chain id changing.
> There are 1 chains and 0 blocks of water and 24 residues with 177 atoms
>
>    chain  #res #atoms
>    1 ' '    24    177
>
> All occupancy fields zero. This is probably not an X-Ray structure
> Opening force field file ./charmm36-modified.ff/atomtypes.atp
> Atomtype 1
> Reading residue database... (charmm36-modified)
> Opening force field file ./charmm36-modified.ff/aminoacids.rtp
> Residue 42
> Sorting it all out...
> Opening force field file ./charmm36-modified.ff/dna.rtp
> Residue 46
> Sorting it all out...
> Opening force field file ./charmm36-modified.ff/lipids.rtp
> Residue 82
> Sorting it all out...
> Opening force field file ./charmm36-modified.ff/rna.rtp
> Residue 86
> Sorting it all out...
> Opening force field file ./charmm36-modified.ff/aminoacids.hdb

... and here's your .hdb file being picked up (and pdb2gmx must have 
found definitions for building hydrogens for FVAL somehow)...

> Opening force field file ./charmm36-modified.ff/dna.hdb
> Opening force field file ./charmm36-modified.ff/lipids.hdb
> Opening force field file ./charmm36-modified.ff/rna.hdb
> Opening force field file ./charmm36-modified.ff/aminoacids.n.tdb
> Opening force field file ./charmm36-modified.ff/dna.n.tdb
> Opening force field file ./charmm36-modified.ff/rna.n.tdb
> Opening force field file ./charmm36-modified.ff/aminoacids.c.tdb
> Opening force field file ./charmm36-modified.ff/dna.c.tdb
> Opening force field file ./charmm36-modified.ff/rna.c.tdb
>
> Back Off! I just backed up topol.top to ./#topol.top.1#
> Processing chain 1 (177 atoms, 24 residues)
> Identified residue FVAL1 as a starting terminus.
> Identified residue GLY24 as a ending terminus.
> 8 out of 8 lines of specbond.dat converted successfully
> Select start terminus type for FVAL-1
>   0: NH3+
>   1: NH2
>   2: None
> 2
> Start terminus FVAL-1: None
> Select end terminus type for GLY-24
>   0: COO-
>   1: COOH
>   2: CT2
>   3: CT3
>   4: None
> 0
> End terminus GLY-24: COO-
> Opening force field file ./charmm36-modified.ff/aminoacids.arn
> Opening force field file ./charmm36-modified.ff/dna.arn
> Opening force field file ./charmm36-modified.ff/rna.arn
> Checking for duplicate atoms....
> Now there are 24 residues with 360 atoms
> Making bonds...
> Warning: Long Bond (1-18 = 0.357049 nm)

...but H1 is built as atom 18 and too far away...

> Number of bonds was 361, now 361
> Generating angles, dihedrals and pairs...
> Before cleaning: 920 pairs
> Before cleaning: 925 dihedrals
> Keeping all generated dihedrals
> Making cmap torsions...There are   22 cmap torsion pairs
> There are  925 dihedrals,   49 impropers,  642 angles
>             911 pairs,      361 bonds and     0 virtual sites
> Total mass 2510.906 a.m.u.
> Total charge 1.000 e
> Writing topology
>
> Back Off! I just backed up posre.itp to ./#posre.itp.1#
>
> Writing coordinate file...
>
>                  --------- PLEASE NOTE ------------
> You have successfully generated a topology from: monomer.pdb.
> The Charmm36-modified force field and the tip3p water model are used.
>                  --------- ETON ESAELP ------------
>
> 6. N-terminal fragment of the grompp -c input
>      1FVAL    CN    1   3.040   1.903   2.415
>      1FVAL    ON    2   3.028   2.024   2.435
>      1FVAL     N    3   2.989   1.816   2.503
>      1FVAL    HN    4   3.030   1.758   2.432
>      1FVAL    CA    5   2.920   1.859   2.622
>      1FVAL    HA    6   2.829   1.820   2.613
>      1FVAL    CB    7   2.990   1.807   2.748
>      1FVAL    HB    8   2.992   1.707   2.746
>      1FVAL   CG1    9   3.133   1.861   2.752
>      1FVAL  HG11   10   3.179   1.827   2.834
>      1FVAL  HG12   11   3.182   1.830   2.671
>      1FVAL  HG13   12   3.131   1.961   2.753
>      1FVAL   CG2   13   2.914   1.855   2.871
>      1FVAL  HG21   14   2.960   1.821   2.953
>      1FVAL  HG22   15   2.912   1.955   2.873
>      1FVAL  HG23   16   2.821   1.820   2.868
>      1FVAL     C   17   2.915   2.011   2.631
>      1FVAL    H1   18   2.825   2.031   2.670

... and this shows atom H1 around 0.1nm from atom C. I can think of no 
legitimate reason how this could occur, given your stated .hdb, but I 
seem to recall an old .hdb version you showed had a line

1    1    H1    C    N    ON

rather than the correct

1    1    H1    CN    N    ON

and this would explain the weird H1 position perfectly (but not how it still gets bonded to atom 1). I'd like you to double check that ./charmm36-modified.ff/aminoacids.hdb has (only) the correct FVAL entry. Also, I'd like to see the first 30 or so lines of the [bonds] section of this [moleculetype] in topol.top, so we can really see what bonds are generated.

I think you may have encountered a new bug in pdb2gmx.

>      1FVAL     O   19   3.019   2.076   2.635
>
> 7.grompp command line
> # grompp -f ions.mdp -c monomer_solv.gro -p topol.top -o ions.tpr
>
>
> 8. grompp output
> Generated 21528 of the 21528 non-bonded parameter combinations
> Generating 1-4 interactions: fudge = 1
> Generated 18355 of the 21528 1-4 parameter combinations
>
> ERROR 1 [file topol.top, line 416]:
>    No default Bond types
>
>
> ERROR 2 [file topol.top, line 1693]:
>    No default U-B types
>
>
> ERROR 3 [file topol.top, line 1694]:
>    No default U-B types
>
>
> ERROR 4 [file topol.top, line 2338]:
>    No default Proper Dih. types
>
>
> ERROR 5 [file topol.top, line 3265]:
>    No default Improper Dih. types

All this probably follows on from the confused state of bonding to H1, 
e.g. if CN-H1 is defined from the .rtp bond, and there's now a C-H1 from 
the possibly erroneous .hdb construction. What is line 416 of topol.top?

Mark

>
> Excluding 3 bonded neighbours molecule type 'Protein'
> Excluding 2 bonded neighbours molecule type 'SOL'
>
> NOTE 1 [file topol.top, line 3365]:
>    System has non-zero total charge: 1.000000
>    Total charge should normally be an integer. See
>    http://www.gromacs.org/Documentation/Floating_Point_Arithmetic
>    for discussion on how close it should be to an integer.
>    
>
>
>
> There was 1 note
>
> -------------------------------------------------------
> Program grompp, VERSION 4.5.5
> Source code file: /home/abuild/rpmbuild/BUILD/gromacs-4.5.5/src/kernel/grompp.c, line: 1372
>
> Fatal error:
> There were 5 errors in input file(s)
> For more information and tips for troubleshooting, please check the GROMACS
> website at http://www.gromacs.org/Documentation/Errors
>
>
>
> I hope there are enough. If anything else, please let me know.
>
>
> Thanks,
> Shima




More information about the gromacs.org_gmx-users mailing list