Parameters calculated for MgB2 and for B-doped diamond.
We present three systematic approaches to use of Density Functional Theory (DFT) for interpretation and prediction of superconductivity in new or existing materials. These approaches do not require estimates of free parameters but utilize standard input values that significantly influence computational resolution of reciprocal space Fermi surfaces and that reduce the meV-scale energy variability of calculated values. Systematic calculations on conventional superconductors show that to attain a level of resolution comparable to the energy gap, two key parameters, Δk and the cut-off energy, must be optimized for a specific compound. The optimal level of resolution is achieved with k-grids smaller than the minimum reciprocal space separation between key parallel Fermi surfaces. These approaches enable estimates of superconducting properties including the transition temperature (Tc) via (i) measurement of the equivalent thermal energy of a phonon anomaly (if present), (ii) the distribution of electrons and effect on Fermi energy (EF) when subjected to a deformation potential and (iii) use of parabolic, or higher order quartic, approximations for key electronic bands implicated in electron–phonon interactions. We demonstrate these approaches for the conventional superconductors MgB2, metal substituted MgB2 and boron-doped diamond.
- Fermi energy
- Fermi level
- Fermi surface
- reciprocal space
- density functional theory
- parabolic equations
- phonon dispersions
- transition temperature
- magnesium diboride
Design and synthesis of new materials requires translation of reciprocal space detail in electronic band structures (EBSs) and phonon dispersions (PDs) to equivalent real space representations [1, 2, 3, 4]. We have shown that for conventional superconductors (SCs), the format and depth of modes in PDs associated with a Kohn anomaly are strongly influenced by the computational resolution of DFT models . Our view is that EBS calculations, performed with appropriate resolution, may also provide critical information on the superconducting gap, which for many SC materials, is in the meV range.
The Fermi level, Fermi energy (EF), and the density of nearly free-electron carriers calculated by DFT are key values that, to date, have been variously reported with differences of many hundreds of meV for the same compound. For example, the value of EF for a well studied compound such as MgB2 has been variously reported as “several eV” , 0.55 eV or 0.122 eV , and more recently, as 0.428 eV . In comparison, many publications on MgB2 and software packages such as CASTEP and ADF consistently report a value for EF of ~8.4 eV [9, 10, 11]. Calculated variations of this magnitude have garnered limited attention  for SCs due to an underlying assumption that EF/kθ > > 1 (where k is Boltzman’s constant and θ is the Debye temperature). Due to this assumption and the deceptive influence of average values for phonon frequencies, a value for EF is not considered in the simplified McMillan version of the Eliashberg model for superconductivity . However, as noted by Malik , equations that explicitly include EF and/or critical current (jo) values may provide clues on how to increase or modify Tc. Malik suggests that regardless of physical attributes, SCs may be distinguished by their values of EF .
Approximations to the Eliashberg model that minimize computational cost require estimates of the electron–phonon interaction, λ, and the Coulomb strength, μ* [13, 14]. For many conventional SCs, these parameters are limited to a narrow range of values and provide reasonable estimates for superconducting properties of known materials [15, 16, 17, 18]. More recently, Sanna
In our search for new SCs, we have evaluated the computational resolution and electron–phonon detail possible with DFT codes readily available in well-known software packages such as CASTEP or Quantum Espresso, to name a couple of examples. Using MgB2 and similar Bardeen-Cooper-Schrieffer (BCS) compounds, we have systematically explored the sensitivity and use of PD and EBS constructs to calculate key superconducting properties without recourse to free parameter estimates or modification of functionals. We initially explored use of a phonon anomaly to estimate Tc [23, 24]; an approach that appears effective for strong phonon mediated superconductivity including for metal substituted MgB2 . In other work, we extended this systematic approach to evaluate PDs and EBSs for a wide range of metal diboride compounds (
More recently, we have examined the link between PDs and EBSs and, in particular, the topology of the Fermi Surface (FS) with pressure  and the change in electron density distributions as MgB2 transitions to the superconducting state . In both approaches, we are able to confirm experimentally determined superconducting properties for a range of conventional (BCS) compounds and then, to predict Tc for new metal substituted analogues of MgB2 [23, 25]. These approaches are, in essence, empirical methods, which systematically identify regular dispersion patterns in calculated PDs and EBSs, based entirely on accepted codification of the DFT [9, 10] and, equally, a clear understanding of input parameter limitations that determine computational resolution. Thus, this use of DFT software, underpinned by elegant formalism and constructs by Kohn and colleagues [29, 30], is a complement to approximations of the Eliashberg model [19, 20].
In this work, we show why computational resolution for DFT models of EBSs and PDs for conventional SCs is critical. In addition, we delineate a third approach to estimate the superconducting gap using parabolic, or higher order quartic, approximations to key bands in the EBS. This approach requires examination of an extended Brillouin zone (BZ) schema and demands lower computational cost compared to equivalent PD calculations for similar outcome(s). When applied to MgB2, at sufficiently fine k-grid and reciprocal space cut-off, this approach directly estimates the superconducting gap and assists identification of valence bands and the origin for EF. In combination, these three approaches provide reliable property predictions for unknown, or theoretical, structures for materials researchers.
2. Calculation methods
EBS calculations are undertaken using DFT as implemented in the Cambridge Serial Total Energy Package (CASTEP) of Materials Studio (MS) 2017 and 2018 [9, 10]. All structures are optimized for geometry, including cell parameters, starting with crystal information files (.cif) available in standard databases. In general, the local density approximation (LDA) and generalized gradient approximation (GGA), with norm-conserving pseudopotentials, are used in DFT calculations. The typical setup for calculations uses a k-grid ranging from 0.06 Å−1 to 0.02 Å−1 or smaller, with a plane wave basis set cut-off of 990 eV, ultra-fine (or better) customized setup to ensure total energy convergence of less than 5 × 10−6 eV/atom, a maximum force of less than 0.01 eV/Å, a maximum stress of less than 0.02 GPa and maximum displacements of less than 5 × 10−4 Å.
All outputs meet convergence criteria at the same fine tolerance level for geometry optimization. The effects of input parameters to DFT calculations described in this work are not related to differences in calculation convergence, but critically, are due to the discreteness, or the finite number of the reciprocal space points, used to select plane waves as basis functions. To illustrate particular points, we also vary specific parameters such as basis set cut-off values, Δk values or software versions as identified in the text.
We also perform numerical interpolations of EBSs to validate higher order trends and to delineate fine structure in computed outcomes. To obtain parabolic, or higher order polynomial, approximations to the electronic bands, the DFT calculated data from MS is exported in csv/excel format. For MgB2, sections of particular bands in the energy range − 14 to 4 eV along the Γ-M (and Γ-K) directions are selected and mirrored across the vertical axis at Γ. Individual parabolic, or higher order quartic, trendline fittings are obtained and used to overlay for comparison with the (periodically repeated) extended BZ scheme of the DFT calculated EBS. Effective masses are also calculated and evaluated for parabolic approximations of different branches of the EBS.
3. Computational resolution
We provide outputs from a series of
3.1 Band structure – variation with k-grid
We have been intrigued by the potential to directly determine the superconducting gap energy for a BCS SC using an appropriate resolution EBS. In this regard, MgB2 offers good opportunity to evaluate this potential due to well defined crystallography and the key role ascribed to σ bands and superconductivity [12, 34]. DFT calculations for a compositional suite, or structure type, produce EBSs that show, in general, similar formats even when two different functional approximations are used . A typical outcome for MgB2 using the LDA approximation in the CASTEP module of Materials Studio using a k-grid value of 0.008 Å−1 is shown in Figure 1. These band structures convey useful information for elucidation of potential for superconductivity. For example, the σ bands appear as approximate inverted parabolas (in red and blue lines; green dotted box) near the Γ center point and cross the Fermi level on either side of Γ (Figure 1). These bands display a strong electron–phonon coupling to the E2g phonon modes and are implicated in superconductivity for MgB2
For MgB2, the relationship of these σ bands near the Fermi level on either side of Γ is shown in the schematic in Figure 2. The reciprocal space projection of degenerate σ bands at the Fermi level correspond to the three-dimensional reciprocal space representations of two FSs at Γ parallel with the
Variations in the calculated energies of specific bands, for example in MgB2 where electron–phonon coupling is predominantly linked to the σ bands [17, 39, 40], are strongly influenced by the k-grid value used in DFT computations. Figure 3 demonstrates the effect of k-grid value, or the sensitivity of DFT calculations, on the EBS for MgB2 using the LDA functional for a series of Δk values 0.02 Å−1, 0.04 Å−1 and 0.06 Å−1. The k-grid value affects calculated energies for bands near the Fermi level particularly those σ bands associated with superconductivity highlighted in Figure 1 for MgB2. The differences in energy at Γ or A between calculations range from tens of meV to hundreds of meV for three different k-grid values as shown in Figure 3.
In Table 1, we show substantial meV shifts in enthalpy and EF for MgB2 calculated at different k-grid values using the same functional and the same ultra-fine tolerance for geometry optimization convergence. For MgB2, Table 1 shows the difference in energy, ΔEv (in eV), between the Fermi level and the vertex of the parabola at Γ for different values of Δk. Differences in lattice parameters (
|Compound||k-grid value [Å−1]||Lattice parameters [Å]||Enthalpy [eV]||Fermi energy [eV]||ΔEv [eV]|
Tc ~ 40 K
Tc ~ 4–7.5 K
We also show calculated values for B-doped diamond using different k-grid values in Table 1. In this case, k-grid intervals are smaller (0.005 Å−1 < Δk < 0.020 Å−1) than those used for MgB2 with corresponding smaller shifts in enthalpy and EF. This lower magnitude impact of the k-grid is in part due to fewer degrees of freedom (
We show a more detailed systematic comparison of calculated enthalpies as a function of the k-grid value for MgB2 in Figure 4. As noted above, all calculations are converged to the same ultra-fine criteria, or tolerance, where self-consistency is achieved. The variability of results shown in Figure 4 is due to the discreteness of functions and values used to derive the full solution of the Schroedinger equation. The variability is not due to lack of convergence which in all cases is defined in Section 2 above.
Figure 4a shows that the value of enthalpy for MgB2 oscillates around a consistent minimal value of −1,748.502 eV as the k-grid value is decreased to <0.015 Å−1. The context for this variation in enthalpy is shown in Figure 4b where the k-grid value is extended to 0.2 Å−1– a value that has been used in some DFT calculations as criterion for machine learning algorithms . At these higher values for Δk, enthalpy calculations do not provide useful information on subtle structural variations due to superlattices or of order/disorder. For example, we have shown using DFT calculations with appropriate k-grid values for (Mg1-xAlx)B2 that ordered motifs with adjacent Al-layers are thermodynamically favored by ~0.15 eV over more complex, disordered configuration(s) . A discrepancy of ~0.2 eV due to incorrect choice of k-grid value (Figure 4) does not enable such distinction to be made with confidence.
3.2 Influence of atom displacements
In a dynamic system, other factors may also influence the position of key electronic bands with respect to the Fermi level. For example, atoms in all solids at temperatures above absolute zero vibrate  and in some cases, the resulting phonons may align with specific crystallographic real space features such as inter-atom bonds. This circumstance occurs for MgB2 in which one of the dominant E2g phonon modes – shown to be intimately involved in electron–phonon coupling at the onset of superconductivity [17, 34] – aligns with B–B bonds in the
Figure 5 shows the effect on the EBS for MgB2 of atom displacement along the B–B bond by ~0.6% (
The σ bands for MgB2 consist of two bands degenerate at the Γ point, but degeneracy is lost when the vector k does not equal 0. The two bands thus have different effective masses (or curvatures) which appear to vary as a deformation potential is applied. In fact, the curvatures become very similar and/or identical along Γ–K at the point where a deformation shifts one band tangential to the Fermi level as shown in Figure 5. Along Γ–M, the effective masses are shown to cross over, where the top band has the curvature of the original light effective mass band, and for k-vectors away from the origin, the curvature remains the same as the original heavy effective mass band.
The response of the effective masses (
The calculated coefficients for the polynomials, fitted to σH, σL and the average trend line for these bands, σM, are dependent on the k-grid values used in DFT calculations as shown for the equilibrium condition (Figure 6a) in Table 2. The terms of these polynomial coefficients (
|Grid Value (Å−1)||Coefficient X4||Coefficient X2||Coefficient X0||Coefficient X0|
|At Equilibrium (Dx = 0.0)|
|Displaced along E2g (Dx = 0.006)|
3.3 Brillouin zone schemes and high order quartic approximations
Conceptually, it is generally accepted that two approaches: (a) the free-electron theory and (b) the tight binding, or linear combination of atomic orbitals (LCAO), offer reasonably good approximations to the conduction and valence bands, respectively, in the electronic structure of materials [1, 2, 47, 48, 49, 50]. With increased computational power, the distinction between these two models becomes negligible. In general, the actual EBS should be similar to an average of respective contributions from these two types of approximations.
The detailed origin of particular bands and that of the zero of EF can be more readily appreciated when we examine an extended Brillouin zone scheme, instead of a reduced zone scheme [1, 2, 48, 51, 52]. For example, Figures 7a–c show periodically repeated reciprocal unit cells with reference to the extended Brillouin zone schemes for the electronic bands of MgB2 along the Γ-M and Γ-K directions, respectively. Figure 7a shows a 2D representation of multiple reciprocal unit cells viewed along c* and identifies reciprocal directions for the calculated EBSs for extended BZs in Figure 7b and c.
The origin of EF determined by the parabolic approximation is identified at (i) the cross-over of two parabolas “1” at M in Figure 7b and (ii) the inflection point of parabola 3 at K + M in Figure 7c. The calculated EF for MgB2 at zero pressure using the LDA functional in CASTEP (and Δk =0.01A−1) is 8.4055 eV. Figure 7b and c show the location for the origin of EF, at a K + M type reciprocal space position or the midpoint between two reciprocal space Γ vectors along Γ–K. This location is difficult to infer from a reduced BZ scheme, particularly for complex structures. Nevertheless, these two locations at M and ± K ± M nodal points, directly relate to the real space B–B hexagonal plane in the MgB2 structure.
For MgB2, sections of the σ bands along the Γ–M and Γ–K directions are approximated by upward facing parabolas, even when inside the valence band region. Deviations from parabolas occur particularly at zone boundaries where the periodic crystal potential primarily influences free-electron like level crossings [49, 52, 53]. Along Γ–M, a parabola with convex inflection at Γ occurs at −12.55 eV (parabola 1, Figure 7b) and its translated homologs reproduce large sections of the light effective mass σ band distant from Γ.
Similarly, a parabola with convex inflection at M and at -M at −2.152 eV (parabola 3, Figure 7b) reproduces the heavy effective mass σ band. Along Γ–K, differentiation of the σ bands is less pronounced in the extended zones but is apparent at Γ (Figure 7c). In Figure 7c along Γ–K, an additional K–M section is shown because a hexagonal boundary edge, equivalent to K–M by symmetry, transects an adjacent reciprocal space point outside the first BZ (node M1 in Figure 7a). Both Figure 7b and c show the origin of EF for this structure. Table 3 summarizes values of key parameters associated with these parabolic approximations to extended BZ schemes for the Γ–M and Γ–K directions, respectively.
|Orientation||Band Type||Energy at Γ|
|Energy at M|
|Energy at K|
3.4 Fermi energy values
The value of EF for a particular DFT calculation is not only sensitive to the k-grid value as shown above but also to other extrinsic conditions such as compositional substitution in a solid-solution and/or changes in applied pressure. Figure 8 displays the calculated Fermi energies for Al and Sc substitution in Mg1-xAlxB2 and Mg1-xScxB2 determined with the LDA and GGA functionals. For each calculated series, the value of k-grid is constant (i.e. Δk = 0.02 Å−1). Figure 8 also shows the calculated Fermi energies for other end-member compositions with AlB2-type structure. Note that NbB2 and ZrB2 are reported superconductors at very low temperatures (Tc = 0.6 K and 5.5 K, respectively) albeit non-stoichiometric or substituted niobium diboride (
Table 4 lists key parameters based on EBS calculations for MgB2 with applied external pressure. For each calculation, the LDA functional and k-grid value is constant and values are computed after geometry optimization. The value for the effective mass, , is determined from the parabolic approximations described below in Section 3.4.
|Unit cell volume|
|Fermi energy [eV]||EF(P)|
Table 4 shows that as pressure is applied, EF largely conforms to the textbook equation:
In general, Eq. (2) provides estimates of EF at ~93% or more of the DFT calculated value without corrections for charge redistribution along bond directions which take place as pressure is applied . Use of parabolic approximations in this manner may provide a useful benchmark (or “rule of thumb”) for models of predicted new compounds.
3.5 Phonon dispersions – variation with k-grid
As noted in earlier publications [17, 23, 25, 34], the k-grid value also influences the form and mode order of phonons in a DFT calculated PD. Figure 6 demonstrates this influence on the MgB2 PD for the range 0.02 Å−1 < Δk < 0.06 Å−1. The more regularly shaped phonon anomaly becomes apparent with smaller k-grid value and is evident for Δk = 0.02 Å−1 (circled; Figure 8a). For values of Δk > 0.05 Å−1, the calculated PD for MgB2 implies that the phase is unstable yet we know from experimental evidence that this is not the case. For SC compounds with lower values of Tc and/or where Fermi surfaces closely intersect the Fermi level with minimal difference in reciprocal space, the sensitivity of the PD to k-grid value will be shifted towards smaller k-grids compared to the effect with MgB2.
Measurement of the parameter, δ, shown in Figure 9a provides a reliable estimate of Tc for MgB2 when the PD is calculated with Δk < 0.02 Å−1 . This approach, which determines the thermal energy, Tδ, of the key E2g phonon mode using an empirical formula , precisely tracks the experimentally determined reduction of Tc for metal-substituted forms of MgB2 such as (Mg1-xAlx)B2 and (Mg1-xScx)B2 [23, 25]. When the value for δ is determined with two sequential DFT calculations using the LDA and the GGA functionals, error estimates (in terms of the amplitude of the spread of the DFT approximations) for the value of Tδ at each level of metal substitution can be obtained.
We have used this approach to estimate the likely value(s) of Tc for other metal-substituted forms of MgB2 that have received limited attention or have not been identified previously in the literature. For example, we determined the PD for (Mg1-xBax)B2, and for (Mg1-xCdx)B2 where x = 0.33, 0.5 or 0.66 [23, 25]. Figure 10a shows the PD for (Mg0.5Ba0.5)B2 calculated using the LDA functional with Δk = 0.02 Å−1. Measurement of the four values for δ in Figure 10a (
The calculated EF and Fermi level allow systematic comparison of EBSs from a structural family or group of materials with varying properties. Moreover, the FS is generally of a well-defined orbital character and topology determined by the value of EF as a result of bands that cross the Fermi level. The Fermi level determined from DFT calculations is defined as at zero energy while the calculated value of EF obtained after an accurate DFT calculation is seldom described in the published literature. By definition, the Fermi level is determined in the ground state by the filling of lower energy electronic states by all the (nearly) free electrons up to a highest possible value of the energy, which in practice should correspond to EF [49, 50, 55, 57, 58].
4.1 Fermi energy
The superconducting gap in many compounds (
We demonstrate this issue using the EBS for MgB2 as shown in Figure 11. In this figure, we have reproduced the EBS for MgB2 as calculated using the LDA functional for Δk = 0.018 Å−1. The Fermi level is set at 0 eV and a notional “Fermi level 2” is also shown as a red dotted line at −250 meV. As noted in Table 1, a change in calculated EF ~ 200 meV may occur with choice of Δk > 0.04 Å−1. The intersection of σ bands with the calculated Fermi level are separated by a distance λ1, which also defines the separation between Fermi surfaces for MgB2.
However, if the value of Δk, or the calculated value for EF, results in a shift of the Fermi level by ~250 meV, the separation of Fermi surfaces, illustrated by λ2, is different. We estimate that this order of Fermi level or EF shift may result in discrepancies between 20% and 35% of the value(s) for λ. As shown in Figure 11, the shape of the σ bands around Γ are asymmetric. The difference in value(s) for λ with variation in EF, will accordingly be dependent on the form, or shape, of these parabolas. Projections of the density of states at the Fermi level will also be affected by this shift of EF as will the outcomes of Eliashberg or McMillan equations for the determination of Tc.
The apparent discrepancy in determination of the value for EF, noted in the Introduction, may be elucidated by examination of Figure 11. For example, we suggest that some researchers define a value for EF as equivalent to the energy shown as d1 in Figure 11 (
4.2 Computational resolution
Table 5 provides a summary of reports on previous DFT calculations for MgB2 and of systematic calculations from this study. This table highlights the diversity of computational methods used to date as well as wide variations in parameters such as k-grid value and the cut-off energy. Systematic evaluation of these two parameters shows that the value for EF may differ by several hundred meV for the same cut-off energy with change in Δk value. For our systematic calculations of these parameters shown in Table 5, the LDA functional is used for consistency. Calculations with the GGA functional show similar trends albeit at different absolute values (by ~0.2 eV) for EF.
|DFT Code||Δk value (Å−1)||Grid||No. of k-points||Energy Cut-off (eV)||Fermi Energy(eV)||ΔEv# (eV)||LDA or GGA||Ref.|
Table 5 shows that a low value for cut-off energy (
Table 5 also lists the variation in energy, ΔEv (in eV), between the Fermi level and the vertex of the parabola at Γ for different values of Δk and for two cut-off energies using the LDA functional for the EBS of MgB2 (in Figure 11, this energy is represented as d1). As we have noted for EF, there are substantial variations (
The calculations by de la Pena-Seaman  on the transformation of Fermi surfaces with substitution of Al and C into MgB2 and recent work by Pesic
4.3 Phonon dispersions and k-grid
We have examined the changes in PD form and mode order for the substitutional series Mg1-xAlxB2  and Mg1-xScxB2  where 0 < x < 1. For PDs, the value of k-grid in a DFT calculation may obscure phenomena that imply superconductivity such as the presence or absence of a phonon anomaly [5, 34]. We have also demonstrated for MgB2 that the change in the E2g phonon anomaly varies with applied pressure and correlates with the experimentally determined change in Tc . For these cases, we show that a temperature, calculated from the extent of the anomaly, Tδ, is a reliable
In an earlier publication , we compare for MgB2 the calculation of Tc (i) using the McMillan formalism of the Eliashberg model  and (ii) using the E2g phonon anomaly energy, Tδ, as noted above [23, 27, 28]. In both cases, with suitable assumptions for the McMillan formalism, the “predictive” fidelity of either method adequately matches experimental data. However, the Eliashberg model requires an estimate for two key parameters, λ and μ*, based on average values of electron–phonon behavior summed over all orientations. In practice, determination of λ and/or μ* by
An and Pickett  estimate that the influence of the E2g mode is at least a factor of 25 times greater than all other phonon modes in MgB2. The E2g mode is predominantly associated with movement within the boron planes of MgB2; that is, along specific orientations . Nevertheless, use of an average value for phonon frequencies integrated over all directions in reciprocal space is a feature of the McMillan formalism that provides a reasonable “
The predictive value of the approaches we advocate to estimate Tc that utilizes calculation of a value for Tδ using a phonon anomaly [23, 24, 25, 74] is evident for Ba-substitution into MgB2 . Our estimates for (Mg1-xBax)B2 at three levels of Ba substitution (x = 0.33, 0.5 and 0.66) and using both LDA and GGA approximations suggest that 62.1 K < Tδ < 64.4 K with an error of ±4.9 K. These estimates are higher by ~15 K than the experimentally determined value of ~45 K by Palnichenko
The presence of multiple phases in the Rb- and Cs- substituted forms of MgB2 synthesized by Palnichenko
Fully converged PDs are a useful indicator of phase stability [26, 74]. The sensitivity of PDs to changes in stoichiometry, composition or Δk is significantly higher than typically encountered in an EBS . The PD calculated at a deliberately large k-grid value 0.06 Å−1 in Figure 9d may be interpreted as a dynamic instability. MgB2 is a well-studied case and we know that this is not correct; however, for unknown or other materials with closer FSs in reciprocal space, we would expect similar phenomena to be manifest at smaller k-grids. Thus, sometimes conclusions about phase transitions may be artifacts of the DFT calculation if k-grids of insufficient resolution are used for materials with approximately parallel FSs in close reciprocal space proximity .
4.4 Fermi surfaces and superconductivity
Electronic bands and FSs of constant energy possess all point symmetries of a crystal as a function of position in reciprocal space [48, 51]. The intersections of the σ bands with the Fermi level, as shown in Figure 1, determine points that, by construction, belong to the FSs. The FS corresponding to σ bands in the reduced BZ become two approximately parallel tubes , as schematically represented in Figure 12 below. As shown for MgB2 in earlier work [23, 24, 27, 28] and by others [65, 75, 76], these σ band FSs are not strictly cylindrical, but form as warped tubes with a narrowing in all directions towards Γ (sketched more accurately in Figure 2).
Since the FS tubes represent hole carrier sections, their interior will be empty in the ground state, while their exterior will be occupied. In a reduced zone schema, this construct creates ambiguous electron/hole character for the inter-tubular region. Ambiguity arises because this inter-tubular region should be, in the ground state, empty (
This notion creates an apparent dilemma, although according to Ziman , “
Given the indeterminate nature of the origin in reciprocal space, specific diameter tubes may be selected interchangeably by the DFT calculation; thus, implying a potential resonating behavior . Analysis of electron–phonon behavior determined by DFT calculations suggests that this inter-tubular region of FSs (or other regions enclosed by parallel surfaces of different topology) is a region in reciprocal space that reveals the extent of superconductivity in typical BCS-type materials [5, 23, 24, 26, 27, 28]. Our calculations for both MgB2 and B-doped diamond show that this inter-tubular region is of meV energy scale from the Fermi energy.
Parallel FSs are common features of superconducting compounds albeit their identification is dependent on crystal symmetry and the choice of k-grid value for DFT calculations [5, 23, 25]. The “resolution” of reciprocal space calculations using DFT (
Thus, the value of k-grid used for DFT calculations is paramount. For PDs, this computational requirement has previously been well documented [14, 17, 20, 24] and, we suggest, is equally requisite for the use of EBS to predict, or design, new superconducting materials.
The recent development of an ML-based scheme to efficiently assimilate the function of the Kohn-Sham equation, and to directly and rapidly, predict the electronic structure of a material or a molecule, given its atomic configuration  is of salient interest with regard to k-grid value. This ML approach maps the atomic environment around a grid-point to the electron density and local density of states at that grid-point. The method clearly demonstrates more than two orders of magnitude improvement in computational time over conventional DFT calculations to generate accurate electronic structure details . Utilization of this methodology at a k-point spacing <0.2 Å−1 to initialise ML-training for charge density  may enable very rapid determination of potential SC materials with many hundreds of atoms in the base structure. Nevertheless, as we have shown in this article, caution in the use of such values for Δk using ML is suggested because “false positives” for superconductivity may emerge and valid “hits” may be missed.
Thermal effects on electronic properties are generally included in DFT calculations as a smearing of electron behavior. However, high structural symmetry, or the lack of it, may impose significant anisotropy and/or preferred directionality of ionic movement that remains active even as temperature is increased. For reference, thermal excitation of the free-electron gas is kBT or about 26 meV at ambient temperatures [51, 57]. As noted above, variations in EF for superconducting phases may be in the meV range depending on the structure. We also note the importance of the smearing parameter in DFT calculations. We suggest that for particular superconducting cases where the Tc and/or phonon energy is low (
Calculated Fermi energies and Fermi levels are essential attributes for determination of materials properties in a range of other applications, such as for the energy band alignment of components in solar cell materials [78, 79], with solid-electrolyte interfaces , as well as for interface induced phenomena such as the substantial increase in Tc of monolayer FeSe on SrTiO3 substrates . Improved interpretation and understanding of electronic behavior in SCs and SC systems can be achieved with reliable calculated output values determined by
The EBS encapsulates a wealth of information for superconductivity that may be misinterpreted due to the quality, or resolution, of DFT computations. A tendency to be satisfied with poor or limited computational resolution is evident in superconductivity literature unlike other fields that compute electronic properties using DFT. Translation of reciprocal space detail to real space periodicity for DFT-based design of new materials in an EBS with appropriate k-grid resolution can provide evidence for structures that may be viable SCs. As we have shown above, the EF value is explicitly determined in DFT computations and, with consistent use of k-grid resolution, can provide comparable estimates of SC properties for proposed structures of a compositional suite. We encourage inclusion of these DFT calculated parameters in reports of SC materials.
We have described three fundamental approaches, based on
We are uncertain whether these approaches to DFT calculations apply to all SCs recognizing that now hundreds of compounds have been identified. Hardware and software limitations may restrict the use of these approaches to small unit cell structures of simple composition and higher symmetry. Nevertheless, in combination, these systematic and simple approaches to use of a well-known theory of electron distribution in solids suggest that prediction of properties for unknown, or hypothesized, SC structures is well within the reach of many materials researchers.
We appreciate discussions with Peter Talbot on molecular orbital calculations and DFT models. We also appreciate access to, and ongoing assistance with, QUT’s HPC facilities, particularly from Hamish Macintosh, Abdul Sharif and Ashley Wright of the e-Research office. This research did not receive any specific grant from funding agencies in the public, commercial or not-for-profit sectors. AA is grateful for generous financial support for higher degree studies from the University of Hafr AlBatin, Saudi Arabia.
All data referred to in and underpinning this publication are available in QUT Research Data Finder and can be found at DOI: 10.25912/5c8b2cc59a2d9.
Conflict of interest
All authors declare that they have no conflict of interest.