PyironAdsorptionGener.py

# This is the Script for adding adsorbates (including H, O or water molecules on the surface)
# To be Continued ...
# Author: Bingxin Li
# Date: 19/08/2024
# Contact: b.li@mpie.de

import numpy as np
import copy as cp
from pyiron.project import Project

class Add_Adsorbate():

    """
    Class for constructing adsorption molecules/radicals/atoms on the Metal/HEAs Surface

    Attributes:
            @Prim_cell(): Find the lattice parameter of the smallest surface unit cell
            @SurfSite_Deter(): Determine the poistion of the specific adsorption site
            @AddAdsorbate(): Add the target adsorbate on the surface and generate a new Pyiron structure                                                         

    """

    def __init__(self, surface, proj='pr', adsorbate='H'):

        '''
        @proj:input the Pyiron project here
        @surface:input the surface structure (from Pyiron) into the class
        @basis:input the structure generated by the Pyiron modules
        '''
        self.project = proj
        self.surface =cp.deepcopy(surface)
        self.adsorbate = adsorbate

        return

    def Prim_cell(self, species='Pt', type='fcc111', a=3.99, vac=12):

        self.prim_surf = self.project.create_surface(species, surface_type=type, a=a, vacuum=vac)
        self.prim_lattice = self.prim_surf.cell[:2,:2]

        return

    def SurfSite_Deter(self, site='ontop', z_bias=0, offsite=np.array([0,0])):
        self.site = site
        self.z_bias = z_bias
        self.offsite= offsite

        if self.site == 'ontop':
            self.z_bias = 2
            self.pos = np.append(self.prim_lattice.T @ (np.array([0, 0])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_bias)
        elif self.site == 'fcc':
            self.pos = np.append(self.prim_lattice.T @ (np.array([1/3, 1/3])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_bias)
        elif self.site == 'hcp':
            self.pos = np.append(self.prim_lattice.T @ (np.array([2/3, 2/3])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_bias)
        elif self.site == 'bridge':
            self.pos = np.append(self.prim_lattice.T @ (np.array([1/2, 0])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_bias)
        elif self.site == 'tetra1':
            self.z_t1_bias = - self.prim_lattice[0][0] * (np.sqrt(6) / 3 - np.sqrt(6) / 4)
            self.pos = np.append(self.prim_lattice.T @ (np.array([2/3, 2/3])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_t1_bias)
        elif self.site == 'tetra2':
            self.z_t2_bias = - self.prim_lattice[0][0] * np.sqrt(6) / 4
            self.pos = np.append(self.prim_lattice.T @ (np.array([0, 0])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_t2_bias)
        elif self.site == 'octa':
            self.z_o_bias = - self.prim_lattice[0][0] / 2
            self.pos = np.append(self.prim_lattice.T @ (np.array([1/3, 1/3])+self.offsite), np.max(self.surface.positions[:, 2]) + self.z_o_bias)
        else:
            raise ValueError(f"Adsorption site type {self.site} not recognized.")

        return(self.pos)

    def AddAdsorbate(self):
        # Consider the PBC in the unit cell
        if self.pos[0] >= self.surface.cell[0][0]:
            self.pos[0] = self.pos[0] - self.surface.cell[0][0]

        if self.pos[1] >= self.surface.cell[1][1]:
            self.pos[1] = self.pos[1] - self.Surface.cell[1][1] 

        self.ads_atoms = self.project.create_atoms(self.adsorbate, positions=[self.pos], cell=self.surface.cell)
        # self.newstruct = self.surface + self.ads_atoms

        return(self.ads_atoms) # Note here we only return a unit cell with the addsorbate, we should use surface + ads_atoms as the final structure