1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
|
import numpy as np
from ase import Atoms
from ase.units import Ha, Bohr
# from ase.io import read
from ase.calculators.openmx.reader import read_openmx
openmx_out_sample = """
System.CurrentDirectory ./
System.Name ch4
Atoms.SpeciesAndCoordinates.Unit Ang
<Definition.of.Atomic.Species
C C5.0-s1p1 C_PBE19
H H5.0-s1 H_PBE19
Definition.of.Atomic.Species>
Atoms.SpeciesAndCoordinates.Unit Ang
<Atoms.SpeciesAndCoordinates
1 C 0.0 0.0 0.1 2.0 2.0
2 H 0.682793 0.682793 0.682793 0.5 0.5
3 H -0.682793 -0.682793 0.68279 0.5 0.5
4 H -0.682793 0.682793 -0.682793 0.5 0.5
5 H 0.682793 -0.682793 -0.682793 0.5 0.5
Atoms.SpeciesAndCoordinates>
<Atoms.UnitVectors
10.0 0.0 0.0
0.0 10.0 0.0
0.0 0.0 10.0
Atoms.UnitVectors>
scf.EigenvalueSolver Band
...
Utot. -8.055096450113
...
Chemical potential (Hartree) -0.156250000000
...
*************************************************************************
*************************************************************************
Decomposed energies in Hartree unit
Utot = Utot(up) + Utot(dn)
= Ukin(up) + Ukin(dn) + Uv(up) + Uv(dn)
+ Ucon(up)+ Ucon(dn) + Ucore+UH0 + Uvdw
Uele = Ukin(up) + Ukin(dn) + Uv(up) + Uv(dn)
Ucon arizes from a constant potential added in the formalism
up: up spin state, dn: down spin state
*************************************************************************
*************************************************************************
Total energy (Hartree) = -8.055096425922011
Decomposed.energies.(Hartree).with.respect.to.atom
Utot
1 C -6.261242355014 ...
2 H -0.445956460556 ...
3 H -0.445956145906 ...
4 H -0.450970732231 ...
5 H -0.450970732215 ...
...
<coordinates.forces
5
1 C 0.00 0.00 0.10 0.00000 0.00000 -0.091659
2 H 0.68 0.68 0.68 0.02700 0.02700 0.029454
3 H -0.68 -0.68 0.68 -0.02700 -0.02700 0.029455
4 H -0.68 0.68 -0.68 0.00894 -0.00894 0.016362
5 H 0.68 -0.68 -0.68 -0.00894 0.00894 0.016362
coordinates.forces>
...
***********************************************************
***********************************************************
Fractional coordinates of the final structure
***********************************************************
***********************************************************
1 C 0.00000 0.00000 0.01000
2 H 0.06827 0.06827 0.06827
3 H 0.93172 0.93172 0.06827
4 H 0.93172 0.06827 0.93172
5 H 0.06827 0.93172 0.93172
"""
def test_openmx_out():
with open('openmx_fio_test.out', 'w') as f:
f.write(openmx_out_sample)
atoms = read_openmx('openmx_fio_test')
tol = 1e-2
# Expected values
energy = -8.0551
energies = np.array([-6.2612, -0.4459, -0.4459, -0.4509, -0.4509])
forces = np.array([[0.00000, 0.00000, -0.091659],
[0.02700, 0.02700, 0.029454],
[-0.02700, -0.02700, 0.029455],
[0.00894, -0.00894, 0.016362],
[-0.00894, 0.00894, 0.016362]])
assert isinstance(atoms, Atoms)
assert np.isclose(atoms.calc.results['energy'], energy * Ha, atol=tol)
assert np.all(np.isclose(atoms.calc.results['energies'],
energies * Ha, atol=tol))
assert np.all(np.isclose(atoms.calc.results['forces'],
forces * Ha / Bohr, atol=tol))
|
