Skip to content

feat: support virial in qe/traj #859

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 5 commits into
base: devel
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
26 changes: 26 additions & 0 deletions dpdata/qe/traj.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@
if TYPE_CHECKING:
from dpdata.utils import FileType

import os

from ..unit import (
EnergyConversion,
ForceConversion,
Expand All @@ -20,6 +22,7 @@

ry2ev = EnergyConversion("rydberg", "eV").value()
kbar2evperang3 = PressureConversion("kbar", "eV/angstrom^3").value()
gpa2evperbohr = PressureConversion("GPa", "eV/bohr^3").value()

length_convert = LengthConversion("bohr", "angstrom").value()
energy_convert = EnergyConversion("hartree", "eV").value()
Expand Down Expand Up @@ -228,6 +231,29 @@ def to_system_data(input_name, prefix, begin=0, step=1):
)
except FileNotFoundError:
data["cells"] = np.tile(cell, (data["coords"].shape[0], 1, 1))

# handle virial
stress_fname = prefix + ".str"
if os.path.exists(stress_fname):
# 1. Read stress tensor (in GPa) for each structure
stress, vsteps = load_data(stress_fname, 3, begin=begin, step=step, convert=1.0)
if csteps != vsteps:
csteps.append(None)
vsteps.append(None)
for int_id in range(len(csteps)):
if csteps[int_id] != vsteps[int_id]:
break
step_id = begin + int_id * step
raise RuntimeError(
f"the step key between files are not consistent. "
f"The difference locates at step: {step_id}, "
f".pos is {csteps[int_id]}, .str is {vsteps[int_id]}"
)
# 2. Calculate volume from cell. revert unit to bohr before taking det
volumes = np.linalg.det(data["cells"] / length_convert).reshape(-1)
# 3. Calculate virials for each structure, shape [nf x 3 x 3]
data["virials"] = gpa2evperbohr * volumes[:, None, None] * stress

return data, csteps


Expand Down
8 changes: 8 additions & 0 deletions tests/qe.traj/si.wrongstr/si.cel
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
20 0.00241888
10.64076206 0.00000000 0.00000000
0.00000000 10.64076206 0.00000000
0.00000000 0.00000000 10.64076206
30 0.00362833
10.64076206 0.00000000 0.00000000
0.00000000 10.64076206 0.00000000
0.00000000 0.00000000 10.64076206
3 changes: 3 additions & 0 deletions tests/qe.traj/si.wrongstr/si.evp
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
# nfi time(ps) ekinc T_cell(K) Tion(K) etot enthal econs econt Volume Pressure(GPa) EXX EVDW
20 2.418884E-03 2.666551E-08 0.000000E+00 1.400063E-01 -45.33322432 -45.33322432 -45.33321900 -45.33321908 1.204809E+03 0.27729
30 3.628326E-03 8.130638E-08 0.000000E+00 3.362418E-01 -45.33323102 -45.33323102 -45.33321824 -45.33321907 1.204809E+03 0.27691
18 changes: 18 additions & 0 deletions tests/qe.traj/si.wrongstr/si.for
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
20 0.00241888
0.79213120739417E-03 0.17846941660539E-02 0.25659395900843E-02
0.12296860773369E-02 -0.17575294070971E-02 -0.12077752684736E-02
0.53449804771536E-03 0.15716728646909E-02 0.11885992380187E-02
-0.20546369136946E-02 -0.31469738318253E-03 -0.13292260295548E-02
0.92996278884107E-03 -0.19904509183874E-02 -0.94000783512887E-03
0.25825132715153E-02 -0.17883694381162E-02 0.11238432666008E-02
0.63315248883314E-03 0.14328752516010E-02 -0.24575206497776E-02
-0.77380112683612E-03 -0.11964478769743E-02 0.14594402556008E-02
30 0.00362833
0.79771922134157E-03 0.17499394041424E-02 0.25186776447507E-02
0.12196024429585E-02 -0.17540520532661E-02 -0.12026567384908E-02
0.55342117414254E-03 0.15344688593791E-02 0.11794075108087E-02
-0.20168463577813E-02 -0.33437869352471E-03 -0.13161886410156E-02
0.92254287169080E-03 -0.19521664753444E-02 -0.92530766852405E-03
0.25425245627026E-02 -0.17646733377723E-02 0.11236784179976E-02
0.62964311964633E-03 0.14348901801160E-02 -0.24175435226922E-02
-0.77656387976694E-03 -0.11736525993363E-02 0.14442073730829E-02
56 changes: 56 additions & 0 deletions tests/qe.traj/si.wrongstr/si.in
Original file line number Diff line number Diff line change
@@ -0,0 +1,56 @@
&control
calculation = 'cp'
restart_mode = 'reset_counters'
prefix = 'si'
outdir = './outdir'
pseudo_dir='./'
nstep=3000
iprint=10
isave=100
dt=5.0
ndr=51, ndw=52 !folder for reading and writing
tstress=.true., tprnfor=.true.
!etot_conv_thr=1.d-9
!ekin_conv_thr=1.d-7
/

&system
ibrav=8,
celldm(1)=10.6407620646
celldm(2) = 1
celldm(3) = 1
nat=8
ntyp=1
ecutwfc=50
ecutrho=400,
nr1b=20, nr2b=20, nr3b=20 !For US PP specify no. of points for FFT of aug. Charge
/

&electrons
electron_dynamics = 'verlet'
electron_velocities = 'zero'
emass=300
orthogonalization='ortho'
ortho_eps=1d-11
/

&ions
ion_dynamics = 'verlet'
ion_velocities='zero'
ion_temperature='nose'
fnosep=25.0
tempw=40.0
/

ATOMIC_SPECIES
Si 28.0855 Si.pbe-n-rrkjus_psl.1.0.0.UPF

ATOMIC_POSITIONS (crystal)
Si 0.0000 0.0000 0.000 1 1 1
Si 0.0000 0.5000 0.5000 1 1 1
Si 0.5000 0.5000 0.0000 1 1 1
Si 0.5000 0.0000 0.5000 1 1 1
Si 0.25000 0.2500 0.2500 1 1 1
Si 0.25000 0.7500 0.7500 1 1 1
Si 0.75000 0.7500 0.2500 1 1 1
Si 0.75000 0.2500 0.7500 1 1 1
18 changes: 18 additions & 0 deletions tests/qe.traj/si.wrongstr/si.pos
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
20 0.00241888
-0.55955613343180E-03 -0.13280468681160E-01 -0.12170559744533E-01
-0.13498425299770E-01 0.53348556364374E+01 0.53265822820792E+01
0.53177169488520E+01 0.53125521646402E+01 -0.13687957707679E-02
0.53300642326909E+01 0.90214941469176E-04 0.53336756691314E+01
0.26530964223874E+01 0.26670146882442E+01 0.26685063320998E+01
0.26474643772774E+01 0.79813721194718E+01 0.79737688869610E+01
0.79764012327704E+01 0.79748937809524E+01 0.26729378003260E+01
0.79873283657577E+01 0.26616463190589E+01 0.79793202591917E+01
30 0.00362833
-0.52107908654115E-03 -0.13021934024520E-01 -0.11855267397130E-01
-0.13404740927597E-01 0.53346716430867E+01 0.53264251309777E+01
0.53177219927670E+01 0.53127844169680E+01 -0.12267029514939E-02
0.53297463825603E+01 0.86918798429361E-04 0.53335032260915E+01
0.26531522876578E+01 0.26668005206457E+01 0.26683824157731E+01
0.26477271882324E+01 0.79811839520919E+01 0.79739025595218E+01
0.79764197425739E+01 0.79751074086469E+01 0.26726239769995E+01
0.79871718245254E+01 0.26615315288521E+01 0.79794965352588E+01
12 changes: 12 additions & 0 deletions tests/qe.traj/si.wrongstr/si.str
Original file line number Diff line number Diff line change
@@ -0,0 +1,12 @@
20 0.00241888
0.27927908 -0.02926875 -0.02277043
-0.02926875 0.27823138 0.03779023
-0.02277043 0.03779026 0.27434825
30 0.00362833
0.27885067 -0.02827892 -0.02066094
-0.02827892 0.27775865 0.03660021
-0.02066094 0.03660025 0.27411244
40 0.00362833
0.27885067 -0.02827892 -0.02066094
-0.02827892 0.27775865 0.03660021
-0.02066094 0.03660025 0.27411244
8 changes: 8 additions & 0 deletions tests/qe.traj/si/si.cel
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
20 0.00241888
10.64076206 0.00000000 0.00000000
0.00000000 10.64076206 0.00000000
0.00000000 0.00000000 10.64076206
30 0.00362833
10.64076206 0.00000000 0.00000000
0.00000000 10.64076206 0.00000000
0.00000000 0.00000000 10.64076206
3 changes: 3 additions & 0 deletions tests/qe.traj/si/si.evp
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
# nfi time(ps) ekinc T_cell(K) Tion(K) etot enthal econs econt Volume Pressure(GPa) EXX EVDW
20 2.418884E-03 2.666551E-08 0.000000E+00 1.400063E-01 -45.33322432 -45.33322432 -45.33321900 -45.33321908 1.204809E+03 0.27729
30 3.628326E-03 8.130638E-08 0.000000E+00 3.362418E-01 -45.33323102 -45.33323102 -45.33321824 -45.33321907 1.204809E+03 0.27691
18 changes: 18 additions & 0 deletions tests/qe.traj/si/si.for
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
20 0.00241888
0.79213120739417E-03 0.17846941660539E-02 0.25659395900843E-02
0.12296860773369E-02 -0.17575294070971E-02 -0.12077752684736E-02
0.53449804771536E-03 0.15716728646909E-02 0.11885992380187E-02
-0.20546369136946E-02 -0.31469738318253E-03 -0.13292260295548E-02
0.92996278884107E-03 -0.19904509183874E-02 -0.94000783512887E-03
0.25825132715153E-02 -0.17883694381162E-02 0.11238432666008E-02
0.63315248883314E-03 0.14328752516010E-02 -0.24575206497776E-02
-0.77380112683612E-03 -0.11964478769743E-02 0.14594402556008E-02
30 0.00362833
0.79771922134157E-03 0.17499394041424E-02 0.25186776447507E-02
0.12196024429585E-02 -0.17540520532661E-02 -0.12026567384908E-02
0.55342117414254E-03 0.15344688593791E-02 0.11794075108087E-02
-0.20168463577813E-02 -0.33437869352471E-03 -0.13161886410156E-02
0.92254287169080E-03 -0.19521664753444E-02 -0.92530766852405E-03
0.25425245627026E-02 -0.17646733377723E-02 0.11236784179976E-02
0.62964311964633E-03 0.14348901801160E-02 -0.24175435226922E-02
-0.77656387976694E-03 -0.11736525993363E-02 0.14442073730829E-02
56 changes: 56 additions & 0 deletions tests/qe.traj/si/si.in
Original file line number Diff line number Diff line change
@@ -0,0 +1,56 @@
&control
calculation = 'cp'
restart_mode = 'reset_counters'
prefix = 'si'
outdir = './outdir'
pseudo_dir='./'
nstep=3000
iprint=10
isave=100
dt=5.0
ndr=51, ndw=52 !folder for reading and writing
tstress=.true., tprnfor=.true.
!etot_conv_thr=1.d-9
!ekin_conv_thr=1.d-7
/

&system
ibrav=8,
celldm(1)=10.6407620646
celldm(2) = 1
celldm(3) = 1
nat=8
ntyp=1
ecutwfc=50
ecutrho=400,
nr1b=20, nr2b=20, nr3b=20 !For US PP specify no. of points for FFT of aug. Charge
/

&electrons
electron_dynamics = 'verlet'
electron_velocities = 'zero'
emass=300
orthogonalization='ortho'
ortho_eps=1d-11
/

&ions
ion_dynamics = 'verlet'
ion_velocities='zero'
ion_temperature='nose'
fnosep=25.0
tempw=40.0
/

ATOMIC_SPECIES
Si 28.0855 Si.pbe-n-rrkjus_psl.1.0.0.UPF

ATOMIC_POSITIONS (crystal)
Si 0.0000 0.0000 0.000 1 1 1
Si 0.0000 0.5000 0.5000 1 1 1
Si 0.5000 0.5000 0.0000 1 1 1
Si 0.5000 0.0000 0.5000 1 1 1
Si 0.25000 0.2500 0.2500 1 1 1
Si 0.25000 0.7500 0.7500 1 1 1
Si 0.75000 0.7500 0.2500 1 1 1
Si 0.75000 0.2500 0.7500 1 1 1
18 changes: 18 additions & 0 deletions tests/qe.traj/si/si.pos
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
20 0.00241888
-0.55955613343180E-03 -0.13280468681160E-01 -0.12170559744533E-01
-0.13498425299770E-01 0.53348556364374E+01 0.53265822820792E+01
0.53177169488520E+01 0.53125521646402E+01 -0.13687957707679E-02
0.53300642326909E+01 0.90214941469176E-04 0.53336756691314E+01
0.26530964223874E+01 0.26670146882442E+01 0.26685063320998E+01
0.26474643772774E+01 0.79813721194718E+01 0.79737688869610E+01
0.79764012327704E+01 0.79748937809524E+01 0.26729378003260E+01
0.79873283657577E+01 0.26616463190589E+01 0.79793202591917E+01
30 0.00362833
-0.52107908654115E-03 -0.13021934024520E-01 -0.11855267397130E-01
-0.13404740927597E-01 0.53346716430867E+01 0.53264251309777E+01
0.53177219927670E+01 0.53127844169680E+01 -0.12267029514939E-02
0.53297463825603E+01 0.86918798429361E-04 0.53335032260915E+01
0.26531522876578E+01 0.26668005206457E+01 0.26683824157731E+01
0.26477271882324E+01 0.79811839520919E+01 0.79739025595218E+01
0.79764197425739E+01 0.79751074086469E+01 0.26726239769995E+01
0.79871718245254E+01 0.26615315288521E+01 0.79794965352588E+01
8 changes: 8 additions & 0 deletions tests/qe.traj/si/si.str
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
20 0.00241888
0.27927908 -0.02926875 -0.02277043
-0.02926875 0.27823138 0.03779023
-0.02277043 0.03779026 0.27434825
30 0.00362833
0.27885067 -0.02827892 -0.02066094
-0.02827892 0.27775865 0.03660021
-0.02066094 0.03660025 0.27411244
35 changes: 35 additions & 0 deletions tests/test_qe_cp_traj.py
Original file line number Diff line number Diff line change
Expand Up @@ -68,5 +68,40 @@ def test_case_null(self):
self.assertAlmostEqual(cell[ii][jj], ref[ii][jj])


class TestVirial(unittest.TestCase):
def test(self):
self.system = dpdata.LabeledSystem("qe.traj/si/si", fmt="qe/cp/traj")
self.assertEqual(self.system["virials"].shape, (2, 3, 3))
np.testing.assert_almost_equal(
self.system["virials"][0],
np.array(
[
[0.31120718, -0.03261485, -0.02537362],
[-0.03261485, 0.3100397, 0.04211053],
[-0.02537362, 0.04211057, 0.30571264],
]
),
)
np.testing.assert_almost_equal(
self.system["virials"][1],
np.array(
[
[0.31072979, -0.03151186, -0.02302297],
[-0.03151186, 0.30951293, 0.04078447],
[-0.02302297, 0.04078451, 0.30544987],
]
),
)

def test_raise(self):
with self.assertRaises(RuntimeError) as c:
self.system = dpdata.LabeledSystem(
"qe.traj/si.wrongstr/si", fmt="qe/cp/traj"
)
self.assertTrue(
"the step key between files are not consistent." in str(c.exception)
)


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