ABINIT Calculations for Silicon

Introduction

ABINIT is a powerful density functional theory (DFT) software package widely used for electronic structure calculations. This article presents five distinct calculation setups for silicon, each demonstrating different capabilities and computational strategies. The examples progress from basic single-dataset calculations to more sophisticated multi-dataset workflows, illustrating best practices for k-point convergence studies, structural optimisation, and band structure computations.

Basic Silicon Structure Definition

All examples in this guide use a common structural definition for silicon in its face-centred cubic (FCC) diamond structure:

acell 3*10.18
rprim  0.0  0.5  0.5
       0.5  0.0  0.5
       0.5  0.5  0.0

ntypat 1
znucl 14
pp_dirpath "$ABI_PSPDIR"
pseudos "Psdj_nc_sr_04_pw_std_psp8/Si.psp8"

natom 2
typat 1 1
xred
   0.0  0.0  0.0
   1/4  1/4  1/4

ecut  12.0

The lattice parameter acell is set to 10.18 Bohr radii (approximately 5.39 Ångströms), which is close to the experimental value for silicon. The primitive lattice vectors rprim define the FCC unit cell, and two silicon atoms are positioned at the conventional FCC positions in reduced coordinates.

Example 1: Basic K-point Grid Calculation

The first example demonstrates a straightforward self-consistent field (SCF) calculation with a coarse k-point grid:

kptopt 1
ngkpt 2 2 2
nshiftk 4
shiftk 0.5 0.5 0.5
       0.5 0.0 0.0
       0.0 0.5 0.0
       0.0 0.0 0.5

nstep 10
toldfe 1.0d-6
diemac 12.0

Key Parameters:

• kptopt 1: Automatic generation of k-points using symmetry
• ngkpt 2 2 2: A 2×2×2 Monkhorst-Pack grid
• nshiftk 4 with four different shifts: Ensures proper sampling of the Brillouin zone for FCC structures
• toldfe 1.0d-6: Convergence tolerance on total energy difference
• diemac 12.0: Model dielectric constant for silicon, improving SCF convergence

Example 2: Denser K-point Mesh

The second example uses an identical structure but with a denser k-point grid:

nkpt 2
kptopt 1
ngkpt 4 4 4
nshiftk 4
shiftk 0.5 0.5 0.5
       0.5 0.0 0.0
       0.0 0.5 0.0
       0.0 0.0 0.5

nstep 10
toldfe 1.0d-6
diemac 12.0

The ngkpt 4 4 4 grid provides 64 k-points in the full Brillouin zone before symmetry reduction, offering better convergence for properties sensitive to k-point density.

Example 3: K-point Convergence Study Using Multi-Dataset

This example demonstrates ABINIT's powerful multi-dataset functionality to perform a k-point convergence study in a single run:

ndtset 4

kptopt 1
nshiftk 4
shiftk  0.5 0.5 0.5
        0.5 0.0 0.0
        0.0 0.5 0.0
        0.0 0.0 0.5

ngkpt1  2 2 2
ngkpt2  4 4 4
ngkpt3  6 6 6
ngkpt4  8 8 8

getwfk -1

nstep 10
toldfe 1.0d-6
diemac 12.0

Multi-Dataset Features:

• ndtset 4: Defines four sequential datasets
• ngkpt1 through ngkpt4: Systematically increases k-point density
• getwfk -1: Each dataset uses the wavefunction from the previous calculation as a starting point, accelerating convergence

This approach efficiently tests k-point convergence by running calculations with 2×2×2, 4×4×4, 6×6×6, and 8×8×8 grids in sequence, with each calculation benefiting from the converged wavefunction of the previous one.

Example 4: Structural Optimisation with Variable K-points

This example performs structural relaxation with two different k-point densities:

ndtset 2

optcell 1
geoopt  "bfgs"
ntime  10
dilatmx 1.05
ecutsm  0.5

kptopt 1
nshiftk 4
shiftk  0.5 0.5 0.5
        0.5 0.0 0.0
        0.0 0.5 0.0
        0.0 0.0 0.5

ngkpt1  2 2 2
ngkpt2  4 4 4

getwfk -1

nstep 10
tolvrs 1.0d-14
diemac 12.0

Structural Optimisation Parameters:

• optcell 1: Optimises unit cell volume whilst maintaining shape
• geoopt "bfgs": Uses the BFGS algorithm for geometry optimisation
• ntime 10: Maximum number of structural relaxation steps
• dilatmx 1.05: Maximum lattice parameter change factor (5%)
• ecutsm 0.5: Smooths energy changes during cell optimisation
• tolvrs 1.0d-14: Stringent convergence on residual potential

The strategy here is to perform a coarse structural optimisation with dataset 1, then refine with the denser k-point grid in dataset 2, starting from the optimised wavefunction.

Example 5: SCF Calculation and Band Structure

The final example demonstrates a complete workflow for band structure calculation:

ndtset 2

acell 3*10.195

ecut 12.0

nstep 20
diemac 12.0

kptopt1 1
nshiftk1 4
shiftk1  0.5 0.5 0.5
         0.5 0.0 0.0
         0.0 0.5 0.0
         0.0 0.0 0.5
ngkpt1  4 4 4
prtden1  1
toldfe1  1.0d-6

iscf2    -2
getden2  -1
kptopt2  -3
nband2   8
ndivsm2   10
kptbounds2  0.5  0.0  0.0
            0.0  0.0  0.0
            0.0  0.5  0.5
            1.0  1.0  1.0
tolwfr2  1.0d-12
enunit2  1

Dataset 1 (SCF Calculation):

• Standard SCF calculation with a 4×4×4 k-point grid
• prtden1 1: Outputs the charge density for use in dataset 2

Dataset 2 (Band Structure):

• iscf2 -2: Non-self-consistent calculation (band structure mode)
• getden2 -1: Reads charge density from dataset 1
• kptopt2 -3: Band structure k-point generation along specified paths
• nband2 8: Calculates 8 bands
• ndivsm2 10: 10 divisions between each pair of special k-points
• kptbounds2: Defines the path L → Γ → X → Γ in the Brillouin zone
• tolwfr2 1.0d-12: Convergence tolerance on wavefunction residuals
• enunit2 1: Output energies in electron volts

This path samples the standard high-symmetry points for the FCC Brillouin zone, providing a complete picture of silicon's electronic band structure.

Code for Plotting Band Structures

Plotting Band Structure

import sys
import matplotlib.pyplot as plt
from abipy.abilab import abiopen

# File path
gsr_path = "tbase3_5o_DS2_GSR.nc"

with abiopen(gsr_path) as gsr:
    ebands = gsr.ebands
    ebands.plot(show=False)

    output_file = "bands_plot.png"
    plt.savefig(output_file)

Plotting Band Structure with Corresponding DOS

import sys
import matplotlib.pyplot as plt
from abipy.abilab import abiopen

# File paths (using corrected calculation results)
bands_file = "tbase3_5o_DS2_GSR.nc"
dos_file = "tbase3_5o_DS1_GSR.nc"

with abiopen(bands_file) as gsr_bands, abiopen(dos_file) as gsr_dos:
    # Get band structure (calculated along k-path)
    ebands_kpath = gsr_bands.ebands
    
    # Get DOS (calculated on k-mesh)
    ebands_mesh = gsr_dos.ebands
    dos = ebands_mesh.get_edos()
    
    print(f"Band structure energy range: {ebands_kpath.enemin():.2f} eV to {ebands_kpath.enemax():.2f} eV")
    print("DEBUG: dos attributes:", dir(dos))
    # print(f"DOS energy range: {dos.energies.min():.2f} eV to {dos.energies.max():.2f} eV")
    print(f"Number of bands: {ebands_kpath.nband}")
    
    # Plot bands and DOS, setting a larger energy range
    fig = ebands_kpath.plot_with_edos(dos, 
                                       e_min=-15,  # Energy lower limit
                                       e_max=30,   # Energy upper limit, showing high-energy part
                                       show=False)
    
    output_file = "bands_dos_plot.png"
    plt.savefig(output_file, dpi=300, bbox_inches='tight')
    print(f"Image saved to {output_file}")
    
    #plt.show()

Computational Strategy and Best Practices

The examples presented illustrate several important computational strategies:

1. Progressive Refinement: Start with coarse parameters and progressively refine, using previous results as starting points (getwfk -1, getden -1)

2. Multi-Dataset Efficiency: Consolidate related calculations into a single input file, reducing file management overhead and leveraging warm-start capabilities

3. Appropriate Convergence Criteria: Use toldfe for standard SCF, tolvrs for structural optimisation, and tolwfr for band structure calculations

4. K-point Sampling: For silicon's FCC structure, using four shifted grids ensures proper Brillouin zone sampling and avoids missing important regions of reciprocal space

5. Structural Optimisation Safeguards: The dilatmx parameter prevents runaway optimisation, whilst ecutsm ensures smooth energy surfaces during cell parameter changes

Conclusion

These ABINIT input examples provide a practical foundation for performing various types of electronic structure calculations on silicon. The progression from basic single calculations to sophisticated multi-dataset workflows demonstrates how ABINIT's features can be leveraged to create efficient and reproducible computational workflows. Researchers can adapt these templates for other materials by modifying the structural parameters, pseudopotentials, and computational settings as appropriate for their specific systems.

The multi-dataset approach is particularly powerful for systematic studies, such as convergence testing or workflows that combine different calculation types. By chaining datasets together and reusing computed quantities, significant computational time can be saved whilst maintaining organised and manageable input files.

References

ABINIT is developed by the ABINIT group and its documentation provides comprehensive guidance on all parameters and methodologies. For detailed information about ABINIT's capabilities and parameter definitions, consult the official ABINIT documentation.

The Monkhorst-Pack k-point sampling scheme used throughout these examples is described in: Monkhorst, H. J., and J. D. Pack. "Special Points for Brillouin-Zone Integrations." Physical Review B 13.12 (1976): 5188-5192.

For theoretical background on density functional theory and pseudopotentials, see: Martin, R. M. Electronic Structure: Basic Theory and Practical Methods. Cambridge University Press, 2004.