From fffabcb727ef8e3cae80068e8c28f69aa6229b10 Mon Sep 17 00:00:00 2001 From: abel Date: Thu, 2 Nov 2017 18:25:36 +0900 Subject: [PATCH] Initial commit --- LICENSE | 21 + plugins/__init__.py | 0 plugins/jobs/__init__.py | 0 plugins/jobs/dynaphopy.py | 284 ++++++++++++ plugins/jobs/lammps/__init__.py | 0 plugins/jobs/lammps/combinate.py | 433 ++++++++++++++++++ plugins/jobs/lammps/force.py | 229 +++++++++ plugins/jobs/lammps/md.py | 260 +++++++++++ plugins/jobs/lammps/optimize.py | 267 +++++++++++ plugins/jobs/lammps/potentials/__init__.py | 37 ++ plugins/jobs/lammps/potentials/eam.py | 17 + .../jobs/lammps/potentials/lennard_jones.py | 16 + plugins/jobs/lammps/potentials/tersoff.py | 16 + plugins/launcher/launch_dynaphopy_si.py | 82 ++++ plugins/launcher/launch_lammps_combinate.py | 110 +++++ plugins/launcher/launch_lammps_force_gan.py | 90 ++++ plugins/launcher/launch_lammps_md_si.py | 91 ++++ .../launcher/launch_lammps_optimization_fe.py | 83 ++++ .../launch_lammps_optimization_gan.py | 96 ++++ .../launcher/launch_lammps_optimization_lj.py | 83 ++++ plugins/launcher/launch_phonopy_gan.py | 67 +++ plugins/launcher/launch_vasp_si.py | 151 ++++++ plugins/parsers/__init__.py | 0 plugins/parsers/dynaphopy.py | 176 +++++++ plugins/parsers/lammps/__init__.py | 0 plugins/parsers/lammps/force.py | 177 +++++++ plugins/parsers/lammps/md.py | 209 +++++++++ plugins/parsers/lammps/optimize.py | 232 ++++++++++ plugins/parsers/phonopy.py | 90 ++++ plugins/parsers/phonopy_dos.py | 0 setup.py | 28 ++ 31 files changed, 3345 insertions(+) create mode 100644 LICENSE create mode 100644 plugins/__init__.py create mode 100644 plugins/jobs/__init__.py create mode 100644 plugins/jobs/dynaphopy.py create mode 100644 plugins/jobs/lammps/__init__.py create mode 100644 plugins/jobs/lammps/combinate.py create mode 100644 plugins/jobs/lammps/force.py create mode 100644 plugins/jobs/lammps/md.py create mode 100644 plugins/jobs/lammps/optimize.py create mode 100644 plugins/jobs/lammps/potentials/__init__.py create mode 100644 plugins/jobs/lammps/potentials/eam.py create mode 100644 plugins/jobs/lammps/potentials/lennard_jones.py create mode 100644 plugins/jobs/lammps/potentials/tersoff.py create mode 100644 plugins/launcher/launch_dynaphopy_si.py create mode 100644 plugins/launcher/launch_lammps_combinate.py create mode 100644 plugins/launcher/launch_lammps_force_gan.py create mode 100644 plugins/launcher/launch_lammps_md_si.py create mode 100644 plugins/launcher/launch_lammps_optimization_fe.py create mode 100644 plugins/launcher/launch_lammps_optimization_gan.py create mode 100644 plugins/launcher/launch_lammps_optimization_lj.py create mode 100644 plugins/launcher/launch_phonopy_gan.py create mode 100644 plugins/launcher/launch_vasp_si.py create mode 100644 plugins/parsers/__init__.py create mode 100644 plugins/parsers/dynaphopy.py create mode 100644 plugins/parsers/lammps/__init__.py create mode 100644 plugins/parsers/lammps/force.py create mode 100644 plugins/parsers/lammps/md.py create mode 100644 plugins/parsers/lammps/optimize.py create mode 100644 plugins/parsers/phonopy.py create mode 100644 plugins/parsers/phonopy_dos.py create mode 100644 setup.py diff --git a/LICENSE b/LICENSE new file mode 100644 index 0000000..2253d69 --- /dev/null +++ b/LICENSE @@ -0,0 +1,21 @@ +The MIT License (MIT) + +Copyright (c) 2017 Abel Carreras Conill + +Permission is hereby granted, free of charge, to any person obtaining a copy +of this software and associated documentation files (the "Software"), to deal +in the Software without restriction, including without limitation the rights +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the Software is +furnished to do so, subject to the following conditions: + +The above copyright notice and this permission notice shall be included in +all copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +THE SOFTWARE. diff --git a/plugins/__init__.py b/plugins/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/plugins/jobs/__init__.py b/plugins/jobs/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/plugins/jobs/dynaphopy.py b/plugins/jobs/dynaphopy.py new file mode 100644 index 0000000..92ff6f8 --- /dev/null +++ b/plugins/jobs/dynaphopy.py @@ -0,0 +1,284 @@ +from aiida.orm.calculation.job import JobCalculation +from aiida.orm.data.parameter import ParameterData +from aiida.orm.data.structure import StructureData +from aiida.orm.data.array.trajectory import TrajectoryData +from aiida.orm.data.array import ArrayData +from aiida.common.exceptions import InputValidationError +from aiida.common.datastructures import CalcInfo, CodeInfo +from aiida.common.utils import classproperty + +import numpy as np + + +def get_FORCE_CONSTANTS_txt(force_constants): + + force_constants = force_constants.get_array('force_constants') + + fc_shape = force_constants.shape + fc_txt = "%4d\n" % (fc_shape[0]) + for i in range(fc_shape[0]): + for j in range(fc_shape[1]): + fc_txt += "%4d%4d\n" % (i+1, j+1) + for vec in force_constants[i][j]: + fc_txt +=("%22.15f"*3 + "\n") % tuple(vec) + + return fc_txt + + +def get_trajectory_txt(trajectory): + + cell = trajectory.get_cells()[0] + + a = np.linalg.norm(cell[0]) + b = np.linalg.norm(cell[1]) + c = np.linalg.norm(cell[2]) + + alpha = np.arccos(np.dot(cell[1], cell[2])/(c*b)) + gamma = np.arccos(np.dot(cell[1], cell[0])/(a*b)) + beta = np.arccos(np.dot(cell[2], cell[0])/(a*c)) + + xhi = a + xy = b * np.cos(gamma) + xz = c * np.cos(beta) + yhi = np.sqrt(pow(b,2)- pow(xy,2)) + yz = (b*c*np.cos(alpha)-xy * xz)/yhi + zhi = np.sqrt(pow(c,2)-pow(xz,2)-pow(yz,2)) + + xhi = xhi + max(0,0, xy, xz, xy+xz) + yhi = yhi + max(0,0, yz) + + xlo_bound = np.min([0.0, xy, xz, xy+xz]) + xhi_bound = xhi + np.max([0.0, xy, xz, xy+xz]) + ylo_bound = np.min([0.0, yz]) + yhi_bound = yhi + np.max([0.0, yz]) + zlo_bound = 0 + zhi_bound = zhi + + ind = trajectory.get_array('steps') + lammps_data_file = '' + for i, position_step in enumerate(trajectory.get_positions()): + lammps_data_file += 'ITEM: TIMESTEP\n' + lammps_data_file += '{}\n'.format(ind[i]) + lammps_data_file += 'ITEM: NUMBER OF ATOMS\n' + lammps_data_file += '{}\n'.format(len(position_step)) + lammps_data_file += 'ITEM: BOX BOUNDS xy xz yz pp pp pp\n' + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f}\n'.format(xlo_bound, xhi_bound, xy) + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f}\n'.format(ylo_bound, yhi_bound, xz) + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f}\n'.format(zlo_bound, zhi_bound, yz) + lammps_data_file += ('ITEM: ATOMS x y z\n') + for position in position_step: + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f}\n'.format(*position) + return lammps_data_file + + +def structure_to_poscar(structure): + + types = [site.kind_name for site in structure.sites] + atom_type_unique = np.unique(types, return_index=True) + sort_index = np.argsort(atom_type_unique[1]) + elements = np.array(atom_type_unique[0])[sort_index] + elements_count= np.diff(np.append(np.array(atom_type_unique[1])[sort_index], [len(types)])) + + poscar = '# VASP POSCAR generated using aiida workflow ' + poscar += '\n1.0\n' + cell = structure.cell + for row in cell: + poscar += '{0: 22.16f} {1: 22.16f} {2: 22.16f}\n'.format(*row) + poscar += ' '.join([str(e) for e in elements]) + '\n' + poscar += ' '.join([str(e) for e in elements_count]) + '\n' + poscar += 'Cartesian\n' + for site in structure.sites: + poscar += '{0: 22.16f} {1: 22.16f} {2: 22.16f}\n'.format(*site.position) + + return poscar + + + +def parameters_to_input_file(parameters_object): + + parameters = parameters_object.get_dict() + input_file = ('STRUCTURE FILE POSCAR\nPOSCAR\n\n') + input_file += ('FORCE CONSTANTS\nFORCE_CONSTANTS\n\n') + input_file += ('PRIMITIVE MATRIX\n') + input_file += ('{} {} {} \n').format(*np.array(parameters['primitive'])[0]) + input_file += ('{} {} {} \n').format(*np.array(parameters['primitive'])[1]) + input_file += ('{} {} {} \n').format(*np.array(parameters['primitive'])[2]) + input_file += ('\n') + input_file += ('SUPERCELL MATRIX PHONOPY\n') + input_file += ('{} {} {} \n').format(*np.array(parameters['supercell'])[0]) + input_file += ('{} {} {} \n').format(*np.array(parameters['supercell'])[1]) + input_file += ('{} {} {} \n').format(*np.array(parameters['supercell'])[2]) + input_file += ('\n') + + return input_file + + +class DynaphopyCalculation(JobCalculation): + """ + A basic plugin for calculating force constants using Phonopy. + + Requirement: the node should be able to import phonopy + """ + + def _init_internal_params(self): + super(DynaphopyCalculation, self)._init_internal_params() + + self._INPUT_FILE_NAME = 'input_dynaphopy' + self._INPUT_TRAJECTORY = 'trajectory' + self._INPUT_CELL = 'POSCAR' + self._INPUT_FORCE_CONSTANTS = 'FORCE_CONSTANTS' + + self._OUTPUT_FORCE_CONSTANTS = 'FORCE_CONSTANTS_OUT' + self._OUTPUT_FILE_NAME = 'OUTPUT' + self._OUTPUT_QUASIPARTICLES = 'quasiparticles_data.yaml' + + self._default_parser = 'dynaphopy' + + + @classproperty + def _use_methods(cls): + """ + Additional use_* methods for the namelists class. + """ + retdict = JobCalculation._use_methods + retdict.update({ + "parameters": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'parameters', + 'docstring': ("Use a node that specifies the dynaphopy input " + "for the namelists"), + }, + "trajectory": { + 'valid_types': TrajectoryData, + 'additional_parameter': None, + 'linkname': 'trajectory', + 'docstring': ("Use a node that specifies the trajectory data " + "for the namelists"), + }, + "force_constants": { + 'valid_types': ArrayData, + 'additional_parameter': None, + 'linkname': 'force_constants', + 'docstring': ("Use a node that specifies the force_constants " + "for the namelists"), + }, + "structure": { + 'valid_types': StructureData, + 'additional_parameter': None, + 'linkname': 'structure', + 'docstring': "Use a node for the structure", + }, + }) + return retdict + + def _prepare_for_submission(self,tempfolder, inputdict): + """ + This is the routine to be called when you want to create + the input files and related stuff with a plugin. + + :param tempfolder: a aiida.common.folders.Folder subclass where + the plugin should put all its files. + :param inputdict: a dictionary with the input nodes, as they would + be returned by get_inputdata_dict (without the Code!) + """ + + try: + parameters_data = inputdict.pop(self.get_linkname('parameters')) + except KeyError: + pass + #raise InputValidationError("No parameters specified for this " + # "calculation") + if not isinstance(parameters_data, ParameterData): + raise InputValidationError("parameters is not of type " + "ParameterData") + + try: + structure = inputdict.pop(self.get_linkname('structure')) + except KeyError: + raise InputValidationError("no structure is specified for this calculation") + + try: + trajectory = inputdict.pop(self.get_linkname('trajectory')) + except KeyError: + raise InputValidationError("trajectory is specified for this calculation") + + try: + force_constants = inputdict.pop(self.get_linkname('force_constants')) + except KeyError: + raise InputValidationError("no force_constants is specified for this calculation") + + try: + code = inputdict.pop(self.get_linkname('code')) + except KeyError: + raise InputValidationError("no code is specified for this calculation") + + + time_step = trajectory.get_times()[1]-trajectory.get_times()[0] + + ############################## + # END OF INITIAL INPUT CHECK # + ############################## + + # =================== prepare the python input files ===================== + + cell_txt = structure_to_poscar(structure) + input_txt = parameters_to_input_file(parameters_data) + force_constants_txt = get_FORCE_CONSTANTS_txt(force_constants) + trajectory_txt = get_trajectory_txt(trajectory) + + # =========================== dump to file ============================= + + input_filename = tempfolder.get_abs_path(self._INPUT_FILE_NAME) + with open(input_filename, 'w') as infile: + infile.write(input_txt) + + cell_filename = tempfolder.get_abs_path(self._INPUT_CELL) + with open(cell_filename, 'w') as infile: + infile.write(cell_txt) + + force_constants_filename = tempfolder.get_abs_path(self._INPUT_FORCE_CONSTANTS) + with open(force_constants_filename, 'w') as infile: + infile.write(force_constants_txt) + + trajectory_filename = tempfolder.get_abs_path(self._INPUT_TRAJECTORY) + with open(trajectory_filename, 'w') as infile: + infile.write(trajectory_txt) + + # ============================ calcinfo ================================ + + local_copy_list = [] + remote_copy_list = [] + # additional_retrieve_list = settings_dict.pop("ADDITIONAL_RETRIEVE_LIST",[]) + + calcinfo = CalcInfo() + + calcinfo.uuid = self.uuid + # Empty command line by default + calcinfo.local_copy_list = local_copy_list + calcinfo.remote_copy_list = remote_copy_list + + # Retrieve files + calcinfo.retrieve_list = [self._OUTPUT_FILE_NAME, + self._OUTPUT_FORCE_CONSTANTS, + self._OUTPUT_QUASIPARTICLES] + + codeinfo = CodeInfo() + codeinfo.cmdline_params = [self._INPUT_FILE_NAME, self._INPUT_TRAJECTORY, + '-ts', '{}'.format(time_step), '--silent', + '-sfc', self._OUTPUT_FORCE_CONSTANTS, '-thm', # '--resolution 0.01', + '-psm','2', '--normalize_dos', '-sdata'] + + if 'temperature' in parameters_data.get_dict(): + codeinfo.cmdline_params.append('--temperature') + codeinfo.cmdline_params.append('{}'.format(parameters_data.dict.temperature)) + + if 'md_commensurate' in parameters_data.get_dict(): + if parameters_data.dict.md_commensurate: + codeinfo.cmdline_params.append('--MD_commensurate') + + codeinfo.stdout_name = self._OUTPUT_FILE_NAME + codeinfo.code_uuid = code.uuid + codeinfo.withmpi = False + calcinfo.codes_info = [codeinfo] + return calcinfo diff --git a/plugins/jobs/lammps/__init__.py b/plugins/jobs/lammps/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/plugins/jobs/lammps/combinate.py b/plugins/jobs/lammps/combinate.py new file mode 100644 index 0000000..2f49819 --- /dev/null +++ b/plugins/jobs/lammps/combinate.py @@ -0,0 +1,433 @@ +from aiida.orm.calculation.job import JobCalculation +from aiida.orm.data.parameter import ParameterData +from aiida.orm.data.structure import StructureData +from aiida.orm.data.array import ArrayData +from aiida.common.exceptions import InputValidationError +from aiida.common.datastructures import CalcInfo, CodeInfo +from aiida.common.utils import classproperty + +from potentials import LammpsPotential +import numpy as np + + +def get_supercell(structure, supercell_shape): + import itertools + + symbols = np.array([site.kind_name for site in structure.sites]) + positions = np.array([site.position for site in structure.sites]) + cell = np.array(structure.cell) + supercell_shape = np.array(supercell_shape.dict.shape) + + supercell_array = np.dot(cell, np.diag(supercell_shape)) + + supercell = StructureData(cell=supercell_array) + for k in range(positions.shape[0]): + for r in itertools.product(*[range(i) for i in supercell_shape[::-1]]): + position = positions[k, :] + np.dot(np.array(r[::-1]), cell) + symbol = symbols[k] + supercell.append_atom(position=position, symbols=symbol) + + return supercell + + +def get_FORCE_CONSTANTS_txt(force_constants): + + force_constants = force_constants.get_array('force_constants') + + fc_shape = force_constants.shape + fc_txt = "%4d\n" % (fc_shape[0]) + for i in range(fc_shape[0]): + for j in range(fc_shape[1]): + fc_txt += "%4d%4d\n" % (i+1, j+1) + for vec in force_constants[i][j]: + fc_txt +=("%22.15f"*3 + "\n") % tuple(vec) + + return fc_txt + + +def structure_to_poscar(structure): + + atom_type_unique = np.unique([site.kind_name for site in structure.sites], return_index=True)[1] + labels = np.diff(np.append(atom_type_unique, [len(structure.sites)])) + + poscar = ' '.join(np.unique([site.kind_name for site in structure.sites])) + poscar += '\n1.0\n' + cell = structure.cell + for row in cell: + poscar += '{0: 22.16f} {1: 22.16f} {2: 22.16f}\n'.format(*row) + poscar += ' '.join(np.unique([site.kind_name for site in structure.sites]))+'\n' + poscar += ' '.join(np.array(labels, dtype=str))+'\n' + poscar += 'Cartesian\n' + for site in structure.sites: + poscar += '{0: 22.16f} {1: 22.16f} {2: 22.16f}\n'.format(*site.position) + + return poscar + + +def parameters_to_input_file(parameters_object): + + parameters = parameters_object.get_dict() + input_file = ('STRUCTURE FILE POSCAR\nPOSCAR\n\n') + input_file += ('FORCE CONSTANTS\nFORCE_CONSTANTS\n\n') + input_file += ('PRIMITIVE MATRIX\n') + input_file += ('{} {} {} \n').format(*np.array(parameters['primitive'])[0]) + input_file += ('{} {} {} \n').format(*np.array(parameters['primitive'])[1]) + input_file += ('{} {} {} \n').format(*np.array(parameters['primitive'])[2]) + input_file += ('\n') + input_file += ('SUPERCELL MATRIX PHONOPY\n') + input_file += ('{} {} {} \n').format(*np.array(parameters['supercell'])[0]) + input_file += ('{} {} {} \n').format(*np.array(parameters['supercell'])[1]) + input_file += ('{} {} {} \n').format(*np.array(parameters['supercell'])[2]) + input_file += ('\n') + + return input_file + + +def generate_LAMMPS_structure(structure): + import numpy as np + + types = [site.kind_name for site in structure.sites] + + type_index_unique = np.unique(types, return_index=True)[1] + count_index_unique = np.diff(np.append(type_index_unique, [len(types)])) + + atom_index = [] + for i, index in enumerate(count_index_unique): + atom_index += [i for j in range(index)] + + masses = [site.mass for site in structure.kinds] + positions = [site.position for site in structure.sites] + + number_of_atoms = len(positions) + + lammps_data_file = 'Generated using dynaphopy\n\n' + lammps_data_file += '{0} atoms\n\n'.format(number_of_atoms) + lammps_data_file += '{0} atom types\n\n'.format(len(masses)) + + cell = np.array(structure.cell) + + a = np.linalg.norm(cell[0]) + b = np.linalg.norm(cell[1]) + c = np.linalg.norm(cell[2]) + + alpha = np.arccos(np.dot(cell[1], cell[2])/(c*b)) + gamma = np.arccos(np.dot(cell[1], cell[0])/(a*b)) + beta = np.arccos(np.dot(cell[2], cell[0])/(a*c)) + + xhi = a + xy = b * np.cos(gamma) + xz = c * np.cos(beta) + yhi = np.sqrt(pow(b,2)- pow(xy,2)) + yz = (b*c*np.cos(alpha)-xy * xz)/yhi + zhi = np.sqrt(pow(c,2)-pow(xz,2)-pow(yz,2)) + + xhi = xhi + max(0,0, xy, xz, xy+xz) + yhi = yhi + max(0,0, yz) + + lammps_data_file += '\n{0:20.10f} {1:20.10f} xlo xhi\n'.format(0, xhi) + lammps_data_file += '{0:20.10f} {1:20.10f} ylo yhi\n'.format(0, yhi) + lammps_data_file += '{0:20.10f} {1:20.10f} zlo zhi\n'.format(0, zhi) + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f} xy xz yz\n\n'.format(xy, xz, yz) + + lammps_data_file += 'Masses\n\n' + + for i, mass in enumerate(masses): + lammps_data_file += '{0} {1:20.10f} \n'.format(i+1, mass) + + + lammps_data_file += '\nAtoms\n\n' + for i, row in enumerate(positions): + lammps_data_file += '{0} {1} {2:20.10f} {3:20.10f} {4:20.10f}\n'.format(i+1, atom_index[i]+1, row[0],row[1],row[2]) + + return lammps_data_file + + +def generate_LAMMPS_input(parameters, + potential_obj, + structure_file='potential.pot', + trajectory_file='trajectory.lammpstr', + command=None): + + random_number = np.random.randint(10000000) + + names_str = ' '.join(potential_obj._names) + + lammps_input_file = 'units metal\n' + lammps_input_file += 'boundary p p p\n' + lammps_input_file += 'box tilt large\n' + lammps_input_file += 'atom_style atomic\n' + lammps_input_file += 'read_data {}\n'.format(structure_file) + + lammps_input_file += potential_obj.get_input_potential_lines() + + lammps_input_file += 'neighbor 0.3 bin\n' + lammps_input_file += 'neigh_modify every 1 delay 0 check no\n' + + lammps_input_file += 'timestep {}\n'.format(parameters.dict.timestep) + + lammps_input_file += 'thermo_style custom step etotal temp vol press\n' + lammps_input_file += 'thermo 100000\n' + + lammps_input_file += 'velocity all create {0} {1} dist gaussian mom yes\n'.format(parameters.dict.temperature, random_number) + lammps_input_file += 'velocity all scale {}\n'.format(parameters.dict.temperature) + + lammps_input_file += 'fix int all nvt temp {0} {0} {1}\n'.format(parameters.dict.temperature, parameters.dict.thermostat_variable) + + lammps_input_file += 'run {}\n'.format(parameters.dict.equilibrium_steps) + lammps_input_file += 'reset_timestep 0\n' + + lammps_input_file += 'dump aiida all custom {0} {1} x y z\n'.format(parameters.dict.dump_rate, trajectory_file) + lammps_input_file += 'dump_modify aiida format "%16.10f %16.10f %16.10f"\n' + lammps_input_file += 'dump_modify aiida sort id\n' + lammps_input_file += 'dump_modify aiida element {}\n'.format(names_str) + + lammps_input_file += 'run {}\n'.format(parameters.dict.total_steps) + + if command: + lammps_input_file += 'shell {}\n'.format(command) + lammps_input_file += 'shell rm {}\n'.format(trajectory_file) + + lammps_input_file += 'print "end of script" \n'.format(1000) + lammps_input_file += 'run {}\n'.format(1000) + + + return lammps_input_file + + +def generate_LAMMPS_potential(pair_style): + + + potential_file = '# Potential file generated by aiida plugin (please check citation in the orignal file)\n' + for key, value in pair_style.dict.data.iteritems(): + potential_file += '{} {}\n'.format(key, value) + + return potential_file + + +class CombinateCalculation(JobCalculation): + """ + A basic plugin for calculating force constants using Lammps. + + Requirement: the node should be able to import phonopy + """ + + def _init_internal_params(self): + super(CombinateCalculation, self)._init_internal_params() + + self._INPUT_FILE_NAME = 'input.in' + self._INPUT_POTENTIAL = 'potential.pot' + self._INPUT_STRUCTURE = 'input.data' + + self._OUTPUT_TRAJECTORY_FILE_NAME = 'trajectory.lampstrj' + + # self._default_parser = "lammps.md" + + self._INPUT_CELL = 'POSCAR' + self._INPUT_FORCE_CONSTANTS = 'FORCE_CONSTANTS' + self._INPUT_FILE_NAME_DYNA = 'input_dynaphopy' + self._OUTPUT_FORCE_CONSTANTS = 'FORCE_CONSTANTS_OUT' + self._OUTPUT_FILE_NAME = 'OUTPUT' + self._OUTPUT_QUASIPARTICLES = 'quasiparticles_data.yaml' + + self._default_parser = 'dynaphopy' + + + @classproperty + def _use_methods(cls): + """ + Additional use_* methods for the namelists class. + """ + retdict = JobCalculation._use_methods + retdict.update({ + "parameters": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'parameters', + 'docstring': ("Use a node that specifies the lammps input data " + "for the namelists"), + }, + "parameters_dynaphopy": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'parameters_phonopy', + 'docstring': ("Use a node that specifies the dynaphopy input data " + "for the namelists"), + }, + "force_constants": { + 'valid_types': ArrayData, + 'additional_parameter': None, + 'linkname': 'force_constants', + 'docstring': ("Use a node that specifies the force_constants " + "for the namelists"), + }, + "potential": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'potential', + 'docstring': ("Use a node that specifies the lammps potential " + "for the namelists"), + }, + "structure": { + 'valid_types': StructureData, + 'additional_parameter': None, + 'linkname': 'structure', + 'docstring': "Use a node for the structure", + }, + "supercell_md": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'supercell_md', + 'docstring': "Use a node for the supercell MD shape", + }, + }) + return retdict + + def _prepare_for_submission(self,tempfolder, inputdict): + """ + This is the routine to be called when you want to create + the input files and related stuff with a plugin. + + :param tempfolder: a aiida.common.folders.Folder subclass where + the plugin should put all its files. + :param inputdict: a dictionary with the input nodes, as they would + be returned by get_inputdata_dict (without the Code!) + """ + + try: + parameters_data = inputdict.pop(self.get_linkname('parameters')) + except KeyError: + raise InputValidationError("No parameters specified for this " + "calculation") + + if not isinstance(parameters_data, ParameterData): + raise InputValidationError("parameters is not of type " + "ParameterData") + + try: + parameters_data_dynaphopy = inputdict.pop(self.get_linkname('parameters_dynaphopy')) + force_constants = inputdict.pop(self.get_linkname('force_constants')) + + except KeyError: + raise InputValidationError("No dynaphopy parameters specified for this " + "calculation") + + try: + potential_data = inputdict.pop(self.get_linkname('potential')) + except KeyError: + raise InputValidationError("No potential specified for this " + "calculation") + + if not isinstance(potential_data, ParameterData): + raise InputValidationError("potential is not of type " + "ParameterData") + + try: + structure = inputdict.pop(self.get_linkname('structure')) + except KeyError: + raise InputValidationError("no structure is specified for this calculation") + + try: + supercell_shape = inputdict.pop(self.get_linkname('supercell_md')) + except KeyError: + raise InputValidationError("no supercell is specified for this calculation") + + + try: + code = inputdict.pop(self.get_linkname('code')) + except KeyError: + raise InputValidationError("no code is specified for this calculation") + + ############################## + # END OF INITIAL INPUT CHECK # + ############################## + + time_step = parameters_data.dict.timestep + + # Dynaphopy command + cmdline_params = ['/usr/bin/python', '/home/abel/BIN/dynaphopy', self._INPUT_FILE_NAME_DYNA, + self._OUTPUT_TRAJECTORY_FILE_NAME, + '-ts', '{}'.format(time_step), '--silent', + '-sfc', self._OUTPUT_FORCE_CONSTANTS, '-thm', # '--resolution 0.05', + '-psm', '2', '--normalize_dos', '-sdata'] # PS algorithm + + if 'temperature' in parameters_data.get_dict(): + cmdline_params.append('--temperature') + cmdline_params.append('{}'.format(parameters_data.dict.temperature)) + + if 'md_commensurate' in parameters_data.get_dict(): + if parameters_data.dict.md_commensurate: + cmdline_params.append('--MD_commensurate') + + cmdline_params.append('> OUTPUT') + + # =================== prepare the python input files ===================== + + structure_md = get_supercell(structure, supercell_shape) + potential_object = LammpsPotential(potential_data, structure_md, potential_filename=self._INPUT_POTENTIAL) + + structure_txt = generate_LAMMPS_structure(structure_md) + input_txt = generate_LAMMPS_input(parameters_data, + potential_object, + structure_file=self._INPUT_STRUCTURE, + trajectory_file=self._OUTPUT_TRAJECTORY_FILE_NAME, + command=' '.join(cmdline_params)) + + potential_txt = potential_object.get_potential_file() + + # =========================== dump to file ============================= + + input_filename = tempfolder.get_abs_path(self._INPUT_FILE_NAME) + with open(input_filename, 'w') as infile: + infile.write(input_txt) + + structure_filename = tempfolder.get_abs_path(self._INPUT_STRUCTURE) + with open(structure_filename, 'w') as infile: + infile.write(structure_txt) + + potential_filename = tempfolder.get_abs_path(self._INPUT_POTENTIAL) + with open(potential_filename, 'w') as infile: + infile.write(potential_txt) + + # =+========================= Dynaphopy =+============================== + + cell_txt = structure_to_poscar(structure) + input_txt = parameters_to_input_file(parameters_data_dynaphopy) + force_constants_txt = get_FORCE_CONSTANTS_txt(force_constants) + + input_filename = tempfolder.get_abs_path(self._INPUT_FILE_NAME_DYNA) + with open(input_filename, 'w') as infile: + infile.write(input_txt) + + cell_filename = tempfolder.get_abs_path(self._INPUT_CELL) + with open(cell_filename, 'w') as infile: + infile.write(cell_txt) + + force_constants_filename = tempfolder.get_abs_path(self._INPUT_FORCE_CONSTANTS) + with open(force_constants_filename, 'w') as infile: + infile.write(force_constants_txt) + + # ============================ calcinfo ================================ + + local_copy_list = [] + remote_copy_list = [] + # additional_retrieve_list = settings_dict.pop("ADDITIONAL_RETRIEVE_LIST",[]) + + calcinfo = CalcInfo() + + calcinfo.uuid = self.uuid + # Empty command line by default + calcinfo.local_copy_list = local_copy_list + calcinfo.remote_copy_list = remote_copy_list + + # Retrieve files + calcinfo.retrieve_list = [self._OUTPUT_FORCE_CONSTANTS, + self._OUTPUT_FILE_NAME, + self._OUTPUT_QUASIPARTICLES] + + codeinfo = CodeInfo() + codeinfo.cmdline_params = ['-in', self._INPUT_FILE_NAME] + + codeinfo.code_uuid = code.uuid + codeinfo.withmpi = True + calcinfo.codes_info = [codeinfo] + return calcinfo diff --git a/plugins/jobs/lammps/force.py b/plugins/jobs/lammps/force.py new file mode 100644 index 0000000..5d72598 --- /dev/null +++ b/plugins/jobs/lammps/force.py @@ -0,0 +1,229 @@ +from aiida.orm.calculation.job import JobCalculation +from aiida.common.exceptions import InputValidationError +from aiida.common.datastructures import CalcInfo, CodeInfo +from aiida.common.utils import classproperty +from aiida.orm import DataFactory + +from potentials import LammpsPotential + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + + +def generate_LAMMPS_structure(structure): + import numpy as np + + types = [site.kind_name for site in structure.sites] + + type_index_unique = np.unique(types, return_index=True)[1] + count_index_unique = np.diff(np.append(type_index_unique, [len(types)])) + + atom_index = [] + for i, index in enumerate(count_index_unique): + atom_index += [i for j in range(index)] + + masses = [site.mass for site in structure.kinds] + positions = [site.position for site in structure.sites] + + number_of_atoms = len(positions) + + lammps_data_file = 'Generated using dynaphopy\n\n' + lammps_data_file += '{0} atoms\n\n'.format(number_of_atoms) + lammps_data_file += '{0} atom types\n\n'.format(len(masses)) + + cell = np.array(structure.cell) + + a = np.linalg.norm(cell[0]) + b = np.linalg.norm(cell[1]) + c = np.linalg.norm(cell[2]) + + alpha = np.arccos(np.dot(cell[1], cell[2])/(c*b)) + gamma = np.arccos(np.dot(cell[1], cell[0])/(a*b)) + beta = np.arccos(np.dot(cell[2], cell[0])/(a*c)) + + xhi = a + xy = b * np.cos(gamma) + xz = c * np.cos(beta) + yhi = np.sqrt(pow(b,2)- pow(xy,2)) + yz = (b*c*np.cos(alpha)-xy * xz)/yhi + zhi = np.sqrt(pow(c,2)-pow(xz,2)-pow(yz,2)) + + xhi = xhi + max(0,0, xy, xz, xy+xz) + yhi = yhi + max(0,0, yz) + + lammps_data_file += '\n{0:20.10f} {1:20.10f} xlo xhi\n'.format(0, xhi) + lammps_data_file += '{0:20.10f} {1:20.10f} ylo yhi\n'.format(0, yhi) + lammps_data_file += '{0:20.10f} {1:20.10f} zlo zhi\n'.format(0, zhi) + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f} xy xz yz\n\n'.format(xy, xz, yz) + + lammps_data_file += 'Masses\n\n' + + for i, mass in enumerate(masses): + lammps_data_file += '{0} {1:20.10f} \n'.format(i+1, mass) + + + lammps_data_file += '\nAtoms\n\n' + for i, row in enumerate(positions): + lammps_data_file += '{0} {1} {2:20.10f} {3:20.10f} {4:20.10f}\n'.format(i+1, atom_index[i]+1, row[0],row[1],row[2]) + + return lammps_data_file + + +def generate_LAMMPS_input(potential_obj, + structure_file='data.gan', + trajectory_file='trajectory.lammpstr'): + + names_str = ' '.join(potential_obj._names) + + lammps_input_file = 'units metal\n' + lammps_input_file += 'boundary p p p\n' + lammps_input_file += 'box tilt large\n' + lammps_input_file += 'atom_style atomic\n' + + lammps_input_file += 'read_data {}\n'.format(structure_file) + + lammps_input_file += potential_obj.get_input_potential_lines() + + lammps_input_file += 'neighbor 0.3 bin\n' + lammps_input_file += 'neigh_modify every 1 delay 0 check no\n' + lammps_input_file += 'dump aiida all custom 1 {0} element fx fy fz\n'.format(trajectory_file) + lammps_input_file += 'dump_modify aiida format "%4s %16.10f %16.10f %16.10f"\n' + lammps_input_file += 'dump_modify aiida sort id\n' + lammps_input_file += 'dump_modify aiida element {}\n'.format(names_str) + + lammps_input_file += 'run 0' + + return lammps_input_file + + +class ForceCalculation(JobCalculation): + """ + A basic plugin for calculating force constants using Lammps. + + Requirement: the node should be able to import phonopy + """ + + def _init_internal_params(self): + super(ForceCalculation, self)._init_internal_params() + + self._INPUT_FILE_NAME = 'input.in' + self._INPUT_POTENTIAL = 'potential.pot' + self._INPUT_STRUCTURE = 'input.data' + + self._OUTPUT_TRAJECTORY_FILE_NAME = 'trajectory.lammpstrj' + + self._OUTPUT_FILE_NAME = 'log.lammps' + self._default_parser = 'lammps.force' + + + @classproperty + def _use_methods(cls): + """ + Additional use_* methods for the namelists class. + """ + retdict = JobCalculation._use_methods + retdict.update({ + "potential": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'potential', + 'docstring': ("Use a node that specifies the lammps potential " + "for the namelists"), + }, + "structure": { + 'valid_types': StructureData, + 'additional_parameter': None, + 'linkname': 'structure', + 'docstring': "Use a node for the structure", + }, + }) + return retdict + + def _prepare_for_submission(self,tempfolder, inputdict): + """ + This is the routine to be called when you want to create + the input files and related stuff with a plugin. + + :param tempfolder: a aiida.common.folders.Folder subclass where + the plugin should put all its files. + :param inputdict: a dictionary with the input nodes, as they would + be returned by get_inputdata_dict (without the Code!) + """ + + try: + potential_data = inputdict.pop(self.get_linkname('potential')) + except KeyError: + raise InputValidationError("No potential specified for this " + "calculation") + + if not isinstance(potential_data, ParameterData): + raise InputValidationError("potential is not of type " + "ParameterData") + + try: + structure = inputdict.pop(self.get_linkname('structure')) + except KeyError: + raise InputValidationError("no structure is specified for this calculation") + + try: + code = inputdict.pop(self.get_linkname('code')) + except KeyError: + raise InputValidationError("no code is specified for this calculation") + + ############################## + # END OF INITIAL INPUT CHECK # + ############################## + + # =================== prepare the python input files ===================== + + potential_object = LammpsPotential(potential_data, structure, potential_filename=self._INPUT_POTENTIAL) + + + structure_txt = generate_LAMMPS_structure(structure) + input_txt = generate_LAMMPS_input(potential_object, + structure_file=self._INPUT_STRUCTURE, + trajectory_file=self._OUTPUT_TRAJECTORY_FILE_NAME,) + + potential_txt = potential_object.get_potential_file() + # =========================== dump to file ============================= + + input_filename = tempfolder.get_abs_path(self._INPUT_FILE_NAME) + with open(input_filename, 'w') as infile: + infile.write(input_txt) + + structure_filename = tempfolder.get_abs_path(self._INPUT_STRUCTURE) + with open(structure_filename, 'w') as infile: + infile.write(structure_txt) + + if potential_txt is not None: + potential_filename = tempfolder.get_abs_path(self._INPUT_POTENTIAL) + with open(potential_filename, 'w') as infile: + infile.write(potential_txt) + + # ============================ calcinfo ================================ + + local_copy_list = [] + remote_copy_list = [] + # additional_retrieve_list = settings_dict.pop("ADDITIONAL_RETRIEVE_LIST",[]) + + calcinfo = CalcInfo() + + calcinfo.uuid = self.uuid + # Empty command line by default + calcinfo.local_copy_list = local_copy_list + calcinfo.remote_copy_list = remote_copy_list + + # Retrieve files + calcinfo.retrieve_list = [] + calcinfo.retrieve_list.append(self._OUTPUT_TRAJECTORY_FILE_NAME) + calcinfo.retrieve_list.append(self._OUTPUT_FILE_NAME) + + codeinfo = CodeInfo() + codeinfo.cmdline_params = ['-in', self._INPUT_FILE_NAME] + codeinfo.code_uuid = code.uuid + codeinfo.withmpi = False + calcinfo.codes_info = [codeinfo] + return calcinfo + + +#$MPI -n $NSLOTS $LAMMPS -sf gpu -pk gpu 2 neigh no -in in.md_data diff --git a/plugins/jobs/lammps/md.py b/plugins/jobs/lammps/md.py new file mode 100644 index 0000000..29e78de --- /dev/null +++ b/plugins/jobs/lammps/md.py @@ -0,0 +1,260 @@ +from aiida.orm.calculation.job import JobCalculation +from aiida.orm.data.parameter import ParameterData +from aiida.orm.data.structure import StructureData +from aiida.common.exceptions import InputValidationError +from aiida.common.datastructures import CalcInfo, CodeInfo +from aiida.common.utils import classproperty + +from potentials import LammpsPotential +import numpy as np + +def generate_LAMMPS_structure(structure): + import numpy as np + + types = [site.kind_name for site in structure.sites] + + type_index_unique = np.unique(types, return_index=True)[1] + count_index_unique = np.diff(np.append(type_index_unique, [len(types)])) + + atom_index = [] + for i, index in enumerate(count_index_unique): + atom_index += [i for j in range(index)] + + masses = [site.mass for site in structure.kinds] + positions = [site.position for site in structure.sites] + + number_of_atoms = len(positions) + + lammps_data_file = 'Generated using dynaphopy\n\n' + lammps_data_file += '{0} atoms\n\n'.format(number_of_atoms) + lammps_data_file += '{0} atom types\n\n'.format(len(masses)) + + cell = np.array(structure.cell) + + a = np.linalg.norm(cell[0]) + b = np.linalg.norm(cell[1]) + c = np.linalg.norm(cell[2]) + + alpha = np.arccos(np.dot(cell[1], cell[2])/(c*b)) + gamma = np.arccos(np.dot(cell[1], cell[0])/(a*b)) + beta = np.arccos(np.dot(cell[2], cell[0])/(a*c)) + + xhi = a + xy = b * np.cos(gamma) + xz = c * np.cos(beta) + yhi = np.sqrt(pow(b,2)- pow(xy,2)) + yz = (b*c*np.cos(alpha)-xy * xz)/yhi + zhi = np.sqrt(pow(c,2)-pow(xz,2)-pow(yz,2)) + + xhi = xhi + max(0,0, xy, xz, xy+xz) + yhi = yhi + max(0,0, yz) + + lammps_data_file += '\n{0:20.10f} {1:20.10f} xlo xhi\n'.format(0, xhi) + lammps_data_file += '{0:20.10f} {1:20.10f} ylo yhi\n'.format(0, yhi) + lammps_data_file += '{0:20.10f} {1:20.10f} zlo zhi\n'.format(0, zhi) + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f} xy xz yz\n\n'.format(xy, xz, yz) + + lammps_data_file += 'Masses\n\n' + + for i, mass in enumerate(masses): + lammps_data_file += '{0} {1:20.10f} \n'.format(i+1, mass) + + + lammps_data_file += '\nAtoms\n\n' + for i, row in enumerate(positions): + lammps_data_file += '{0} {1} {2:20.10f} {3:20.10f} {4:20.10f}\n'.format(i+1, atom_index[i]+1, row[0],row[1],row[2]) + + return lammps_data_file + + +def generate_LAMMPS_input(parameters, + potential_obj, + structure_file='potential.pot', + trajectory_file='trajectory.lammpstr'): + + random_number = np.random.randint(10000000) + + names_str = ' '.join(potential_obj._names) + + lammps_input_file = 'units metal\n' + lammps_input_file += 'boundary p p p\n' + lammps_input_file += 'box tilt large\n' + lammps_input_file += 'atom_style atomic\n' + lammps_input_file += 'read_data {}\n'.format(structure_file) + + lammps_input_file += potential_obj.get_input_potential_lines() + + lammps_input_file += 'neighbor 0.3 bin\n' + lammps_input_file += 'neigh_modify every 1 delay 0 check no\n' + + lammps_input_file += 'timestep {}\n'.format(parameters.dict.timestep) + + lammps_input_file += 'thermo_style custom step etotal temp vol press\n' + lammps_input_file += 'thermo 1000\n' + + lammps_input_file += 'velocity all create {0} {1} dist gaussian mom yes\n'.format(parameters.dict.temperature, random_number) + lammps_input_file += 'velocity all scale {}\n'.format(parameters.dict.temperature) + + lammps_input_file += 'fix int all nvt temp {0} {0} {1}\n'.format(parameters.dict.temperature, parameters.dict.thermostat_variable) + + lammps_input_file += 'run {}\n'.format(parameters.dict.equilibrium_steps) + lammps_input_file += 'reset_timestep 0\n' + + lammps_input_file += 'dump aiida all custom {0} {1} element x y z\n'.format(parameters.dict.dump_rate, trajectory_file) + lammps_input_file += 'dump_modify aiida format "%4s %16.10f %16.10f %16.10f"\n' + lammps_input_file += 'dump_modify aiida sort id\n' + lammps_input_file += 'dump_modify aiida element {}\n'.format(names_str) + + lammps_input_file += 'run {}\n'.format(parameters.dict.total_steps) + + return lammps_input_file + + +class MdCalculation(JobCalculation): + """ + A basic plugin for calculating force constants using Lammps. + + Requirement: the node should be able to import phonopy + """ + + def _init_internal_params(self): + super(MdCalculation, self)._init_internal_params() + + self._INPUT_FILE_NAME = 'input.in' + self._INPUT_POTENTIAL = 'potential.pot' + self._INPUT_STRUCTURE = 'input.data' + + self._OUTPUT_TRAJECTORY_FILE_NAME = 'trajectory.lammpstrj' + + self._OUTPUT_FILE_NAME = 'log.lammps' + self._default_parser = 'lammps.md' + + @classproperty + def _use_methods(cls): + """ + Additional use_* methods for the namelists class. + """ + retdict = JobCalculation._use_methods + retdict.update({ + "parameters": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'parameters', + 'docstring': ("Use a node that specifies the lammps input data " + "for the namelists"), + }, + "potential": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'potential', + 'docstring': ("Use a node that specifies the lammps potential " + "for the namelists"), + }, + "structure": { + 'valid_types': StructureData, + 'additional_parameter': None, + 'linkname': 'structure', + 'docstring': "Use a node for the structure", + }, + }) + return retdict + + def _prepare_for_submission(self,tempfolder, inputdict): + """ + This is the routine to be called when you want to create + the input files and related stuff with a plugin. + + :param tempfolder: a aiida.common.folders.Folder subclass where + the plugin should put all its files. + :param inputdict: a dictionary with the input nodes, as they would + be returned by get_inputdata_dict (without the Code!) + """ + + try: + parameters_data = inputdict.pop(self.get_linkname('parameters')) + except KeyError: + raise InputValidationError("No parameters specified for this " + "calculation") + + if not isinstance(parameters_data, ParameterData): + raise InputValidationError("parameters is not of type " + "ParameterData") + + try: + potential_data = inputdict.pop(self.get_linkname('potential')) + except KeyError: + raise InputValidationError("No potential specified for this " + "calculation") + + if not isinstance(potential_data, ParameterData): + raise InputValidationError("potential is not of type " + "ParameterData") + + try: + structure = inputdict.pop(self.get_linkname('structure')) + except KeyError: + raise InputValidationError("no structure is specified for this calculation") + + try: + code = inputdict.pop(self.get_linkname('code')) + except KeyError: + raise InputValidationError("no code is specified for this calculation") + + ############################## + # END OF INITIAL INPUT CHECK # + ############################## + + # =================== prepare the python input files ===================== + + potential_object = LammpsPotential(potential_data, structure, potential_filename=self._INPUT_POTENTIAL) + + structure_txt = generate_LAMMPS_structure(structure) + input_txt = generate_LAMMPS_input(parameters_data, + potential_object, + structure_file=self._INPUT_STRUCTURE, + trajectory_file=self._OUTPUT_TRAJECTORY_FILE_NAME) + + potential_txt = potential_object.get_potential_file() + + # =========================== dump to file ============================= + + input_filename = tempfolder.get_abs_path(self._INPUT_FILE_NAME) + with open(input_filename, 'w') as infile: + infile.write(input_txt) + + structure_filename = tempfolder.get_abs_path(self._INPUT_STRUCTURE) + with open(structure_filename, 'w') as infile: + infile.write(structure_txt) + + if potential_txt is not None: + potential_filename = tempfolder.get_abs_path(self._INPUT_POTENTIAL) + with open(potential_filename, 'w') as infile: + infile.write(potential_txt) + + # ============================ calcinfo ================================ + + local_copy_list = [] + remote_copy_list = [] + # additional_retrieve_list = settings_dict.pop("ADDITIONAL_RETRIEVE_LIST",[]) + + calcinfo = CalcInfo() + + calcinfo.uuid = self.uuid + # Empty command line by default + calcinfo.local_copy_list = local_copy_list + calcinfo.remote_copy_list = remote_copy_list + + # Retrieve files + calcinfo.retrieve_list = [] + calcinfo.retrieve_list.append(self._OUTPUT_TRAJECTORY_FILE_NAME) + calcinfo.retrieve_list.append(self._OUTPUT_FILE_NAME) + + codeinfo = CodeInfo() + codeinfo.cmdline_params = ['-in', self._INPUT_FILE_NAME] + codeinfo.code_uuid = code.uuid + codeinfo.withmpi = True + calcinfo.codes_info = [codeinfo] + return calcinfo + + +#$MPI -n $NSLOTS $LAMMPS -sf gpu -pk gpu 2 neigh no -in in.md_data diff --git a/plugins/jobs/lammps/optimize.py b/plugins/jobs/lammps/optimize.py new file mode 100644 index 0000000..1b0f8fc --- /dev/null +++ b/plugins/jobs/lammps/optimize.py @@ -0,0 +1,267 @@ +from aiida.orm.calculation.job import JobCalculation +from aiida.orm.data.parameter import ParameterData +from aiida.orm.data.structure import StructureData +from aiida.common.exceptions import InputValidationError +from aiida.common.datastructures import CalcInfo, CodeInfo +from aiida.common.utils import classproperty + +from potentials import LammpsPotential + + +def generate_LAMMPS_structure(structure): + import numpy as np + + types = [site.kind_name for site in structure.sites] + + type_index_unique = np.unique(types, return_index=True)[1] + count_index_unique = np.diff(np.append(type_index_unique, [len(types)])) + + atom_index = [] + for i, index in enumerate(count_index_unique): + atom_index += [i for j in range(index)] + + masses = [site.mass for site in structure.kinds] + positions = [site.position for site in structure.sites] + + number_of_atoms = len(positions) + + lammps_data_file = 'Generated using dynaphopy\n\n' + lammps_data_file += '{0} atoms\n\n'.format(number_of_atoms) + lammps_data_file += '{0} atom types\n\n'.format(len(masses)) + + cell = np.array(structure.cell) + + a = np.linalg.norm(cell[0]) + b = np.linalg.norm(cell[1]) + c = np.linalg.norm(cell[2]) + + alpha = np.arccos(np.dot(cell[1], cell[2])/(c*b)) + gamma = np.arccos(np.dot(cell[1], cell[0])/(a*b)) + beta = np.arccos(np.dot(cell[2], cell[0])/(a*c)) + + xhi = a + xy = b * np.cos(gamma) + xz = c * np.cos(beta) + yhi = np.sqrt(pow(b,2)- pow(xy,2)) + yz = (b*c*np.cos(alpha)-xy * xz)/yhi + zhi = np.sqrt(pow(c,2)-pow(xz,2)-pow(yz,2)) + + xhi = xhi + max(0,0, xy, xz, xy+xz) + yhi = yhi + max(0,0, yz) + + lammps_data_file += '\n{0:20.10f} {1:20.10f} xlo xhi\n'.format(0, xhi) + lammps_data_file += '{0:20.10f} {1:20.10f} ylo yhi\n'.format(0, yhi) + lammps_data_file += '{0:20.10f} {1:20.10f} zlo zhi\n'.format(0, zhi) + lammps_data_file += '{0:20.10f} {1:20.10f} {2:20.10f} xy xz yz\n\n'.format(xy, xz, yz) + + lammps_data_file += 'Masses\n\n' + + for i, mass in enumerate(masses): + lammps_data_file += '{0} {1:20.10f} \n'.format(i+1, mass) + + + lammps_data_file += '\nAtoms\n\n' + for i, row in enumerate(positions): + lammps_data_file += '{0} {1} {2:20.10f} {3:20.10f} {4:20.10f}\n'.format(i+1, atom_index[i]+1, row[0],row[1],row[2]) + + return lammps_data_file + + +def generate_LAMMPS_input(potential_obj, + parameters_data, + structure_file='data.gan', + optimize_path_file='path.lammpstrj'): + + names_str = ' '.join(potential_obj._names) + + + + parameters = parameters_data.get_dict() + + lammps_input_file = 'units metal\n' + lammps_input_file += 'boundary p p p\n' + lammps_input_file += 'box tilt large\n' + lammps_input_file += 'atom_style atomic\n' + lammps_input_file += 'read_data {}\n'.format(structure_file) + + lammps_input_file += potential_obj.get_input_potential_lines() + + lammps_input_file += 'fix int all box/relax {} {} vmax {}\n'.format(parameters['relaxation'], + parameters['pressure'] * 1000, # pressure kb -> bar + parameters['vmax']) + + lammps_input_file += 'compute stpa all stress/atom NULL\n' + # xx, yy, zz, xy, xz, yz + lammps_input_file += 'compute stgb all reduce sum c_stpa[1] c_stpa[2] c_stpa[3] c_stpa[4] c_stpa[5] c_stpa[6]\n' + lammps_input_file += 'variable pr equal -(c_stgb[1]+c_stgb[2]+c_stgb[3])/(3*vol)\n' + lammps_input_file += 'thermo_style custom step temp press v_pr etotal c_stgb[1] c_stgb[2] c_stgb[3] c_stgb[4] c_stgb[5] c_stgb[6]\n' + + lammps_input_file += 'dump aiida all custom 1 {0} element x y z fx fy fz\n'.format(optimize_path_file) + lammps_input_file += 'dump_modify aiida format "%4s %16.10f %16.10f %16.10f %16.10f %16.10f %16.10f"\n' + lammps_input_file += 'dump_modify aiida sort id\n' + lammps_input_file += 'dump_modify aiida element {}\n'.format(names_str) + lammps_input_file += 'min_style cg\n' + lammps_input_file += 'minimize {} {} {} {}\n'.format(parameters['energy_tolerance'], + parameters['force_tolerance'], + parameters['max_iterations'], + parameters['max_evaluations']) + # lammps_input_file += 'print "$(xlo - xhi) $(xy) $(xz)"\n' + # lammps_input_file += 'print "0.000 $(yhi - ylo) $(yz)"\n' + # lammps_input_file += 'print "0.000 0.000 $(zhi-zlo)"\n' + lammps_input_file += 'print "$(xlo) $(xhi) $(xy)"\n' + lammps_input_file += 'print "$(ylo) $(yhi) $(xz)"\n' + lammps_input_file += 'print "$(zlo) $(zhi) $(yz)"\n' + + return lammps_input_file + + +class OptimizeCalculation(JobCalculation): + """ + A basic plugin for calculating force constants using Lammps. + + Requirement: the node should be able to import phonopy + """ + + def _init_internal_params(self): + super(OptimizeCalculation, self)._init_internal_params() + + self._INPUT_FILE_NAME = 'input.in' + self._INPUT_POTENTIAL = 'potential.pot' + self._INPUT_STRUCTURE = 'input.data' + + self._OUTPUT_TRAJECTORY_FILE_NAME = 'path.lammpstrj' + + self._OUTPUT_FILE_NAME = 'log.lammps' + self._default_parser = 'lammps.optimize' + + @classproperty + def _use_methods(cls): + """ + Additional use_* methods for the namelists class. + """ + retdict = JobCalculation._use_methods + retdict.update({ + "potential": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'potential', + 'docstring': ("Use a node that specifies the lammps potential " + "for the namelists"), + }, + "parameters": { + 'valid_types': ParameterData, + 'additional_parameter': None, + 'linkname': 'parameters', + 'docstring': ("Use a node that specifies the lammps parameters " + "for the namelists"), + }, + "structure": { + 'valid_types': StructureData, + 'additional_parameter': None, + 'linkname': 'structure', + 'docstring': "Use a node for the structure", + }, + }) + return retdict + + def _prepare_for_submission(self,tempfolder, inputdict): + """ + This is the routine to be called when you want to create + the input files and related stuff with a plugin. + + :param tempfolder: a aiida.common.folders.Folder subclass where + the plugin should put all its files. + :param inputdict: a dictionary with the input nodes, as they would + be returned by get_inputdata_dict (without the Code!) + """ + + try: + potential_data = inputdict.pop(self.get_linkname('potential')) + except KeyError: + raise InputValidationError("No potential specified for this " + "calculation") + + if not isinstance(potential_data, ParameterData): + raise InputValidationError("potential is not of type " + "ParameterData") + + try: + parameters_data = inputdict.pop(self.get_linkname('parameters')) + except KeyError: + raise InputValidationError("No parameters specified for this " + "calculation") + + if not isinstance(parameters_data, ParameterData): + raise InputValidationError("parameters is not of type " + "ParameterData") + + + try: + structure = inputdict.pop(self.get_linkname('structure')) + except KeyError: + raise InputValidationError("no structure is specified for this calculation") + + try: + code = inputdict.pop(self.get_linkname('code')) + except KeyError: + raise InputValidationError("no code is specified for this calculation") + + ############################## + # END OF INITIAL INPUT CHECK # + ############################## + + # =================== prepare the python input files ===================== + + potential_object = LammpsPotential(potential_data, structure, potential_filename=self._INPUT_POTENTIAL) + + + structure_txt = generate_LAMMPS_structure(structure) + input_txt = generate_LAMMPS_input(potential_object, + parameters_data, + structure_file=self._INPUT_STRUCTURE, + optimize_path_file=self._OUTPUT_TRAJECTORY_FILE_NAME) + + potential_txt = potential_object.get_potential_file() + + # =========================== dump to file ============================= + + input_filename = tempfolder.get_abs_path(self._INPUT_FILE_NAME) + with open(input_filename, 'w') as infile: + infile.write(input_txt) + + structure_filename = tempfolder.get_abs_path(self._INPUT_STRUCTURE) + with open(structure_filename, 'w') as infile: + infile.write(structure_txt) + + if potential_txt is not None: + potential_filename = tempfolder.get_abs_path(self._INPUT_POTENTIAL) + with open(potential_filename, 'w') as infile: + infile.write(potential_txt) + + # ============================ calcinfo ================================ + + local_copy_list = [] + remote_copy_list = [] + # additional_retrieve_list = settings_dict.pop("ADDITIONAL_RETRIEVE_LIST",[]) + + calcinfo = CalcInfo() + + calcinfo.uuid = self.uuid + # Empty command line by default + calcinfo.local_copy_list = local_copy_list + calcinfo.remote_copy_list = remote_copy_list + + # Retrieve files + calcinfo.retrieve_list = [] + calcinfo.retrieve_list.append(self._OUTPUT_TRAJECTORY_FILE_NAME) + calcinfo.retrieve_list.append(self._OUTPUT_FILE_NAME) + + codeinfo = CodeInfo() + codeinfo.cmdline_params = ['-in', self._INPUT_FILE_NAME] + codeinfo.code_uuid = code.uuid + codeinfo.withmpi = False + calcinfo.codes_info = [codeinfo] + return calcinfo + + +#$MPI -n $NSLOTS $LAMMPS -sf gpu -pk gpu 2 neigh no -in in.md_data diff --git a/plugins/jobs/lammps/potentials/__init__.py b/plugins/jobs/lammps/potentials/__init__.py new file mode 100644 index 0000000..a648daa --- /dev/null +++ b/plugins/jobs/lammps/potentials/__init__.py @@ -0,0 +1,37 @@ +import importlib + + +class LammpsPotential: + def __init__(self, + potential_data, + structure, + potential_filename='None'): + + self._names = [site.name for site in structure.kinds] + self._type = potential_data.dict.pair_style + self._data = potential_data.dict.data + + self._potential_filename = potential_filename + + try: + self._potential_module = importlib.import_module('.{}'.format(self._type), __name__) + except ImportError: + raise ImportError('This lammps potential is not implemented') + +# if self._type == 'tersoff': +# from tersoff import _generate_LAMMPS_potential, _get_input_potential_lines +# elif self._type == 'lennard_jones': +# from lennard_jones import _generate_LAMMPS_potential, _get_input_potential_lines +# else: +# raise ValueError('This lammps potential is not implemented') + +# self._generate_LAMMPS_potential = _generate_LAMMPS_potential +# self._get_input_potential_lines = _get_input_potential_lines + + def get_potential_file(self): + return self._potential_module.generate_LAMMPS_potential(self._data) + + def get_input_potential_lines(self): + return self._potential_module.get_input_potential_lines(self._data, + potential_filename=self._potential_filename, + names=self._names) diff --git a/plugins/jobs/lammps/potentials/eam.py b/plugins/jobs/lammps/potentials/eam.py new file mode 100644 index 0000000..6146c6c --- /dev/null +++ b/plugins/jobs/lammps/potentials/eam.py @@ -0,0 +1,17 @@ + +def generate_LAMMPS_potential(data): + #potential_file = '# Potential file generated by aiida plugin (please check citation in the orignal file)\n' + + potential_file = '' + for line in data['file_contents']: + potential_file += '{}'.format(line) + + return potential_file + + +def get_input_potential_lines(data, names=None, potential_filename='potential.pot'): + + lammps_input_text = 'pair_style eam/{}\n'.format(data['type']) + lammps_input_text += 'pair_coeff * * {} {}\n'.format(potential_filename, ' '.join(names)) + + return lammps_input_text \ No newline at end of file diff --git a/plugins/jobs/lammps/potentials/lennard_jones.py b/plugins/jobs/lammps/potentials/lennard_jones.py new file mode 100644 index 0000000..02b2ba5 --- /dev/null +++ b/plugins/jobs/lammps/potentials/lennard_jones.py @@ -0,0 +1,16 @@ +import numpy as np + + +def generate_LAMMPS_potential(data): + return None + + +def get_input_potential_lines(data, names=None, potential_filename='potential.pot'): + + cut = np.max([float(i.split()[2]) for i in data.values()]) + + lammps_input_text = 'pair_style lj/cut {}\n'.format(cut) + + for key, value in data.iteritems(): + lammps_input_text += 'pair_coeff {} {}\n'.format(key, value) + return lammps_input_text diff --git a/plugins/jobs/lammps/potentials/tersoff.py b/plugins/jobs/lammps/potentials/tersoff.py new file mode 100644 index 0000000..8b8af19 --- /dev/null +++ b/plugins/jobs/lammps/potentials/tersoff.py @@ -0,0 +1,16 @@ + + +def generate_LAMMPS_potential(data): + potential_file = '# Potential file generated by aiida plugin (please check citation in the orignal file)\n' + for key, value in data.iteritems(): + potential_file += '{} {}\n'.format(key, value) + + return potential_file + + +def get_input_potential_lines(data, names=None, potential_filename='potential.pot'): + + lammps_input_text = 'pair_style tersoff\n' + lammps_input_text += 'pair_coeff * * {} {}\n'.format(potential_filename, ' '.join(names)) + + return lammps_input_text \ No newline at end of file diff --git a/plugins/launcher/launch_dynaphopy_si.py b/plugins/launcher/launch_dynaphopy_si.py new file mode 100644 index 0000000..f89cbcf --- /dev/null +++ b/plugins/launcher/launch_dynaphopy_si.py @@ -0,0 +1,82 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'dynaphopy@stern' + +############################ +# Define input parameters # +############################ + +a = 5.404 +cell = [[a, 0, 0], + [0, a, 0], + [0, 0, a]] + +symbols=['Si'] * 8 +scaled_positions = [(0.875, 0.875, 0.875), + (0.875, 0.375, 0.375), + (0.375, 0.875, 0.375), + (0.375, 0.375, 0.875), + (0.125, 0.125, 0.125), + (0.125, 0.625, 0.625), + (0.625, 0.125, 0.625), + (0.625, 0.625, 0.125)] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +dynaphopy_parameters ={'supercell': [[2, 0, 0], + [0, 2, 0], + [0, 0, 2]], + 'primitive': [[1.0, 0.0, 0.0], + [0.0, 1.0, 0.0], + [0.0, 0.0, 1.0]], + 'mesh': [40, 40, 40], + 'md_commensurate': True, + 'temperature': 300} # Temperature can be omitted (If ommited calculated from Max.-Boltz.) + + +dynaphopy_machine = { + 'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + + +from aiida.orm import load_node +force_constants = load_node(20569) # Loads node that contains the harmonic force constants (Array data) +trajectory = load_node(20528) # Loads node that constains the MD trajectory (TrajectoryData) + + +codename = codename +code = Code.get_from_string(codename) +calc = code.new_calc(max_wallclock_seconds=3600, + resources=dynaphopy_machine) + +calc.label = "test dynaphopy calculation" +calc.description = "A much longer description" + +calc.use_code(code) + +calc.use_structure(structure) +calc.use_parameters(ParameterData(dict=dynaphopy_parameters)) +calc.use_force_constants(force_constants) +calc.use_trajectory(trajectory) + +calc.store_all() + + +calc.submit() +print "submitted calculation with PK={}".format(calc.dbnode.pk) diff --git a/plugins/launcher/launch_lammps_combinate.py b/plugins/launcher/launch_lammps_combinate.py new file mode 100644 index 0000000..e47a179 --- /dev/null +++ b/plugins/launcher/launch_lammps_combinate.py @@ -0,0 +1,110 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'lammps_comb@boston' + +############################ +# Define input parameters # +############################ + +a = 5.404 +cell = [[a, 0, 0], + [0, a, 0], + [0, 0, a]] + +symbols=['Si'] * 8 +scaled_positions = [(0.875, 0.875, 0.875), + (0.875, 0.375, 0.375), + (0.375, 0.875, 0.375), + (0.375, 0.375, 0.875), + (0.125, 0.125, 0.125), + (0.125, 0.625, 0.625), + (0.625, 0.125, 0.625), + (0.625, 0.625, 0.125)] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +# Silicon(C) Tersoff +tersoff_si = {'Si Si Si ': '3.0 1.0 1.7322 1.0039e5 16.218 -0.59826 0.78734 1.0999e-6 1.7322 471.18 2.85 0.15 2.4799 1830.8'} + + +potential ={'pair_style': 'tersoff', + 'data': tersoff_si} + + +dynaphopy_parameters ={'supercell': [[2, 0, 0], + [0, 2, 0], + [0, 0, 2]], + 'primitive': [[1.0, 0.0, 0.0], + [0.0, 1.0, 0.0], + [0.0, 0.0, 1.0]], + 'mesh': [40, 40, 40], + 'md_commensurate': True, + 'temperature': 300} # Temperature can be omitted (If ommited calculated from Max.-Boltz.) + + + +from aiida.orm import load_node +force_constants = load_node(20569) # Loads node that contains the harmonic force constants (Array data) + +machine = { + 'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + + +parameters_md = {'timestep': 0.001, + 'temperature': 300, + 'thermostat_variable': 0.5, + 'equilibrium_steps': 100, + 'total_steps': 2000, + 'dump_rate': 1} + + +code = Code.get_from_string(codename) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources=machine) + +calc.label = "test lammps calculation" +calc.description = "A much longer description" +calc.use_code(code) +calc.use_structure(structure) +calc.use_potential(ParameterData(dict=potential)) +calc.use_force_constants(force_constants) +calc.use_parameters_dynaphopy(ParameterData(dict=dynaphopy_parameters)) +calc.use_supercell_md(ParameterData(dict={'shape': [2, 2, 2]})) + +calc.use_parameters(ParameterData(dict=parameters_md)) + + +test_only = False + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/launcher/launch_lammps_force_gan.py b/plugins/launcher/launch_lammps_force_gan.py new file mode 100644 index 0000000..c9a3ebf --- /dev/null +++ b/plugins/launcher/launch_lammps_force_gan.py @@ -0,0 +1,90 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'lammps_force@boston' + +############################ +# Define input parameters # +############################ + +# GaN +cell = [[ 3.1900000572, 0, 0], + [-1.5950000286, 2.762621076, 0], + [ 0.0, 0, 5.1890001297]] + + +scaled_positions=[(0.6666669, 0.3333334, 0.0000000), + (0.3333331, 0.6666663, 0.5000000), + (0.6666669, 0.3333334, 0.3750000), + (0.3333331, 0.6666663, 0.8750000)] + +symbols=['Ga', 'Ga', 'N', 'N'] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +# GaN Tersoff +tersoff_gan = {'Ga Ga Ga': '1.0 0.007874 1.846 1.918000 0.75000 -0.301300 1.0 1.0 1.44970 410.132 2.87 0.15 1.60916 535.199', + 'N N N' : '1.0 0.766120 0.000 0.178493 0.20172 -0.045238 1.0 1.0 2.38426 423.769 2.20 0.20 3.55779 1044.77', + 'Ga Ga N' : '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 0.0 0.00000 0.00000 2.90 0.20 0.00000 0.00000', + 'Ga N N' : '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 1.0 2.63906 3864.27 2.90 0.20 2.93516 6136.44', + 'N Ga Ga': '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 1.0 2.63906 3864.27 2.90 0.20 2.93516 6136.44', + 'N Ga N ': '1.0 0.766120 0.000 0.178493 0.20172 -0.045238 1.0 0.0 0.00000 0.00000 2.20 0.20 0.00000 0.00000', + 'N N Ga': '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 0.0 0.00000 0.00000 2.90 0.20 0.00000 0.00000', + 'Ga N Ga': '1.0 0.007874 1.846 1.918000 0.75000 -0.301300 1.0 0.0 0.00000 0.00000 2.87 0.15 0.00000 0.00000'} + +# Silicon(C) Tersoff +# tersoff_si = {'Si Si Si ': '3.0 1.0 1.7322 1.0039e5 16.218 -0.59826 0.78734 1.0999e-6 1.7322 471.18 2.85 0.15 2.4799 1830.8'} + + +potential ={'pair_style': 'tersoff', + 'data': tersoff_gan} + +lammps_machine = { + 'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + + +code = Code.get_from_string(codename) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources=lammps_machine) + +calc.label = "test lammps calculation" +calc.description = "A much longer description" +calc.use_code(code) +calc.use_structure(structure) +calc.use_potential(ParameterData(dict=potential)) + +print "submitted calculation with PK={}".format(calc.dbnode.pk) + +test_only = False + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/launcher/launch_lammps_md_si.py b/plugins/launcher/launch_lammps_md_si.py new file mode 100644 index 0000000..af009a0 --- /dev/null +++ b/plugins/launcher/launch_lammps_md_si.py @@ -0,0 +1,91 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'lammps_md@boston' + +############################ +# Define input parameters # +############################ + +a = 5.404 +cell = [[a, 0, 0], + [0, a, 0], + [0, 0, a]] + +symbols=['Si'] * 8 +scaled_positions = [(0.875, 0.875, 0.875), + (0.875, 0.375, 0.375), + (0.375, 0.875, 0.375), + (0.375, 0.375, 0.875), + (0.125, 0.125, 0.125), + (0.125, 0.625, 0.625), + (0.625, 0.125, 0.625), + (0.625, 0.625, 0.125)] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +# Silicon(C) Tersoff +tersoff_si = {'Si Si Si ': '3.0 1.0 1.7322 1.0039e5 16.218 -0.59826 0.78734 1.0999e-6 1.7322 471.18 2.85 0.15 2.4799 1830.8'} + + +potential ={'pair_style': 'tersoff', + 'data': tersoff_si} + +lammps_machine = { + 'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + + +parameters_md = {'timestep': 0.001, + 'temperature' : 300, + 'thermostat_variable': 0.5, + 'equilibrium_steps': 100, + 'total_steps': 2000, + 'dump_rate': 1} + + +code = Code.get_from_string(codename) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources=lammps_machine) + +calc.label = "test lammps calculation" +calc.description = "A much longer description" +calc.use_code(code) +calc.use_structure(structure) +calc.use_potential(ParameterData(dict=potential)) + +calc.use_parameters(ParameterData(dict=parameters_md)) + + +test_only = False + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/launcher/launch_lammps_optimization_fe.py b/plugins/launcher/launch_lammps_optimization_fe.py new file mode 100644 index 0000000..7bfde1c --- /dev/null +++ b/plugins/launcher/launch_lammps_optimization_fe.py @@ -0,0 +1,83 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'lammps_optimize@boston' + +############################ +# Define input parameters # +############################ + +# Fe BCC +cell = [[2.848116, 0.000000, 0.000000], + [0.000000, 2.848116, 0.000000], + [0.000000, 0.000000, 2.848116]] + +scaled_positions=[(0.0000000, 0.0000000, 0.0000000), + (0.5000000, 0.5000000, 0.5000000)] + +symbols=['Fe', 'Fe'] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + + +eam_data = {'type': 'fs', + 'file_contents': open('Fe_mm.eam.fs').readlines()} + +potential ={'pair_style': 'eam', + 'data': eam_data} + +lammps_machine = {'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + +parameters_opt = {'relaxation': 'tri', # iso/aniso/tri + 'pressure': 0.0, # kbars + 'vmax': 0.000001, # Angstrom^3 + 'energy_tolerance': 1.0e-25, # eV + 'force_tolerance': 1.0e-25, # eV angstrom + 'max_evaluations': 1000000, + 'max_iterations': 500000} + +code = Code.get_from_string(codename) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources=lammps_machine) + +calc.label = "test lammps calculation" +calc.description = "A much longer description" +calc.use_code(code) +calc.use_structure(structure) +calc.use_potential(ParameterData(dict=potential)) + +calc.use_parameters(ParameterData(dict=parameters_opt)) + +test_only = False + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/launcher/launch_lammps_optimization_gan.py b/plugins/launcher/launch_lammps_optimization_gan.py new file mode 100644 index 0000000..c966a8e --- /dev/null +++ b/plugins/launcher/launch_lammps_optimization_gan.py @@ -0,0 +1,96 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'lammps_optimize@boston' + +############################ +# Define input parameters # +############################ + +# GaN +cell = [[ 3.1900000572, 0, 0], + [-1.5950000286, 2.762621076, 0], + [ 0.0, 0, 5.1890001297]] + + +scaled_positions=[(0.6666669, 0.3333334, 0.0000000), + (0.3333331, 0.6666663, 0.5000000), + (0.6666669, 0.3333334, 0.3750000), + (0.3333331, 0.6666663, 0.8750000)] + +symbols=['Ga', 'Ga', 'N', 'N'] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +# GaN Tersoff +tersoff_gan = {'Ga Ga Ga': '1.0 0.007874 1.846 1.918000 0.75000 -0.301300 1.0 1.0 1.44970 410.132 2.87 0.15 1.60916 535.199', + 'N N N' : '1.0 0.766120 0.000 0.178493 0.20172 -0.045238 1.0 1.0 2.38426 423.769 2.20 0.20 3.55779 1044.77', + 'Ga Ga N' : '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 0.0 0.00000 0.00000 2.90 0.20 0.00000 0.00000', + 'Ga N N' : '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 1.0 2.63906 3864.27 2.90 0.20 2.93516 6136.44', + 'N Ga Ga': '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 1.0 2.63906 3864.27 2.90 0.20 2.93516 6136.44', + 'N Ga N ': '1.0 0.766120 0.000 0.178493 0.20172 -0.045238 1.0 0.0 0.00000 0.00000 2.20 0.20 0.00000 0.00000', + 'N N Ga': '1.0 0.001632 0.000 65.20700 2.82100 -0.518000 1.0 0.0 0.00000 0.00000 2.90 0.20 0.00000 0.00000', + 'Ga N Ga': '1.0 0.007874 1.846 1.918000 0.75000 -0.301300 1.0 0.0 0.00000 0.00000 2.87 0.15 0.00000 0.00000'} + +# Silicon(C) Tersoff +# tersoff_si = {'Si Si Si ': '3.0 1.0 1.7322 1.0039e5 16.218 -0.59826 0.78734 1.0999e-6 1.7322 471.18 2.85 0.15 2.4799 1830.8'} + + +potential ={'pair_style': 'tersoff', + 'data': tersoff_gan} + +lammps_machine = {'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + +parameters_opt = {'relaxation': 'tri', # iso/aniso/tri + 'pressure': 0.0, # kbars + 'vmax': 0.000001, # Angstrom^3 + 'energy_tolerance': 1.0e-25, # eV + 'force_tolerance': 1.0e-25, # eV angstrom + 'max_evaluations': 1000000, + 'max_iterations': 500000} + +code = Code.get_from_string(codename) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources=lammps_machine) + +calc.label = "test lammps calculation" +calc.description = "A much longer description" +calc.use_code(code) +calc.use_structure(structure) +calc.use_potential(ParameterData(dict=potential)) + +calc.use_parameters(ParameterData(dict=parameters_opt)) + +test_only = False + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/launcher/launch_lammps_optimization_lj.py b/plugins/launcher/launch_lammps_optimization_lj.py new file mode 100644 index 0000000..57ee1ea --- /dev/null +++ b/plugins/launcher/launch_lammps_optimization_lj.py @@ -0,0 +1,83 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +import numpy as np + + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +codename = 'lammps_md@boston' + +############################ +# Define input parameters # +############################ + +cell = [[ 3.987594, 0.000000, 0.000000], + [-1.993797, 3.453358, 0.000000], + [ 0.000000, 0.000000, 6.538394]] + +symbols=['Ar'] * 2 +scaled_positions = [(0.33333, 0.66666, 0.25000), + (0.66667, 0.33333, 0.75000)] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +# Example LJ parameters for Argon. These may not be accurate at all +potential ={'pair_style': 'lennard_jones', + # epsilon, sigma, cutoff + 'data': {'1 1': '0.01029 3.4 2.5', + #'2 2': '1.0 1.0 2.5', + #'1 2': '1.0 1.0 2.5' + }} + +lammps_machine = { + 'num_machines': 1, + 'parallel_env': 'mpi*', + 'tot_num_mpiprocs': 16} + +parameters_opt = {'relaxation': 'tri', # iso/aniso/tri + 'pressure': 0.0, # kbars + 'vmax': 0.000001, # Angstrom^3 + 'energy_tolerance': 1.0e-25, # eV + 'force_tolerance': 1.0e-25, # eV angstrom + 'max_evaluations': 1000000, + 'max_iterations': 500000} + +code = Code.get_from_string(codename) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources=lammps_machine) + +calc.label = "test lammps calculation" +calc.description = "A much longer description" +calc.use_code(code) +calc.use_structure(structure) +calc.use_potential(ParameterData(dict=potential)) + +calc.use_parameters(ParameterData(dict=parameters_opt)) + +test_only = False + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/launcher/launch_phonopy_gan.py b/plugins/launcher/launch_phonopy_gan.py new file mode 100644 index 0000000..b97353f --- /dev/null +++ b/plugins/launcher/launch_phonopy_gan.py @@ -0,0 +1,67 @@ +from aiida import load_dbenv +load_dbenv() + +from aiida.orm import Code, DataFactory, load_node +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +import numpy as np +import os + + + +codename = 'phonopy@stern' +code = Code.get_from_string(codename) + +cell = [[ 3.1900000572, 0, 0], + [-1.5950000286, 2.762621076, 0], + [ 0.0, 0, 5.1890001297]] + +scaled_positions=[(0.6666669, 0.3333334, 0.0000000), + (0.3333331, 0.6666663, 0.5000000), + (0.6666669, 0.3333334, 0.3750000), + (0.3333331, 0.6666663, 0.8750000)] + +symbols=['Ga', 'Ga', 'N', 'N'] + +s = StructureData(cell=cell) + +for i, scaled_position in enumerate(scaled_positions): + s.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +parameters = ParameterData(dict={'supercell': [[3,0,0], + [0,3,0], + [0,0,3]], + 'primitive': [[1.0, 0.0, 0.0], + [0.0, 1.0, 0.0], + [0.0, 0.0, 1.0]], + 'distance': 0.01, + 'mesh': [40, 40, 40], + 'symmetry_precision': 1e-5} + ) + +calc = code.new_calc(max_wallclock_seconds=3600, + resources={"num_machines": 1, + "parallel_env":"localmpi", + "tot_num_mpiprocs": 6}) + + +calc.label = "test phonopy calculation" +calc.description = "A much longer description" + +calc.use_structure(s) +calc.use_code(code) +calc.use_parameters(parameters) +calc.use_data_sets(load_node(23913)) + +if False: + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation with PK={}".format(calc.pk) + calc.submit() \ No newline at end of file diff --git a/plugins/launcher/launch_vasp_si.py b/plugins/launcher/launch_vasp_si.py new file mode 100644 index 0000000..d947f39 --- /dev/null +++ b/plugins/launcher/launch_vasp_si.py @@ -0,0 +1,151 @@ +from aiida import load_dbenv +load_dbenv() +from aiida.orm import Code, DataFactory + +from pymatgen.io import vasp as vaspio + +StructureData = DataFactory('structure') +ParameterData = DataFactory('parameter') + +import numpy as np + +codename = 'vasp541mpi@stern' + +############################ +# Define input parameters # +############################ + +a = 5.40400123456789 +cell = [[a, 0, 0], + [0, a, 0], + [0, 0, a]] + +symbols=['Si'] * 8 +scaled_positions = [(0.875, 0.875, 0.875), + (0.875, 0.375, 0.375), + (0.375, 0.875, 0.375), + (0.375, 0.375, 0.875), + (0.125, 0.125, 0.125), + (0.125, 0.625, 0.625), + (0.625, 0.125, 0.625), + (0.625, 0.625, 0.125)] + +structure = StructureData(cell=cell) +positions = np.dot(scaled_positions, cell) + +for i, scaled_position in enumerate(scaled_positions): + structure.append_atom(position=np.dot(scaled_position, cell).tolist(), + symbols=symbols[i]) + +structure.store() + +# VASP input parameters +incar_dict = { + 'NELMIN' : 5, + 'NELM' : 100, + 'ENCUT' : 400, + 'ALGO' : 38, + 'ISMEAR' : 0, + 'SIGMA' : 0.01, + 'GGA' : 'PS' +} + +pseudo_dict = {'functional': 'PBE', + 'symbols': np.unique(symbols).tolist()} # Can be a list of strings with the potwpaw pseudopotentials folder name. Ex: ['Si'] + +# supported_style_modes: "Gamma", "Monkhorst", "Automatic", "Line_mode", "Cartesian" & "Reciprocal" (pymatgen) +kpoints_dict = {'style' : 'Monkhorst', + 'points': [2, 2, 2], + 'shift' : [0.0, 0.0, 0.0]} + +# Cluster information +machine_dict = { + 'num_machines': 1, + 'parallel_env':'mpi*', + 'tot_num_mpiprocs' : 16} + +test_only = True + + +################### +# Set calculation # +################### + +code = Code.get_from_string(codename) +calc = code.new_calc( + max_wallclock_seconds=3600, + resources=machine_dict +) +calc.set_withmpi(True) + +calc.label = 'VASP' +calc.label = 'Silicon VASP' +calc.description = "This is an example calculation of VASP" + +# POSCAR +calc.use_structure(structure) + +incar = vaspio.Incar(incar_dict) +calc.use_incar(ParameterData(dict=incar.as_dict())) + +# KPOINTS +kpoints = kpoints_dict + +kpoints = vaspio.Kpoints(comment='aiida generated', + style=kpoints['style'], + kpts=(kpoints['points'],), kpts_shift=kpoints['shift']) + +calc.use_kpoints(ParameterData(dict=kpoints.as_dict())) + +# POTCAR +pseudo = pseudo_dict +potcar = vaspio.Potcar(symbols=pseudo['symbols'], + functional=pseudo['functional']) +calc.use_potcar(ParameterData(dict=potcar.as_dict())) + +# Define parsers to use +settings = {'PARSER_INSTRUCTIONS': []} +pinstr = settings['PARSER_INSTRUCTIONS'] +pinstr += [{ + 'instr': 'array_data_parser', + 'type': 'data', + 'params': {}}, + { + 'instr': 'output_parameters', + 'type': 'data', + 'params': {}}, + { + 'instr': 'dummy_error_parser', + 'type': 'error', + 'params': {}}, + { + 'instr': 'default_structure_parser', + 'type': 'structure', + 'params': {}} +] + +# additional files to return +settings.setdefault( + 'ADDITIONAL_RETRIEVE_LIST', [ + 'OSZICAR', + 'CONTCAR', + 'OUTCAR', +# 'vasprun.xml' + ] +) +calc.use_settings(ParameterData(dict=settings)) + +if test_only: # It will not be submitted + import os + subfolder, script_filename = calc.submit_test() + print "Test_submit for calculation (uuid='{}')".format(calc.uuid) + print "Submit file in {}".format(os.path.join( + os.path.relpath(subfolder.abspath), + script_filename)) +else: + calc.store_all() + print "created calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) + calc.submit() + print "submitted calculation; calc=Calculation(uuid='{}') # ID={}".format( + calc.uuid, calc.dbnode.pk) diff --git a/plugins/parsers/__init__.py b/plugins/parsers/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/plugins/parsers/dynaphopy.py b/plugins/parsers/dynaphopy.py new file mode 100644 index 0000000..3b9bac4 --- /dev/null +++ b/plugins/parsers/dynaphopy.py @@ -0,0 +1,176 @@ +from aiida.parsers.parser import Parser +from aiida.parsers.exceptions import OutputParsingError +from aiida.orm.data.array import ArrayData +from aiida.orm.data.parameter import ParameterData + +import numpy as np + + +def parse_FORCE_CONSTANTS(filename): + + fcfile = open(filename) + num = int((fcfile.readline().strip().split())[0]) + force_constants = np.zeros((num, num, 3, 3), dtype=float) + for i in range(num): + for j in range(num): + fcfile.readline() + tensor = [] + for k in range(3): + tensor.append([float(x) for x in fcfile.readline().strip().split()]) + force_constants[i, j] = np.array(tensor) + fcfile.close() + return force_constants + + +def parse_quasiparticle_data(qp_file): + import yaml + + f = open(qp_file, "r") + quasiparticle_data = yaml.load(f) + f.close() + return quasiparticle_data + + +def parse_dynaphopy_output(file): + + thermal_properties = None + f = open(file, 'r') + data_lines = f.readlines() + + indices = [] + q_points = [] + for i, line in enumerate(data_lines): + if 'Q-point' in line: + # print i, np.array(line.replace(']', '').replace('[', '').split()[4:8], dtype=float) + indices.append(i) + q_points.append(np.array(line.replace(']', '').replace('[', '').split()[4:8],dtype=float)) + + indices.append(len(data_lines)) + + phonons = {} + for i, index in enumerate(indices[:-1]): + + fragment = data_lines[indices[i]: indices[i+1]] + if 'kipped' in fragment: + continue + # print q_points[i], i + phonon_modes = {} + for j, line in enumerate(fragment): + if 'Peak' in line: + number = line.split()[2] + phonon_mode = {'width': float(fragment[j+2].split()[1]), + 'positions': float(fragment[j+3].split()[1]), + 'shift': float(fragment[j+12].split()[2])} + phonon_modes.update({number: phonon_mode}) + + if 'Thermal' in line: + free_energy = float(fragment[j+4].split()[4]) + entropy = float(fragment[j+5].split()[3]) + cv = float(fragment[j+6].split()[3]) + total_energy = float(fragment[j+7].split()[4]) + + temperature = float(fragment[j].split()[5].replace('(','')) + + thermal_properties = {'temperature': temperature, + 'free_energy': free_energy, + 'entropy': entropy, + 'cv': cv, + 'total_energy': total_energy} + + + phonon_modes.update({'q_point': q_points[i].tolist()}) + + phonons.update({'wave_vector_'+str(i): phonon_modes}) + + f.close() + + return phonons, thermal_properties + + +class DynaphopyParser(Parser): + """ + Simple Parser for LAMMPS. + """ + + def __init__(self, calc): + """ + Initialize the instance of LammpsParser + """ + super(DynaphopyParser, self).__init__(calc) + + def parse_with_retrieved(self, retrieved): + """ + Parses the datafolder, stores results. + """ + + # suppose at the start that the job is successful + successful = True + + # select the folder object + # Check that the retrieved folder is there + try: + out_folder = retrieved[self._calc._get_linkname_retrieved()] + except KeyError: + self.logger.error("No retrieved folder found") + return False, () + + # check what is inside the folder + list_of_files = out_folder.get_folder_list() + + # OUTPUT file should exist + if not self._calc._OUTPUT_FILE_NAME in list_of_files: + successful = False + self.logger.error("Output file not found") + return successful, () + + # Get file and do the parsing + outfile = out_folder.get_abs_path( self._calc._OUTPUT_FILE_NAME) + force_constants_file = out_folder.get_abs_path(self._calc._OUTPUT_FORCE_CONSTANTS) + qp_file = out_folder.get_abs_path(self._calc._OUTPUT_QUASIPARTICLES) + + try: + quasiparticle_data, thermal_properties = parse_dynaphopy_output(outfile) + quasiparticle_data = parse_quasiparticle_data(qp_file) + except ValueError: + pass + + try: + force_constants = parse_FORCE_CONSTANTS(force_constants_file) + except: + pass + + # look at warnings + warnings = [] + with open(out_folder.get_abs_path( self._calc._SCHED_ERROR_FILE )) as f: + errors = f.read() + if errors: + warnings = [errors] + + # ====================== prepare the output node ====================== + + # save the outputs + new_nodes_list = [] + + # save phonon data into node + try: + new_nodes_list.append(('quasiparticle_data', ParameterData(dict=quasiparticle_data))) + except KeyError: # keys not + pass + + try: + new_nodes_list.append(('thermal_properties', ParameterData(dict=thermal_properties))) + except KeyError: # keys not + pass + + try: + array_data = ArrayData() + array_data.set_array('force_constants', force_constants) + new_nodes_list.append(('array_data', array_data)) + except KeyError: # keys not + pass + + + # add the dictionary with warnings + new_nodes_list.append((self.get_linkname_outparams(), ParameterData(dict={'warnings': warnings}))) + + return successful, new_nodes_list diff --git a/plugins/parsers/lammps/__init__.py b/plugins/parsers/lammps/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/plugins/parsers/lammps/force.py b/plugins/parsers/lammps/force.py new file mode 100644 index 0000000..3a293cf --- /dev/null +++ b/plugins/parsers/lammps/force.py @@ -0,0 +1,177 @@ +from aiida.parsers.parser import Parser +from aiida.parsers.exceptions import OutputParsingError + +from aiida.orm import DataFactory + +ArrayData = DataFactory('array') +ParameterData = DataFactory('parameter') + + +import numpy as np + +def read_log_file(logfile): + + f = open(logfile, 'r') + data = f.readlines() + + data_dict = {} + for i, line in enumerate(data): + if 'Loop time' in line: + energy = float(data[i-1].split()[4]) + data_dict['energy'] = energy + + return data_dict + +def read_lammps_forces(file_name): + + import mmap + # Time in picoseconds + # Coordinates in Angstroms + + # Starting reading + + # Dimensionality of LAMMP calculation + number_of_dimensions = 3 + + cells = [] + + with open(file_name, "r+") as f: + + file_map = mmap.mmap(f.fileno(), 0) + + # Read time steps + position_number=file_map.find('TIMESTEP') + + file_map.seek(position_number) + file_map.readline() + + + #Read number of atoms + position_number=file_map.find('NUMBER OF ATOMS') + file_map.seek(position_number) + file_map.readline() + number_of_atoms = int(file_map.readline()) + + + #Read cell + position_number=file_map.find('ITEM: BOX') + file_map.seek(position_number) + file_map.readline() + + + bounds = [] + for i in range(3): + bounds.append(file_map.readline().split()) + + bounds = np.array(bounds, dtype=float) + if bounds.shape[1] == 2: + bounds = np.append(bounds, np.array([0, 0, 0])[None].T ,axis=1) + + + xy = bounds[0, 2] + xz = bounds[1, 2] + yz = bounds[2, 2] + + xlo = bounds[0, 0] - np.min([0.0, xy, xz, xy+xz]) + xhi = bounds[0, 1] - np.max([0.0, xy, xz, xy+xz]) + ylo = bounds[1, 0] - np.min([0.0, yz]) + yhi = bounds[1, 1] - np.max([0.0, yz]) + zlo = bounds[2, 0] + zhi = bounds[2, 1] + + super_cell = np.array([[xhi-xlo, xy, xz], + [0, yhi-ylo, yz], + [0, 0, zhi-zlo]]) + + cells.append(super_cell.T) + + position_number = file_map.find('ITEM: ATOMS') + file_map.seek(position_number) + file_map.readline() + + #Reading forces + forces = [] + read_elements = [] + for i in range (number_of_atoms): + line = file_map.readline().split()[0:number_of_dimensions+1] + forces.append(line[1:number_of_dimensions+1]) + read_elements.append(line[0]) + + file_map.close() + + forces = np.array([forces], dtype=float) + + return forces + + +class ForceParser(Parser): + """ + Simple Parser for LAMMPS. + """ + + def __init__(self, calc): + """ + Initialize the instance of LammpsParser + """ + super(ForceParser, self).__init__(calc) + + def parse_with_retrieved(self, retrieved): + """ + Parses the datafolder, stores results. + """ + + # suppose at the start that the job is successful + successful = True + + # select the folder object + # Check that the retrieved folder is there + try: + out_folder = retrieved[self._calc._get_linkname_retrieved()] + except KeyError: + self.logger.error("No retrieved folder found") + return False, () + + # check what is inside the folder + list_of_files = out_folder.get_folder_list() + + # OUTPUT file should exist + if not self._calc._OUTPUT_FILE_NAME in list_of_files: + successful = False + self.logger.error("Output file not found") + return successful, () + + # Get file and do the parsing + outfile = out_folder.get_abs_path( self._calc._OUTPUT_FILE_NAME) + ouput_trajectory = out_folder.get_abs_path( self._calc._OUTPUT_TRAJECTORY_FILE_NAME) + + outputa_data = read_log_file(outfile) + forces = read_lammps_forces(ouput_trajectory) + + # look at warnings + warnings = [] + with open(out_folder.get_abs_path( self._calc._SCHED_ERROR_FILE )) as f: + errors = f.read() + if errors: + warnings = [errors] + + # ====================== prepare the output node ====================== + + # save the outputs + new_nodes_list = [] + + # save trajectory into node + + array_data = ArrayData() + array_data.set_array('forces', forces) + new_nodes_list.append(('output_array', array_data)) + + # add the dictionary with warnings + outputa_data.update({'warnings': warnings}) + + parameters_data = ParameterData(dict=outputa_data) + new_nodes_list.append((self.get_linkname_outparams(), parameters_data)) + + # add the dictionary with warnings + # new_nodes_list.append((self.get_linkname_outparams(), ParameterData(dict={'warnings': warnings}))) + + return successful, new_nodes_list diff --git a/plugins/parsers/lammps/md.py b/plugins/parsers/lammps/md.py new file mode 100644 index 0000000..0c47eb8 --- /dev/null +++ b/plugins/parsers/lammps/md.py @@ -0,0 +1,209 @@ +from aiida.parsers.parser import Parser +from aiida.parsers.exceptions import OutputParsingError +from aiida.orm.data.array import ArrayData +from aiida.orm.data.parameter import ParameterData +from aiida.orm.data.array.trajectory import TrajectoryData + +import numpy as np + +def read_lammps_trajectory(file_name, + limit_number_steps=100000000, + initial_cut=1, + end_cut=None, + timestep=1): + + import mmap + #Time in picoseconds + #Coordinates in Angstroms + + #Starting reading + print("Reading LAMMPS trajectory") + print("This could take long, please wait..") + + #Dimensionality of LAMMP calculation + number_of_dimensions = 3 + + step_ids = [] + data = [] + cells = [] + counter = 0 + bounds = None + number_of_atoms = None + + lammps_labels = False + + with open(file_name, "r+") as f: + + file_map = mmap.mmap(f.fileno(), 0) + + while True: + + counter += 1 + + #Read time steps + position_number=file_map.find('TIMESTEP') + if position_number < 0: break + + file_map.seek(position_number) + file_map.readline() + step_ids.append(float(file_map.readline())) + + if number_of_atoms is None: + #Read number of atoms + position_number=file_map.find('NUMBER OF ATOMS') + file_map.seek(position_number) + file_map.readline() + number_of_atoms = int(file_map.readline()) + + if True: + #Read cell + position_number=file_map.find('ITEM: BOX') + file_map.seek(position_number) + file_map.readline() + + + bounds = [] + for i in range(3): + bounds.append(file_map.readline().split()) + + bounds = np.array(bounds, dtype=float) + if bounds.shape[1] == 2: + bounds = np.append(bounds, np.array([0, 0, 0])[None].T ,axis=1) + + xy = bounds[0, 2] + xz = bounds[1, 2] + yz = bounds[2, 2] + + xlo = bounds[0, 0] - np.min([0.0, xy, xz, xy+xz]) + xhi = bounds[0, 1] - np.max([0.0, xy, xz, xy+xz]) + ylo = bounds[1, 0] - np.min([0.0, yz]) + yhi = bounds[1, 1] - np.max([0.0, yz]) + zlo = bounds[2, 0] + zhi = bounds[2, 1] + + super_cell = np.array([[xhi-xlo, xy, xz], + [0, yhi-ylo, yz], + [0, 0, zhi-zlo]]) + + cells.append(super_cell.T) + + + position_number = file_map.find('ITEM: ATOMS') + file_map.seek(position_number) + lammps_labels=file_map.readline() + + + #Initial cut control + if initial_cut > counter: + continue + + #Reading coordinates + read_coordinates = [] + read_elements = [] + for i in range (number_of_atoms): + line = file_map.readline().split()[0:number_of_dimensions+1] + read_coordinates.append(line[1:number_of_dimensions+1]) + read_elements.append(line[0]) + try: + data.append(np.array(read_coordinates, dtype=float)) #in angstroms + # print read_coordinates + except ValueError: + print("Error reading step {0}".format(counter)) + break + # print(read_coordinates) + #security routine to limit maximum of steps to read and put in memory + if limit_number_steps+initial_cut < counter: + print("Warning! maximum number of steps reached! No more steps will be read") + break + + if end_cut is not None and end_cut <= counter: + break + + + file_map.close() + + data = np.array(data) + step_ids = np.array(step_ids, dtype=int) + cells = np.array(cells) + elements = np.array(read_elements) + time = np.array(step_ids)*timestep + return data, step_ids, cells, elements, time + + +class MdParser(Parser): + """ + Simple Parser for LAMMPS. + """ + + def __init__(self, calc): + """ + Initialize the instance of LammpsParser + """ + super(MdParser, self).__init__(calc) + + def parse_with_retrieved(self, retrieved): + """ + Parses the datafolder, stores results. + """ + + # suppose at the start that the job is successful + successful = True + + # select the folder object + # Check that the retrieved folder is there + try: + out_folder = retrieved[self._calc._get_linkname_retrieved()] + except KeyError: + self.logger.error("No retrieved folder found") + return False, () + + # check what is inside the folder + list_of_files = out_folder.get_folder_list() + + # OUTPUT file should exist + if not self._calc._OUTPUT_FILE_NAME in list_of_files: + successful = False + self.logger.error("Output file not found") + return successful, () + + # Get file and do the parsing + outfile = out_folder.get_abs_path( self._calc._OUTPUT_FILE_NAME) + ouput_trajectory = out_folder.get_abs_path( self._calc._OUTPUT_TRAJECTORY_FILE_NAME) + + timestep = self._calc.inp.parameters.dict.timestep + positions, step_ids, cells, symbols, time = read_lammps_trajectory(ouput_trajectory, timestep=timestep) + + # Delete trajectory once parsed + try: + import os + os.remove(ouput_trajectory) + except: + pass + +# force_constants = parse_FORCE_CONSTANTS(outfile) + + # look at warnings + warnings = [] + with open(out_folder.get_abs_path( self._calc._SCHED_ERROR_FILE )) as f: + errors = f.read() + if errors: + warnings = [errors] + + # ====================== prepare the output node ====================== + + # save the outputs + new_nodes_list = [] + + # save trajectory into node + try: + trajectory_data = TrajectoryData() + trajectory_data.set_trajectory(step_ids, cells, symbols, positions, times=time) + new_nodes_list.append(('trajectory_data', trajectory_data)) + except KeyError: # keys not found in json + pass + + # add the dictionary with warnings + new_nodes_list.append((self.get_linkname_outparams(), ParameterData(dict={'warnings': warnings}))) + + return successful, new_nodes_list + diff --git a/plugins/parsers/lammps/optimize.py b/plugins/parsers/lammps/optimize.py new file mode 100644 index 0000000..7cbbf82 --- /dev/null +++ b/plugins/parsers/lammps/optimize.py @@ -0,0 +1,232 @@ +from aiida.parsers.parser import Parser +from aiida.parsers.exceptions import OutputParsingError +from aiida.orm.data.array import ArrayData +from aiida.orm.data.structure import StructureData +from aiida.orm.data.parameter import ParameterData + +import numpy as np + +def read_log_file(logfile): + + f = open(logfile, 'r') + data = f.readlines() + + data_dict = {} + for i, line in enumerate(data): + if 'Loop time' in line: + energy = float(data[i-1].split()[4]) + data_dict['energy'] = energy + xx, yy, zz, xy, xz, yz = data[i-1].split()[5:11] + stress = np.array([[xx, xy, xz], + [xy, yy, yz], + [xz, yz, zz]], dtype=float) + + + if '$(xlo)' in line: + a = data[i+1].split() + if '$(ylo)' in line: + b = data[i+1].split() + if '$(zlo)' in line: + c = data[i+1].split() + + bounds = np.array([a, b, c], dtype=float) + + # lammps_input_file += 'print "$(xlo) $(xhi) $(xy)"\n' + # lammps_input_file += 'print "$(ylo) $(yhi) $(xz)"\n' + # lammps_input_file += 'print "$(zlo) $(zhi) $(yz)"\n' + + + xy = bounds[0, 2] + xz = bounds[1, 2] + yz = bounds[2, 2] + + xlo = bounds[0, 0] + xhi = bounds[0, 1] + ylo = bounds[1, 0] + yhi = bounds[1, 1] + zlo = bounds[2, 0] + zhi = bounds[2, 1] + + super_cell = np.array([[xhi-xlo, xy, xz], + [0, yhi-ylo, yz], + [0, 0, zhi-zlo]]) + + cell=super_cell.T + + if np.linalg.det(cell) < 0: + cell = -1.0*cell + + volume = np.linalg.det(cell) + stress = -stress/volume * 1.e-3 # bar*A^3 -> kbar + + return data_dict, cell, stress + +def read_lammps_positions_and_forces(file_name): + + import mmap + # Time in picoseconds + # Coordinates in Angstroms + + # Starting reading + + # Dimensionality of LAMMP calculation + number_of_dimensions = 3 + + + with open(file_name, "r+") as f: + + file_map = mmap.mmap(f.fileno(), 0) + + # Read time steps + while True: + position_number=file_map.find('TIMESTEP') + try: + file_map.seek(position_number) + file_map.readline() + except ValueError: + break + + #Read number of atoms + position_number=file_map.find('NUMBER OF ATOMS') + file_map.seek(position_number) + file_map.readline() + number_of_atoms = int(file_map.readline()) + + + #Read cell + position_number=file_map.find('ITEM: BOX') + file_map.seek(position_number) + file_map.readline() + + bounds = [] + for i in range(3): + bounds.append(file_map.readline().split()) + + bounds = np.array(bounds, dtype=float) + if bounds.shape[1] == 2: + bounds = np.append(bounds, np.array([0, 0, 0])[None].T ,axis=1) + + xy = bounds[0, 2] + xz = bounds[1, 2] + yz = bounds[2, 2] + + xlo = bounds[0, 0] - np.min([0.0, xy, xz, xy+xz]) + xhi = bounds[0, 1] - np.max([0.0, xy, xz, xy+xz]) + ylo = bounds[1, 0] - np.min([0.0, yz]) + yhi = bounds[1, 1] - np.max([0.0, yz]) + zlo = bounds[2, 0] + zhi = bounds[2, 1] + + super_cell = np.array([[xhi-xlo, xy, xz], + [0, yhi-ylo, yz], + [0, 0, zhi-zlo]]) + + cell=super_cell.T + + + + + + position_number = file_map.find('ITEM: ATOMS') + file_map.seek(position_number) + file_map.readline() + + #Reading positions + positions = [] + forces = [] + read_elements = [] + for i in range (number_of_atoms): + line = file_map.readline().split()[0:number_of_dimensions*2+1] + positions.append(line[1:number_of_dimensions+1]) + forces.append(line[1+number_of_dimensions:number_of_dimensions*2+1]) + read_elements.append(line[0]) + + file_map.close() + + positions = np.array([positions]) + forces = np.array([forces], dtype=float) + + return positions, forces, read_elements, cell + + +class OptimizeParser(Parser): + """ + Simple Parser for LAMMPS. + """ + + def __init__(self, calc): + """ + Initialize the instance of LammpsParser + """ + super(OptimizeParser, self).__init__(calc) + + def parse_with_retrieved(self, retrieved): + """ + Parses the datafolder, stores results. + """ + + # suppose at the start that the job is successful + successful = True + + # select the folder object + # Check that the retrieved folder is there + try: + out_folder = retrieved[self._calc._get_linkname_retrieved()] + except KeyError: + self.logger.error("No retrieved folder found") + return False, () + + # check what is inside the folder + list_of_files = out_folder.get_folder_list() + + # OUTPUT file should exist + if not self._calc._OUTPUT_FILE_NAME in list_of_files: + successful = False + self.logger.error("Output file not found") + return successful, () + + # Get file and do the parsing + outfile = out_folder.get_abs_path( self._calc._OUTPUT_FILE_NAME) + ouput_trajectory = out_folder.get_abs_path( self._calc._OUTPUT_TRAJECTORY_FILE_NAME) + + output_data, cell, stress_tensor = read_log_file(outfile) + positions, forces, symbols, cell2 = read_lammps_positions_and_forces(ouput_trajectory) + + # look at warnings + warnings = [] + with open(out_folder.get_abs_path( self._calc._SCHED_ERROR_FILE )) as f: + errors = f.read() + if errors: + warnings = [errors] + + # ====================== prepare the output node ====================== + + # save the outputs + new_nodes_list = [] + + # save optimized structure into node + structure = StructureData(cell=cell) + + for i, position in enumerate(positions[-1]): + structure.append_atom(position=position.tolist(), + symbols=symbols[i]) + + new_nodes_list.append(('output_structure', structure)) + + # save forces into node + array_data = ArrayData() + array_data.set_array('forces', forces) + array_data.set_array('stress', stress_tensor) + + new_nodes_list.append(('output_array', array_data)) + + # add the dictionary with warnings + output_data.update({'warnings': warnings}) + + parameters_data = ParameterData(dict=output_data) + new_nodes_list.append((self.get_linkname_outparams(), parameters_data)) + + # add the dictionary with warnings + # new_nodes_list.append((self.get_linkname_outparams(), ParameterData(dict={'warnings': warnings}))) + + return successful, new_nodes_list diff --git a/plugins/parsers/phonopy.py b/plugins/parsers/phonopy.py new file mode 100644 index 0000000..41ac62c --- /dev/null +++ b/plugins/parsers/phonopy.py @@ -0,0 +1,90 @@ +from aiida.parsers.parser import Parser +from aiida.parsers.exceptions import OutputParsingError +from aiida.orm.data.array import ArrayData +from aiida.orm.data.force_constants import ForceConstants + +from aiida.orm.data.parameter import ParameterData + +import numpy as np + + +def parse_FORCE_CONSTANTS(filename): + fcfile = open(filename) + num = int((fcfile.readline().strip().split())[0]) + force_constants = np.zeros((num, num, 3, 3), dtype=float) + for i in range(num): + for j in range(num): + fcfile.readline() + tensor = [] + for k in range(3): + tensor.append([float(x) for x in fcfile.readline().strip().split()]) + force_constants[i, j] = np.array(tensor) + return force_constants + + +class PhonopyParser(Parser): + """ + Parser the FORCE_CONSTANTS file of phonopy. + """ + + def __init__(self, calc): + """ + Initialize the instance of PhonopyParser + """ + super(PhonopyParser, self).__init__(calc) + + def parse_with_retrieved(self, retrieved): + """ + Parses the datafolder, stores results. + """ + + # suppose at the start that the job is successful + successful = True + + # select the folder object + # Check that the retrieved folder is there + try: + out_folder = retrieved[self._calc._get_linkname_retrieved()] + except KeyError: + self.logger.error("No retrieved folder found") + return False, () + + # check what is inside the folder + list_of_files = out_folder.get_folder_list() + + # OUTPUT file should exist + if not self._calc._OUTPUT_FILE_NAME in list_of_files: + successful = False + self.logger.error("Output file not found") + return successful, () + + # Get file and do the parsing + outfile = out_folder.get_abs_path(self._calc._OUTPUT_FILE_NAME) + + force_constants = parse_FORCE_CONSTANTS(outfile) + + # look at warnings + warnings = [] + with open(out_folder.get_abs_path(self._calc._SCHED_ERROR_FILE)) as f: + errors = f.read() + if errors: + warnings = [errors] + + # ====================== prepare the output node ====================== + + # save the outputs + new_nodes_list = [] + + # save force constants into node + try: + #array_data = ArrayData() + #array_data.set_array('force_constants', force_constants) + + new_nodes_list.append(('force_constants', ForceConstants(array=force_constants))) + except KeyError: # keys not found in json + pass + + # add the dictionary with warnings + new_nodes_list.append((self.get_linkname_outparams(), ParameterData(dict={'warnings': warnings}))) + + return successful, new_nodes_list diff --git a/plugins/parsers/phonopy_dos.py b/plugins/parsers/phonopy_dos.py new file mode 100644 index 0000000..e69de29 diff --git a/setup.py b/setup.py new file mode 100644 index 0000000..a785c1f --- /dev/null +++ b/setup.py @@ -0,0 +1,28 @@ +from setuptools import setup, find_packages + +setup( + name='aiida-lammps', + version='0.1', + description='AiiDA plugin for LAMMPS', + url='https://github.com/abelcarreras/aiida_extensions', + author='Abel Carreras', + author_email='abelcarreras83@gmail.com', + license='MIT license', + packages=find_packages(exclude=['aiida']), + requires=['phonopy', 'numpy', 'dynaphopy'], + setup_requires=['reentry'], + reentry_register=True, + entry_points={ + 'aiida.calculations': [ + 'lammps.combinate = plugins.jobs.lammps.combinate:CombinateCalculation', + 'lammps.force = plugins.jobs.lammps.force:ForceCalculation', + 'lammps.md = plugins.jobs.lammps.md:MdCalculation', + 'lammps.optimize = plugins.jobs.lammps.optimize:OptimizeCalculation', + 'dynaphopy = plugins.jobs.dynaphopy: DynaphopyCalculation'], + 'aiida.parsers': [ + 'lammps.force = plugins.parsers.lammps.force:ForceParser', + 'lammps.md = plugins.parsers.lammps.md:MdParser', + 'lammps.optimize = plugins.parsers.lammps.optimize:OptimizeParser', + 'dynaphopy = plugins.parsers.dynaphopy: DynaphopyParser'] + } + ) \ No newline at end of file