Hello,
I’m using ASE for the diagonalization of the hessian matrix. My code is:
import os
from ase.units import Ry
from ase.vibrations import Vibrations
from ase.io.trajectory import Trajectory
from ase.calculators.siesta import Siesta
from ase.calculators.siesta.parameters import Species
from basis import basisSet
from arguments import singlePointArguments
traj = Trajectory('./PBEstructures.traj')
oxygenIndex = [atom.index for atom in traj[0] if atom.symbol == 'O'][0]
hydrogen1Index = [atom.index for atom in traj[0] if atom.symbol == 'H'][0]
hydrogen2Index = [atom.index for atom in traj[0] if atom.symbol == 'H'][1]
for j in range(len(traj)):
system = traj[j]
for k in range(0, oxygenIndex):
if k <= 31:
system[k].tag = 1
elif k >= 32:
system[k].tag = 2
siesta = Siesta(
label=f'structure-{j}',
mesh_cutoff=700*Ry,
kpts=[1, 1, 1],
xc="PBE",
species=[
Species(symbol="H", basis_set=basisSet("H"), pseudopotential="H.psf"),
Species(symbol="O", basis_set=basisSet("O"), pseudopotential="O.psf"),
Species(symbol='Au', basis_set=basisSet('Au2'), pseudopotential='Au.gga.psf', tag=1),
Species(symbol='Au', basis_set=basisSet('Au3'), pseudopotential='Au.gga.psf', tag=1)
],
fdf_arguments=singlePointArguments(),
)
system.calc = siesta
if not os.path.exists(f'structure-{j}'):
os.makedirs(f'structure-{j}')
os.chdir(f'image-{j}')
vib = Vibrations(
atoms=system,
indices=[ oxygenIndex, hydrogen1Index, hydrogen2Index ],
name=f'structure-{j}',
delta=0.02,
nfree=2
)
vib.run()
print(vib.get_vibrations().get_energies_and_modes())
vib.summary()
os.chdir('../')
system = None
siesta = None
PBEStructures.traj is a file containing a box with 64 gold atoms and a water a molecular on top of it.
My goal here is to understand the units of measure (mainly of the eigenvectors) used when I print vib.get_vibrations().get_energies_and_modes(). The output I receive is:
(array([0. +0.02172963j, 0. +0.00898803j,
0.0047047 +0.j , 0.00515926+0.j ,
0.01071083+0.j , 0.01819307+0.j ,
0.20329078+0.j , 0.46236624+0.j ,
0.4736801 +0.j ]), array([[[-4.42204141e-02, 8.21599997e-04, -8.20717964e-04],
[ 1.88182925e-01, 4.34003246e-02, 6.77329557e-01],
[ 1.87967471e-01, -2.82722309e-02, -6.54830052e-01]],
[[-1.30904884e-03, -1.17128050e-01, 8.54004482e-02],
[-1.48444828e-02, 1.62819355e-01, 5.50239486e-01],
[-2.37987090e-02, 1.16354127e-01, 5.61187346e-01]],
[[ 2.30667250e-01, 1.59479522e-02, 1.83724813e-02],
[ 2.57152279e-01, 2.91414458e-03, 6.41077081e-02],
[ 2.57173126e-01, 2.47372183e-02, -3.41282499e-02]],
[[-2.12276397e-02, 2.00084623e-01, 1.18940029e-01],
[-3.12744446e-02, 2.17942020e-01, 1.27328885e-01],
[-3.23472198e-02, 1.96342333e-01, 1.43769534e-01]],
[[ 9.45816964e-03, 5.59756452e-02, -1.97466996e-01],
[-8.46615478e-03, 3.37616114e-01, 2.52215042e-01],
[-1.63751228e-02, 2.77428558e-01, 2.59051203e-01]],
[[ 4.77412315e-02, -1.63801442e-03, 4.83822863e-03],
[-3.29851589e-01, 5.74999313e-01, 4.96722335e-02],
[-3.28998731e-01, -6.24454755e-01, -1.24080348e-01]],
[[ 1.21337881e-05, 5.97117696e-02, -3.33461626e-02],
[ 4.09683326e-01, -4.74611830e-01, 2.56556369e-01],
[-4.09675251e-01, -4.74696304e-01, 2.56387055e-01]],
[[-2.03379881e-04, -4.23625906e-02, 2.44361192e-02],
[ 5.74482847e-01, 3.37161265e-01, -1.90362467e-01],
[-5.71214032e-01, 3.34886317e-01, -1.89000171e-01]],
[[-6.72311778e-02, 1.32489072e-04, -6.62638940e-05],
[ 5.31582749e-01, 3.64169723e-01, -2.05093709e-01],
[ 5.35119980e-01, -3.66289334e-01, 2.06244020e-01]]]))
---------------------
# meV cm^-1
---------------------
0 21.7i 175.3i
1 9.0i 72.5i
2 4.7 37.9
3 5.2 41.6
4 10.7 86.4
5 18.2 146.7
6 203.3 1639.7
7 462.4 3729.2
8 473.7 3820.5
---------------------
Zero-point energy: 0.589 eV
Thank you in advance,
Regards,
João Rheinheimer