Source code for pyiron.vasp.metadyn

# coding: utf-8
# Copyright (c) Max-Planck-Institut für Eisenforschung GmbH - Computational Materials Design (CM) Department
# Distributed under the terms of "New BSD License", see the LICENSE file.

import numpy as np
from pyiron_base import GenericParameters
from pyiron.vasp.vasp import Vasp
from pyiron.vasp.base import Input, Output
from pyiron.vasp.report import Report
import os
import posixpath

__author__ = "Sudarsan Surendralal"
__copyright__ = (
    "Copyright 2020, Max-Planck-Institut für Eisenforschung GmbH - "
    "Computational Materials Design (CM) Department"
)
__version__ = "1.0"
__maintainer__ = "Sudarsan Surendralal"
__email__ = "surendralal@mpie.de"
__status__ = "testing"
__date__ = "March 1, 2020"


[docs]class VaspMetadyn(Vasp): """ Class to setup and run and analyze VASP and VASP metadynamics simulations. For more details see the appropriate `VASP documentation`_ This class is a derivative of pyiron.objects.job.generic.GenericJob. The functions in these modules are written in such the function names and attributes are very generic (get_structure(), molecular_dynamics(), version) but the functions are written to handle VASP specific input/output. Args: project (pyiron.project.Project instance): Specifies the project path among other attributes job_name (str): Name of the job Examples: Let's say you need to run a vasp simulation where you would like to control the input parameters manually. To set up a static dft run with Gaussian smearing and a k-point MP mesh of [6, 6, 6]. You would have to set it up as shown below: >>> ham = Vasp(job_name="trial_job") >>> ham.input.incar[IBRION] = -1 >>> ham.input.incar[ISMEAR] = 0 >>> ham.input.kpoints.set(size_of_mesh=[6, 6, 6]) However, the according to pyiron's philosophy, it is recommended to avoid using code specific tags like IBRION, ISMEAR etc. Therefore the recommended way to set this calculation is as follows: >>> ham = Vasp(job_name="trial_job") >>> ham.calc_static() >>> ham.set_occupancy_smearing(smearing="gaussian") >>> ham.set_kpoints(mesh=[6, 6, 6]) The exact same tags as in the first examples are set automatically. .. _VASP documentation: https://cms.mpi.univie.ac.at/vasp/vasp/Advanced_MD_techniques.html """ def __init__(self, project, job_name): super(VaspMetadyn, self).__init__(project, job_name) self.__name__ = "VaspMetadyn" self.input = MetadynInput() self._output_parser = MetadynOutput() self.supported_primitive_constraints = ["bond", "angle", "torsion", "x_pos", "y_pos", "z_pos"] self.supported_complex_constraints = ["linear_combination", "norm", "coordination_number"] self._constraint_dict = dict() self._complex_constraints = dict() # This is necessary to write the constrained dynamics output to the REPORT file self.input.incar["LBLUEOUT"] = True
[docs] def set_primitive_constraint(self, name, constraint_type, atom_indices, biased=False, increment=0.0): """ Function to set primitive geometric constraints in VASP. Args: name (str): Name of the constraint constraint_type (str): Type of the constraint (has to be part of `supported_primitive_constraints` atom_indices (int/list/numpy.ndarray): Indices of the atoms for which the constraint should be applied biased (bool): True if potential bias is to be applied (biased MD and metadynamics calculations) increment (float): Increment in the constraint variable at every time step """ if self.structure is None: raise ValueError("The structure has to be set before a dynamic constraint is assigned") self._constraint_dict[name] = dict() if constraint_type in self.supported_primitive_constraints: self._constraint_dict[name]["constraint_type"] = constraint_type else: raise ValueError("The constraint type '{}' is not supported".format(constraint_type)) self._constraint_dict[name]["atom_indices"] = atom_indices self._constraint_dict[name]["biased"] = biased self._constraint_dict[name]["increment"] = increment
def _set_primitive_constraint_iconst(self, constraint_type, atom_indices, biased=False): if self.structure is None: raise ValueError("The structure has to be set before a dynamic constraint is assigned") constraint_dict = {"bond": "R", "angle": "A", "torsion": "T", "x_pos": "X", "y_pos": "Y", "z_pos": "Z"} if constraint_type not in list(constraint_dict.keys()): raise ValueError("Use a compatible constraint type") line = len(self.input.iconst._dataset["Value"]) status = 0 if biased: status = 7 if constraint_type in ["x_pos", "y_pos", "z_pos"]: if isinstance(atom_indices, (list, np.ndarray)): a_ind = str(self.sorted_indices[atom_indices[0]] + 1) else: a_ind = str(self.sorted_indices[atom_indices] + 1) elif constraint_type in ["bond"]: if len(atom_indices) != 2: raise ValueError("For this constraint the atom_indices must be a list or numpy array with 2 values") a_ind = ' '.join(map(str, self.sorted_indices[atom_indices] + 1)) elif constraint_type in ["angle"]: if len(atom_indices) != 3: raise ValueError("For this constraint the atom_indices must be a list or numpy array with 3 values") a_ind = ' '.join(map(str, self.sorted_indices[atom_indices] + 1)) else: raise ValueError("The constraint {} is not implemented!".format(constraint_type)) constraint_string = "{} {} {}".format(constraint_dict[constraint_type], a_ind, status) self.input.iconst.set_value(line, constraint_string)
[docs] def set_complex_constraint(self, name, constraint_type, coefficient_dict, biased=False, increment=0.0): """ Set complex constraints based on defined primitive constraints. Args: name (str): Name of the complex constraint constraint_type (str): Type of the constraint (has to be part of `supported_complex_constraints` coefficient_dict (dict): Dictionary containing the primitive constraint name as the key and the constraint coefficient as the values biased (bool): True if potential bias is to be applied (biased MD and metadynamics calculations) increment (float): Increment in the constraint variable at every time step """ if self.structure is None: raise ValueError("The structure has to be set before a dynamic constraint is assigned") if constraint_type not in self.supported_complex_constraints: raise ValueError("Use a compatible complex constraint type") self._complex_constraints[name] = dict() self._complex_constraints[name]["constraint_type"] = constraint_type if len(list(coefficient_dict.keys())) != len(list(self._constraint_dict.keys())): raise ValueError("The number of coefficients should match the number of primitive contstraints") self._complex_constraints[name]["coefficient_dict"] = coefficient_dict self._complex_constraints[name]["biased"] = biased self._complex_constraints[name]["increment"] = increment
def _set_complex_constraint(self, constraint_type, coefficients, biased=False): constraint_dict = {"linear_combination": "S", "norm": "C", "coordination_number": "D"} if constraint_type not in list(constraint_dict.keys()): raise ValueError("Use a compatible constraint type") status = 0 if biased: status = 5 coeffs = " ".join(map(str, coefficients)) constraint_string = "{} {} {}".format(constraint_dict[constraint_type], coeffs, status) line = len(self.input.iconst._dataset["Value"]) self.input.iconst.set_value(line, constraint_string)
[docs] def write_constraints(self): """ Function to write constraint parameters in the INCAR and ICONST files """ increment_list = list() linear_constraint_order = list() # Clear the existing ICONST file self.input.iconst.clear_all() for key, constraint in self._constraint_dict.items(): linear_constraint_order.append(key) self._set_primitive_constraint_iconst(constraint_type=constraint["constraint_type"], atom_indices=constraint["atom_indices"], biased=constraint["biased"]) increment_list.append(constraint["increment"]) for constraint in self._complex_constraints.values(): coefficients = [constraint["coefficient_dict"][val] for val in linear_constraint_order] self._set_complex_constraint(constraint_type=constraint["constraint_type"], coefficients=coefficients, biased=constraint["biased"]) increment_list.append(constraint["increment"]) self.input.incar["INCREM"] = " ".join(map(str, increment_list))
[docs] def write_input(self): """ Call routines that generate the INCAR, POTCAR, KPOINTS and POSCAR input files """ self.write_constraints() self.to_hdf() super(VaspMetadyn, self).write_input()
[docs] def to_hdf(self, hdf=None, group_name=None): self.input._constraint_dict = self._constraint_dict self.input._complex_constraints = self._complex_constraints super(VaspMetadyn, self).to_hdf(hdf=hdf, group_name=group_name)
[docs] def from_hdf(self, hdf=None, group_name=None): super(VaspMetadyn, self).from_hdf(hdf=hdf, group_name=group_name) self._constraint_dict = self.input._constraint_dict self._complex_constraints = self.input._complex_constraints
[docs]class MetadynInput(Input): def __init__(self): super(MetadynInput, self).__init__() self.iconst = GenericParameters(input_file_name=None, table_name="iconst", val_only=True, comment_char="!") self.penaltypot = GenericParameters(input_file_name=None, table_name="penaltypot", val_only=True, comment_char="!") self._constraint_dict = dict() self._complex_constraints = dict()
[docs] def write(self, structure, modified_elements, directory=None): """ Writes all the input files to a specified directory Args: structure (atomistics.structure.atoms.Atoms instance): Structure to be written directory (str): The working directory for the VASP run """ # Writing the constraints, increments, and penalty potentials super(MetadynInput, self).write(structure, modified_elements, directory) self.iconst.write_file(file_name="ICONST", cwd=directory) self.penaltypot.write_file(file_name="PENALTYPOT", cwd=directory)
[docs] def to_hdf(self, hdf): super(MetadynInput, self).to_hdf(hdf) with hdf.open("input") as hdf5_input: self.iconst.to_hdf(hdf5_input) self.penaltypot.to_hdf(hdf5_input) hdf5_input["constraint_dict"] = self._constraint_dict hdf5_input["complex_constraint_dict"] = self._complex_constraints
[docs] def from_hdf(self, hdf): super(MetadynInput, self).from_hdf(hdf) with hdf.open("input") as hdf5_input: self.iconst.from_hdf(hdf5_input) self.penaltypot.from_hdf(hdf5_input) if "constraint_dict" in hdf5_input.list_nodes(): self._constraint_dict = hdf5_input["constraint_dict"] if "complex_constraint_dict" in hdf5_input.list_nodes(): self._complex_constraints = hdf5_input["complex_constraint_dict"]
[docs]class MetadynOutput(Output): def __init__(self): super(MetadynOutput, self).__init__() self.report_file = Report()
[docs] def collect(self, directory=os.getcwd(), sorted_indices=None): super(MetadynOutput, self).collect(directory=directory, sorted_indices=sorted_indices) files_present = os.listdir(directory) if "REPORT" in files_present: self.report_file.from_file(filename=posixpath.join(directory, "REPORT"))
[docs] def to_hdf(self, hdf): """ Save the object in a HDF5 file Args: hdf (pyiron_base.generic.hdfio.ProjectHDFio): HDF path to which the object is to be saved """ super(MetadynOutput, self).to_hdf(hdf) with hdf.open("output") as hdf5_output: with hdf5_output.open("constrained_md") as hdf5_const: for key, value in self.report_file.parse_dict.items(): hdf5_const[key] = value