AMBER Archive (2008)

Subject: RE: AMBER: Using AMBER forces

From: Ross Walker (ross_at_rosswalker.co.uk)
Date: Tue Jun 24 2008 - 11:22:25 CDT

Hi Steven,

> The only modification I have made to sander is the following, placed in
> the force subroutine in force.f, between "call trace_exit( 'force' )"
> and "return":

> so sander never gets further than the first call to force(). Is there
> any reason why this shouldn't give the expected results, even when using
> QM/MM? In particular, I see in the runmd() subroutine that the first MD
> step (where this first force() call occurs) is treated differently to
> the others.

As far as I can tell this should work fine even with QM/MM. The reason the
first QM call is slightly different from subsequent calls stems from a
number of reasons. Firstly there is a whole set of memory allocation that is
done within the qm_mm routine, it is not necessary to deallocate and
reallocate this on every call and so it is just done the once on the first
call. The second reason is that for speed the QM/MM SCF reuses the previous
density matrix on every call. Of course on the first call there is no
previous matrix so it has to create a guess. There are also some other
issues with loading in the semi-empirical parameters etc. None of this
should effect you though so what you are doing should be fine.

I would suggest checking it manually though. You could run sander with the
debugf option turned on which will compare analytical and numerical forces.
This will print you out both the analytical and numerical forces for the
initial structure and you can compare this with what is in your force file.
You may need to do some kind of unit conversion but perhaps not.

E.g. for a simple gas phase system

&cntrl
  imin=0, nstlim=1,
  ntf=2, ntc=2, dt=0.002,
  ntb=0, tempi=0.0,
  temp0=0.0, ntt=3, gamma_ln=1.0,
  cut=999.0,
  ifqnt=1,
 /
 &qmmm
  qmmask='@1-10',
  qmtheory=1,
  qmcharge=0,
  tight_p_conv=1,
  scfconv=1.0d-10,
 /
 &debugf
  do_debugf=1,
  atomn=1,2,3,4,5,6,7,8,9,10,11,12,13,14,
 /

This will print out the forces for atoms 1 to 14 for you for comparisson.

Note, if you are specifically interested in the forces from QM/MM then it is
probably advisable to tighten the scf convergence criteria as in the example
above (tight_p_conv=1,scfconv=1.0d-10).

Good luck,
Ross

/\
\/
|\oss Walker

| Assistant Research Professor |
| San Diego Supercomputer Center |
| Tel: +1 858 822 0854 | EMail:- ross_at_rosswalker.co.uk |
| http://www.rosswalker.co.uk | PGP Key available on request |

Note: Electronic Mail is not secure, has no guarantee of delivery, may not
be read every day, and should not be used for urgent or sensitive issues.

  

-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
      to majordomo_at_scripps.edu