[ase-users] emt giving nonzero interaction energies at large separation

[Jens Jørgen Mortensen]

man, 12 09 2011 kl. 21:23 +0200, skrev Tom Haxton:
> Hello,
> I want to use ase to calculate the interaction energy between a metal
> slab and a substituted benzene molecule as a function of relative
> position.  I expect that the interaction energy should go to zero as
> the separation between the molecules increases.  However, using the
> emt calculator, this is not the case.  At large separations, the sum
> of the energies of the isolated slab and the isolated molecule is not
> equal to the energy of the combined system.  I use the same cell size,
> shape, and periodic boundary conditions in all three calculations.
> For large separations and large cell sizes, the result is independent
> of separation, cell size, and boundary conditions.
> 
> I am using ase-3.4.  I also tried ase-3.5 and got the same problem,
> though a different number for the energy discrepancy.  I am working on
> a x86_64 machine running Linux 2.6, python 2.6.5, and numpy 1.3.0.
> 
> Below is the script that produces the error.  It is the more minimal
> case of a single isolated Au atom interacting with a Benzene molecule,
> separated by 50 angstroms, in a 100 x 100 x 100 angstrom cell.  The
> energies are:
> 
> Benzene    -3.47792756546
> Au                3.8
> Combined    -1.46565762747
> Difference   -1.78773006201
> 
> Can anyone tell me if there is something I am doing wrong or if this
> is a bug?

I get only a 2 meV error, but still - I would call it a bug.  I think I
have fixed it in the development version of ASE.  Thanks for finding
this.

See here for a simple way to build a benzene molecule:

https://trac.fysik.dtu.dk/projects/ase/browser/trunk/ase/test/emt2.py?rev=2295

Please keep in mind that the parameters for H and C in ASE's EMT
potential should not be used for real work.  They are there only for
testing and demonstration purposes.

Jens Jørgen

>   Also, is there any way to get out any useful intermediate
> information from the emt calculation?  Does ASAP have the same
> problem?  (I don't have ASAP installed yet.)  Thank you very much.
> 
> Best,
> Tom
> 
> 
> 
> #/usr/bin/python
> 
> import sys
> import math
> import copy
> from ase import Atoms
> from ase.visualize import view
> from ase.calculators.emt import EMT
> from ase.constraints import FixAtoms
> from ase.optimize import QuasiNewton
> from ase.lattice.surface import fcc111,add_adsorbate
> 
> mycalculator = EMT()
> mybigcell = (100., 100., 100.)
> separation = 50.
> 
> r_CC = 1.374
> r_HC = 1.084
> 
> benzene = Atoms('CCCCCCHHHHHH', positions=[(r_CC, 0., 0.), (-r_CC, 0.,
> 0), (0.5 * r_CC, 0.5 * math.sqrt(3.0) * r_CC, 0.), (0.5 * r_CC, -0.5 *
> math.sqrt(3.0) * r_CC, 0.), (-0.5 * r_CC, 0.5 * math.sqrt(3.0) * r_CC,
> 0.), (-0.5 * r_CC, -0.5 * math.sqrt(3.0) * r_CC, 0.), (r_CC + r_HC,
> 0., 0.), (-(r_CC + r_HC), 0., 0.), (0.5 * (r_CC + r_HC), 0.5 *
> math.sqrt(3.0) * (r_CC + r_HC), 0), (0.5 * (r_CC + r_HC), -0.5 *
> math.sqrt(3.0) * (r_CC + r_HC), 0), (-0.5 * (r_CC + r_HC), 0.5 *
> math.sqrt(3.0) * (r_CC + r_HC), 0), (-0.5 * (r_CC + r_HC), -0.5 *
> math.sqrt(3.0) * (r_CC + r_HC), 0)])
> 
> gold = Atoms('Au', positions=[(0., 0., 0.)])
> 
> molecule1 = copy.deepcopy(benzene)
> molecule2 = copy.deepcopy(gold)
> 
> molecule1.set_cell(mybigcell)
> molecule1.set_pbc(False)
> molecule1.set_calculator(mycalculator)
> molecule2.set_cell(mybigcell)
> molecule2.set_pbc(False)
> molecule2.set_calculator(mycalculator)
> molecule2.translate((0., 0., separation))
> moleculecombined = molecule1 + molecule2
> molecule1.center()
> molecule2.center()
> moleculecombined.center()
> view(molecule1)
> view(molecule2)
> 
> e1 = molecule1.get_potential_energy()
> e2 = molecule2.get_potential_energy()
> 
> moleculecombined.set_cell(mybigcell)
> moleculecombined.set_pbc(False)
> moleculecombined.set_calculator(mycalculator)
> 
> ecombined = moleculecombined.get_potential_energy()
> 
> print ecombined - e1 - e2, e1, e2, ecombined
> 
> #view(moleculecombined)
> 
> 
>