Re: [AMBER] QM/MM using Amber and Gaussian: Discrepancies between the QM energy from the gaussian output file and the energy reported in the mdout file as EXTERNESCF

From: Goetz, Andreas <agoetz.sdsc.edu>
Date: Wed, 5 Dec 2018 19:38:02 +0000

Hi Antoine,

Yes, you are right. Gaussian of course does not know the topology, thus does not exclude 1-2 and 1-3 pairs or scale charges of 1-4 pairs, which means that even if all MM charges are included as external electrostatic potential in the QM calculation, the self-energy that Gaussian computes will differ from EEL. Thank you for pointing this out.

All the best,
Andy


Dr. Andreas W. Goetz
Assistant Research Scientist
San Diego Supercomputer Center
Tel: +1-858-822-4771
Email: agoetz.sdsc.edu
Web: www.awgoetz.de

> On Dec 5, 2018, at 2:16 AM, Antoine MARION <amarion.metu.edu.tr> wrote:
>
> Hi Andy,
>
> Just out of curiosity: shouldn’t the self-energy and EEL differ even if conditions a) and b) hold due to the exclude list (bonds and angles) and the scaling of 1-4 interactions?
>
> Best,
> Antoine
>
>
>> On 4 Dec 2018, at 22:48, Goetz, Andreas <agoetz.sdsc.edu> wrote:
>>
>> Hi Eric,
>>
>>> Does the self-energy of the charges only corresponds to the energy for the interaction *between* the point charges or does it also include the interaction between point charges and the atom of the QM region?
>>
>> Yes, this is the electrostatic energy between the point charges. See also the comments in the source code.
>>
>>> So then the energy between MM point charges is included in the EEL energy, is that correct?
>>
>> Yes, this is part of the force field, calculated by Amber, and included in EEL.
>>
>>> Are they calculated in the same manner by Amber than they are by Gaussian? The values in kcal/mol are very different so I wouldn't think so.
>>
>> EEL and the self-energy are identical if two conditions hold: a) no periodic boundary conditions (for instance a molecule in gas phase or in a droplet) and b) a QM/MM cutoff that is sufficiently large to include all point charges. Otherwise there is a truncation of electrostatic interactions between QM and MM region: Point charges beyond the cutoff are not passed to the QM code and long-range electrostatic interactions between QM and MM region are not accounted. Of course Amber handles long-range interactions between MM point charges correctly using PME.
>>
>>> Also, I had a look at other qm2_extern modules and it looks like for example for Orca this subtraction is not performed. Is it because Orca remove this form the total energy on its own? Or perhaps I missed it in the code.
>>
>> None of the other supported QM codes calculates the self energy between the external point charges. It seems somewhat pointless (no pun intended) but Gaussian does it so we have to subtract it out from the total energy that Gaussian reports.
>>
>> All the best,
>> Andy
>>
>> —
>> Dr. Andreas W. Goetz
>> Assistant Research Scientist
>> San Diego Supercomputer Center
>> Tel: +1-858-822-4771
>> Email: agoetz.sdsc.edu
>> Web: www.awgoetz.de
>>
>>> On Dec 4, 2018, at 8:10 AM, Eric Lang <eric.lang.bristol.ac.uk> wrote:
>>>
>>> Dear Antoine,
>>>
>>>
>>> Thank you very much for your reply.
>>>
>>> The results make sense now. I had a look at qm2_extern_gau_module.F90 and I can see where it comes from and I can recalculate the value in the mdout file from the gaussian file.
>>>
>>>
>>> However, I wouldn't mind some clarifications on why this is done.
>>>
>>> Does the self-energy of the charges only corresponds to the energy for the interaction *between* the point charges or does it also include the interaction between point charges and the atom of the QM region?
>>>
>>> I guess that it is only between charges otherwise what would be the point of doing electrostatic embedding.
>>>
>>> So then the energy between MM point charges is included in the EEL energy, is that correct?
>>>
>>> Are they calculated in the same manner by Amber than they are by Gaussian? The values in kcal/mol are very different so I wouldn't think so.
>>>
>>> Then I guess if this self energy is only between the point charges, it correspond purely to the MM part so it make sense to have it calculated from the MM program based on the force field
>>>
>>> Also, I had a look at other qm2_extern modules and it looks like for example for Orca this subtraction is not performed. Is it because Orca remove this form the total energy on its own? Or perhaps I missed it in the code.
>>>
>>>
>>> Many thanks in advance for your help,
>>>
>>>
>>> Eric
>>>
>>>
>>>
>>> --
>>>
>>> Eric Lang
>>>
>>>
>>> Centre for Computational Chemistry
>>> School of Chemistry - University of Bristol
>>> Bristol BS8 1TS - United Kingdom
>>>
>>>
>>> ________________________________
>>> From: Antoine MARION <amarion.metu.edu.tr>
>>> Sent: 04 December 2018 14:22:56
>>> To: AMBER Mailing List
>>> Subject: Re: [AMBER] QM/MM using Amber and Gaussian: Discrepancies between the QM energy from the gaussian output file and the energy reported in the mdout file as EXTERNESCF
>>>
>>> Dear Eric,
>>>
>>> As far as I see from the code, you should also subtract the “Self energy of the charges” (as found in the Gaussian output) to match EXTERNESCF.
>>> Probably to avoid double counting. (?)
>>>
>>> EXTERNESCF = E(RwB97XD) - self_energy
>>>
>>> Best,
>>> Antoine
>>>
>>> ____________________________
>>> -- Dr. Antoine MARION
>>> Assistant Professor
>>> Doktor Öğretim Üyesi
>>>
>>> Department of Chemistry
>>> Middle East Technical University
>>> 06800, Ankara, Turkey
>>>
>>> Office: O-310
>>> Email: amarion.metu.edu.tr <mailto:amarion.metu.edu.tr>
>>> Phone: +90 312 210 51 44
>>>
>>>> On 4 Dec 2018, at 14:36, Eric Lang <eric.lang.bristol.ac.uk> wrote:
>>>>
>>>> Dear Amber Developers and Users,
>>>>
>>>>
>>>>
>>>> I am running some single point QM/MM in Sander (compiled as Amber16+AmberTools17) using Gaussian 09.
>>>>
>>>> The calculations run ok, but I think something is not quite right with the extraction of the QM energy from the gaussian output file: the EXTERNESCF energy in the mdout file in kcal/mol doesn’t match the final energy (last ‘SCF Done’ line or ‘HF’ at the very end of the Gaussian file) when converted to kcal/mol using the conversion used in Amber (AU_TO_KCAL = 6.2750946943E+02).
>>>>
>>>> I would think they should be the same.
>>>>
>>>> Am I missing something obvious in the way the calculations are done or is there a problem with the way the gaussian output files are parsed?
>>>>
>>>>
>>>>
>>>> Here is the output form Amber:
>>>>
>>>>
>>>>
>>>> FINAL RESULTS
>>>>
>>>> NSTEP ENERGY RMS GMAX NAME NUMBER
>>>>
>>>> 1 -5.0519E+05 1.3953E+01 1.1652E+02 CA 1579
>>>>
>>>>
>>>>
>>>> BOND = 456.4977 ANGLE = 1327.2559 DIHED = 2177.6750
>>>>
>>>> VDWAALS = 8267.2236 EEL = -89612.1756 HBOND = 0.0000
>>>>
>>>> 1-4 VDW = 625.9136 1-4 EEL = 7948.6914 RESTRAINT = 0.0000
>>>>
>>>> EXTERNESCF = -436379.7510
>>>>
>>>>
>>>>
>>>> And here is the final energy from Gaussian
>>>>
>>>>
>>>>
>>>> SCF Done: E(RwB97XD) = -3178.95512525 A.U.
>>>>
>>>>
>>>>
>>>> Which gives -1994824.4 kcal/mol and not -436379.7510 kcal/mol as found in EXTERNESCF.
>>>>
>>>>
>>>>
>>>> Many thanks in advance for your help,
>>>>
>>>>
>>>>
>>>> Eric
>>>>
>>>>
>>>>
>>>> --
>>>>
>>>> Eric Lang
>>>>
>>>> Centre for Computational Chemistry
>>>> School of Chemistry - University of Bristol
>>>> Bristol BS8 1TS - United Kingdom
>>>>
>>>>
>>>>
>>>> _______________________________________________
>>>> AMBER mailing list
>>>> AMBER.ambermd.org
>>>> http://lists.ambermd.org/mailman/listinfo/amber
>>>
>>> _______________________________________________
>>> AMBER mailing list
>>> AMBER.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber
>>> _______________________________________________
>>> AMBER mailing list
>>> AMBER.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber
>>
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
>
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber

_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Dec 05 2018 - 12:00:03 PST
Custom Search