Pymatgen sub-package for lifting a surface layer in z (c) direction

Hi everyone,

I wonder if there is any subpackage in pymatgen can help me lifting an surface layer of a slab structure with 0.2Å increment (from 2Å - 5Å) and write it into POSCAR format to perform crack initiation test.

Here is an example
POSCAR file
1.00000000000000
0.0000000000000000 0.0000000000000000 -4.9935999999999998
-4.3245844599999996 0.0000000001000000 2.4967999999999999
0.0000000000000000 41.3633999999999986 -0.0000000004000000
Cr Fe O
10 24 15
Selective dynamics
Direct
0.3514277164107921 0.6907049529237036 0.3360449354706682 T T T
0.6723760709429016 0.4075823983139003 0.3736943951592519 T T T
0.9479207188265740 0.9724474079673007 0.3706400244759484 T T T
0.2617487971845416 0.6907609695464666 0.4233634401920597 T T T
0.9881900877704727 -0.0116999314849779 0.4280450420952487 T T T
0.2893505895297185 0.6931712775084026 0.4802452474270836 T T T
0.5851658341304642 0.2906045538818027 0.4762435910074906 T T T
0.9272200336548112 0.9653393177405936 0.5309163533682427 T T T
0.6339373444534714 0.3820419944223448 0.5332686281224046 T T T
0.3140539881041046 0.6721136696657830 0.5683855293813306 T T T
0.1666666650000010 0.3333333349999990 0.0559285849999966 F F F
0.1666666650000010 0.8333333349999990 0.0559285849999966 F F F
0.6666666650000010 0.3333333349999990 0.0559285849999966 F F F
0.6666666650000010 0.8333333349999990 0.0559285849999966 F F F
0.3369065168727638 0.1854341988110731 0.2938774571544731 T T T
0.3374054300480825 0.6665161054082195 0.2904954615644242 T T T
0.8244321165308205 0.1681676604025589 0.2937505101724555 T T T
0.8481311017587185 0.6731526193764762 0.2937362747893189 T T T
0.0030563489535511 0.0067540173480085 0.2484276577972868 T T T
0.9847826829734412 0.4931786586376498 0.2463107060608104 T T T
0.5034911770606328 0.0212029737920214 0.2466004471029293 T T T
0.5174067960293125 0.5005212639591873 0.2462724475614931 T T T

I would like to lift the BOLD section of the poscar in Z-direction by 0.2Å increment up to 5Å.

Does anyone has study similar problem? I have been doing it manually and would like to speed up the process since I have 30+ POSCAR structures that I need to perform crack test.

I tried mpinterface.interface.Interface but it gave me this error

---------------------------------------------------------------------------

ImportError Traceback (most recent call last)
in
----> 1 from mpinterfaces.calibrate import CalibrateSlab
2 from mpinterfaces.interface import Interface
3 from mpinterfaces.transformations import *
4 from mpinterfaces.utils import *
5

~/.local/lib/python3.8/site-packages/mpinterfaces/calibrate.py in
36 import numpy as np
37
—> 38 from pymatgen import Lattice
39 from pymatgen.core.structure import Structure
40 from pymatgen.core.surface import SlabGenerator

ImportError: cannot import name ‘Lattice’ from ‘pymatgen’ (unknown location)

Try from pymatgen.core.lattice import Lattice – it looks like MPInterfaces needs to be updated for the latest version of pymatgen (v2022 onwards). Alternatively, you can install a slightly older version of pymatgen (one of the v2021 versions).

Hi,

as @mkhorton said, MPInterfaces is not really updated a lot any more. But for your problem there should be a rather simple solution if your slabs are structured in a similar manner.

I think the hardest part in the general problem is to find the topmost layer, but if you slab sites are ordered like in you POSCAR example (the first n sites are make up the first layer in all your systems, or you know the site indices if they are non consecutive), you could do something like this:

def shift_sites(slab, sites_to_shift, shift):
    shifted_coords = []
    for i, s in enumerate(slab.cart_coords):
        if i in sites_to_shift:
            shifted_coords.append(s+shift)
        else:
            shifted_coords.append(s)
    shifted_slab = Slab(lattice=slab.lattice,
                        species=slab.species,
                        coords=shifted_coords,
                        miller_index=slab.miller_index,
                        oriented_unit_cell=slab.oriented_unit_cell,
                        shift=slab.shift,
                        scale_factor=slab.scale_factor,
                        reconstruction=slab.reconstruction,
                        coords_are_cartesian=True,
                        site_properties=slab.site_properties)
    return shifted_slab

where slab is your pymatgen slab object without shifts, sites_to_shift is a list of site indices (in your example list(range(10)) should work), and shift is a list of shifts, e.g. [delta_x, delta_y, delta_z], in your case [0,0,0.2].

If you do not know the site indices of the first layer, you might be able to find them using some functions in the pymatgen.core.surface module, e.g. get_d and/or get_slab_regions.

Good luck, Michael

3 Likes