Source code for calorine.gpumd.io

from typing import List

import numpy as np
from ase import Atoms
from ase.io import read, write
from pandas import DataFrame


[docs]def read_kappa(filename: str) -> DataFrame: """Parses a file in ``kappa.out`` format from GPUMD and returns the content as a data frame. More information concerning file format, content and units can be found `here <https://gpumd.org/gpumd/output_files/kappa_out.html>`__. Parameters ---------- filename input file name """ data = np.loadtxt(filename) if isinstance(data[0], np.float64): # If only a single row in kappa.out, append a dimension data = data.reshape(1, -1) tags = 'kx_in kx_out ky_in ky_out kz_tot'.split() if len(data[0]) != len(tags): raise ValueError(f'Input file contains {len(data[0])} data columns.' f' Expected {len(tags)} columns.') df = DataFrame(data=data, columns=tags) df['kx_tot'] = df.kx_in + df.kx_out df['ky_tot'] = df.ky_in + df.ky_out return df
[docs]def read_hac(filename: str) -> DataFrame: """Parses a file in ``hac.out`` format from GPUMD and returns the content as a data frame. More information concerning file format, content and units can be found `here <https://gpumd.org/gpumd/output_files/hac_out.html>`__. Parameters ---------- filename input file name """ data = np.loadtxt(filename) if isinstance(data[0], np.float64): # If only a single row in hac.out, append a dimension data = data.reshape(1, -1) tags = 'time' tags += ' jin_jtot_x jout_jtot_x jin_jtot_y jout_jtot_y jtot_jtot_z' tags += ' kx_in kx_out ky_in ky_out kz_tot' tags = tags.split() if len(data[0]) != len(tags): raise ValueError(f'Input file contains {len(data[0])} data columns.' f' Expected {len(tags)} columns.') df = DataFrame(data=data, columns=tags) df['kx_tot'] = df.kx_in + df.kx_out df['ky_tot'] = df.ky_in + df.ky_out # remove columns with less relevant data to save space for col in df: if 'jtot' in col or '_in' in col: del df[col] return df
[docs]def read_thermo(filename: str, natoms: int = 1) -> DataFrame: """Parses a file in ``thermo.out`` format from GPUMD and returns the content as a data frame. More information concerning file format, content and units can be found `here <https://gpumd.org/gpumd/output_files/thermo_out.html>`__. Parameters ---------- filename input file name natoms number of atoms; used to normalize energies """ data = np.loadtxt(filename) if isinstance(data[0], np.float64): # If only a single row in loss.out, append a dimension data = data.reshape(1, -1) if len(data[0]) == 9: # orthorhombic box tags = 'temperature kinetic_energy potential_energy' tags += ' stress_xx stress_yy stress_zz' tags += ' cell_xx cell_yy cell_zz' elif len(data[0]) == 12: # orthorhombic box with stresses in Voigt notation (v3.3.1+) tags = 'temperature kinetic_energy potential_energy' tags += ' stress_xx stress_yy stress_zz stress_yz stress_xz stress_xy' tags += ' cell_xx cell_yy cell_zz' elif len(data[0]) == 15: # triclinic box tags = 'temperature kinetic_energy potential_energy' tags += ' stress_xx stress_yy stress_zz' tags += ' cell_xx cell_xy cell_xz cell_yx cell_yy cell_yz cell_zx cell_zy cell_zz' elif len(data[0]) == 18: # triclinic box with stresses in Voigt notation (v3.3.1+) tags = 'temperature kinetic_energy potential_energy' tags += ' stress_xx stress_yy stress_zz stress_yz stress_xz stress_xy' tags += ' cell_xx cell_xy cell_xz cell_yx cell_yy cell_yz cell_zx cell_zy cell_zz' else: raise ValueError(f'Input file contains {len(data[0])} data columns.' ' Expected 9, 12, 15 or 18 columns.') df = DataFrame(data=data, columns=tags.split()) assert natoms > 0, 'natoms must be positive' df.kinetic_energy /= natoms df.potential_energy /= natoms return df
[docs]def read_xyz(filename: str) -> Atoms: """ Read the structure input file (`model.xyz`) for GPUMD and return the structure along with run input parameters from the file. This is a wrapper function around :func:`ase.io.read_xyz` since the ASE implementation does not read velocities properly. Parameters ---------- filename Name of file from which to read the structure Returns ------- structure as ASE Atoms object with additional per-atom arrays representing atomic masses, velocities etc. """ structure = read(filename, format='extxyz') if structure.has('vel'): structure.set_velocities(structure.get_array('vel')) return structure
[docs]def write_xyz(filename: str, structure: Atoms, groupings: List[List[List[int]]] = None): """ Writes a structure into GPUMD input format (`model.xyz`). Parameters ---------- filename Name of file to which the structure should be written structure Input structure groupings Groups into which the individual atoms should be divided in the form of a list of list of lists. Specifically, the outer list corresponds to the grouping methods, of which there can be three at the most, which contains a list of groups in the form of lists of site indices. The sum of the lengths of the latter must be the same as the total number of atoms. Raises ------ ValueError Raised if parameters are incompatible """ # Make a local copy of the atoms object _structure = structure.copy() # Check velocties parameter velocities = _structure.get_velocities() if velocities is None or np.max(np.abs(velocities)) < 1e-6: has_velocity = 0 else: has_velocity = 1 # Check groupings parameter if groupings is None: number_of_grouping_methods = 0 else: number_of_grouping_methods = len(groupings) if number_of_grouping_methods > 3: raise ValueError('There can be no more than 3 grouping methods!') for g, grouping in enumerate(groupings): all_indices = [i for group in grouping for i in group] if len(all_indices) != len(_structure) or \ set(all_indices) != set(range(len(_structure))): raise ValueError(f'The indices listed in grouping method {g} are' ' not compatible with the input' ' structure!') # Allowed keyword=value pairs. Use ASEs extyz write functionality.: # pbc="pbc_a pbc_b pbc_c" # lattice="ax ay az bx by bz cx cy cz" # properties=property_name:data_type:number_of_columns # species:S:1 # pos:R:3 # mass:R:1 # vel:R:3 # group:I:number_of_grouping_methods _structure.new_array('mass', _structure.get_masses()) if has_velocity: _structure.new_array('vel', _structure.get_velocities()) if groupings is not None: group_indices = np.array([ [ [group_index for group_index, group in enumerate( grouping) if structure_idx in group] for grouping in groupings ] for structure_idx in range(len(_structure))]).squeeze() # pythoniccc _structure.new_array('group', group_indices) write(filename=filename, images=_structure, write_info=True, format='extxyz')