Commit cbed8fdb authored by Lauri Himanen's avatar Lauri Himanen

More general information parsed from geometry optimization.

parent 52e6a24d
import ase.io
import logging
logger = logging.getLogger("nomad")
#===============================================================================
def iread(filename, index=None, file_format=None):
"""
"""
try:
generator = ase.io.iread(filename, index, file_format)
except ValueError:
logger.error("ASE could not read the file '{}'.".format(filename))
raise
# Return the positions in a numpy array instead of the atoms object
for atoms in generator:
pos = atoms.positions
yield pos
......@@ -423,3 +423,15 @@ class CommonMatcher(object):
for attr, callback in extractOnCloseTriggers(self).items():
onClose[attr] = [callback]
return onClose
#===========================================================================
def onclosecatcher(func):
"""Used to catch exceptions in onClose
"""
def dec():
try:
func()
except Exception, e:
print 'Decorator handled exception %s' % e
raise e
return dec
from nomadcore.simple_parser import SimpleMatcher as SM
from nomadcore.baseclasses import MainHierarchicalParser
from commonmatcher import CommonMatcher
from cp2kparser.generic.configurationreading import iread
from nomadcore.caching_backend import CachingLevel
import logging
import ase.io
......@@ -28,21 +29,26 @@ class CP2KGeoOptParser(MainHierarchicalParser):
# Cache levels
self.caching_level_for_metaname.update({
'x_cp2k_optimization_energy': CachingLevel.ForwardAndCache,
'x_cp2k_optimization_step_size_convergence_limit': CachingLevel.ForwardAndCache,
'x_cp2k_section_geometry_optimization_information': CachingLevel.ForwardAndCache,
})
#=======================================================================
# SimpleMatchers
self.geo_opt = SM(
" *** STARTING GEOMETRY OPTIMIZATION ***".replace("*", "\*"),
sections=["section_frame_sequence", "section_sampling_method"],
sections=["section_frame_sequence"],
subMatchers=[
SM( " -------- Informations at step =\s+{}\s+------------".format(self.cm.regex_i),
SM( " *** CONJUGATE GRADIENTS ***".replace("*", "\*"),
adHoc=self.adHoc_conjugate_gradient()
),
SM( " -------- Informations at step",
forwardMatch=True,
name="geooptstep",
repeats=True,
sections=["section_single_configuration_calculation", "section_system"],
subMatchers=[
SM( " -------- Informations at step =\s+{}\s+------------".format(self.cm.regex_i),
SM( " -------- Informations at step",
sections=["x_cp2k_section_geometry_optimization_information"],
subMatchers=[
SM( " Optimization Method =\s+(?P<x_cp2k_optimization_method>{})".format(self.cm.regex_word)),
......@@ -51,10 +57,12 @@ class CP2KGeoOptParser(MainHierarchicalParser):
SM( " Decrease in energy =\s+(?P<x_cp2k_optimization_energy_decrease>{})".format(self.cm.regex_word)),
SM( " Used time =\s+(?P<x_cp2k_optimization_used_time>{})".format(self.cm.regex_f)),
SM( " Max. step size =\s+(?P<x_cp2k_optimization_max_step_size__bohr>{})".format(self.cm.regex_f)),
SM( " Conv. limit for step size =\s+(?P<x_cp2k_optimization_step_size_convergence_limit__bohr>{})".format(self.cm.regex_f)),
SM( " Convergence in step size =\s+(?P<x_cp2k_optimization_step_size_convergence>{})".format(self.cm.regex_word)),
SM( " RMS step size =\s+(?P<x_cp2k_optimization_rms_step_size__bohr>{})".format(self.cm.regex_f)),
SM( " Convergence in RMS step =\s+(?P<x_cp2k_optimization_rms_step_size_convergence>{})".format(self.cm.regex_word)),
SM( " Max. gradient =\s+(?P<x_cp2k_optimization_max_gradient__bohr_1hartree>{})".format(self.cm.regex_f)),
SM( " Conv. limit for gradients =\s+(?P<x_cp2k_optimization_gradient_convergence_limit__bohr_1hartree>{})".format(self.cm.regex_f)),
SM( " Conv. for gradients =\s+(?P<x_cp2k_optimization_max_gradient_convergence>{})".format(self.cm.regex_word)),
SM( " RMS gradient =\s+(?P<x_cp2k_optimization_rms_gradient__bohr_1hartree>{})".format(self.cm.regex_f)),
SM( " Conv. in RMS gradients =\s+(?P<x_cp2k_optimization_rms_gradient_convergence>{})".format(self.cm.regex_word)),
......@@ -73,7 +81,7 @@ class CP2KGeoOptParser(MainHierarchicalParser):
# computational time.
self.root_matcher = SM("",
forwardMatch=True,
sections=["section_run"],
sections=["section_run", "section_sampling_method"],
subMatchers=[
SM( "",
forwardMatch=True,
......@@ -93,18 +101,31 @@ class CP2KGeoOptParser(MainHierarchicalParser):
self.cache_service.push_value("number_of_frames_in_sequence")
self.cache_service.push_array_values("frame_sequence_potential_energy")
opt_section = section["x_cp2k_section_geometry_optimization_information"]
if opt_section is not None:
opt_section = opt_section[-1]
geo_limit = opt_section["x_cp2k_optimization_step_size_convergence_limit"]
if geo_limit is not None:
self.backend.addValue("geometry_optimization_geometry_change", geo_limit[0])
force_limit = opt_section["x_cp2k_optimization_gradient_convergence_limit"]
if force_limit is not None:
self.backend.addValue("geometry_optimization_threshold_force", force_limit[0])
def onClose_section_sampling_method(self, backend, gIndex, section):
self.backend.addValue("sampling_method", "geometry_optimization")
def onClose_x_cp2k_section_geometry_optimization_information(self, backend, gIndex, section):
energy = section["x_cp2k_optimization_energy"][0]
self.cache_service["frame_sequence_potential_energy"].append(energy)
def onClose_section_method(self, backend, gIndex, section):
traj_file = self.file_service.get_file_by_id("trajectory")
try:
if traj_file is not None:
self.traj_iterator = ase.io.iread(traj_file)
except ValueError:
# The format was not supported by ase
pass
if traj_file is not None:
try:
self.traj_iterator = iread(traj_file)
except ValueError:
# The format was not supported by ase
pass
#===========================================================================
# adHoc functions
......@@ -122,6 +143,13 @@ class CP2KGeoOptParser(MainHierarchicalParser):
parser.backend.addValue("geometry_optimization_converged", False)
return wrapper
def adHoc_conjugate_gradient(self):
"""Called when conjugate gradient method is used.
"""
def wrapper(parser):
parser.backend.addValue("geometry_optimization_method", "conjugate_gradient")
return wrapper
def adHoc_step(self):
"""Called when all the step information has been retrieved from the
output file. Here further information is gathered from external files.
......@@ -131,8 +159,7 @@ class CP2KGeoOptParser(MainHierarchicalParser):
# Get the next position from the trajectory file
if self.traj_iterator is not None:
atoms = next(self.traj_iterator)
pos = atoms.positions
pos = next(self.traj_iterator)
self.cache_service["atom_positions"] = pos
return wrapper
......
3
i = 1, E = -17.1534159246
O 12.2353220000 1.3766420000 10.8698800000
H 12.4175775999 2.2362362573 11.2616216864
H 11.9271436933 1.5723516602 10.0115134757
3
i = 2, E = -17.1941015290
O 12.2353220000 1.3766420000 10.8698800000
H 12.4618486015 2.2314871691 11.3335607388
H 11.9990227122 1.5776813026 10.0384213366
3
i = 3, E = -17.2092321965
O 12.2353220000 1.3766420000 10.8698800000
H 12.4884579626 2.2312983427 11.3362172362
H 11.9955329604 1.5759482262 10.0026671854
3
i = 4, E = -17.2097667733
O 12.2353220000 1.3766420000 10.8698800000
H 12.4962705528 2.2308411983 11.3355758433
H 11.9975151486 1.5746309898 10.0054430868
3
i = 5, E = -17.2097743028
O 12.2353220000 1.3766420000 10.8698800000
H 12.4957759240 2.2307207009 11.3354305795
H 11.9975745036 1.5747800582 10.0062718875
3
i = 6, E = -17.2097743229
O 12.2353220000 1.3766420000 10.8698800000
H 12.4958168364 2.2307249171 11.3354322532
H 11.9975556812 1.5748088251 10.0062793864
3
i = 7, E = -17.2097782066
O 12.2353220000 1.3766420000 10.8698800000
H 12.4958164689 2.2307248873 11.3354322515
H 11.9975558616 1.5748085240 10.0062792262
&GLOBAL
PROJECT H2O
RUN_TYPE GEO_OPT
PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
METHOD QS
&SUBSYS
&CELL
ABC 12.4138 12.4138 12.4138
&END CELL
&COORD
O 12.235322 1.376642 10.869880
H 12.415139 2.233125 11.257611
H 11.922476 1.573799 9.986994
&END COORD
&KIND H
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q1
&END KIND
&KIND O
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q6
&END KIND
&END SUBSYS
&DFT
BASIS_SET_FILE_NAME ../../BASIS_SET
POTENTIAL_FILE_NAME ../../GTH_POTENTIALS
&QS
EPS_DEFAULT 1.0E-7
&END QS
&MGRID
CUTOFF 100
NGRIDS 4
REL_CUTOFF 30
&END MGRID
&SCF
SCF_GUESS ATOMIC
EPS_SCF 1.0E-05
MAX_SCF 200
&DIAGONALIZATION T
ALGORITHM STANDARD
&END DIAGONALIZATION
&MIXING T
ALPHA 0.5
METHOD PULAY_MIXING
NPULAY 5
&END MIXING
&PRINT
&RESTART OFF
&END RESTART
&END PRINT
&END SCF
&XC
&XC_FUNCTIONAL PADE
&END XC_FUNCTIONAL
&END XC
&END DFT
&END FORCE_EVAL
&MOTION
&GEO_OPT
TYPE MINIMIZATION
MAX_DR 1.0E-03
MAX_FORCE 1.0E-03
RMS_DR 1.0E-03
RMS_FORCE 1.0E-03
MAX_ITER 200
OPTIMIZER CG
&CG
MAX_STEEP_STEPS 0
RESTART_LIMIT 9.0E-01
&END CG
&END GEO_OPT
&CONSTRAINT
&FIXED_ATOMS
COMPONENTS_TO_FIX XYZ
LIST 1
&END FIXED_ATOMS
&END CONSTRAINT
&END MOTION
This diff is collapsed.
3
i = 1, E = -17.1643447508
O 12.2353220000 1.3766420000 10.8698800000
H 12.4189404412 2.2491702452 11.2652097702
H 11.9149490419 1.5748472981 9.9693615826
3
i = 2, E = -17.1643578234
O 12.2353220000 1.3766420000 10.8698800000
H 12.4200101767 2.2468980462 11.2643381000
H 11.9115895109 1.5703371433 9.9671556911
3
i = 3, E = -17.1643766833
O 12.2353220000 1.3766420000 10.8698800000
H 12.4233662348 2.2477209582 11.2631307672
H 11.9069167231 1.5604096790 9.9689943947
3
i = 4, E = -17.1644565615
O 12.2353220000 1.3766420000 10.8698800000
H 12.4423251455 2.2644837735 11.2309777118
H 11.8922856584 1.5514786478 9.9754905025
3
i = 5, E = -17.1645203350
O 12.2353220000 1.3766420000 10.8698800000
H 12.4431348352 2.2577511412 11.2267313861
H 11.8919233389 1.5518052357 9.9721549140
3
i = 6, E = -17.1646190312
O 12.2353220000 1.3766420000 10.8698800000
H 12.4736086753 2.2567157451 11.2132013565
H 11.9047281262 1.5469224568 9.9655629396
3
i = 7, E = -17.1646198446
O 12.2353220000 1.3766420000 10.8698800000
H 12.4742379375 2.2563137117 11.2125606428
H 11.9051180271 1.5469372555 9.9662564829
3
i = 8, E = -17.1646203294
O 12.2353220000 1.3766420000 10.8698800000
H 12.4749724176 2.2565225180 11.2119388324
H 11.9039582865 1.5481084051 9.9667841313
3
i = 9, E = -17.1646204547
O 12.2353220000 1.3766420000 10.8698800000
H 12.4752917518 2.2561254118 11.2121163217
H 11.9031806595 1.5481019598 9.9671120050
3
i = 10, E = -17.1646204237
O 12.2353220000 1.3766420000 10.8698800000
H 12.4755864059 2.2562726297 11.2125118440
H 11.9026358443 1.5480115451 9.9673498890
3
i = 11, E = -17.1646204766
O 12.2353220000 1.3766420000 10.8698800000
H 12.4754906682 2.2560912723 11.2123985687
H 11.9026556670 1.5480524898 9.9673305122
3
i = 12, E = -17.1646347711
O 12.2353220000 1.3766420000 10.8698800000
H 12.4754916182 2.2560930719 11.2123996927
H 11.9026554703 1.5480520836 9.9673307045
&GLOBAL
PROJECT H2O
RUN_TYPE GEO_OPT
PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
METHOD QS
&SUBSYS
&CELL
ABC 12.4138 12.4138 12.4138
&END CELL
&COORD
O 12.235322 1.376642 10.869880
H 12.415139 2.233125 11.257611
H 11.922476 1.573799 9.986994
&END COORD
&KIND H
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q1
&END KIND
&KIND O
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q6
&END KIND
&END SUBSYS
&DFT
BASIS_SET_FILE_NAME ../../BASIS_SET
POTENTIAL_FILE_NAME ../../GTH_POTENTIALS
&QS
EPS_DEFAULT 1.0E-7
&END QS
&MGRID
CUTOFF 200
NGRIDS 4
REL_CUTOFF 30
&END MGRID
&SCF
SCF_GUESS ATOMIC
EPS_SCF 1.0E-05
MAX_SCF 200
&DIAGONALIZATION T
ALGORITHM STANDARD
&END DIAGONALIZATION
&MIXING T
ALPHA 0.5
METHOD PULAY_MIXING
NPULAY 5
&END MIXING
&PRINT
&RESTART OFF
&END RESTART
&END PRINT
&END SCF
&XC
&XC_FUNCTIONAL PADE
&END XC_FUNCTIONAL
&END XC
&END DFT
&END FORCE_EVAL
&MOTION
&GEO_OPT
TYPE MINIMIZATION
MAX_DR 1.0E-03
MAX_FORCE 1.0E-03
RMS_DR 1.0E-03
RMS_FORCE 1.0E-03
MAX_ITER 200
OPTIMIZER CG
&CG
MAX_STEEP_STEPS 0
RESTART_LIMIT 9.0E-01
&END CG
&END GEO_OPT
&CONSTRAINT
&FIXED_ATOMS
COMPONENTS_TO_FIX XYZ
LIST 1
&END FIXED_ATOMS
&END CONSTRAINT
&PRINT
&TRAJECTORY
FILENAME customname
&END TRAJECTORY
&END PRINT
&END MOTION
This diff is collapsed.
......@@ -458,6 +458,46 @@ class TestGeoOpt(unittest.TestCase):
result = self.results["number_of_frames_in_sequence"]
self.assertEqual(result, 7)
def test_sampling_method(self):
result = self.results["sampling_method"]
self.assertEqual(result, "geometry_optimization")
def test_geometry_optimization_method(self):
result = self.results["geometry_optimization_method"]
self.assertEqual(result, "conjugate_gradient")
def test_geometry_optimization_geometry_change(self):
result = self.results["geometry_optimization_geometry_change"]
expected_result = convert_unit(
0.0010000000,
"bohr"
)
self.assertEqual(result, expected_result)
def test_geometry_optimization_threshold_force(self):
result = self.results["geometry_optimization_threshold_force"]
expected_result = convert_unit(
0.0010000000,
"bohr^-1*hartree"
)
self.assertEqual(result, expected_result)
def test_frame_sequence_potential_energy(self):
result = self.results["frame_sequence_potential_energy"]
expected_result = convert_unit(
np.array([
-17.1488237635,
-17.1534159246,
-17.1941015290,
-17.2092321965,
-17.2097667733,
-17.2097743028,
-17.2097743229,
]),
"hartree"
)
self.assertTrue(np.array_equal(result, expected_result))
def test_atom_positions(self):
result = self.results["atom_positions"]
expected_start = convert_unit(
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment