Skip to content

Commit

Permalink
Make improvements to the format for OpenMX (#619)
Browse files Browse the repository at this point in the history
Dear developers,

I would like to make some improvements on the format for OpenMX.

1.  Change the format name from `openmx/out` to `openmx/md`
2. Add the if statement to read geometry optimization trajectories
3. Add the if statement to check if the SCF calculation has converged

Thank you in advance.
  • Loading branch information
shigeandtomo authored Mar 20, 2024
1 parent 0b8014f commit b048b37
Show file tree
Hide file tree
Showing 7 changed files with 1,591 additions and 1,407 deletions.
25 changes: 22 additions & 3 deletions dpdata/openmx/omx.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
energy_convert = EnergyConversion("hartree", "eV").value()
force_convert = ForceConversion("hartree/bohr", "eV/angstrom").value()

import warnings
from collections import OrderedDict

### iterout.c from OpenMX soure code: column numbers and physical quantities ###
Expand All @@ -28,6 +29,12 @@
# /* 12: magnetic moment (muB) */
# /* 13,14: angles of spin */

# 15: scf_convergence_flag (optional)
#
# 1. Move the declaration of `scf_convergence_flag` in `DFT.c` to `openmx_common.h`.
# 2. Add `scf_convergence_flag` output to the end of `iterout.c` where `*.md` is written.
# 3. Recompile OpenMX.


def load_atom(lines):
atom_names = []
Expand Down Expand Up @@ -70,9 +77,18 @@ def load_cells(lines):
for index, line in enumerate(lines):
if "Cell_Vectors=" in line:
parts = line.split()
cell.append([float(parts[12]), float(parts[13]), float(parts[14])])
cell.append([float(parts[15]), float(parts[16]), float(parts[17])])
cell.append([float(parts[18]), float(parts[19]), float(parts[20])])
if len(parts) == 21: # MD.Type is NVT_NH
cell.append([float(parts[12]), float(parts[13]), float(parts[14])])
cell.append([float(parts[15]), float(parts[16]), float(parts[17])])
cell.append([float(parts[18]), float(parts[19]), float(parts[20])])
elif len(parts) == 16: # MD.Type is Opt
cell.append([float(parts[7]), float(parts[8]), float(parts[9])])
cell.append([float(parts[10]), float(parts[11]), float(parts[12])])
cell.append([float(parts[13]), float(parts[14]), float(parts[15])])
else:
raise RuntimeError(
"Does the file System.Name.md contain unsupported calculation results?"
)
cells.append(cell)
cell = []
cells = np.array(cells)
Expand Down Expand Up @@ -104,6 +120,9 @@ def load_coords(lines, atom_names, natoms):
parts = line.split()
for_line = [float(parts[1]), float(parts[2]), float(parts[3])]
coord.append(for_line)
# It may be necessary to recompile OpenMX to make scf convergence determination.
if len(parts) == 15 and parts[14] == "0":
warnings.warn("SCF in System.Name.md has not converged!")
if cnt == natoms:
coords.append(coord)
cnt = 0
Expand Down
4 changes: 2 additions & 2 deletions dpdata/plugins/openmx.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,15 +3,15 @@
from dpdata.format import Format


@Format.register("openmx/out")
@Format.register("openmx/md")
class OPENMXFormat(Format):
"""Format for the `OpenMX <https://www.openmx-square.org/>`.
OpenMX (Open source package for Material eXplorer) is a nano-scale material simulation package based on DFT, norm-conserving pseudopotentials, and pseudo-atomic localized basis functions.
Note that two output files, System.Name.dat and System.Name.md, are required.
Use the `openmx/out` keyword argument to supply this format.
Use the `openmx/md` keyword argument to supply this format.
"""

@Format.post("rot_lower_triangular")
Expand Down
2,800 changes: 1,400 additions & 1,400 deletions tests/openmx/Methane.md

Large diffs are not rendered by default.

68 changes: 68 additions & 0 deletions tests/openmx/Methane2.dat
Original file line number Diff line number Diff line change
@@ -0,0 +1,68 @@
#
# File Name
#

System.CurrrentDirectory ./ # default=./
System.Name Methane2
level.of.stdout 1 # default=1 (1-3)
level.of.fileout 1 # default=1 (0-2)

#
# Definition of Atomic Species
#

Species.Number 2
<Definition.of.Atomic.Species
H H5.0-s2 H_PBE19
C C5.0-s2p2 C_PBE19
Definition.of.Atomic.Species>

#
# Atoms
#

Atoms.Number 5
Atoms.SpeciesAndCoordinates.Unit Ang # Ang|AU
<Atoms.SpeciesAndCoordinates
1 C 0.300000 0.000000 0.000000 2.0 2.0
2 H -0.889981 -0.629312 0.000000 0.5 0.5
3 H 0.000000 0.629312 -0.889981 0.5 0.5
4 H 0.000000 0.629312 0.889981 0.5 0.5
5 H 0.889981 -0.629312 0.000000 0.5 0.5
Atoms.SpeciesAndCoordinates>
Atoms.UnitVectors.Unit Ang # Ang|AU
#<Atoms.UnitVectors
# 10.0 0.0 0.0
# 0.0 10.0 0.0
# 0.0 0.0 10.0
#Atoms.UnitVectors>

#
# SCF or Electronic System
#

scf.XcType GGA-PBE # LDA|LSDA-CA|LSDA-PW|GGA-PBE
scf.SpinPolarization off # On|Off|NC
scf.ElectronicTemperature 100.0 # default=300 (K)
scf.energycutoff 200.0 # default=150 (Ry)
scf.maxIter 1 # default=40
scf.EigenvalueSolver cluster # DC|GDC|Cluster|Band
scf.Kgrid 1 1 1 # means n1 x n2 x n3
scf.Mixing.Type rmm-diis # Simple|Rmm-Diis|Gr-Pulay|Kerker|Rmm-Diisk
scf.Init.Mixing.Weight 0.30 # default=0.30
scf.Min.Mixing.Weight 0.001 # default=0.001
scf.Max.Mixing.Weight 0.400 # default=0.40
scf.Mixing.History 15 # default=5
scf.Mixing.StartPulay 4 # default=6
scf.criterion 1.0e-8 # default=1.0e-6 (Hartree)

#
# MD or Geometry Optimization
#

MD.Type opt # Nomd|Opt|DIIS|NVE|NVT_VS|NVT_NH
MD.Opt.DIIS.History 7 # default=7
MD.Opt.StartDIIS 5 # default=5
MD.maxIter 5 # default=1
MD.TimeStep 1.0 # default=0.5 (fs)
MD.Opt.criterion 1.0e-4 # default=1.0e-4 (Hartree/bohr)
35 changes: 35 additions & 0 deletions tests/openmx/Methane2.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,35 @@
5
time= 0.000 (fs) Energy= -8.15440 (Hartree) Cell_Vectors= 10.00000 0.00000 0.00000 0.00000 10.00000 0.00000 0.00000 0.00000 10.00000
C 0.30000 0.00000 0.00000 -0.36382 0.22843 -0.00000 0.00000 0.00000 0.00000 0.17513 0.00000 0.00000 0.00000 0
H -0.88998 -0.62931 0.00000 0.04918 0.01544 0.00000 0.00000 0.00000 0.00000 -0.07837 0.00000 0.00000 0.00000 0
H 0.00000 0.62931 -0.88998 0.02120 -0.00206 -0.00338 0.00000 0.00000 0.00000 -0.01532 0.00000 0.00000 0.00000 0
H 0.00000 0.62931 0.88998 0.02120 -0.00206 0.00338 0.00000 0.00000 0.00000 -0.01532 0.00000 0.00000 0.00000 0
H 0.88998 -0.62931 0.00000 0.23151 -0.22432 -0.00000 0.00000 0.00000 0.00000 -0.06612 0.00000 0.00000 0.00000 0
5
time= 1.000 (fs) Energy= -8.22382 (Hartree) Cell_Vectors= 10.00000 0.00000 0.00000 0.00000 10.00000 0.00000 0.00000 0.00000 10.00000
C 0.21037 0.05628 -0.00000 -0.12208 -0.03279 0.00000 0.00000 0.00000 0.00000 0.18366 0.00000 0.00000 0.00000 0
H -0.87786 -0.62551 0.00000 0.03791 0.01857 0.00000 0.00000 0.00000 0.00000 -0.08412 0.00000 0.00000 0.00000 0
H 0.00522 0.62881 -0.89081 0.00960 0.01866 -0.02578 0.00000 0.00000 0.00000 -0.03522 0.00000 0.00000 0.00000 0
H 0.00522 0.62881 0.89081 0.00960 0.01866 0.02578 0.00000 0.00000 0.00000 -0.03522 0.00000 0.00000 0.00000 0
H 0.94702 -0.68458 -0.00000 0.04642 -0.02939 -0.00000 0.00000 0.00000 0.00000 -0.02910 0.00000 0.00000 0.00000 0
5
time= 2.000 (fs) Energy= -8.24265 (Hartree) Cell_Vectors= 10.00000 0.00000 0.00000 0.00000 10.00000 0.00000 0.00000 0.00000 10.00000
C 0.12898 0.03442 0.00000 -0.04237 -0.04096 0.00000 0.00000 0.00000 0.00000 0.18611 0.00000 0.00000 0.00000 0
H -0.85259 -0.61313 0.00000 0.01505 0.00643 0.00000 0.00000 0.00000 0.00000 -0.06398 0.00000 0.00000 0.00000 0
H 0.01162 0.64125 -0.90800 0.00781 0.01155 -0.01439 0.00000 0.00000 0.00000 -0.04343 0.00000 0.00000 0.00000 0
H 0.01162 0.64125 0.90800 0.00781 0.01155 0.01439 0.00000 0.00000 0.00000 -0.04343 0.00000 0.00000 0.00000 0
H 0.97796 -0.70417 -0.00000 0.00515 0.00555 -0.00000 0.00000 0.00000 0.00000 -0.03526 0.00000 0.00000 0.00000 0
5
time= 3.000 (fs) Energy= -8.24626 (Hartree) Cell_Vectors= 10.00000 0.00000 0.00000 0.00000 10.00000 0.00000 0.00000 0.00000 10.00000
C 0.10073 0.00711 0.00000 -0.01636 -0.00597 0.00000 0.00000 0.00000 0.00000 0.18796 0.00000 0.00000 0.00000 0
H -0.84256 -0.60884 0.00000 -0.00199 -0.00521 -0.00000 0.00000 0.00000 0.00000 -0.05693 0.00000 0.00000 0.00000 0
H 0.01683 0.64895 -0.91760 0.00694 0.00262 -0.00333 0.00000 0.00000 0.00000 -0.04635 0.00000 0.00000 0.00000 0
H 0.01683 0.64895 0.91760 0.00694 0.00262 0.00333 0.00000 0.00000 0.00000 -0.04635 0.00000 0.00000 0.00000 0
H 0.98139 -0.70048 -0.00000 0.00184 0.00495 -0.00000 0.00000 0.00000 0.00000 -0.03833 0.00000 0.00000 0.00000 0
5
time= 4.000 (fs) Energy= -8.24682 (Hartree) Cell_Vectors= 10.00000 0.00000 0.00000 0.00000 10.00000 0.00000 0.00000 0.00000 10.00000
C 0.08982 0.00313 0.00000 -0.00704 -0.00134 0.00000 0.00000 0.00000 0.00000 0.19038 0.00000 0.00000 0.00000 0
H -0.84389 -0.61232 0.00000 -0.00554 -0.00705 0.00000 0.00000 0.00000 0.00000 -0.05558 0.00000 0.00000 0.00000 0
H 0.02145 0.65069 -0.91982 0.00591 0.00121 -0.00163 0.00000 0.00000 0.00000 -0.04727 0.00000 0.00000 0.00000 0
H 0.02145 0.65069 0.91982 0.00591 0.00121 0.00163 0.00000 0.00000 0.00000 -0.04727 0.00000 0.00000 0.00000 0
H 0.98262 -0.69717 -0.00000 -0.00051 0.00564 -0.00000 0.00000 0.00000 0.00000 -0.04026 0.00000 0.00000 0.00000 0
4 changes: 2 additions & 2 deletions tests/test_openmx.py
Original file line number Diff line number Diff line change
Expand Up @@ -50,12 +50,12 @@ def test_coord(self):

class TestOPENMXTraj(unittest.TestCase, TestOPENMXTRAJProps):
def setUp(self):
self.system = dpdata.System("openmx/Methane", fmt="openmx/out")
self.system = dpdata.System("openmx/Methane", fmt="openmx/md")


class TestOPENMXLabeledTraj(unittest.TestCase, TestOPENMXTRAJProps):
def setUp(self):
self.system = dpdata.LabeledSystem("openmx/Methane", fmt="openmx/out")
self.system = dpdata.LabeledSystem("openmx/Methane", fmt="openmx/md")


if __name__ == "__main__":
Expand Down
62 changes: 62 additions & 0 deletions tests/test_openmx_check_convergence.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,62 @@
import unittest

import numpy as np
from context import dpdata


class TestOPENMXTRAJProps:
def test_atom_names(self):
self.assertEqual(self.system.data["atom_names"], ["C", "H"])

def test_atom_numbs(self):
self.assertEqual(self.system.data["atom_numbs"], [1, 4])

def test_atom_types(self):
for ii in range(0, 1):
self.assertEqual(self.system.data["atom_types"][ii], 0)
for ii in range(1, 5):
self.assertEqual(self.system.data["atom_types"][ii], 1)

def test_cell(self):
ref = 10.0 * np.eye(3)
self.assertEqual(self.system.get_nframes(), 5)
for ff in range(self.system.get_nframes()):
for ii in range(3):
for jj in range(3):
self.assertEqual(self.system["cells"][ff][ii][jj], ref[ii][jj])

def test_coord(self):
with open("openmx/Methane2.md") as md_file:
lines = md_file.readlines()
lines = lines[-5:]
coords = []
for line in lines:
parts = line.split()
for_line = [float(parts[1]), float(parts[2]), float(parts[3])]
coords.append(for_line)
coords = np.array(coords)
celll = 10.0
## Applying PBC ##
for ii in range(5):
for jj in range(3):
if coords[ii][jj] < 0:
coords[ii][jj] += celll
elif coords[ii][jj] >= celll:
coords[ii][jj] -= celll
self.assertAlmostEqual(
self.system["coords"][-1][ii][jj], coords[ii][jj]
)


class TestOPENMXTraj(unittest.TestCase, TestOPENMXTRAJProps):
def setUp(self):
self.system = dpdata.System("openmx/Methane2", fmt="openmx/md")


class TestOPENMXLabeledTraj(unittest.TestCase, TestOPENMXTRAJProps):
def setUp(self):
self.system = dpdata.LabeledSystem("openmx/Methane2", fmt="openmx/md")


if __name__ == "__main__":
unittest.main()

0 comments on commit b048b37

Please sign in to comment.