Search for adsorption structure using Minima Hopping method

Using the Minima Hopping method implemented in ASE (Atomic Simulation Environment), we will explore the adsorption structure of H 2 </ sub> molecules on the Pt surface. This article was created with reference to the Official Documentation Tutorial.

The reason why EMT is used is because of the calculation speed. Especially for organic molecules (C, H, O, N), as mentioned in the Official Documents, It's just a play model. It is not reliable. Select an appropriate theory when using it in practice.

Create a slab model

Create a model in which H 2 </ sub> molecules are adsorbed on the surface of Pt (211).

from ase.build import surface, add_adsorbate
from ase.data import atomic_numbers, covalent_radii
from ase.visualize import view

#Pt bulk
Pt_bulk= bulk('Pt', 'fcc', a=3.9, cubic=True) 

# Pt(211)Slab
Pt_211 = surface(Pt_bulk, (2,1,1), 3, vacuum=10, periodic=True) 
d_PtH = (covalent_radii[atomic_numbers['H']] + covalent_radii[atomic_numbers['Pt']])
z_max =  np.max(Pt_211.positions[:,2]) #Maximum z-coordinate value of Pt atom
xy_center = ((Pt_211.cell[0]+Pt_211.cell[1])/2)[0:2]

#Add H2 molecule
add_adsorbate(Pt_211, molecule('H2'), d_PtH, xy_center)

#Check the structure
view(Pt_211)

fig1.png

Set constraint conditions

Use the FixAtoms and Hookean classes to set constraints. Under the constraint condition of Hookean, a spring-like force is applied between atoms.

--Pt slab is completely fixed --H 2 </ sub> Molecules should not break and be not too far from the surface

from ase.calculators.emt import EMT
from ase.constraints import FixAtoms, Hookean
from ase.optimize.minimahopping import MinimaHopping

#Restraint condition
const_Pt = [FixAtoms(indices=[atom.index for atom in Pt_211 if atom.symbol=='Pt'])]
d_HH = covalent_radii[atomic_numbers['H']] * 2
H_indices = [atom.index for atom in Pt_211 if atom.symbol=='H']
const_H2 = [Hookean(a1=H_indices[0], a2=H_indices[1], rt=d_HH*1.3, k=15.), #rt is the distance threshold[A],k is the spring constant
            Hookean(a1=H_indices[0], a2=(0., 0., 1., -(z_max+d_PtH*5)), k=5.),
            Hookean(a1=H_indices[1], a2=(0., 0., 1., -(z_max+d_PtH*5)), k=5.)]
surf.set_constraint(const_Pt + const_H2)
surf.set_calculator(EMT())

Run Minima Hopping

Execute Minima Hopping. By default, the log is output to the hop.log file.

hop = MinimaHopping(surf, Ediff0=0.5, T0=2000)
hop(totalsteps=10)

Visualization of results

The ʻase.optimize.minimahoppingmodule implements theMHPlot` function that visualizes the result of MinimaHopping. It is a wrapper for matplotlib. A black check mark indicates the newly found tiny structure.

import matplotlib.pyplot as plt
from ase.optimize.minimahopping import MHPlot

plt.rcParams["font.size"] = 18 #Adjust font size
MHPlot()

fig2.png

Let's also check the tiny structure found by Minima Hopping.

fig3.png

reference

Recommended Posts

Search for adsorption structure using Minima Hopping method
Search for profitable brands using COTOHA
Saddle point search using the gradient method
Python pandas: Search for DataFrame using regular expressions
Refined search for Pokemon race values using Python
Directory structure for test-driven development using pytest in python