Hi Dave,
Yes I know exactly what is wrong here and I have been 'suggesting' for a
long time how we might go about fixing it but it has never been put in.
Essentially when you run a periodic PME calculation with NDDO based QM/MM
inside the cutoff you compute the full QM/MM interaction. In the 'PME' part
of the calculation though you use a mulliken charge approximation. Therefore
there is a discontinuity in the electrostatics as you cross the QM/MM
boundary. You can check the QM/MM is okay by running a solvent cap with no
cutoff. Things get better for QM/MM when you run a larger system (MM atom
count), use qmcut=12.0 and also if you have a buried QM/MM region inside a
protein where there is less effective crossing of the QM/MM boundary.
The reason DFTB is fine is that it uses a mulliken charge interaction
(q[i]q[j]/r[ij]) for both the direct and reciprocal space calculations. One
possible solution is to add a force switching function to the direct space
sum in QM/MM which would then conserve energy but for the wrong reasons. The
other that I have been encouraging 'certain' people to implement for ages is
to add an option (supposed to be qmmm_int=3) to use a mulliken charge
interaction in the direct space for the QM/MM interaction with NDDO. I
believe this is what Jorgensen has been using for years. The potential would
then be smooth and energy would be conserved. It would be interesting to
compare the difference between the current approach and a Mulliken charge
QM/MM interaction. My bet is the difference would not be huge and the
improvement in energy conservation would make up for the more approximate
QM/MM interaction.
If Insuk would like to try implementing this I would be happy to help out.
All the best
Ross
> -----Original Message-----
> From: David A Case [mailto:case.biomaps.rutgers.edu]
> Sent: Tuesday, April 05, 2011 7:22 AM
> To: amber-developers.ambermd.org
> Subject: [AMBER-Developers] problems with QM/MM energy conservation
>
> InSuk and I are having problems seeing energy conservation in simple
> QM/MM calculations with the NDDO Hamiltonians. Out simplest example is
> NMA in water, and here is a typical input file:
>
>
> production
> &cntrl
> imin=0, ntx=7, ntpr=100,
> ntf=2, ntc=2, tol=1.d-7,
> ntb=1, ntp=0,
> ntt=0,
> nstlim=50000, dt=0.001, ntave=25000,
> cut=9.0,
> irest=1,
> iwrap=1,
> ifqnt=1,
> /
> &ewald
> dsum_tol=1.d-6,
> /
> &qmmm
> qmmask=':1', qmcharge=0, qm_theory='AM1',
> scfconv=1.0d-11, tight_p_conv=1, peptide_corr=0,
> /
>
> This shows severe lack of energy conservation (change in total energy
> of
> several dozen kcal/mol in just 50000 steps. If I set ifqnt=0, or
> qm_theory="DFTB", then I get essentially no drift in total energy (i.e.
> less
> than 0.1 kcal/mol over 50000 steps. I've attached the prmtop and
> inpcrd
> files.
>
> Does anyone see what is going wrong here? Our executables pass all the
> tests,
> and the fact that DFTB seems to work but PM3 and AM1 fail is
> intriguing.
>
> ...thanks for any help...dac
_______________________________________________
AMBER-Developers mailing list
AMBER-Developers.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber-developers
Received on Tue Apr 05 2011 - 08:30:03 PDT
