Dear Ross,
Thanks for your reply. I just thought I'd point out that a while ago I
tried obtaining the forces just by using the debugf options, but they
don't play nicely with QM/MM:
http://amber.ch.ic.ac.uk/archive/200702/0386.html
Steve.
Ross Walker wrote:
> 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.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.scripps.edu
> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
> to majordomo.scripps.edu
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber.scripps.edu
To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
to majordomo.scripps.edu
Received on Wed Jun 25 2008 - 06:07:49 PDT
