Back to the main Index

Structure object

The AbiPy structure inherits from the pymatgen structure. One has therefore access to all the methods and tools already available in pymatgen. In this notebook, we mainly focus on the extensions added by AbiPy. For the features provided by pymatgen, please consult the official pymatgen documentation

Table of Contents

[back to top]

In [1]:
from __future__ import division, print_function, unicode_literals

import warnings
warnings.filterwarnings("ignore") # to get rid of deprecation warnings

# Import abipy modules
from abipy import abilab
from abipy.abilab import Structure
import abipy.data as abidata

# Useful tools we'll need later on.
from pprint import pprint
import numpy as np

# This line configures matplotlib to show figures embedded in the notebook.
# Replace `inline` with `notebook` in classic notebook
%matplotlib inline   

# Option available in jupyterlab. See https://github.com/matplotlib/jupyter-matplotlib
#%matplotlib widget  

Reading a structure from file

[back to top]

It is possible to initialize a crystalline structure from different file formats:

  • CIF
  • POSCAR/CONTCAR
  • CHGCAR
  • LOCPOT,
  • vasprun.xml
  • CSSR
  • ABINIT Netcdf files
  • pymatgen's JSON serialized structures

Note, in particular, that one can initialize the structure from the netcdf files
produced by Abinit (GSR.nc, WFK.nc, etc) as well as output files in text format such as the Abinit input/output files or even the DDB file.

To initialize the structure from a CIF file use the from_file method:

In [2]:
# abidata.cif_file returns one of the CIF files shipped with AbiPy.
structure = Structure.from_file(abidata.cif_file("si.cif"))
print(structure)
Full Formula (Si2)
Reduced Formula: Si
abc   :   3.866975   3.866975   3.866975
angles:  60.000000  60.000000  60.000000
Sites (2)
  #  SP       a     b     c
---  ----  ----  ----  ----
  0  Si    0     0     0
  1  Si    0.25  0.25  0.25

To read the structure from a netcdf file:

In [3]:
structure = Structure.from_file(abidata.ref_file("si_nscf_GSR.nc"))

# Use to_string with verbose > 0 to get more info 
print(structure.to_string(verbose=1))
Full Formula (Si2)
Reduced Formula: Si
abc   :   3.866975   3.866975   3.866975
angles:  60.000000  60.000000  60.000000

Spglib space group info (magnetic symmetries not taken into account).
Spacegroup: Fd-3m (227), Hall: F 4d 2 3 -1d, Abinit spg_number: 227
Crystal_system: cubic, Lattice_type: cubic, Point_group: m-3m

  Idx  Symbol    Reduced_Coords              Wyckoff      EqIdx
-----  --------  --------------------------  ---------  -------
    0  Si        +0.00000 +0.00000 +0.00000  (2a)             0
    1  Si        +0.25000 +0.25000 +0.25000  (2a)             0

Abinit Spacegroup: spgid: 227, num_spatial_symmetries: 48, has_timerev: True, symmorphic: True

Use to_abivars to get the list of Abinit variables in a python dictionary:

In [4]:
structure.to_abivars()
Out[4]:
{'natom': 2,
 'ntypat': 1,
 'typat': array([1, 1]),
 'znucl': [14],
 'xred': array([[0.  , 0.  , 0.  ],
        [0.25, 0.25, 0.25]]),
 'acell': [1.0, 1.0, 1.0],
 'rprim': array([[6.32850055, 0.        , 3.6537615 ],
        [2.10950018, 5.96656754, 3.6537615 ],
        [0.        , 0.        , 7.30752299]])}

and abi_string to get a string that can be used directly in the input file:

In [5]:
print(structure.abi_string)
 natom 2
 ntypat 1
 typat 1 1
 znucl 14
 xred
    0.0000000000    0.0000000000    0.0000000000
    0.2500000000    0.2500000000    0.2500000000
 acell    1.0    1.0    1.0
 rprim
    6.3285005521    0.0000000000    3.6537614973
    2.1095001840    5.9665675402    3.6537614973
    0.0000000000    0.0000000000    7.3075229946

To visualize the structure with matplotlib, use:

In [6]:
structure.plot();

The matplotlib version is minimalistic but it plays well with jupyter notebooks. For a more advanced visualization we suggest using a specialized graphical applications. Fortunately, one can invoke (already installed) external applications directly from AbiPy with e.g.

In [7]:
# structure.visualize("vesta")

To get a structure from the materials project database (https://www.materialsproject.org ), use:

In [8]:
# You can pass the api_key or set the env variable PMG_MAPI_KEY in your ~/.pmgrc.yaml files.
si2_mp = Structure.from_mpid("mp-149", api_key=None)
print(si2_mp)
Full Formula (Si2)
Reduced Formula: Si
abc   :   3.866975   3.866975   3.866975
angles:  60.000000  60.000000  60.000000
Sites (2)
  #  SP        a      b      c    magmom
---  ----  -----  -----  -----  --------
  0  Si    0.875  0.875  0.875         0
  1  Si    0.125  0.125  0.125         0

In some cases, we have multiple structures and we need to compare the lattice parameters. Use dataframes_from_structures to build a pandas DataFrame:

In [9]:
dfs = abilab.dataframes_from_structures([structure, si2_mp], index=["CIF", "MP"])

then we can compare the lattice parameters with:

In [10]:
dfs.lattice
Out[10]:
formula natom alpha beta gamma a b c volume abispg_num spglib_symb spglib_num spglib_lattice_type
CIF Si2 2 60.0 60.0 60.0 3.866975 3.866975 3.866975 40.888292 227 Fd-3m 227 cubic
MP Si2 2 60.0 60.0 60.0 3.866975 3.866975 3.866975 40.888292 None Fd-3m 227 cubic

Note that all AbiPy robots have this feature built-in. Sometimes it is much easier to build a robot directly from files and then compare the structures with e.g. robot.get_lattice_dataframe().

Converting to other formats

[back to top]

Use structure.convert(format) to get the string representation in the new format:

In [11]:
for fmt in ["cif", "POSCAR", "qe"]:
    print((" Abinit --> %s " % fmt).center(80, "*"))
    print(structure.convert(fmt=fmt))
******************************** Abinit --> cif ********************************
# generated using pymatgen
data_Si
_symmetry_space_group_name_H-M   'P 1'
_cell_length_a   3.86697464
_cell_length_b   3.86697464
_cell_length_c   3.86697464
_cell_angle_alpha   60.00000000
_cell_angle_beta   60.00000000
_cell_angle_gamma   60.00000000
_symmetry_Int_Tables_number   1
_chemical_formula_structural   Si
_chemical_formula_sum   Si2
_cell_volume   40.88829228
_cell_formula_units_Z   2
loop_
 _symmetry_equiv_pos_site_id
 _symmetry_equiv_pos_as_xyz
  1  'x, y, z'
loop_
 _atom_site_type_symbol
 _atom_site_label
 _atom_site_symmetry_multiplicity
 _atom_site_fract_x
 _atom_site_fract_y
 _atom_site_fract_z
 _atom_site_occupancy
  Si  Si0  1  0.000000  0.000000  0.000000  1
  Si  Si1  1  0.250000  0.250000  0.250000  1

****************************** Abinit --> POSCAR *******************************
Si2
1.0
3.348898 0.000000 1.933487
1.116299 3.157372 1.933487
0.000000 0.000000 3.866975
Si
2
direct
0.000000 0.000000 0.000000 Si
0.250000 0.250000 0.250000 Si

******************************** Abinit --> qe *********************************
&CONTROL
  calculation = 'scf',
/
&SYSTEM
  ibrav = 0,
  nat = 2,
  ntyp = 1,
/
&ELECTRONS
/
&IONS
/
&CELL
/
ATOMIC_SPECIES
  Si  28.0855 Si.pseudo
ATOMIC_POSITIONS crystal
  Si 0.000000 0.000000 0.000000
  Si 0.250000 0.250000 0.250000
K_POINTS automatic
  1 1 1 0 0 0
CELL_PARAMETERS angstrom
  3.348898 0.000000 1.933487
  1.116299 3.157372 1.933487
  0.000000 0.000000 3.866975

Getting information on the structure

[back to top]

In [12]:
print(structure.reciprocal_lattice)
1.876195 -0.663335 0.000000
0.000000 1.990005 0.000000
-0.938097 -0.663335 1.624832
In [13]:
structure.reciprocal_lattice.matrix.T @ structure.lattice.matrix / (2 * np.pi)
Out[13]:
array([[ 1.00000000e+00,  0.00000000e+00, -2.49862819e-17],
       [-2.24284041e-19,  1.00000000e+00,  4.44094187e-17],
       [ 0.00000000e+00,  0.00000000e+00,  1.00000000e+00]])
In [14]:
# List of high-symmetry k-points.
print(structure.hsym_kpoints)
0) [+0.000, +0.000, +0.000], name: $\Gamma$, weight: 0.000
1) [+0.500, +0.000, +0.500], name: X, weight: 0.000
2) [+0.500, +0.250, +0.750], name: W, weight: 0.000
3) [+0.375, +0.375, +0.750], name: K, weight: 0.000
4) [+0.000, +0.000, +0.000], name: $\Gamma$, weight: 0.000
5) [+0.500, +0.500, +0.500], name: L, weight: 0.000
6) [+0.625, +0.250, +0.625], name: U, weight: 0.000
7) [+0.500, +0.250, +0.750], name: W, weight: 0.000
8) [+0.500, +0.500, +0.500], name: L, weight: 0.000
9) [+0.375, +0.375, +0.750], name: K, weight: 0.000
10) [+0.625, +0.250, +0.625], name: U, weight: 0.000
11) [+0.500, +0.000, +0.500], name: X, weight: 0.000

The method calc_ksampling allows one to get an efficient sampling of the Brillouin zone by just specifying the number of divisions to be used for the smallest lattice vector of the reciprocal lattice:

In [15]:
pprint(structure.calc_ksampling(nksmall=10))
{'ngkpt': array([10, 10, 10]),
 'shiftk': array([[0.5, 0.5, 0.5],
       [0.5, 0. , 0. ],
       [0. , 0.5, 0. ],
       [0. , 0. , 0.5]])}

To get the recommended high symmetry $k$-path in reduced coordinates:

In [16]:
structure.calc_kptbounds()
Out[16]:
array([[0.   , 0.   , 0.   ],
       [0.5  , 0.   , 0.5  ],
       [0.5  , 0.25 , 0.75 ],
       [0.375, 0.375, 0.75 ],
       [0.   , 0.   , 0.   ],
       [0.5  , 0.5  , 0.5  ],
       [0.625, 0.25 , 0.625],
       [0.5  , 0.25 , 0.75 ],
       [0.5  , 0.5  , 0.5  ],
       [0.375, 0.375, 0.75 ],
       [0.625, 0.25 , 0.625],
       [0.5  , 0.   , 0.5  ]])

The high-symmetry q-path is automatically selected assuming the structure fulfills the convention described in Setyawan2010

To visualize the Brillouin zone with matplotlib:

In [17]:
structure.plot_bz();

To get the number of valence electrons for a given set of pseudopotentials:

In [18]:
structure.num_valence_electrons(pseudos=abidata.pseudos("14si.pspnc"))
Out[18]:
8

To visualize the X-ray diffraction plot with pymatgen XRDCalculator

In [19]:
structure.plot_xrd();

abistruct.py

[back to top]

abistruct.py provides a handy command line interface to operate on structure objects constructed from external files. There are several options available as well an interface to the materials project and the COD database.

In [20]:
!abistruct.py --help
usage: abistruct.py [-h] [-V]
                    {spglib,abispg,convert,supercell,abisanitize,irefine,conventional,proto,wyckoff,tensor_site,neighbors,interpolate,xrd,oxistate,ipython,notebook,kpath,bz,ngkpt,ktables,abikmesh,lgk,kstar,keq,visualize,mp_id,mp_match,mp_search,mp_pd,mp_ebands,cod_search,cod_id,animate}
                    ...

optional arguments:
  -h, --help            show this help message and exit
  -V, --version         show program's version number and exit

subcommands:
  Valid subcommands, use command --help for help

  {spglib,abispg,convert,supercell,abisanitize,irefine,conventional,proto,wyckoff,tensor_site,neighbors,interpolate,xrd,oxistate,ipython,notebook,kpath,bz,ngkpt,ktables,abikmesh,lgk,kstar,keq,visualize,mp_id,mp_match,mp_search,mp_pd,mp_ebands,cod_search,cod_id,animate}
                        sub-command help
    spglib              Analyze structure with spglib.
    abispg              Extract/Compute Abinit space group from file with
                        structure.
    convert             Convert structure to the specified format.
    supercell           Generate supercell.
    abisanitize         Sanitize structure with abi_sanitize, compare
                        structures and save result to file.
    irefine             Refine structure with abi_sanitize iteratively, stop
                        if target space group is obtained.
    conventional        Gives a structure with a conventional cell according
                        to certain standards. The standards are defined in
                        doi:10.1016/j.commatsci.2010.05.010
    proto               Find prototype in the AFLOW LIBRARY OF
                        CRYSTALLOGRAPHIC PROTOTYPES.
                        http://doi.org/10.1016/j.commatsci.2017.01.017
    wyckoff             Print wyckoff positions. WARNING: still under
                        development!
    tensor_site         Print symmetry properties of tensors due to site-
                        symmetries. WARNING: still under development!
    neighbors           Get neighbors for each atom in the unit cell, out to a
                        distance radius.
    interpolate         Interpolate between two structures. Useful for the
                        construction of NEB inputs.
    xrd                 X-ray diffraction plot.
    oxistate            Estimate oxidation states with pymatgen bond valence
                        analysis.
    ipython             Open IPython shell for advanced operations on
                        structure object.
    notebook            Read structure from file and generate jupyter
                        notebook.
    kpath               Read structure from file, generate k-path for band-
                        structure calculations.
    bz                  Read structure from file, plot Brillouin zone with
                        matplotlib.
    ngkpt               Return the Abinit k-point sampling variables from the
                        number of divisions used to sample the smallest
                        lattice vector of the reciprocal lattice.
    ktables             Read structure from filepath, call spglib to sample
                        the BZ, print k-points in the IBZ with weights.
    abikmesh            Read structure from file, call Abinit to sample the BZ
                        with ngkpt, shiftk, and kptopt. Print k-points in the
                        IBZ with weights.
    lgk                 Read structure from file, find little group of
                        k-point, print Bilbao character table.
    kstar               Read structure from file, print star of k-point.
    keq                 Read structure from file, check whether two k-points
                        are equivalent by symmetry.
    visualize           Visualize the structure with the specified
                        application. Requires external app or optional python
                        modules (mayavi, vtk).
    mp_id               Get structure from the pymatgen database. Export to
                        format. Requires internet connection and PMG_MAPI_KEY.
    mp_match            Get structure from the pymatgen database. Requires
                        internet connection and PMG_MAPI_KEY.
    mp_search           Get structure from the pymatgen database. Requires
                        internet connection and PMG_MAPI_KEY
    mp_pd               Generate phase diagram with entries from the Materials
                        Project. Requires internet connection and PMG_MAPI_KEY
    mp_ebands           Get structure from the pymatgen database. Export to
                        format. Requires internet connection and PMG_MAPI_KEY.
    cod_search          Get structure from COD database. Requires internet
                        connection and mysql
    cod_id              Get structure from COD database. Requires internet
                        connection and mysql
    animate             Read structures from HIST.nc or XDATCAR. Print
                        structures in Xcrysden AXSF format to stdout.

Usage example:

###################
# Space group tools
###################

  abistruct.py spglib FILE                 => Read structure from FILE and analyze it with spglib.
  abistruct.py abispg FILE                 => Read structure from FILE, and compute ABINIT space group.
  abistruct.py abisanitize FILE            => Read structure from FILE, call abisanitize, compare structures
                                              and save "abisanitized" structure to file.
  abistruct.py conventional FILE           => Read structure from FILE, generate conventional structure
                                              following doi:10.1016/j.commatsci.2010.05.010
  abistruct.py proto FILE                 => Read structure from FILE, find possible crystallographic prototypes:
                                             in the AFLOW library of crystallographic prototypes.
                                             http://doi.org/10.1016/j.commatsci.2017.01.017

##################
# Conversion tools
##################

  abistruct.py convert FILE                => Print the ABINIT variables defining the structure.
  abistruct.py convert FILE -f cif         => Read structure from FILE and output CIF file
                                              (Use convert --help to get list of formats supported)
  abistruct.py convert out_HIST.nc         => Read FINAL structure from the HIST file and
                                              print the corresponding ABINIT variables.
  abistruct.py supercell FILE -s 2 2 1     => Read structure from FILE and build [2, 2, 1] supercell,
                                              print new structure using --format (default abivars).
################
# K-points tools
################

  abistruct.py kpath FILE                  => Read structure from FILE and print ABINIT variables for k-path.
  abistruct.py bz FILE                     => Read structure from FILE, plot BZ with matplotlib.
  abistruct.py ngkpt FILE -n 4             => Compute `ngkpt` and `shiftk` from the number of divisions used to sample
                                              the smallest reciprocal lattice vector.
  abistruct.py abikmesh FILE --ngkpt 2 2 2 => Read structure from FILE, call Abinit to sample the BZ
                                              with a 2, 2, 2 k-mesh, print points in IBZ with weights.
  abistruct.py ktables FILE -m 2 2 2       => Read structure from FILE, call spglib to sample the BZ
                                              with a 2, 2, 2 k-mesh, print points in IBZ with weights.
  abistruct.py lgk FILE -k 0.25 0 0        => Read structure from FILE, find little group of k-point,
                                              print Bilbao character table.
  abistruct.py kstar FILE -k 0.25 0 0      => Read structure from FILE, print star of k-point.
  abistruct.py keq FILE -k 0.5 0 0 0 0.5 0  => Read structure from FILE, test whether k1 and k2 are
                                               symmetry-equivalent k-points.

###############
# Miscelleanous
###############

  abistruct.py neighbors FILE              => Get neighbors for each atom in the unit cell, out to a distance radius.
  abistruct.py interpolate FILE1 FILE2     => Interpolate between two structures. Useful for the construction of NEB inputs.
  abistruct.py xrd FILE                    => X-ray diffraction plot.
  abistruct.py oxistate FILE               => Estimate oxidation states with pymatgen bond valence analysis.

###############
# Visualization
###############

  abistruct.py visualize FILE -a vesta     => Visualize the structure with vesta (default if -a is not given)
                                              Supports also ovito, xcrysden, vtk, mayavi, matplotlib See --help
  abistruct.py ipython FILE                => Read structure from FILE and open it in the Ipython terminal.
  abistruct.py notebook FILE               => Read structure from FILE and generate jupyter notebook.

###########
# Databases
###########

  abistruct.py cod_id 1526507              => Get structure from COD database (http://www.crystallography.net/cod).
  abistruct.py cod_search MgB2             => Search for structures in the COD database.
  abistruct.py mp_id mp-149                => Get structure from materials project database and print
                                              JSON representation. Use e.g. `-f abivars` to change format.
  abistruct.py mp_match FILE               => Read structure from FILE and find matching structures on the
                                              Materials Project site. Use e.g. `-f cif` to change output format.
  abistruct.py mp_search LiF               => Connect to the materials project database. Get structures corresponding
                                              to a chemical system or formula e.g. `Fe2O3` or `Li-Fe-O` or
                                              `Ir-O-*` for wildcard pattern matching.
                                              Print info and Abinit input files. Use e.g. `-f POSCAR`
                                              to change output format. `-f None` to disable structure output.
  abistruct.py mp_pd FILE-or-elements      => Generate phase diagram with entries from the Materials Project.
  abistruct.py mp_ebands FILE             => Fetch electron band structure from MP database. Print gaps.
                                              Accept FILE with structure if ebands fro structure is wanted
                                              or mp id e.g. "mp-149 or list of elements e.g `Li-Fe-O` or chemical formula.

`FILE` is any file supported by abipy/pymatgen e.g Netcdf files, Abinit input/output, POSCAR, xsf ...
Use `abistruct.py --help` for help and `abistruct.py COMMAND --help` to get the documentation for `COMMAND`.
Use `-v` to increase verbosity level (can be supplied multiple times e.g -vv).

Back to the main Index