Commit eaa6fd34 authored by Lauri Himanen's avatar Lauri Himanen
Browse files

NPT, NVT parsing.

parent 70e01898
......@@ -5,7 +5,9 @@ from commonmatcher import CommonMatcher
import cp2kparser.generic.configurationreading
import cp2kparser.generic.csvparsing
from nomadcore.caching_backend import CachingLevel
from nomadcore.unit_conversion.unit_conversion import convert_unit
import logging
import math
logger = logging.getLogger("nomad")
......@@ -53,16 +55,19 @@ class CP2KMDParser(MainHierarchicalParser):
subMatchers=[
SM( " MD\| Molecular Dynamics Protocol",
forwardMatch=True,
sections=["x_cp2k_section_md_settings"],
sections=["section_sampling_method", "x_cp2k_section_md_settings"],
subMatchers=[
SM( " MD\| Ensemble Type\s+(?P<x_cp2k_md_ensemble_type>{})".format(self.cm.regex_word)),
SM( " MD\| Number of Time Steps\s+(?P<x_cp2k_md_number_of_time_steps>{})".format(self.cm.regex_i)),
SM( " MD\| Time Step \[fs\]\s+(?P<x_cp2k_md_time_step__fs>{})".format(self.cm.regex_f)),
SM( " MD\| Temperature \[K\]\s+(?P<x_cp2k_md_temperature>{})".format(self.cm.regex_f)),
SM( " MD\| Temperature tolerance \[K\]\s+(?P<x_cp2k_md_temperature_tolerance>{})".format(self.cm.regex_f)),
SM( " MD\| Temperature \[K\]\s+(?P<x_cp2k_md_target_temperature>{})".format(self.cm.regex_f)),
SM( " MD\| Temperature tolerance \[K\]\s+(?P<x_cp2k_md_target_temperature_tolerance>{})".format(self.cm.regex_f)),
SM( " MD\| Pressure \[Bar\]\s+(?P<x_cp2k_md_target_pressure>{})".format(self.cm.regex_f)),
SM( " MD\| Barostat time constant \[ fs\]\s+(?P<x_cp2k_md_barostat_time_constant>{})".format(self.cm.regex_f)),
SM( " MD\| Print MD information every\s+(?P<x_cp2k_md_print_frequency>{}) step\(s\)".format(self.cm.regex_i)),
SM( " MD\| File type Print frequency\[steps\] File names"),
SM( " MD\| Coordinates\s+(?P<x_cp2k_md_coordinates_print_frequency>{})\s+(?P<x_cp2k_md_coordinates_filename>{})".format(self.cm.regex_i, self.cm.regex_word)),
SM( " MD\| Simulation Cel\s+(?P<x_cp2k_md_simulation_cell_print_frequency>{})\s+(?P<x_cp2k_md_simulation_cell_filename>{})".format(self.cm.regex_i, self.cm.regex_word)),
SM( " MD\| Velocities\s+(?P<x_cp2k_md_velocities_print_frequency>{})\s+(?P<x_cp2k_md_velocities_filename>{})".format(self.cm.regex_i, self.cm.regex_word)),
SM( " MD\| Energies\s+(?P<x_cp2k_md_energies_print_frequency>{})\s+(?P<x_cp2k_md_energies_filename>{})".format(self.cm.regex_i, self.cm.regex_word)),
SM( " MD\| Dump\s+(?P<x_cp2k_md_dump_print_frequency>{})\s+(?P<x_cp2k_md_dump_filename>{})".format(self.cm.regex_i, self.cm.regex_word)),
......@@ -78,7 +83,11 @@ class CP2KMDParser(MainHierarchicalParser):
SM( " INITIAL POTENTIAL ENERGY\[hartree\] =\s+(?P<x_cp2k_md_potential_energy_instantaneous__hartree>{})".format(self.cm.regex_f)),
SM( " INITIAL KINETIC ENERGY\[hartree\] =\s+(?P<x_cp2k_md_kinetic_energy_instantaneous__hartree>{})".format(self.cm.regex_f)),
SM( " INITIAL TEMPERATURE\[K\] =\s+(?P<x_cp2k_md_temperature_instantaneous>{})".format(self.cm.regex_f)),
SM( " INITIAL VOLUME\[bohr\^3\] =\s+(?P<x_cp2k_md_volume__bohr3>{})".format(self.cm.regex_f)),
SM( " INITIAL BAROSTAT TEMP\[K\] =\s+(?P<x_cp2k_md_barostat_temperature_instantaneous>{})".format(self.cm.regex_f)),
SM( " INITIAL PRESSURE\[bar\] =\s+(?P<x_cp2k_md_pressure_instantaneous__bar>{})".format(self.cm.regex_f)),
SM( " INITIAL VOLUME\[bohr\^3\] =\s+(?P<x_cp2k_md_volume_instantaneous__bohr3>{})".format(self.cm.regex_f)),
SM( " INITIAL CELL LNTHS\[bohr\] =\s+(?P<x_cp2k_md_cell_length_a_instantaneous__bohr>{0})\s+(?P<x_cp2k_md_cell_length_b_instantaneous__bohr>{0})\s+(?P<x_cp2k_md_cell_length_c_instantaneous__bohr>{0})".format(self.cm.regex_f)),
SM( " INITIAL CELL ANGLS\[deg\] =\s+(?P<x_cp2k_md_cell_angle_a_instantaneous__deg>{0})\s+(?P<x_cp2k_md_cell_angle_b_instantaneous__deg>{0})\s+(?P<x_cp2k_md_cell_angle_c_instantaneous__deg>{0})".format(self.cm.regex_f)),
],
),
SM( " SCF WAVEFUNCTION OPTIMIZATION",
......@@ -97,6 +106,13 @@ class CP2KMDParser(MainHierarchicalParser):
SM( " POTENTIAL ENERGY\[hartree\] =\s+(?P<x_cp2k_md_potential_energy_instantaneous__hartree>{})\s+(?P<x_cp2k_md_potential_energy_average__hartree>{})".format(self.cm.regex_f, self.cm.regex_f)),
SM( " KINETIC ENERGY \[hartree\] =\s+(?P<x_cp2k_md_kinetic_energy_instantaneous__hartree>{})\s+(?P<x_cp2k_md_kinetic_energy_average__hartree>{})".format(self.cm.regex_f, self.cm.regex_f)),
SM( " TEMPERATURE \[K\] =\s+(?P<x_cp2k_md_temperature_instantaneous>{})\s+(?P<x_cp2k_md_temperature_average>{})".format(self.cm.regex_f, self.cm.regex_f)),
SM( " PRESSURE \[bar\] =\s+(?P<x_cp2k_md_pressure_instantaneous__bar>{0})\s+(?P<x_cp2k_md_pressure_average__bar>{0})".format(self.cm.regex_f)),
SM( " BAROSTAT TEMP\[K\] =\s+(?P<x_cp2k_md_barostat_temperature_instantaneous>{0})\s+(?P<x_cp2k_md_barostat_temperature_average>{0})".format(self.cm.regex_f)),
SM( " VOLUME\[bohr\^3\] =\s+(?P<x_cp2k_md_volume_instantaneous__bohr3>{0})\s+(?P<x_cp2k_md_volume_average__bohr3>{0})".format(self.cm.regex_f)),
SM( " CELL LNTHS\[bohr\] =\s+(?P<x_cp2k_md_cell_length_a_instantaneous__bohr>{0})\s+(?P<x_cp2k_md_cell_length_b_instantaneous__bohr>{0})\s+(?P<x_cp2k_md_cell_length_c_instantaneous__bohr>{0})".format(self.cm.regex_f)),
SM( " AVE. CELL LNTHS\[bohr\] =\s+(?P<x_cp2k_md_cell_length_a_average__bohr>{0})\s+(?P<x_cp2k_md_cell_length_b_average__bohr>{0})\s+(?P<x_cp2k_md_cell_length_c_average__bohr>{0})".format(self.cm.regex_f)),
SM( " CELL ANGLS\[deg\] =\s+(?P<x_cp2k_md_cell_angle_a_instantaneous__deg>{0})\s+(?P<x_cp2k_md_cell_angle_b_instantaneous__deg>{0})\s+(?P<x_cp2k_md_cell_angle_c_instantaneous__deg>{0})".format(self.cm.regex_f)),
SM( " AVE. CELL ANGLS\[deg\] =\s+(?P<x_cp2k_md_cell_angle_a_average__deg>{0})\s+(?P<x_cp2k_md_cell_angle_b_average__deg>{0})\s+(?P<x_cp2k_md_cell_angle_c_average__deg>{0})".format(self.cm.regex_f)),
],
),
]
......@@ -107,7 +123,7 @@ class CP2KMDParser(MainHierarchicalParser):
# computational time.
self.root_matcher = SM("",
forwardMatch=True,
sections=["section_run", "section_sampling_method"],
sections=["section_run"],
subMatchers=[
SM( "",
forwardMatch=True,
......@@ -131,7 +147,8 @@ class CP2KMDParser(MainHierarchicalParser):
sampling_map = {
"NVE": "NVE",
"NVT": "NVT",
"NPT": "NPT",
"NPT_F": "NPT",
"NPT_I": "NPT",
}
sampling = sampling_map.get(sampling)
if sampling is not None:
......@@ -220,12 +237,15 @@ class CP2KMDParser(MainHierarchicalParser):
frame_sequence_time = []
frame_sequence_kinetic_energy = []
frame_sequence_conserved_quantity = []
frame_sequence_pressure = []
single_conf_gids = []
i_md_step = 0
for i_step in range(self.n_steps + 1):
sectionGID = backend.openSection("section_single_configuration_calculation")
systemGID = backend.openSection("section_system")
single_conf_gids.append(sectionGID)
# Trajectory
if self.traj_iterator is not None:
......@@ -272,13 +292,65 @@ class CP2KMDParser(MainHierarchicalParser):
quickstep.add_latest_value(["x_cp2k_section_scf_iteration", "x_cp2k_energy_XC_scf_iteration"], "energy_XC_scf_iteration")
backend.closeSection("section_scf_iteration", scfGID)
i_md_step += 1
pressure = md_step.get_latest_value("x_cp2k_md_pressure_instantaneous")
if pressure is not None:
frame_sequence_pressure.append(pressure)
backend.closeSection("section_system", systemGID)
backend.closeSection("section_single_configuration_calculation", sectionGID)
# Add the summaries to frame sequence
backend.addArrayValues("frame_sequence_potential_energy", np.array(frame_sequence_potential_energy), unit="hartree")
backend.addArrayValues("frame_sequence_kinetic_energy", np.array(frame_sequence_kinetic_energy), unit="hartree")
backend.addArrayValues("frame_sequence_conserved_quantity", np.array(frame_sequence_conserved_quantity), unit="hartree")
backend.addArrayValues("frame_sequence_temperature", np.array(frame_sequence_temperature))
backend.addArrayValues("frame_sequence_time", np.array(frame_sequence_time), unit="fs")
frame_sequence_potential_energy = convert_unit(np.array(frame_sequence_potential_energy), "hartree")
frame_sequence_kinetic_energy = convert_unit(np.array(frame_sequence_kinetic_energy), "hartree")
frame_sequence_conserved_quantity = convert_unit(np.array(frame_sequence_conserved_quantity), "hartree")
frame_sequence_time = np.array(frame_sequence_time)
frame_sequence_temperature = np.array(frame_sequence_temperature)
frame_sequence_pressure = np.array(frame_sequence_pressure)
backend.addArrayValues("frame_sequence_potential_energy", frame_sequence_potential_energy)
backend.addArrayValues("frame_sequence_kinetic_energy", frame_sequence_kinetic_energy)
backend.addArrayValues("frame_sequence_conserved_quantity", frame_sequence_conserved_quantity)
backend.addArrayValues("frame_sequence_temperature", frame_sequence_temperature)
backend.addArrayValues("frame_sequence_time", frame_sequence_time, unit="fs")
if frame_sequence_pressure.size != 0:
backend.addArrayValues("frame_sequence_pressure", frame_sequence_pressure)
# Number of frames
backend.addValue("number_of_frames_in_sequence", self.n_steps + 1)
# Reference to sampling method
backend.addValue("frame_sequence_to_sampling_ref", 0)
# References to local frames
backend.addArrayValues("frame_sequence_local_frames_ref", np.array(single_conf_gids))
# Temperature stats
mean_temp = frame_sequence_temperature.mean()
c = frame_sequence_temperature - mean_temp
std_temp = math.sqrt(np.dot(c, c)/frame_sequence_temperature.size)
backend.addArrayValues("frame_sequence_temperature_stats", np.array([mean_temp, std_temp]))
# Potential energy stats
mean_pot = frame_sequence_potential_energy.mean()
c = frame_sequence_potential_energy - mean_pot
std_pot = math.sqrt(np.dot(c, c)/frame_sequence_potential_energy.size)
backend.addArrayValues("frame_sequence_potential_energy_stats", np.array([mean_pot, std_pot]))
# Kinetic energy stats
mean_kin = frame_sequence_kinetic_energy.mean()
c = frame_sequence_kinetic_energy - mean_kin
std_kin = math.sqrt(np.dot(c, c)/frame_sequence_kinetic_energy.size)
backend.addArrayValues("frame_sequence_kinetic_energy_stats", np.array([mean_kin, std_kin]))
# Conserved quantity stats
mean_cons = frame_sequence_conserved_quantity.mean()
c = frame_sequence_conserved_quantity - mean_cons
std_cons = math.sqrt(np.dot(c, c)/frame_sequence_conserved_quantity.size)
backend.addArrayValues("frame_sequence_conserved_quantity_stats", np.array([mean_cons, std_cons]))
# Pressure stats
if frame_sequence_pressure.size != 0:
mean_pressure = frame_sequence_pressure.mean()
c = frame_sequence_pressure - mean_pressure
std_pressure = math.sqrt(np.dot(c, c)/frame_sequence_pressure.size)
backend.addArrayValues("frame_sequence_pressure_stats", np.array([mean_pressure, std_pressure]))
# Step Nr. Time[fs] Kin.[a.u.] Temp[K] Pot.[a.u.] Cons Qty[a.u.] UsedTime[s]
0 0.000000 0.007125335 300.000000000 -34.336403530 -34.324047996 0.000000000
1 0.500000 0.007358038 309.797597156 -34.336638312 -34.324054403 2.979168899
2 1.000000 0.007414878 312.190722217 -34.336698674 -34.324056733 0.842197949
3 1.500000 0.007245416 305.055835462 -34.336539787 -34.324059082 1.418614861
4 2.000000 0.006832946 287.689475065 -34.336141815 -34.324057040 1.415143902
5 2.500000 0.006241331 262.780549510 -34.335562912 -34.324043992 1.398936698
6 3.000000 0.005685824 239.391887129 -34.335021181 -34.324022684 0.867581070
7 3.500000 0.005480859 230.762197167 -34.334844452 -34.324007616 1.268830910
8 4.000000 0.005835031 245.673981276 -34.335254239 -34.324014235 1.262329495
9 4.500000 0.006642601 279.675336061 -34.336142063 -34.324043595 1.263161745
10 5.000000 0.007478540 314.871094720 -34.337060810 -34.324078054 1.387408595
6
i = 0, time = 0.000, E = -34.3364035304
O -4.5832976890 5.3339524760 1.5600533860
O -4.5832976890 8.3339524760 1.5600533860
H -3.7776655230 5.3315557480 0.9431994470
H -5.0815598010 4.5898899650 1.1769412080
H -3.7776655230 8.3315557480 0.9431994470
H -5.0815598010 7.5898899650 1.1769412080
6
i = 1, time = 0.500, E = -34.3366383116
O -4.5817445718 5.3355868921 1.5589463670
O -4.5843710351 8.3367641088 1.5602036855
H -3.7808196412 5.3299000381 0.9495250366
H -5.0946313631 4.5838911317 1.1707137947
H -3.7897032638 8.3431761391 0.9355342799
H -5.0853477221 7.6085130826 1.1707220814
6
i = 2, time = 1.000, E = -34.3366986739
O -4.5799918657 5.3371286726 1.5575079911
O -4.5851343385 8.3396223014 1.5600600622
H -3.7861200969 5.3285392590 0.9580740075
H -5.1071763088 4.5781786556 1.1657849252
H -3.8041487990 8.3552228149 0.9305319372
H -5.0901482179 7.6252630035 1.1648052814
6
i = 3, time = 1.500, E = -34.3365397866
O -4.5781525332 5.3385086913 1.5557554785
O -4.5856067912 8.3425610025 1.5596663005
H -3.7927172952 5.3274048616 0.9680705618
H -5.1183861690 4.5738665384 1.1627189257
H -3.8204148813 8.3677541909 0.9278283032
H -5.0964224937 7.6394940377 1.1589505954
6
i = 4, time = 2.000, E = -34.3361418152
O -4.5763474117 5.3396767410 1.5537289640
O -4.5858370053 8.3455949542 1.5590739711
H -3.7994932085 5.3264231795 0.9785678628
H -5.1276395448 4.5717746569 1.1619227294
H -3.8375962604 8.3808222456 0.9268079947
H -5.1044431438 7.6508796083 1.1530569525
6
i = 5, time = 2.500, E = -34.3355629121
O -4.5746973646 5.3406032468 1.5514872810
O -4.5859023088 8.3487166732 1.5583381725
H -3.8051398314 5.3255151479 0.9885418342
H -5.1345026391 4.5724201887 1.1636284032
H -3.8545273685 8.3944656500 0.9266321439
H -5.1142414556 7.6594624260 1.1471963631
6
i = 6, time = 3.000, E = -34.3350211809
O -4.5733108862 5.3412793526 1.5490990859
O -4.5859010716 8.3518981492 1.5575126636
H -3.8083270379 5.3245982175 0.9970179837
H -5.1387347438 4.5760274540 1.1678900873
H -3.8699633888 8.4087054796 0.9263165384
H -5.1256013031 7.6656128586 1.1416025367
6
i = 7, time = 3.500, E = -34.3348444523
O -4.5722672947 5.3417167646 1.5466308712
O -4.5859384471 8.3550968032 1.5566437253
H -3.8079480022 5.3235909353 1.0032208412
H -5.1403047092 4.5825311319 1.1745837856
H -3.8828053198 8.4235445448 0.9248681293
H -5.1380900787 7.6699325915 1.1366207240
6
i = 8, time = 4.000, E = -34.3352542387
O -4.5715987287 5.3419478216 1.5441349254
O -4.5861081221 8.3582640927 1.5557626871
H -3.8033803816 5.3224195671 1.0067267084
H -5.1394103594 4.5915665781 1.1834103145
H -3.8922840023 8.4389692827 0.9214663221
H -5.1511094332 7.6731378661 1.1326376756
6
i = 9, time = 4.500, E = -34.3361420625
O -4.5712774870 5.3420256604 1.5416418597
O -4.5864750634 8.3613547086 1.5548789667
H -3.7946535888 5.3210242633 1.0075618904
H -5.1364921448 4.6024497398 1.1939105245
H -3.8980875126 8.4549504626 0.9156404974
H -5.1639606530 7.6759486146 1.1300132301
6
i = 10, time = 5.000, E = -34.3370608096
O -4.5712148711 5.3420237890 1.5391606731
O -4.5870647368 8.3643343914 1.5539766792
H -3.7824389243 5.3193620987 1.0061894957
H -5.1322314658 4.6141604170 1.2054997435
H -3.9004203933 8.4714379402 0.9073777219
H -5.1759289758 7.6790007799 1.1290280189
&FORCE_EVAL
METHOD Quickstep
STRESS_TENSOR ANALYTICAL
&DFT
BASIS_SET_FILE_NAME ../../BASIS_SET
POTENTIAL_FILE_NAME ../../GTH_POTENTIALS
&MGRID
CUTOFF 200
NGRIDS 4
REL_CUTOFF 20
&END MGRID
&QS
EPS_DEFAULT 1.0E-12
WF_INTERPOLATION PS
EXTRAPOLATION_ORDER 3
&END QS
&SCF
&PRINT
&RESTART OFF
&END
&END
EPS_SCF 1.0E-6
MAX_SCF 150
SCF_GUESS atomic
&END SCF
&XC
&XC_FUNCTIONAL Pade
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 6.000 6.000 6.000
&END CELL
&COORD
O -4.583297689 5.333952476 1.560053386
O -4.583297689 8.333952476 1.560053386
H -3.777665523 5.331555748 0.9431994470
H -5.081559801 4.589889965 1.176941208
H -3.777665523 8.331555748 0.9431994470
H -5.081559801 7.589889965 1.176941208
&END COORD
&KIND O
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q6
&END KIND
&KIND H
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q1
&END KIND
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT H2O-md
RUN_TYPE MD
PRINT_LEVEL MEDIUM
&END GLOBAL
&MOTION
&MD
ENSEMBLE NPT_F
STEPS 10
TIMESTEP 0.5
TEMPERATURE 300.0
&THERMOSTAT
&NOSE
LENGTH 3
YOSHIDA 3
TIMECON 1000.
MTS 2
&END NOSE
&END
&BAROSTAT
PRESSURE 1
&END BAROSTAT
&END MD
&END MOTION
This diff is collapsed.
# Step Nr. Time[fs] Kin.[a.u.] Temp[K] Pot.[a.u.] Cons Qty[a.u.] UsedTime[s]
0 0.000000 0.007125335 300.000000000 -34.336403530 -34.328803174 0.000000000
1 0.500000 0.007392508 311.248871016 -34.336672017 -34.328802635 2.148559709
2 1.000000 0.008067310 339.660255576 -34.337355043 -34.328808886 0.864046893
3 1.500000 0.009049638 381.019502834 -34.338350140 -34.328819470 0.723231513
4 2.000000 0.010159126 427.732574903 -34.339475566 -34.328832955 0.704930791
5 2.500000 0.011176982 470.587657681 -34.340509751 -34.328846556 0.463860581
6 3.000000 0.011910577 501.474453250 -34.341253436 -34.328853690 0.707140392
7 3.500000 0.012291989 517.533117883 -34.341630388 -34.328846129 0.704500439
8 4.000000 0.012459614 524.590704368 -34.341779925 -34.328824864 0.699688757
9 4.500000 0.012709582 535.115189702 -34.342014148 -34.328805883 0.783292716
10 5.000000 0.013296346 559.819878564 -34.342602508 -34.328804131 0.779202085
6
i = 0, time = 0.000, E = -34.3364035304
O -4.5832976890 5.3339524760 1.5600533860
O -4.5832976890 8.3339524760 1.5600533860
H -3.7776655230 5.3315557480 0.9431994470
H -5.0815598010 4.5898899650 1.1769412080
H -3.7776655230 8.3315557480 0.9431994470
H -5.0815598010 7.5898899650 1.1769412080
6
i = 1, time = 0.500, E = -34.3366720170
O -4.5832141237 5.3348615025 1.5594895361
O -4.5851051829 8.3343406371 1.5584697202
H -3.7641052192 5.3220856563 0.9511634220
H -5.0679756011 4.5853534777 1.1803360421
H -3.7750119931 8.3280492899 0.9563835095
H -5.0827396824 7.5857887257 1.1849577344
6
i = 2, time = 1.000, E = -34.3373550429
O -4.5828049383 5.3357691352 1.5586089968
O -4.5866562402 8.3347154046 1.5566020367
H -3.7535977354 5.3130127821 0.9621191549
H -5.0552348681 4.5793407851 1.1841770228
H -3.7747808015 8.3248037224 0.9719224244
H -5.0843282279 7.5806009883 1.1934920789
6
i = 3, time = 1.500, E = -34.3383501397
O -4.5820489290 5.3366838261 1.5574164613
O -4.5879973096 8.3350769589 1.5544875646
H -3.7462220341 5.3044063128 0.9760413735
H -5.0436371645 4.5715894016 1.1883549041
H -3.7763914112 8.3217228178 0.9892247854
H -5.0862018301 7.5743945505 1.2024939792
6
i = 4, time = 2.000, E = -34.3394755665
O -4.5809522660 5.3376017740 1.5559245524
O -4.5891965182 8.3354207910 1.5521718699
H -3.7417785789 5.2963500176 0.9927183401
H -5.0333412437 4.5620220298 1.1928647337
H -3.7789875266 8.3187003624 1.0076153202
H -5.0881783004 7.5673142580 1.2119797771
6
i = 5, time = 2.500, E = -34.3405097510
O -4.5795491273 5.3385062527 1.5541528845
O -4.5903398830 8.3357392133 1.5497065976
H -3.7397632032 5.2889278617 1.0118086580
H -5.0243367934 4.5507617109 1.1978087788
H -3.7814809155 8.3156165145 1.0263624254
H -5.0900309156 7.5595583610 1.2220188093
6
i = 6, time = 3.000, E = -34.3412534359
O -4.5779007434 5.3393694326 1.5521281411
O -4.5915236270 8.3360233929 1.5471452643
H -3.7393559470 5.2822040584 1.0328330950
H -5.0164596950 4.5381206968 1.2033782222
H -3.7826498710 8.3123359540 1.0446600978
H -5.0915076171 7.5513497901 1.2327082964
6
i = 7, time = 3.500, E = -34.3416303882
O -4.5760921479 5.3401563070 1.5498842299
O -4.5928417265 8.3362654016 1.5445366175
H -3.7394416888 5.2762030099 1.0550905498
H -5.0094342338 4.5245773555 1.2098200186
H -3.7813106527 8.3087103779 1.0616871265
H -5.0923522905 7.5429091919 1.2441397137
6
i = 8, time = 4.000, E = -34.3417799250
O -4.5742250712 5.3408297161 1.5474606674
O -4.5943696311 8.3364597947 1.5419166175
H -3.7386894421 5.2708966763 1.0776631045
H -5.0029238013 4.5107413737 1.2173979021
H -3.7765497830 8.3045878080 1.0768004989
H -5.0923255202 7.5344332562 1.2563619268
6
i = 9, time = 4.500, E = -34.3420141482
O -4.5724056254 5.3413551533 1.5448981480
O -4.5961485492 8.3366043602 1.5393019255
H -3.7357153878 5.2662061862 1.0996567950
H -4.9965782504 4.4972892910 1.2263556009
H -3.7679497896 8.2998290242 1.0897340295
H -5.0912224598 7.5260787626 1.2693508542
6
i = 10, time = 5.000, E = -34.3426025085
O -4.5707272065 5.3417043346 1.5422317470
O -4.5981764323 8.3367000926 1.5366881064
H -3.7293273128 5.2620203316 1.1204819383
H -4.9900727878 4.4848628743 1.2368847643
H -3.7557182679 8.2943276832 1.1006576433
H -5.0888863567 7.5179510705 1.2830009137
&FORCE_EVAL
METHOD Quickstep
&DFT
BASIS_SET_FILE_NAME ../../BASIS_SET
POTENTIAL_FILE_NAME ../../GTH_POTENTIALS
&MGRID
CUTOFF 200
NGRIDS 4
REL_CUTOFF 20
&END MGRID
&QS
EPS_DEFAULT 1.0E-12
WF_INTERPOLATION PS
EXTRAPOLATION_ORDER 3
&END QS
&SCF
&PRINT
&RESTART OFF
&END
&END
EPS_SCF 1.0E-6
MAX_SCF 150
SCF_GUESS atomic
&END SCF
&XC
&XC_FUNCTIONAL Pade
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 6.000 6.000 6.000
&END CELL
&COORD
O -4.583297689 5.333952476 1.560053386
O -4.583297689 8.333952476 1.560053386
H -3.777665523 5.331555748 0.9431994470
H -5.081559801 4.589889965 1.176941208
H -3.777665523 8.331555748 0.9431994470
H -5.081559801 7.589889965 1.176941208
&END COORD
&KIND O
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q6
&END KIND
&KIND H
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q1
&END KIND
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT H2O-md
RUN_TYPE MD
PRINT_LEVEL MEDIUM
&END GLOBAL
&MOTION
&MD
ENSEMBLE NVT
STEPS 10
TIMESTEP 0.5
TEMPERATURE 300.0
&THERMOSTAT
&NOSE
LENGTH 3
YOSHIDA 3
TIMECON 1000.
MTS 2
&END NOSE
&END
&END MD
&END MOTION
This diff is collapsed.
......@@ -685,6 +685,70 @@ class TestMD(unittest.TestCase):
@classmethod
def setUpClass(cls):
cls.results = get_results("md/nve", "section_run")
cls.temp = convert_unit(
np.array([
300.000000000,
275.075405378,
235.091633019,
202.752506973,
192.266488819,
201.629598676,
218.299664775,
230.324748558,
232.691881533,
226.146979313,
213.165337396,
]),
"K"
)
cls.cons = convert_unit(
np.array([
-34.323271136,
-34.323245645,
-34.323206964,
-34.323183380,
-34.323187747,
-34.323208962,
-34.323227533,
-34.323233583,
-34.323230715,
-34.323227013,
-34.323224123,
]),
"hartree"
)
cls.pot = convert_unit(
np.array([
-34.330396471,
-34.329778993,
-34.328790653,
-34.327998978,
-34.327754290,
-34.327997890,
-34.328412394,
-34.328704052,
-34.328757407,
-34.328598255,
-34.328287038,
]),
"hartree"
)
cls.kin = convert_unit(
np.array([
0.007125335,
0.006533348,
0.005583688,
0.004815598,
0.004566544,
0.004788928,
0.005184860,
0.005470470,
0.005526692,
0.005371243,
0.005062914,
]),
"hartree"
)
def test_ensemble_type(self):
result = self.results["ensemble_type"]
......@@ -694,6 +758,10 @@ class TestMD(unittest.TestCase):
result = self.results["sampling_method"]
self.assertEqual(result, "molecular_dynamics")
def test_number_of_frames_in_sequence(self):
result = self.results["number_of_frames_in_sequence"]
self.assertEqual(result, 11)
def test_atom_positions(self):
result = self.results["atom_positions"]
expected_start = convert_unit(
......@@ -723,43 +791,99 @@ class TestMD(unittest.TestCase):
def test_frame_sequence_potential_energy(self):
result = self.results["frame_sequence_potential_energy"]