Skip to content

vasp calculator: small fix for xdat2traj to map forces correctly to atoms

Prateek Mehta requested to merge prtkm/ase:xdat2traj into master

This MR contains a small fix to ase.calculators.vasp.xdat2traj to ensure that the mapping of the atoms and the forces are consistently maintained.

I have illustrated the problem with a simple vasp calculation below (performed using the jasp wrapper - https://github.com/jkitchin/jasp).

#+BEGIN_SRC python
from jasp import *
from ase import Atoms, Atom
from ase.visualize import view
from ase.structure import molecule

from ase.io.trajectory import Trajectory


# make an unsorted water molecule
atoms = Atoms([Atom('H', [ 0., 0.76, -0.5]),
               Atom('O', [0., 0., 0.12]),
               Atom('H', [0., -0.76, -0.5])],
              cell = [10, 9, 11])
atoms.center()

JASPRC['queue.nprocs'] = 1

# Cheap calculation to relax atoms
with jasp('/scratch365/pmehta1/calculations/xdat-test4',
          encut = 200,
          xc = 'PBE',
          isif = 2,
          ibrion = 2,
          nsw=100,
          atoms=atoms) as calc:
    calc.calculate()
    
    print "Forces for relaxed configuration:"
    for sym, force in zip(atoms.get_chemical_symbols(), atoms.get_forces()):
        print sym, force

    # Read the trajectory
    xd = xdat2traj(calc=calc)
    xd.convert()
    traj = Trajectory('out.traj')
    
    print "Forces for final traj image:"
    image = traj[-1]
    for sym, force in zip(image.get_chemical_symbols(), image.get_forces()):
        print sym, force
#+END_SRC

#+RESULTS:

With the original code, I get the following output. Note that the atoms for the trajectory image are sorted, but the forces are not mapped correctly.

: Forces for relaxed configuration:
: H [ 0.    -0.025 -0.036]
: O [ 0.     0.     0.072]
: H [ 0.     0.025 -0.036]
: Forces for final traj image:
: H [ 0.    -0.025 -0.036]
: H [ 0.     0.     0.072]
: O [ 0.     0.025 -0.036]

With the modified code, I get this. In this case the atoms are not sorted, and the forces are consistent.

: Forces for relaxed configuration:
: H [ 0.    -0.025 -0.036]
: O [ 0.     0.     0.072]
: H [ 0.     0.025 -0.036]
: Forces for final traj image:
: H [ 0.    -0.025 -0.036]
: O [ 0.     0.     0.072]
: H [ 0.     0.025 -0.036]

Note: To get the sorted trajectory, with we need to simply run xdat2traj without the calculator argument, i.e., xd = xdat2traj().

Merge request reports