diff --git a/TCutility/results/adf.py b/TCutility/results/adf.py index b3cd7ee9..c687df2c 100644 --- a/TCutility/results/adf.py +++ b/TCutility/results/adf.py @@ -22,14 +22,21 @@ def get_calc_settings(info: Result) -> Result: ''' assert info.engine == 'adf', f'This function reads ADF data, not {info.engine} data' - assert 'adf.rkf' in info.files, f'Missing adf.rkf file, [{", ".join([": ".join(item) for item in info.files.items()])}]' + # assert 'adf.rkf' in info.files, f'Missing adf.rkf file, [{", ".join([": ".join(item) for item in info.files.items()])}]' - reader_adf = cache.get(info.files['adf.rkf']) ret = Result() # set the calculation task at a higher level ret.task = info.input.task + # the VibrationalAnalysis task does not produce adf.rkf files + # in that case we end the reading here, since the following + # variables do not apply + if ret.task.lower() == 'vibrationalanalysis': + return ret + + reader_adf = cache.get(info.files['adf.rkf']) + # determine if calculation used relativistic corrections # if it did, variable 'escale' will be present in 'SFOs' # if it didnt, only variable 'energy' will be present @@ -79,12 +86,32 @@ def get_properties(info: Result) -> Result: - **vdd.charges (list[float])** - list of Voronoi Deformation Denisty (VDD) charges in [electrons], being the difference between the final (SCF) and initial VDD charges. ''' + def read_vibrations(reader: cache.TrackKFReader) -> Result: + ret = Result() + ret.number_of_modes = reader.read('Vibrations', 'nNormalModes') + ret.frequencies = ensure_list(reader.read('Vibrations', 'Frequencies[cm-1]')) + if ('Vibrations', 'Intensities[km/mol]') in reader: + ret.intensities = ensure_list(reader.read('Vibrations', 'Intensities[km/mol]')) + ret.number_of_imag_modes = len([freq for freq in ret.frequencies if freq < 0]) + ret.character = 'minimum' if ret.number_of_imag_modes == 0 else 'transitionstate' + ret.modes = [] + for i in range(ret.number_of_modes): + ret.modes.append(reader.read('Vibrations', f'NoWeightNormalMode({i+1})')) + return ret + + assert info.engine == 'adf', f'This function reads ADF data, not {info.engine} data' - assert 'adf.rkf' in info.files, f'Missing adf.rkf file, [{", ".join([": ".join(item) for item in info.files.items()])}]' - reader_adf = cache.get(info.files['adf.rkf']) + ret = Result() + if info.adf.task.lower() == 'vibrationalanalysis': + reader_ams = cache.get(info.files['ams.rkf']) + ret.vibrations = read_vibrations(reader_ams) + return ret + + reader_adf = cache.get(info.files['adf.rkf']) + # read energies (given in Ha in rkf files) ret.energy.bond = reader_adf.read('Energy', 'Bond Energy') * constants.HA2KCALMOL ret.energy.elstat = reader_adf.read('Energy', 'elstat') * constants.HA2KCALMOL @@ -97,15 +124,7 @@ def get_properties(info: Result) -> Result: # vibrational information if ('Vibrations', 'nNormalModes') in reader_adf: - ret.vibrations.number_of_modes = reader_adf.read('Vibrations', 'nNormalModes') - freqs = reader_adf.read('Vibrations', 'Frequencies[cm-1]') - ints = reader_adf.read('Vibrations', 'Intensities[km/mol]') - ret.vibrations.frequencies = freqs if isinstance(freqs, list) else [freqs] - ret.vibrations.intensities = ints if isinstance(ints, list) else [ints] - ret.vibrations.number_of_imag_modes = len([freq for freq in freqs if freq < 0]) - ret.vibrations.modes = [] - for i in range(ret.vibrations.number_of_modes): - ret.vibrations.modes.append(reader_adf.read('Vibrations', f'NoWeightNormalMode({i+1})')) + ret.vibrations = read_vibrations(reader_adf) # read the Voronoi Deformation Charges Deformation (VDD) before and after SCF convergence (being "inital" and "SCF") vdd_scf: List[float] = ensure_list(reader_adf.read('Properties', 'AtomCharge_SCF Voronoi')) # type: ignore since plams does not include typing for KFReader. List[float] is returned diff --git a/TCutility/results/dftb.py b/TCutility/results/dftb.py index 1b189932..a2f2bdc6 100644 --- a/TCutility/results/dftb.py +++ b/TCutility/results/dftb.py @@ -1,4 +1,5 @@ from TCutility.results import cache, Result +from TCutility import ensure_list def get_calc_settings(info: Result) -> Result: @@ -9,19 +10,11 @@ def get_calc_settings(info: Result) -> Result: assert info.engine == 'dftb', f'This function reads DFTB data, not {info.engine} data' assert 'dftb.rkf' in info.files, f'Missing dftb.rkf file, [{", ".join([": ".join(item) for item in info.files.items()])}]' - reader_dftb = cache.get(info.files['dftb.rkf']) ret = Result() # set the calculation task at a higher level ret.task = info.input.task - # read properties of the calculation - number_of_properties = int(reader_dftb.read('Properties', 'nEntries')) # type: ignore plams does not include type hints. Returns int - for i in range(1, number_of_properties + 1): - prop_type = str(reader_dftb.read('Properties', f'Type({i})')).strip() - prop_subtype = str(reader_dftb.read('Properties', f'Subtype({i})')).strip() - ret.properties[prop_type][prop_subtype] = reader_dftb.read('Properties', f'Value({i})') - return ret @@ -46,13 +39,30 @@ def get_properties(info: Result) -> Result: # properties are not given as in ADF. Instead you must first know which index each key is at. For example we might have: # (Properties, SubType(k)) == Coulomb # therefore, to get the Coulomb energy we first have to know the index k - # we therefore first loop through every property entry and store the properties in order in a list + # we therefore first loop through every subtype entry and store the substypes in order in a list + # we also store the main type (e.g. 'Energy', 'Gradient', ...) and the value number_of_properties = reader_dftb.read('Properties', 'nEntries') - properties = [] + subtypes = [] + types = [] + values = [] for i in range(1, number_of_properties+1): - properties.append(reader_dftb.read('Properties', f'Subtype({i})').strip()) - print(properties) + subtypes.append(reader_dftb.read('Properties', f'Subtype({i})').strip()) + types.append(reader_dftb.read('Properties', f'Type({i})').strip()) + values.append(reader_dftb.read('Properties', f'Value({i})')) + + # then simply add the properties to ret + for typ, subtyp, value in zip(types, subtypes, values): + ret[typ.replace(' ', '_')][subtyp] = value + + # we also read vibrations + ret.vibrations.number_of_modes = reader_dftb.read('Vibrations', 'nNormalModes') + ret.vibrations.frequencies = ensure_list(reader_dftb.read('Vibrations', 'Frequencies[cm-1]')) + if ('Vibrations', 'Intensities[km/mol]') in reader_dftb: + ret.vibrations.intensities = ensure_list(reader_dftb.read('Vibrations', 'Intensities[km/mol]')) + ret.vibrations.number_of_imag_modes = len([freq for freq in ret.vibrations.frequencies if freq < 0]) + ret.vibrations.character = 'minimum' if ret.vibrations.number_of_imag_modes == 0 else 'transitionstate' + ret.vibrations.modes = [] + for i in range(ret.vibrations.number_of_modes): + ret.vibrations.modes.append(reader_dftb.read('Vibrations', f'NoWeightNormalMode({i+1})')) - # read energies (given in Ha i return ret -