[gmx-users] stepsize too small ... but potential energy negative!

Anna Marabotti anna.marabotti at isa.cnr.it
Mon May 24 12:26:58 CEST 2010


Dear Luca, dear all,
thank you for your hints. I made some trials with my systems and these are
my answers to your questions:

- my system is a protein (with or w/o ligand) in solvent (water SPC).
Following your suggestions, I tried to perform an EM on the protein w/o
ligand after the editconf step (i.e. I created the topology with pdb2gmx,
created a cubic box with editconf, then I used grompp+mdrun to perform EM).
I received the same error: the minimization starts, performs 37 steps, then
stops with the "stepsize too small...." communication. The energy in this
case is -3.06e+4 (one order of magnitude higher than the solvated system)
but still negative, and the maximum force is still "inf on atom 1". So
minimizing in vacuo does not help to solve the problem. This atom 1 is the
Nter of the protein, it is belonging to a Ser and I did not charge it
explicitly with pdb2gmx (i.e. I did not use the flag -ter) but it is bound
to 3 H in the .gro file.

- If I look at the protein with VMD or other visualizing tools, it seems to
me that no major problems are present on the structure. In particular, atom
1 is not "broken" or something similar, I can't see no "aberrations". I
don't know how to check if something goes wrong apart from visualization, do
you know some tools that could automatically check the file? To the best of
my knowledge, gmxcheck performs checks only on trajectories, or can I use it
also on single structures?

- I also tried to continue anyway after the minimization step with the PR MD
in NVT, but it did not start because of a lot of LINCS error. It suggests me
that some distortions are present in my structure, but if I cannot minimize
it, how to relieve this problem?

I copy&paste below the .log file for this new minimization in vacuo. The
parameters I used are the same for the minimization in solvent.
Any help will be appreciated.
Thank you very much and best regards

Input Parameters:
   integrator           = steep
   nsteps               = 50000
   init_step            = 0
   ns_type              = Grid
   nstlist              = 5
   ndelta               = 2
   nstcomm              = 1
   comm_mode            = Linear
   nstlog               = 100
   nstxout              = 100
   nstvout              = 100
   nstfout              = 0
   nstenergy            = 100
   nstxtcout            = 0
   init_t               = 0
   delta_t              = 0.001
   xtcprec              = 1000
   nkx                  = 70
   nky                  = 70
   nkz                  = 70
   pme_order            = 4
   ewald_rtol           = 1e-05
   ewald_geometry       = 0
   epsilon_surface      = 0
   optimize_fft         = TRUE
   ePBC                 = xyz
   bPeriodicMols        = FALSE
   bContinuation        = FALSE
   bShakeSOR            = FALSE
   etc                  = No
   epc                  = No
   epctype              = Isotropic
   tau_p                = 1
   ref_p (3x3):
      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}
   compress (3x3):
      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}
   refcoord_scaling     = No
   posres_com (3):
      posres_com[0]= 0.00000e+00
      posres_com[1]= 0.00000e+00
      posres_com[2]= 0.00000e+00
   posres_comB (3):
      posres_comB[0]= 0.00000e+00
      posres_comB[1]= 0.00000e+00
      posres_comB[2]= 0.00000e+00
   andersen_seed        = 815131
   rlist                = 1
   rtpi                 = 0.05
   coulombtype          = PME
   rcoulomb_switch      = 0
   rcoulomb             = 1
   vdwtype              = Cut-off
   rvdw_switch          = 0
   rvdw                 = 1.2
   epsilon_r            = 1
   epsilon_rf           = 1
   tabext               = 1
   implicit_solvent     = No
   gb_algorithm         = Still
   gb_epsilon_solvent   = 80
   nstgbradii           = 1
   rgbradii             = 2
   gb_saltconc          = 0
   gb_obc_alpha         = 1
   gb_obc_beta          = 0.8
   gb_obc_gamma         = 4.85
   sa_surface_tension   = 2.092
DispCorr             = No
   free_energy          = no
   init_lambda          = 0
   sc_alpha             = 0
   sc_power             = 0
   sc_sigma             = 0.3
   delta_lambda         = 0
   nwall                = 0
   wall_type            = 9-3
   wall_atomtype[0]     = -1
   wall_atomtype[1]     = -1
   wall_density[0]      = 0
   wall_density[1]      = 0
   wall_ewald_zfac      = 3
   pull                 = no
   disre                = No
   disre_weighting      = Conservative
   disre_mixed          = FALSE
   dr_fc                = 1000
   dr_tau               = 0
   nstdisreout          = 100
   orires_fc            = 0
   orires_tau           = 0
   nstorireout          = 100
   dihre-fc             = 1000
   em_stepsize          = 0.1
   em_tol               = 1000
   niter                = 20
   fc_stepsize          = 0
   nstcgsteep           = 1000
   nbfgscorr            = 10
   ConstAlg             = Lincs
shake_tol            = 0.0001
   lincs_order          = 4
   lincs_warnangle      = 30
   lincs_iter           = 1
   bd_fric              = 0
   ld_seed              = 1993
   cos_accel            = 0
   deform (3x3):
      deform[    0]={ 0.00000e+00,  0.00000e+00,  0.00000e+00}
      deform[    1]={ 0.00000e+00,  0.00000e+00,  0.00000e+00}
      deform[    2]={ 0.00000e+00,  0.00000e+00,  0.00000e+00}
   userint1             = 0
   userint2             = 0
   userint3             = 0
   userint4             = 0
   userreal1            = 0
   userreal2            = 0
   userreal3            = 0
   userreal4            = 0
grpopts:
   nrdf:        6039
   ref_t:           0
   tau_t:           0
anneal:          No
ann_npoints:           0
   acc:            0           0           0
   nfreeze:           N           N           N
   energygrp_flags[  0]: 0
   efield-x:
      n = 0
   efield-xt:
      n = 0
efield-y:
      n = 0
   efield-yt:
      n = 0
   efield-z:
      n = 0
   efield-zt:
      n = 0
   bQMMM                = FALSE
   QMconstraints        = 0
   QMMMscheme           = 0
   scalefactor          = 1
qm_opts:
   ngQM                 = 0
Table routines are used for coulomb: TRUE
Table routines are used for vdw:     FALSE
Will do PME sum in reciprocal space.

++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
U. Essman, L. Perela, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen
A smooth particle mesh Ewald method
J. Chem. Phys. 103 (1995) pp. 8577-8592
-------- -------- --- Thank You --- -------- --------

Using a Gaussian width (1/beta) of 0.320163 nm for Ewald
Cut-off's:   NS: 1   Coulomb: 1   LJ: 1.2
System total charge: -3.000
Generated table with 4400 data points for Ewald.
Tabscale = 2000 points/nm
Generated table with 4400 data points for LJ6.
Tabscale = 2000 points/nm
Generated table with 4400 data points for LJ12.
Tabscale = 2000 points/nm
Generated table with 4400 data points for 1-4 COUL.
Tabscale = 2000 points/nm
Generated table with 4400 data points for 1-4 LJ6.
Tabscale = 2000 points/nm
Generated table with 4400 data points for 1-4 LJ12.
Tabscale = 2000 points/nm
Configuring nonbonded kernels...
Testing x86_64 SSE2 support... present.


Removing pbc first time
Initiating Steepest Descents
Max number of connections per atom is 27
Total number of connections is 27956
Max number of graph edges per atom is 4
Total number of graph edges is 4086
Started Steepest Descents on node 0 Mon May 24 12:03:01 2010
Steepest Descents:
   Tolerance (Fmax)   =  1.00000e+03
   Number of steps    =        50000
Grid: 7 x 7 x 7 cells
           Step           Time         Lambda
              0        0.00000        0.00000

   Energies (kJ/mol)
        G96Bond       G96Angle    Proper Dih.  Improper Dih.          LJ-14
    2.40716e+03    1.23569e+03    6.85831e+02    8.70819e+01   2.11992e-314
     Coulomb-14        LJ (SR)        LJ (LR)   Coulomb (SR)   Coul. recip.
   2.11992e-314   -2.19860e+03   -1.52663e+02   -6.70307e+03   -2.59368e+04
      Potential Pressure (bar)
   -3.05754e+04            nan
           Step           Time         Lambda
              1        1.00000        0.00000

           Step           Time         Lambda
              2        2.00000        0.00000

           Step           Time         Lambda
              3        3.00000        0.00000

           Step           Time         Lambda
              4        4.00000        0.00000

           Step           Time         Lambda
              5        5.00000        0.00000

           Step           Time         Lambda
              6        6.00000        0.00000

           Step           Time         Lambda
              7        7.00000        0.00000

           Step           Time         Lambda
              8        8.00000        0.00000

           Step           Time         Lambda
              9        9.00000        0.00000

           Step           Time         Lambda
             10       10.00000        0.00000

..........(same thing until)

          Step           Time         Lambda
             37       37.00000        0.00000


Stepsize too small, or no change in energy.
Converged to machine precision,
but not to the requested precision Fmax < 1000

Steepest Descents converged to machine precision in 38 steps,
but did not reach the requested Fmax < 1000.
Potential Energy  = -3.05753603909389e+04
Maximum force     =                   inf on atom 1
Norm of force     =                   inf

Anna


Message: 2
Date: Fri, 21 May 2010 14:44:04 +0200
From: Luca Mollica <luca.mollica at ibs.fr>
Subject: Re: [gmx-users] stepsize too small ... but potential energy
	negative!
To: Discussion list for GROMACS users <gmx-users at gromacs.org>
Message-ID: <4BF68014.1060706 at ibs.fr>
Content-Type: text/plain; charset="iso-8859-1"

Dear Anna,
are you talking about a system "in vacuo" or in solvent ?
If you have placed the protein in water w/o minimizing it before 
solvation, probabily an "in vacuo" minimization could be useful for your 
system before moving into the solvated case. Moreover, are you sure that 
the protein does not have "broken" residue or something like that ? 
Sometimes, the completion of topology creation step goes fine but 
something wrong (on the gromacs side, apart from visualization) with 
sidechains or portion of the proteins that are generated by other 
softwares/servers.
I do not think it's a ligand problem, BTW.

Cheers

Luca

Message: 3
Date: Fri, 21 May 2010 14:46:28 +0200
From: Luca Mollica <luca.mollica at ibs.fr>
Subject: Re: [gmx-users] stepsize too small ... but potential energy
	negative!
To: gmx-users at gromacs.org
Message-ID: <4BF680A4.3060702 at ibs.fr>
Content-Type: text/plain; charset="iso-8859-1"

Moreover,
the problem is not the potential energy, but the force that must converge.
The "inf" there about force and about atom 1 tells me that the computed 
force has problems, indeed.
Which is atom 1 ? Is the protein Nter neutral or charged ? Or is atom 1 
an atom from the ligand ?

Cheers

L







More information about the gromacs.org_gmx-users mailing list