2021

[1,2]Kamal Choudhary

[1] Material Measurement Laboratory, Maryland 20899, USA , National Institute
of Standards and Technology, Gaithersburg, 20899, MD, USA

2]Theiss Research, La Jolla, 92037, CA, USA

# Designing High-$T_{C}$ Superconductors with BCS-inspired Screening, Density
Functional Theory and Deep-learning

<EMAIL_ADDRESS>Kevin Garrity  *  [

###### Abstract

We develop a multi-step workflow for the discovery of conventional
superconductors, starting with a Bardeen–Cooper–Schrieffer inspired pre-
screening of 1736 materials with high Debye temperature and electronic density
of states. Next, we perform electron-phonon coupling calculations for 1058 of
them to establish a large and systematic database of BCS superconducting
properties. Using the McMillan-Allen-Dynes formula, we identify 105
dynamically stable materials with transition temperatures, $T_{C}$ $\geq$ 5 K.
Additionally, we analyze trends in our dataset and individual materials
including MoN, VC, VTe, KB6, Ru3NbC, V3Pt, ScN, LaN2, RuO2, and TaC. We
demonstrate that deep-learning(DL) models can predict superconductor
properties faster than direct first principles computations. Notably, we find
that by predicting the Eliashberg function as an intermediate quantity, we can
improve model performance versus a direct DL prediction of $T_{C}$. We apply
the trained models on the crystallographic open database and pre-screen
candidates for further DFT calculations.

###### keywords:

Superconductors, Deep-learning, JARVIS-DFT, Electron-phonon coupling

##  1 Introduction

Since the discovery of superconductivity in 1911 by Onnes
kamerlingh1911resistance , the identification of novel superconducting
materials, especially those with high transition temperatures ($T_{C}$), has
been an active area of research in condensed matter physics
poole2013superconductivity ; rogalla2011100 . The highest temperature
conventional superconductor in ambient conditions, MgB2 ($T_{C}$ = 39 K), was
discovered relatively recentlynagamatsu2001superconductivity , and progress in
high pressure hydrogen-based superconductorspickard2020superconducting  and
topological superconductorssato2017topological  further motivates the search
for additional ambient or near-ambient condition superconductors with high
$T_{C}$ that may be suitable for scientific or industrial applications.

There have been several previous efforts to computationally and/or
experimentally discover superconducting materials falling into certain
materials classes, such as transition metals roberts1976survey , A15 , B1
kihlstrom1985tunneling ; stewart2015superconductivity , AB2 compounds
ivanovskii2003band ; buzea2001review ; nagamatsu2001superconductivity ,
cuprates plakida2010high , iron-based compounds hosono2015iron , hydrides
shipley2020stability ; liu2017potential ; drozdov2015conventional  and many
other material-classes poole2013superconductivity ; shipley2021high ;
yuan2019recent ; nagamatsu2001superconductivity ; rodriguez1990optical ;
subedi2013electron ; duan2019ab ; kolmogorov2010new ; gao2010high ;
duan2014pressure . However, a more systematic data-driven search should help
expedite the discovery of potentially high-$T_{C}$ superconductors. Moreover,
machine learning has become popular in the search for superconductors. There
have been several reports of machine-learning applications for finding
superconductorshutcheon2020predicting ; stanev2018machine ; yuan2019recent ;
xie2022machine , but thus far they are mostly based on chemical formulas, and
lack detailed atomic structure information that can be critical for
superconducting behavior.

Two key ingredients required to computationally identify
Bardeen–Cooper–Schrieffer (BCS) conventional superconductors cooper2010bcs ;
giustino2017electron  with high-$T_{C}$ are: 1) a robust computational
workflow, and 2) a database of curated materials with prior knowledge such as
elastic constants and electronic density of states. Using density functional
theory perturbation theory (DFT-PT), electron-phonon coupling (EPC) can be
calculated and used to predict $T_{C}$ with reasonable accuracy for arbitrary
materialsgiustino2017electron ; PhysRevB.101.134511 . However, the
computational expense of these calculations is very high when compared to a
single DFT self-consistent energy calculation, especially when fully
converging the relevant $k$-point and $q$-point Brillioun zone samplings for
electrons and phonons. Hence, a fast and reliable computational workflow for
identifying BCS type conventional superconductors that balances computational
cost, accuracy, and scope is needed.

In this work, we develop such a computational approach to discover BCS
superconductors, combining several methods at various levels of computational
expense and accuracy. We start with a BCS-inspired pre-screening based for
materials with high Debye temperature ($\theta_{D}$) and high electron density
of states (DOS) at Fermi-level ($N(0)$), using the existing JARVIS-DFT
database choudhary2020joint . We then develop and apply a DFT-PT workflow to
compute $T_{C}$ using electron-phonon coupling and the McMillan-Allen-Dynes
formula mcmillan1968transition , with initially low convergence settings. We
benchmark the DFT workflow on known superconductors and apply it to materials
from our pre-screening step. For the best candidates, we perform additional
convergence tests to validate our predictions. In addition, we use the dataset
to develop deep-learning models using the atomistic line-graph graph neural
network (ALIGNN) choudhary2021atomistic ; choudhary2021recent  to predict the
Debye temperature, electronic DOS, $T_{C}$, and electron-phonon coupling
parameters for arbitrary crystal structures.

We utilize the publicly available JARVIS choudhary2020joint  infrastructure to
achieve the goals mentioned above. JARVIS (Joint Automated Repository for
Various Integrated Simulations, https://jarvis.nist.gov/) is a collection of
databases and tools to automate materials design using classical force-field,
density functional theory, machine learning calculations and experiments. In
particular, we obtain elastic tensor and DOS data from JARVIS-DFT database,
establish the DFT workflow with JARVIS-Tools package and train the deep-
learning model using ALIGNN. JARVIS-DFT is a density functional theory based
database of 55645 materials with several material properties such as formation
energy, bandgap with different level of theories choudhary2018computational ,
solar-cell efficiency choudhary2019accelerated , topological spin-orbit
coupling spillage choudhary2021high ; choudhary2019high ;
choudhary2020computational , elastic tensors choudhary2018elastic , dielectric
tensors, piezoelectric tensors, infrared and Raman spectrum choudhary2020high
, electric field gradients choudhary2020density , exfolation energies
choudhary2017high  etc. with stringent DFT-convergence setup
choudhary2019convergence . Including $T_{C}$ in JARVIS-DFT would greatly
enrich the applicability of the database and aid in guided materials design.
In this work, we add the $T_{C}$ data in JARVIS-DFT.

##  2 Results and discussion

###  2.1 First principles workflow

A flow chart for designing conventional superconductors using BCS inspired
screening (iBCS), JARVIS-Tools based DFT-Screening (J-Scr), and ALIGNN in an
integrated way is shown in Fig. 1. The combination of computationally
expensive DFT-based screening with empirical rules-based screening and deep
learning methods allows us to apply our combined workflow to a wider group of
materials than a brute-force computational approach.

Figure 1: Schematic showing the steps involved in identifying high-$T_{C}$
superconductors. a) statistical distribution of Debye temperature (K) and b)
statistical distribution of electronic density of states (states/eV/total
number of electrons) at Fermi level from the JARVIS-DFT database, c)
probability that compounds containing a given element have $\theta_{D}>$300 K.
The flow chart shows the application of BCS-inspired screening, density
functional theory calculations and deep-learning training.

First, we screen likely superconductors using a criteria inspired by the BCS
equation for $T_{C}$ (Eq.1). According to this equation, high-$T_{C}$
materials usually have high Debye temperature ($\theta_{D}$), and high
electronic DOS (states/eV/total number of electrons) at the Fermi level
($N(0)$). The JARVIS-DFT database provides DFT-based electronic DOS at the
Fermi level choudhary2018computational  and the Debye temperature, as derived
from the elastic tensor databasechoudhary2018elastic .

Currently, there are electronic DOS database available for 55723 materials and
elastic tensors for 17419 materials (using the v08.18.2021 version of JARVIS-
DFT database at the time of writing). We analyze the statistical distribution
of $\theta_{D}$ and $N(0)$ in Fig. 1a and Fig. 1b, respectively. We observe
that the Debye temperature can range up to 1500 K with high peak near 200 K
(Fig. 1a). The DOS at the Fermi level can range up to 10 states/eV/Nelect with
peak value around 1 (Fig. 1b). Out of 17419 materials, we find 5618 of them
with $\theta_{D}$ greater than 300 K. Furthermore, selecting materials with
electronic DOS at the Fermi level greater than 1 states/eV/Nelect, we find
1736 materials. In Fig. 1c, we present the probability that compounds
containing a given element have $\theta_{D}>$ 300 K. We observe that light
elements that form strong covalent bonds, including the $1s$, $2p$, and $3p$
elements, as well as many lighter $3d$ transition metals, tend to have high
$\theta_{D}$, while denser transition metals and elements that form weaker
bonds tend to have low $\theta_{D}$. There are multiple methods to calculate
Debye temperature anderson1963simplified ; garrity2016first . In this work, we
primarily use the elastic tensor data, as it was available before our phonon
calculations and thus available for pre-screening, but we compare to an
approach based on the phonon DOS in Supplementary Figure 1.

Figure 2: Effect of different k-point and q-point selection on EPC parameters
for VC. a) $\lambda$, b) $\omega_{log}$, c) $T_{C}$. Broadening dependent
scaled d) $\lambda$, e) $\omega_{log}$ and f) $T_{C}$ for the all the systems
in the database. The scaling was done with respect to the corresponding values
at broadening of 0.05 Ry. The legend in panel a corresponds to the number of
k-points and q-points, e.g. 10x10x10 k-points and 10x10x10 q-points (denoted
by black star) is same as 1 k-point per q-point; and 20x20x20 k-point and
2x2x2 q-point correspond to 1000 k-point per q-point (denoted by left purple
triangles).  Figure 3: Broadening convergence with respect to experimental
$T_{C}$ data for 14 superconductors.

Performing fully converged electron-phonon coupling calculations using DFT-PT
is very computationally expensive, as the calculations in general require both
high density $k$-point grids to sample to electronic states and high $q$-point
sampling for the phononswierzbowska2005origins , and the number of modes to
calculate at each q-point increases with the number of atoms in the primitive
cell. Therefore, in the supplementary information, we perform a large number
of convergence checks to understand the minimal set of convergence criteria
that provide useful information in identifying possible high-$T_{C}$
conventional superconductors (see Supplementary Table 1-3). We show the
convergence of $\lambda$, $\omega_{log}$ and $T_{C}$ with respect to q-points
in Supplementary Figure 2, Supplementary Figure 3 and Supplementary Figure 4
respectively. The comparison for the same for PBE perdew1996generalized  vs
PBEsol are shown in Supplementary Figure 5, Supplementary Figure 6,
Supplementary Figure 7 and Supplementary Figure 8 respectively.

In Fig. 2, we show the effect of broadening parameter on EPC and $T_{C}$
values. As an example of convergence results, in Fig. 2a, Fig. 2b and Fig. 2c,
we show convergence with several k-point and q-point settings for the VC
compound. From Fig. 2a and 2b, it is clear that denser k-points are necessary
for reasonable $\lambda$ and $\omega_{log}$ values. We find that larger
broadenings are necessary for lower q-point grids. Similar behavior was also
found by Shipley et al. for rare-earth hydride systemshipley2020stability .

In Fig. 2a, we find that though for 20x20x20 k-points and 2x2x2 q-points case,
$\lambda$ and $\omega_{log}$ values (panel a and panel b) are different from
other selections, the $T_{C}$ value is still closer to 10x10x10 and 2x2x2
k-points and q-points (panel c) respectively. We find that $T_{C}$ for
10x10x10 k-points and 10x10x10 q-points (black plus sign) are similar to that
of 10x10x10 k-points and 5x5x5 q-points (orange upper triangle), suggesting
that increasing q-points do not affect much for the 10x10x10 k-points in terms
of predicting $T_{C}$. Similarly, for 10x10x10 k-points, 2x2x2 q-points (right
green triangle); 20x20x20 k-points , 2x2x2 q-points (left purple triangle),
the $T_{C}$ values are closer suggesting in this case increasing k-points do
not change the $T_{C}$ much for 20x20x20 k-points. Interestingly, for the
20x20x20 k-points and 4x4x4 q-points case, the $T_{C}$ lies in-between the
above two extremes of $T_{C}$. We note that though $\omega_{log}$ seem to
converge and overlap for right green and downward red triangles in Fig. 2b,
the $\lambda$ value is still slightly different (as from Fig. 2a) leading to
large difference in $T_{C}$ values indicating that it harder to converge
$\lambda$ that $\omega_{log}$.

Next, we compare the effect of broadening parameters on the EPC parameters and
$T_{C}$ for the entire database in Fig. 2d, Fig. 2e and Fig.2f. We find
relatively consistent results for broadenings of 0.03-0.05 with our relatively
sparse q-point grids, while smaller broadenings lead to erratic unconverged
results. In Fig. 3 we compare the effects of broadening on our prediction of
$T_{C}$ for 14 known superconductor (see table 1), and we find that higher
broadening values give better agreement.

We investigate the effect of $\mu^{*}$ parameter lee1995first ;
giustino2017electron ; marques2005ab ; morel1962calculation ;
zarifi2018crystal ; wierzbowska2005origins  in Eq. 7 on $T_{C}$ for several
materials in Fig. 4. In Fig. 4a, we show the effect of different $\mu^{*}$ on
5 well-known conventional superconductors (Al, Pb, MgB2, NbN, V3Si). We find
that as we increase the $\mu^{*}$ values, the $T_{C}$ decreases, as expected.
The rate of decay somewhat is dependent on the material, but the materials
remain in largely the same order as $\mu^{*}$ is varied. Furthermore, we
compute the mean absolute error in predicting $T_{C}$ of 14 known
superconductors in table 1 using different $\mu^{*}$ values. We find that
minimum MAE was obtained for $\mu^{*}=0.09$. Similarly, we compare the $T_{C}$
for the entire database for different $\mu^{*}$ in Fig. 4c. We compare these
values to $\mu^{*}=0.09$ which is shown as a solid black line. We find that
increasing $\mu^{*}$ decreases the $T_{C}$ similar to Fig. 4a. Now, we
quantify, how many materials in our dataset have $T_{C}\geq 5K$ as we vary
$\mu^{*}$. Relative to $\mu^{*}=0.09$, (denoted by dotted green line), we find
that we gain 85.5 % and 37.4 % materials for $\mu^{*}$ 0.03 and 0.06
respectively. As we increase the $\mu^{*}$ to 0.12, 0.15 and 0.18 we lose 29
%, 43.5 % and 57.3 % materials respectively.

Based on these tests, we find that $q$-point grids as small as $2\times
2\times 2$, combined with k-point grids similar to the typical grids used for
self-consistent DFT total energy calculations are already useful in
identifying candidate materials. Therefore, we adopt the strategy of
performing an initial screening step with low convergence parameters applied
to a high number of materials, which will be further refined later. During
this step, we use the same $k$-points grid found during the JARVIS-DFT total
energy convergence, as well as a $q$-point grid of at least $2\times 2\times
2$ with a broadening of 0.05 Rydberg ($\approx 0.68$ eV).

Figure 4: Effect of $\mu^{*}$ on transition temperatures. a) Change in $T_{C}$
as we increase $\mu^{*}$ for a few representative superconductors, b) change
in $T_{C}$ with respect to $\mu^{*}$ 0.09 for the entire dataset, c) change in
percenatge materials with $T_{c}$ greater equal to 5 K. Here $\mu^{*}$ 0.09 is
shown as 100 %. As we decrease $\mu^{*}$, we add materials and vice-versa.

We compare the present DFT-screening workflow (J-Scr) with other methods such
as superconducting density functional theory (SCDFT) sanna2020combining ,
Lüders-Marques (LM) luders2005ab ; marques2005ab  and experiments
allen1975transition ; roberts1976survey ; webb2015superconductivity  in Table
I. Note that we ignore spin-orbit coupling, spin-polarization, and spin-
fluctuation contributions heid2010effect ; gibson1989evidence ;
mandal2014strong ; gerber2017femtosecond  during the calculations because of
the additional computational cost, but these effects should be considered in
follow-up investigations. Importantly, we find that for the top seven systems
shown in Table I, the $T_{C}$ from the McMillan-Allen-Dynes formula, using
$\mu^{*}=0.09$, from our J-Scr workflow has an excellent agreement with
respect to both SCDFT and the experiments, justifying our approximations. We
find that the mean absolute error (MAE) for J-Scr and experiment for all the
top 14 near ambient condition superconductors is 3.3 K. Now, comparing the top
7 materials in Table I with respect to SCDFT and LM, the MAE of experiment vs
J-Scr, experiment vs SCDFT, experiment vs LM, SCDFT vs J-Scr are: 1.9 K, 1.4
K, 4.9 K and 1.1 K respectively.

Table 1: Comparison of JARVIS-Tools screening workflow based (J-Scr) critical transition temperatures in K with experimental data (Exp) allen1975transition ; roberts1976survey ; webb2015superconductivity ; drzazga2018characteristics ; dew1979formation ; tsindlekht2008linear ; RevModPhysMatthias ; hulm1972superconductivity ; pechen1994tunneling , SCDFT sanna2020combining  and LM methods kruglov2020superconductivity ; liu2017potential . Form. | Spg | JID | Exp | SCDFT | LM | J-Scr
---|---|---|---|---|---|---
Al | 225 | 816 | 1.2 | 1.55 | 0.3 | 1.6
Ta | 229 | 1014 | 4.5 | 5.17 | 2.45 | 7.6
Pb | 225 | 961 | 7.2 | 6.06 | 4.95 | 5.4
Nb | 229 | 934 | 9.3 | 10.29 | 7.0 | 10.7
ZrN | 225 | 19679 | 10.0 | 11.6 | 6.12 | 10.0
V3Si | 223 | 14960 | 17 | 18.1 | 13.1 | 17.6
MgB2 | 191 | 1151 | 39 | 35.4 | 20.04 | 33.0
V | 229 | 14837 | 5.3 |  |  | 18.3
Nb3Si | 223 | 15938 | 18.0 |  |  | 16.5
NbO | 221 | 14492 | 1.38 |  |  | 3.6
NbC | 225 | 19889 | 12 |  |  | 17.1
NbN | 221 | 36335 | 16.0 |  |  | 17.6
YB6 | 221 | 20620 | 7.2 |  |  | 5.1
Nb3Al | 223 | 11981 | 16.8 |  |  | 9.0
YH10(250GPa) | 225 |  | 260 |  |  | 213.5
LaH10(250GPa) | 225 |  | 211 |  |  | 190.0
Table 2: JARVIS-Tools screening workflow based $T_{C}$ for some of the potential candidate superconductors ($T_{C}$), their chemical formula (Form.), spacegroup number (Spg), JARVIS-DFT ID (JID), Inorganic Crystal Structure Database ID (ICSD) belsky2002new  wherever available, JARVIS-DFT based formation energy ($E_{form}$ (eV/atom)) and energy above convex hull ($E_{hull}$ (eV)). For the complete list, refer to the supplementary information  Form. | Spg | JID | ICSD | $E_{form}$ | $E_{hull}$ |  $T_{C}$(K)
---|---|---|---|---|---|---
MoN | 187 | 16897 | 187185 | -0.47 | 0.09 | 33.4
CaB2 | 191 | 36379 | 237011 | -0.25 | 0.09 | 31.0
ZrN | 194 | 13861 | 161885 | -1.76 | 0.18 | 30.0
VC | 225 | 19657 | 619079 | -0.48 | 0.06 | 28.1
V2CN | 123 | 105356 | - | -0.82 | 0.11 | 26.2
Mn | 225 | 25344 | 41509 | 0.08 | 0.08 | 23.0
NbFeB | 187 | 4546 | - | -0.15 | 0.39 | 22.1
NbVC2 | 5 | 102190 | - | -0.46 | 0.08 | 21.9
ScN | 225 | 15086 | 290470 | -2.15 | 0.0 | 20.8
LaN2 | 2 | 118592 |  | -1.05 | 0.0 | 20.4
VRu | 221 | 19694 | 106010 | -0.22 | 0.01 | 20.3
TiReN3 | 161 | 36745 | - | -0.68 | 0.10 | 20.0
B2CN | 51 | 91700 | 183794 | -0.53 | 0.19 | 19.4
KB6 | 221 | 20067 | 98987 | -0.09 | 0.0 | 19.0
ZrMoC2 | 166 | 99893 | - | -0.49 | 0.08 | 17.9
TaB2 | 191 | 20082 | 30420 | -0.60 | 0.0 | 17.2
NbS | 194 | 18923 | 44992 | -0.98 | 0.05 | 17.0
TaVC2 | 166 | 101106 | - | -0.54 | 0.05 | 16.3
TaC | 187 | 36405 |  | -0.24 | 0.40 | 16.1
MgBH | 11 | 120827 | - | -0.03 | 0.11 | 15.5
CoN | 216 | 14724 | 236792 | -0.02 | 0.0 | 15.0
NbRu3C | 221 | 8528 | 77216 | -0.02 | 0.19 | 15.0
TiB4Mo | 191 | 105371 | - | -0.67 | 0.06 | 14.8
Al4CO | 6 | 53086 | - | -0.60 | 0.41 | 14.4
ScC | 194 | 37100 | - | -0.02 | 0.20 | 14.3
TaS | 187 | 4699 | 52114 | -0.74 | 0.11 | 14.2
VOs | 221 | 122961 |  | -0.26 | 0.0 | 14.1
RuO2 | 136 | 19852 | 236962 | -1.23 | 0.0 | 14.0
NbVCN | 160 | 99944 | - | -0.78 | 0.15 | 13.6
NbB2 | 191 | 14726 | 614908 | -0.65 | 0.0 | 13.8
CaRu2N2 | 164 | 108175 | - | -0.32 | 0.10 | 13.3
ZrNbCN | 160 | 102642 | - | -1.18 | 0.06 | 13.1
V3Al | 221 | 36199 |  | -0.15 | 0.14 | 13.0
TiO | 221 | 50092 |  | -2.65 | 0.08 | 12.6
LiBe2Ru | 123 | 71476 |  | -0.14 | 0.08 | 12.5
V2TcRu | 225 | 41603 |  | -0.3 | 0.0 | 12.3
ScSi2 | 191 | 4406 | 651822 | -0.31 | 0.23 | 12.2
VTe | 221 | 122995 |  | 0.22 | 0.50 | 11.3
AsP2 | 191 | 120382 |  | 0.45 | 0.45 | 11.2
Cr3GaN | 221 | 15120 | 626036 | -0.30 | 0.04 | 10.9
CrH | 225 | 18397 | 626118 | 0.06 | 0.06 | 10.7
TiReTc2 | 225 | 39667 |  | -0.27 | 0.0 | 10.7
ZrSc | 187 | 100313 |  | -0.002 | 0.01 | 10.3
V3Pt | 223 | 20434 | 649824 | -0.43 | 0.0 | 10.3
La3SnC | 221 | 101642 | - | -0.60 | 0.0 | 10.2
NbS | 187 | 20636 | 645320 | -1.0 | 0.07 | 10.1
FeN | 216 | 15718 | 236789 | -0.30 | 0.0 | 10.0
InP | 225 | 17716 | 188692 | -0.01 | 0.20 | 10.5
V3Ni | 223 | 15777 | 647028 | -0.16 | 0.0 | 10.0
NbRu3 | 221 | 16588 | 105224 | -0.12 | 0.0 | 10.0
MnN | 216 | 78284 | 236788 | -0.44 | 0.0 | 10.0
CaZrO2 | 99 | 114723 |  | -2.76 | 0.30 | 10.0
ScSe | 225 | 19974 | 44972 | -1.68 | 0.0 | 9.1

This suggests that there is a close agreement between the computational
methods. We also evaluate the workflow for the well-known high-pressure
hydrides shipley2020stability ; liu2017potential ; drozdov2015conventional
such as LaH10 and YH10 at 250 GPa as shown in Table I. We find the MAE to be
25 K for these two systems, which is reasonable for high-throughput type
screeningshipley2020stability ; liu2017potential ; drozdov2015conventional .
Therefore, the workflow can be applied to study high-pressure
superconductivity as well. However, for the present work, we analyze materials
without external pressure.

Next, we apply the J-Scr workflow on materials with number of atoms less than
equal to 5, $\theta_{D}>$300K and $N(0)>$1 as discussed in the previous
section. As of now, we have applied J-Scr to 1058 materials. Out of these,
only 626 of them are dynamically stable (i.e. no imaginary phonon modes) and
105 of them have $T_{C}>$5 K, which is very promising. Table II provides the
chemical formula, spacegroup, JARVIS-DFT ID, Inorganic Crystal Structure
Database (ICSD) ID belsky2002new , and the OptB88vdW-based formation energy
from JARVIS-DFT and predicted $T_{C}$ for some of the best candidates for
experimental synthesis. The JARVIS-IDs can be used to obtain more detailed
information, e.g., https://www.ctcms.nist.gov/ knc6/static/JARVIS-
DFT/JVASP-16897.xml for JARVIS-ID: 16897

We find that several potential candidates are based on nitrides, borides and
Vanadium containing compounds. All of these systems have high symmetry, with
spacegroup numbers ranging from 99 to 225, although this may be an artifact of
pre-selecting systems with few atoms. We visualize some of the crystal
structures in Supplementary Figure 9. Many of these systems have hexagonal and
cubic symmetry. The JARVIS-DFT identifiers for corresponding systems can be
used to obtain further properties. Importantly, most of these systems have
negative formation energies and energy above convex hull less than 0.5
eV/atom. We also show Fermi-surfaces of such example compounds in Fig. S9. The
shape of the Fermi-surface depicts the electron motion inside a material
yamase2021fermi . In earlier works an2021superconductivity  related to MgB2
type compounds, the shape of the Fermi surface helped understand the mechanism
of superconductivity in such compounds.

As shown in table2, some of our highest $T_{C}$ candidate superconductors
include MoN, VC, Mn, MnN, LaN2, KB6, TaC, etc. MoN in a rock-salt structure
has been previously been reported with a $T_{C}$ as high as 30 K
inumaru2008high ; wang2015hardest ; however, here we propose a hexagonal form
which has not been observed experimentally. Similarly, ZrN in rock-salt form
is known to have $T_{C}$ of 10 K; however, here we propose a hexagonal form of
ZrN to have a $T_{C}$ around 30 K. CaB2 was previously theoretically predicted
be a high-$T_{C}$ superconductorchoi2009prediction  . Previous reports of
synthesizing VC in 1:1 ratio has been found to be challenging by previous
researchers zbasnik1973electronic . Lanthanum nitride
vaitheeswaran2002structural  in 1:1 ratio compound has been found to have
superconducting properties. Lanthanum nitride in 1:2 ratio is not
characterized for the superconducting properties to the best of our knowledge.
Vanadium Ruthenium alloy with a different ratio than 1:1 has been reported
earlier as superconducting with $T_{C}$ close to 5 K chu1971electronic . KB6
has been proposed to be a high-$T_{C}$ superconductor as well
katsura2010possibility . Tantalum carbide has been shown to be a potential
topological superconductor with a $T_{C}$ of 10.3 K yan2021superconductivity .
Li et. al. li2020superconducting  recently reported the CoN system as a
potential superconductor with $T_{C}$ of 16 K. Nb-Ru-C alloy system has been
synthesized experimentally before hong2019structural , but no superconducting
data is available to the best of our knowledge. Ta-S system with a 1:2 ratio
is a superconductor sernetz1974superconductivity  with $T_{C}$ of 5.4 K,
however, a 1:1 ratio material i.e Tantalum rich system has not been
characterized for superconducting properties yet. Nb-B system has been
reported earlier with a $T_{C}$ of around 9.8 K ren2008structural . V-Al alloy
system can attain $T_{C}$ of around 11.15 K as reported by Kodess et al.
Kodess71 . The Sc-Si system has been realized experimentally before but their
superconducting properties has not been reportedgokhale1986sc . 2D B2O system
was recently proposed with a $T_{C}$ of 10.3 K yan2020theoretical . Cr3GaN
belongs to the cubic anti-perovskite family of superconductors and has been
predicted to be superconductor by Tutuncu et al. tutuncu2013theoretical . V3Pt
belongs to the A15 family and has been investigated for superconducting
properties earlier ramakrishnan1986resistivity . Fe-N system with 1:2 ratio
was shown with a $T_{C}$ of 8 K chen2018novel .

Next, we show the relationship between EPC parameters for the dynamically
stable compounds and the Eliashberg spectral functions of some of the
candidate materials in Fig. 5. In Fig. 5a, we find that $\lambda$ and
$\omega_{log}$ to have an inverse relationship while in Fig. 5 we see that
$\lambda$ and $T_{C}$ tend to follow a linear relationship,  both of which are
typical for BCS superconductors. For high-$T_{C}$, high $\omega_{log}$ as well
as high $\lambda$ is favorable as evident from the colormap. The EPC
Eliashberg spectral function expresses the electron-phonon interaction in the
form of a spectral density. The weighted area under the EPC function
determines the $\lambda$ as well as $\omega_{log}$ parameters.

Figure 5: Relation between Electron-phonon coupling parameters and EPC
function of some of the potential candidate superconductors. a) $\omega_{log}$
vs $\lambda$, b) $T_{C}$ vs $\lambda$, c) MoN, d) VC, e) KB6, and f) VTe.

###  2.2 Machine Learning

Next, we develop deep-learning models to accelerate both our initial BCS-
inspired screening and our calculation of the electron-phonon coupling
parameters.

Figure 6: Atomistic line graph neural network based deep-learning (DL)
regression model performance on 5 % test set for a) Debye temperature and b)
DOS. Classical force-field descriptor (CFID) (c,d,e) and DL (f,g,h) based
regression model performance on 5 % test set for DFT calculated transition
temperature ($T_{C}$), EPC parameter $\omega_{log}$, and EPC parameter
$\lambda$. In (f), we show performance with direct $T_{C}$ prediction (red
color), $T_{C}$ prediction with direct prediction of wlog and lambda and then
using eq. 7 (green color) and $T_{C}$ prediction with Eliashberg function and
then using eq. 5-8 (black color).

While the BCS pre-screening step is much less expensive than the full EPC
calculation, it still requires the DOS and $\theta_{D}$, which require
substantial computation. Therefore, we develop regression models to predict
these properties directly from an arbitrary crystal structure, using the large
datasets available in the JARVIS-DFT database. There are multiple methods to
establish the structure-property relation for crystal structures, but in this
work we use atomistic line graph neural network (ALIGNN) model
choudhary2021atomistic , which has been shown to outperform many well-known
benchmarks for solids and explicitly capture chemical and many-body physical
interactions in materials. Our results of these models on 5 % held test sets
are shown in Fig. 6 and Table. III. The baseline model MAE was computed by
using mean of the target values in the training dataset and using it as
predictions for all the materials in the 5% test data. We observe that the
mean absolute error for the Debye temperature is 49 K while that for DOS is
1.5 states/eV/Nelect. The baseline model MAE for the Debye temperature and DOS
are 145.5 K and 3.62 states/eV/Nelect respectively. Using these trained
models, we applied it to Crystallography Open Database (COD)
gravzulis2009crystallography  with 431778 materials and with number of atoms
less than 100, and pre-screened 8293 materials with high Debye temperature and
DOS values.

Figure 7: Prediction of Eliashberg function with ALIGNN for the 5 % test set.
ALIGNN can capture peak-positions and heights reasonably well.

Next, we develop machine-learning models to directly predict EPC parameters
and $T_{C}$ using our database of 626 dynamically stable data-points from the
j-Scr developed above. We use two methods: hand-craft descriptors (Classical
force-field descriptor (CFID) choudhary2018machine ) and a deep-learning
approach (atomistic line graph neural network (ALIGNN)) choudhary2021atomistic
. In particular, we train models for the McMillan-Allen-Dynes transition
temperature ($T_{C}$) and the electron phonon coupling parameters-
$\omega_{log}$ and $\lambda$. We note that ML models usually require large
datasets with sample sizes in the range of thousands, but we show preliminary
but useful results already with our current dataset, which will continue to
grow.

Classical force-field descriptor (CFID) choudhary2018machine  based
performances are shown in Fig. 6c, Fig. 6d and Fig. 6e for $T_{C}$,
$\omega_{log}$ and $\lambda$ respectively. Similarly, ALIGNN based
performances are shown in Fig. 6f, Fig. 6g and Fig. 6h. For a perfect
agreement, the data-points should lie on the $y=x$ line. We observe that the
MAE using the CFID approach for $T_{C}$, $\omega_{log}$ and $\lambda$ are 1.84
K, 53.48 K and 0.19. Similarly, the MAEs for ALIGNN approach are 1.84 K, 37.43
K and 0.14 respectively. ALIGNN outperforms CFID in predicting $\omega_{log}$,
but the other performances are similar. Moreover, we notice that it is
significantly easier to learn $\omega_{log}$ than $T_{C}$ and $\lambda$ as
evident from the parity plots and table 3. In particular, the model for
$\lambda$ is only slightly better than the baseline model. Using the ALIGNN
models for $\omega_{log}$ and $\lambda$, and eq. 7, we predict the TC with an
MAE of 1.77 K. We apply the ALIGNN TC predictor model on the COD database and
find that 39595 candidate materials can have $T_{C}\geq 10K$ and 2161 of them
have $T_{C}\geq 15K$.

Table 3:  DL model performance for direct property predictions. We compare the mean absolute error (MAE) of a baseline model, the MAE for the DL model on the 5% test data, and the root mean squared error (RMSE) for the 5 % test data. The baseline model always predicts the mean of the target values in the training dataset. Property | BaselineMAE | TestMAE | TestRMSE
---|---|---|---
$\theta_{D}$ (K) | 145.5 | 49.0 | 124.62
$N(0)$ (states/eV/Nelect) | 3.62 | 1.5 | 4.40
$\lambda$ | 0.20 | 0.14 | 0.24
$\omega_{log}$ (K) | 76.0 | 37.43 | 50.73
$T_{C}$ (K) | 2.72 | 1.84 | 3.31

As predicting $\lambda$ directly using ML models is evidently challenging, we
attempt an alternate method to directly predict Eliashberg function using the
ALIGNN model. We choose an energy range of 0 to 100 meV with 1 meV binsize and
predict the Eliashberg functions. We show the DFT and ALIGNN based Eliashberg
functions for samples in the test set in Fig. 7. We find that the ALIGNN model
does an good job of capturing most of the peaks. We calculate the $T_{C}$
using the ALIGNN based Eliashberg function predictions, and find the MAE is
1.39 K, which improves on our direct prediction method above by 24%.
Interestingly, we have observed a similar behavior for predicting phonon
properties using phonon density of states as an intermediate quantity in
Ref.gurunathan2022rapid . This suggests that learning more fundamental and
information-rich quantities such as Eliashberg functions can be useful for ML
approaches with limited data, as compared to more direct predictions of
integrated quantities. Unfortunately, a drawback of deep-learning models is
that it is difficult to extract physical insight from their internal
parameters, but we hope to investigate these ideas further in future works.

The above ML-exercise indicates that while $\omega_{log}$ is relatively
straightforward to predict using the crystal structure, $\lambda$ and
consequently $T_{C}$ are much more challenging. We note that deep learning
models usually require larger datasets of at least ten thousand entires, which
remains well beyond the data available in this work. In addition to larger
datasets, these results should motivate future works on designing better
descriptors and models.

In summary, we have developed a combined high-throughput DFT and ML approach
to study conventional superconductors, finding over one-hundred candidate
materials with predicted $T_{C}\geq 5K$. We provide data from this work as
well as machine learning models to help accelerate the discovery of
superconductors. Because we employ a high-throughput approach to screen large
databases, we employ several assumptions, but we perform significant benchmark
testing as well as detailed convergence checks on particular materials found
to be promising in our initial screening to verify our results. We have made
our datasets and tools publicly available to enhance the reproducibility and
transparency of our work. We believe that our work can be of great help to
guide future computational as well as experimental efforts to discover and
characterize BCS-superconductors.

##  3 Methods

In this section, we discuss details of various computational methods.

###  3.1 BCS inspired screening

According to BCS-theorybardeen1957theory  , the attractive electron-electron
interaction mediated by phonons gives rise to Cooper pairs, i.e. bound states
are formed by two electrons with opposite spins and momenta. BCS-theory
provides the relation between the superconducting transition temperature
($T_{C}$), Debye temperature ($\theta_{D}$), electronic DOS at Fermi level
$N(0)$, and electron-phonon interaction ($V$) as follows:

$T_{c}=1.14\theta_{D}\exp\bigg{(}-\frac{1}{N(0)V}\bigg{)}$ (1)

$\theta_{D}$ is defined as anderson1963simplified :

$\theta_{D}=\frac{h}{k_{B}}\bigg{[}\frac{3nN_{a}\rho}{4\pi
M}\bigg{]}^{\frac{1}{3}}v_{m}$ (2)

where $h$ is Planck’s constant, $k_{B}$ is the Boltzmann constant, $n$ is the
number of atoms per formula unit, $N_{A}$ is Avogadro constant, $\rho$ is the
crystal structure’s density, $M$ is the molar mass, and $v_{m}$ is the average
sound velocity obtained from the elastic tensor anderson1963simplified .

We use the $\theta_{D}$ obtained from finite-difference calculations of
elastic tensorschoudhary2018elastic , and $N(0)$ is obtained from the DOS at
the Fermi-level as available in the JARVIS-DFT database. Currently, there are
electronic DOS database available for 55723 materials and elastic tensors for
17419 materials (v08.18.2021) and the database is continuously expanding. We
normalize the electronic density of states by the total number valence
electrons in the DFT calculation. JARVIS-DFT is primarily based on Vienna Ab
initio Simulation Package (VASP) kresse1996efficient ; kresse1996efficiency .
software and OptB88vdW klimevs2009chemical  functional but also contains data
with different functionals and methods. In JARVIS-DFT k-points are converged
with respect to total energy choudhary2019convergence . These converged
k-points are also used in subsequent electron-phonon calculations using
Quantum Espresso package giannozzi2009quantum ; giannozzi2020quantum  (see
below).

###  3.2 Density functional theory (DFT) calculations

The superconducting properties of conventional superconductors depend on EPC.
There are several methods to perform EPC calculations such as
interpolated/Gaussian broadening method wierzbowska2005origins , tetrahedron
method kawamura2014improved  and other Wannier-based electron-phonon method
ponce2016epw , etc. In this work, we primarily report the interpolated method,
but we also compare with the tetrahedron method (see supplementary
information). We perform EPC calculations using DFT-PTbaroni1987green ;
gonze1995perturbation  with the Quantum Espresso software package
giannozzi2009quantum  and the GBRV garrity2014pseudopotentials
pseudopotentials. We report results with the PBEsolperdew2008restoring
functional in the main text, and we compare to the PBEperdew1996generalized
functional in the supplementary information. We begin with structures from the
JARVIS-DFT database, and perform full relaxation using our Quantum Espresso
settings. The EPC parameter is derived from spectral function
${\alpha}^{2}F(\omega)$ which is calculated as follows:

${\alpha}^{2}F(\omega)=\frac{1}{2{\pi}N({\epsilon_{F}})}\sum_{qj}\frac{\gamma_{qj}}{\omega_{qj}}\delta(\omega-\omega_{qj})w(q)$
(3)

where $\omega_{qj}$ is the mode frequency, $N({\epsilon_{F}})$ is the DOS at
the Fermi level ${\epsilon_{F}}$, $\delta$ is the Dirac-delta function, $w(q)$
is the weight of the $q$ point, $\gamma_{qj}$ is the linewidth of a phonon
mode $j$ at wave vector $q$ and is given by:

$\gamma_{qj}=2\pi\omega_{qj}\sum_{nm}\int\frac{d^{3}k}{\Omega_{BZ}}\mbox{{}\rm\sf\leavevmode\hbox{}\hbox{}g\\_{kn,k+qm}^j\/}^2
\delta(\epsilon_{kn}-\epsilon_F) \delta(\epsilon_{k+qm}-\epsilon_F)
\end{equation}Here,theintegralisoverthefirstBrillouinzone,\epsilon_{kn}and\epsilon_{k+qm}aretheDFTeigenvalueswithwavevectorkandk+qwithinthenthandmthbandsrespectively,g_{kn,k+qm}^{j}istheelectron-
phononmatrixelement.\gamma_{qj}isrelatedtothemodeEPCparameter\lambda_{qj}by:\par\begin{equation}\lambda_{qj}=\frac{\gamma_{qj}}{\pi
hN(\epsilon_{F})\omega_{qj}^{2}}\end{equation}\par\par
Now,theEPCparameterisgivenby:\par\begin{equation}\lambda=2\int\frac{\alpha^{2}F(\omega)}{\omega}d\omega=\sum_{qj}\lambda_{qj}w(q)\end{equation}withw(q)astheweightofaqpoint.\par\par\par
Thesuperconductingtransitiontemperature,T_{C}canthenbeapproximatedusingMcMillan-
Allen-Dynes\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{mcmillan1968transition}{\@@citephrase{(}}{\@@citephrase{)}}}equationasfollows:\par\begin{equation}T_{c}=\frac{\omega_{log}}{1.2}\exp\bigg{[}-\frac{1.04(1+\lambda)}{\lambda-\mu^{*}(1+0.62\lambda)}\bigg{]}\end{equation}where\begin{equation}\omega_{log}=\exp\bigg{[}\frac{\int
d\omega\frac{\alpha^{2}F(\omega)}{\omega}\ln\omega}{\int
d\omega\frac{\alpha^{2}F(\omega)}{\omega}}\bigg{]}\end{equation}\par\par{\color[rgb]{0,0,0}WenotethatseveralvariantsofMcMillan-
Allen-Dynesformulahavebeenproposedpreviouslytobetteraccountforhigh-\lambda
materials\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{allen1975transition,allen1983theory,margine2013anisotropic,marsiglio1988iterative}{\@@citephrase{(}}{\@@citephrase{)}}},however,inthiswork,weoptforthesimpleroriginalMcMillan-
Allen-
Dynesformulaasourscreeningtool.Furthermore,wenotethatmaterialswithnarrowelectronicbandwidthmaynotbewell-
reproducedbytheformula.}InEq.~{}\ref{eq:mad},theparameter\mu^{*}istheeffectiveCoulombpotentialparameter.{\color[rgb]{0,0,0}Whiletherearewaystocalculatethisparameterfromfirstprinciples\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{lee1995first}{\@@citephrase{(}}{\@@citephrase{)}}},thisparametergenerallyvariesoverarelativelysmallrange(suchas0.09to0.18).Followingpreviousworks\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{lee1995first,giustino2017electron,marques2005ab,morel1962calculation,zarifi2018crystal,wierzbowska2005origins}{\@@citephrase{(}}{\@@citephrase{)}}},wetakeafixedvalue(here\mu^{*}=0.09)whenreportingourmainresults,andwepresentadiscussionontheeffectofvarying\mu^{*}inthelatersection.}\par\par\par\par\par\par\par$
(8)

###  3.3 Deep-learning training

Deep-learning is one of the fastest developing methods especially for
materials science applications choudhary2021recent . Graph neural networks
(GNN) are based on deep-learning framework because they can work with
unstructured and non-Euclidean data i.e. non-grid data such as atomic
structures. In this work, we use the recently developed atomistic line graph
neural network (ALIGNN) choudhary2021atomistic , which is publicly available
at https://github.com/usnistgov/alignn. ALIGNN has been used to train fast and
accurate models for more than 60 properties of solids and molecules with high
accuracy. We use ALIGNN to train models for the Debye temperature as well as
the electronic density of states at the Fermi-level using large datasets from
the JARVIS-DFT database. We also train models to predict the DFT-based
transition temperatures ($T_{C}$) and the underlying electron-phonon coupling
parameters $\omega_{log}$ and $\lambda$, using the smaller datasets computed
in this work. We use both ALIGNN and the hand-crafted (Classical force-field
descriptor (CFID) choudhary2018machine )-based approaches for comparison. CFID
approach utilizes the LightGBM package ke2017lightgbm .

In ALIGNN, a crystal structure is represented as a graph using atomic elements
as nodes and atomic bonds as edges. Each node in the atomistic graph is
assigned 9 input node features based on its atomic species: electronegativity,
group number, covalent radius, valence electrons, first ionization energy,
electron affinity, block and atomic volume. The inter-atomic bond distances
are used as edge features with radial basis function up to 8 Å cut-off. We use
a periodic 12-nearest-neighbor ($N$) graph construction. This atomistic graph
is then used for constructing the corresponding line graph using interatomic
bond-distances as nodes and bond-angles as edge features. ALIGNN uses edge-
gated graph convolution for updating nodes as well as edge features using a
propagation function ($f$) for layer ($l$), atom features ($h$), and node
($i$), details of which can be found in Ref. choudhary2021atomistic :

$h_{i}^{(l+1)}=f(h_{i}^{l}{\\{h_{j}^{l}}\\}_{{}_{i}})$ (9)

Unlike many other conventional GNNs (Xie2018CGCNN ; Schutt2018SchNet ;
Chen2019MEGNet ), ALIGNN uses bond-distances as well as bond-angles to
distinguish atomic structures. One ALIGNN layer composes an edge-gated graph
convolution on the bond graph with an edge-gated graph convolution on the line
graph. The line graph convolution produces bond messages that are propagated
to the atomistic graph, which further updates the bond features in combination
with atom features. The ALIGNN model is implemented in PyTorch
paszke2019pytorch  and deep graph library (DGL) wang2019deep . The
hyperparameters for ALIGNN are kept same as the original paper. We use the
available 17419 Debye temperature data and 55645 electronic DOS at Fermi level
data for model training. In addition, we use the 626 electron-phonon coupling
parameters for dynamically stable materials, their transition temperatures and
Eliashberg functions. For the Debye temperature and electronic density of
states at the Fermi-level, we use a batch size of 64 for 500 epochs and
80:10:10 training-validation-testing data split, while for predicting the EPC
parameters and Eliashberg functions, we use a batch size of 16, 90:5:5 split
and training for 300 epochs. The test set was never used during the training
procedure. We train use Tesla V100 SXM2 32 gigabyte Graphics processing unit
(GPU), with 8 Intel Xeon E5-2698 v4 CPU cores for concurrently fetching and
preprocessing batches of data during training. Please note commercial software
and hardware are identified to specify procedures. Such identification does
not imply recommendation by National Institute of Standards and Technology
(NIST).

##  4 Data Availability Statement

The dataset is available at the Figshare repository:
https://doi.org/10.6084/m9.figshare.21370572.

##  5 Code Availability Statement

Software packages mentioned in the article can be found at
https://github.com/usnistgov/jarvis and https://github.com/usnistgov/alignn.

##  6 Acknowledgments

K.C. and K.F.G. thank the National Institute of Standards and Technology for
funding, computational, and data-management resources. K.C. thanks the
computational support from XSEDE computational resources under allocation
number TG-DMR 190095. Contributions from K.C. were supported by the financial
assistance award 70NANB19H117 from the U.S. Department of Commerce, National
Institute of Standards and Technology.

##  7 Author Contributions

K.C. and K.F.G. jointly developed the workflow. K.C. carried out the high-
throughput DFT calculations and developed machine learning workflow. K.C. and
K.F.G. jointly analyzed the DFT data and contributed in writing the
manuscript.

##  8 Competing interests

The authors declare no competing interests.

## References

  *   * (1) Kamerlingh Onnes, H.  The resistance of pure mercury at helium temperatures.  _Commun. Phys. Lab. Univ. Leiden, b_ 120 (1911) .
  * (2) Poole, C. P., Farach, H. A. & Creswick, R. J.  _Superconductivity_ (Academic press, 2013).
  * (3) Rogalla, H. & Kes, P. H.  _100 years of superconductivity_ (Taylor & Francis, 2011).
  * (4) Nagamatsu, J., Nakagawa, N., Muranaka, T., Zenitani, Y. & Akimitsu, J.  Superconductivity at 39 k in magnesium diboride.  _Nature_ 410 (6824), 63–64 (2001) .
  * (5) Pickard, C. J., Errea, I. & Eremets, M. I.  Superconducting hydrides under pressure.  _Annu. Rev. Condens. Matter Phys_ 11, 57–76 (2020) .
  * (6) Sato, M. & Ando, Y.  Topological superconductors: a review.  _Rep. Prog. Phys._ 80 (7), 076501 (2017) .
  * (7) Roberts, B. W.  Survey of superconductive materials and critical evaluation of selected properties.  _J. Phys. Chem. Ref. Data_ 5 (3), 581–822 (1976) .
  * (8) Kihlstrom, K., Simon, R. & Wolf, S.  Tunneling $\alpha$2f ($\omega$) on high tc a15 and b1 compounds.  _Physica B+ C_ 135 (1-3), 198–202 (1985) .
  * (9) Stewart, G. R.  Superconductivity in the a15 structure.  _Phys. C: Supercond. Appl_ 514, 28–35 (2015) .
  * (10) Ivanovskii, A.  Band structure and properties of superconducting mgb2 and related compounds (a review).  _Phys. Solid State_ 45 (10), 1829–1859 (2003) .
  * (11) Buzea, C. & Yamashita, T.  Review of the superconducting properties of mgb2.  _Supercond. Sci. Technol._ 14 (11), R115 (2001) .
  * (12) Plakida, N.  _High-temperature cuprate superconductors: Experiment, theory, and applications_ Vol. 166 (Springer Science & Business Media, 2010).
  * (13) Hosono, H. & Kuroki, K.  Iron-based superconductors: Current status of materials and pairing mechanism.  _Phys. C: Supercond. Appl_ 514, 399–422 (2015) .
  * (14) Shipley, A. M., Hutcheon, M. J., Johnson, M. S., Needs, R. J. & Pickard, C. J.  Stability and superconductivity of lanthanum and yttrium decahydrides.  _Phys. Rev. B_ 101 (22), 224511 (2020) .
  * (15) Liu, H., Naumov, I. I., Hoffmann, R., Ashcroft, N. & Hemley, R. J.  Potential high-tc superconducting lanthanum and yttrium hydrides at high pressure.  _Proc. Natl. Acad. Sci._ 114 (27), 6990–6995 (2017) .
  * (16) Drozdov, A., Eremets, M., Troyan, I., Ksenofontov, V. & Shylin, S. I.  Conventional superconductivity at 203 kelvin at high pressures in the sulfur hydride system.  _Nature_ 525 (7567), 73–76 (2015) .
  * (17) Shipley, A. M., Hutcheon, M. J., Needs, R. J. & Pickard, C. J.  High-throughput discovery of high-temperature conventional superconductors.  _Phys. Rev. B_ 104 (5), 054501 (2021) .
  * (18) Yuan, J., Stanev, V., Gao, C., Takeuchi, I. & Jin, K.  Recent advances in high-throughput superconductivity research.  _Supercond. Sci. Technol._ 32 (12), 123001 (2019) .
  * (19) Rodriguez, C. _et al._ Optical near-zone-center phonons and their interaction with electrons in yba 2 cu 3 o 7: results of the local-density approximation.  _Phys. Rev. B_ 42 (4), 2692 (1990) .
  * (20) Subedi, A., Ortenzi, L. & Boeri, L.  Electron-phonon superconductivity in a pt 3 p (a= sr, ca, la) compounds: From weak to strong coupling.  _Phys. Rev. B_ 87 (14), 144504 (2013) .
  * (21) Duan, D., Yu, H., Xie, H. & Cui, T.  Ab initio approach and its impact on superconductivity.  _J. Supercond. Nov. Magn._ 32 (1), 53–60 (2019) .
  * (22) Kolmogorov, A. _et al._ New superconducting and semiconducting fe-b compounds predicted with an ab initio evolutionary search.  _Phys. Rev. Lett._ 105 (21), 217003 (2010) .
  * (23) Gao, G. _et al._ High-pressure crystal structures and superconductivity of stannane (snh4).  _Proc. Natl. Acad. Sci._ 107 (4), 1317–1320 (2010) .
  * (24) Duan, D. _et al._ Pressure-induced metallization of dense (h2s) 2h2 with high-tc superconductivity.  _Sci. Rep._ 4 (1), 1–6 (2014) .
  * (25) Hutcheon, M. J., Shipley, A. M. & Needs, R. J.  Predicting novel superconducting hydrides using machine learning approaches.  _Phys. Rev. B_ 101 (14), 144505 (2020) .
  * (26) Stanev, V. _et al._ Machine learning modeling of superconducting critical temperature.  _NPJ Comput. Mater._ 4 (1), 1–14 (2018) .
  * (27) Xie, S. _et al._ Machine learning of superconducting critical temperature from eliashberg theory.  _NPJ Comput. Mater._ 8 (1), 1–8 (2022) .
  * (28) Cooper, L. N. & Feldman, D.  _BCS: 50 years_ (World scientific, 2010).
  * (29) Giustino, F.  Electron-phonon interactions from first principles.  _Rev. Mod. Phys._ 89 (1), 015003 (2017) .
  * (30) Kawamura, M., Hizume, Y. & Ozaki, T.  Benchmark of density functional theory for superconductors in elemental materials.  _Phys. Rev. B_ 101, 134511 (2020) .
  * (31) Choudhary, K. _et al._ The joint automated repository for various integrated simulations (jarvis) for data-driven materials design.  _NPJ Comput. Mater._ 6 (1), 1–13 (2020) .
  * (32) McMillan, W.  Transition temperature of strong-coupled superconductors.  _Phys. Rev._ 167 (2), 331 (1968) .
  * (33) Choudhary, K. & DeCost, B.  Atomistic line graph neural network for improved materials property predictions.  _NPJ Comput. Mater._ 7 (1), 1–8 (2021) .
  * (34) Choudhary, K. _et al._ Recent advances and applications of deep learning methods in materials science.  _NPJ Comput. Mater._ 8 (1), 1–26 (2022) .
  * (35) Choudhary, K. _et al._ Computational screening of high-performance optoelectronic materials using optb88vdw and tb-mbj formalisms.  _Sci. Data_ 5 (1), 1–12 (2018) .
  * (36) Choudhary, K. _et al._ Accelerated discovery of efficient solar cell materials using quantum and machine-learning methods.  _Chem. Mater._ 31 (15), 5900–5908 (2019) .
  * (37) Choudhary, K., Garrity, K. F., Ghimire, N. J., Anand, N. & Tavazza, F.  High-throughput search for magnetic topological materials using spin-orbit spillage, machine learning, and experiments.  _Phys. Rev. B_ 103 (15), 155131 (2021) .
  * (38) Choudhary, K., Garrity, K. F. & Tavazza, F.  High-throughput discovery of topologically non-trivial materials using spin-orbit spillage.  _Sci. Rep._ 9 (1), 1–8 (2019) .
  * (39) Choudhary, K., Garrity, K. F., Jiang, J., Pachter, R. & Tavazza, F.  Computational search for magnetic and non-magnetic 2d topological materials using unified spin–orbit spillage screening.  _NPJ Comput. Mater._ 6 (1), 1–8 (2020) .
  * (40) Choudhary, K., Cheon, G., Reed, E. & Tavazza, F.  Elastic properties of bulk and low-dimensional materials using van der waals density functional.  _Phys. Rev. B_ 98 (1), 014107 (2018) .
  * (41) Choudhary, K. _et al._ High-throughput density functional perturbation theory and machine learning predictions of infrared, piezoelectric, and dielectric responses.  _NPJ Comput. Mater._ 6 (1), 1–13 (2020) .
  * (42) Choudhary, K., Ansari, J. N., Mazin, I. I. & Sauer, K. L.  Density functional theory-based electric field gradient database.  _Sci. Data_ 7 (1), 1–10 (2020) .
  * (43) Choudhary, K., Kalish, I., Beams, R. & Tavazza, F.  High-throughput identification and characterization of two-dimensional materials using density functional theory.  _Sci. Rep._ 7 (1), 1–16 (2017) .
  * (44) Choudhary, K. & Tavazza, F.  Convergence and machine learning predictions of monkhorst-pack k-points and plane-wave cut-off in high-throughput dft calculations.  _Comput. Mater. Sci._ 161, 300–308 (2019) .
  * (45) Anderson, O. L.  A simplified method for calculating the debye temperature from elastic constants.  _JJ Phys Chem Solids_ 24 (7), 909–917 (1963) .
  * (46) Garrity, K. F.  First-principles search for n-type oxide, nitride, and sulfide thermoelectrics.  _Phys. Rev. B_ 94 (4), 045122 (2016) .
  * (47) Wierzbowska, M., de Gironcoli, S. & Giannozzi, P.  Origins of low-and high-pressure discontinuities of $t\\_${$c$} in niobium.  _Preprint at https://arxiv.org/abs/cond-mat/0504077v2_ (2005) .
  * (48) Perdew, J. P., Burke, K. & Ernzerhof, M.  Generalized gradient approximation made simple.  _Phys. Rev. Lett._ 77 (18), 3865 (1996) .
  * (49) Lee, K.-H., Chang, K.-J. & Cohen, M. L.  First-principles calculations of the coulomb pseudopotential $\mu$*: Application to al.  _Phys. Rev. B_ 52 (3), 1425 (1995) .
  * (50) Marques, M. _et al._ Ab initio theory of superconductivity. ii. application to elemental metals.  _Phys. Rev. B_ 72 (2), 024546 (2005) .
  * (51) Morel, P. & Anderson, P.  Calculation of the superconducting state parameters with retarded electron-phonon interaction.  _Phys. Rev._ 125 (4), 1263 (1962) .
  * (52) Zarifi, N., Bi, T., Liu, H. & Zurek, E.  Crystal structures and properties of iron hydrides at high pressure.  _J. Phys. Chem. C._ 122 (42), 24262–24269 (2018) .
  * (53) Sanna, A., Pellegrini, C. & Gross, E.  Combining eliashberg theory with density functional theory for the accurate prediction of superconducting transition temperatures and gap functions.  _Phys. Rev. Lett._ 125 (5), 057001 (2020) .
  * (54) Lüders, M. _et al._ Ab initio theory of superconductivity. i. density functional formalism and approximate functionals.  _Phys. Rev. B_ 72 (2), 024545 (2005) .
  * (55) Allen, P. B. & Dynes, R.  Transition temperature of strong-coupled superconductors reanalyzed.  _Phys. Rev. B_ 12 (3), 905 (1975) .
  * (56) Webb, G., Marsiglio, F. & Hirsch, J.  Superconductivity in the elements, alloys and simple compounds.  _Phys. C: Supercond. Appl_ 514, 17–27 (2015) .
  * (57) Heid, R., Bohnen, K.-P., Sklyadneva, I. Y. & Chulkov, E.  Effect of spin-orbit coupling on the electron-phonon interaction of the superconductors pb and tl.  _Phys. Rev. B_ 81 (17), 174527 (2010) .
  * (58) Gibson, G. & Meservey, R.  Evidence for spin fluctuations in vanadium from a tunneling study of fermi-liquid effects.  _Phys. Rev. B_ 40 (13), 8705 (1989) .
  * (59) Mandal, S., Cohen, R. E. & Haule, K.  Strong pressure-dependent electron-phonon coupling in fese.  _Phys. Rev. B_ 89 (22), 220502 (2014) .
  * (60) Gerber, S. _et al._ Femtosecond electron-phonon lock-in by photoemission and x-ray free-electron laser.  _Science_ 357 (6346), 71–75 (2017) .
  * (61) Drzazga, E., Domagalska, I., Jarosik, M., Szczeniak, R. & Kalaga, J.  Characteristics of superconducting state in vanadium: The eliashberg equations and semi-analytical formulas.  _J. Supercond. Nov. Magn._ 31 (4), 1029–1034 (2018) .
  * (62) Dew-Hughes, D. & Linse, V.  Formation of superconducting nb3si by explosive compression.  _J. Appl. Phys._ 50 (5), 3500–3504 (1979) .
  * (63) Tsindlekht, M. _et al._ Linear and nonlinear low-frequency electrodynamics of surface superconducting states in an yttrium hexaboride single crystal.  _Phys. Rev. B_ 78 (2), 024522 (2008) .
  * (64) Matthias, B. T., Geballe, T. H. & Compton, V. B.  Superconductivity.  _Rev. Mod. Phys._ 35, 1–22 (1963) .
  * (65) Hulm, J., Jones, C., Hein, R. & Gibson, J.  Superconductivity in the tio and nbo systems.  _J. Low Temp. Phys._ 7 (3), 291–307 (1972) .
  * (66) Pechen, E. _et al._ Tunneling and critical-magnetic-field study of superconducting nbc thin films.  _Physica C: Superconductivity_ 235, 2511–2512 (1994) .
  * (67) Kruglov, I. A. _et al._ Superconductivity of lah 10 and lah 16 polyhydrides.  _Phys. Rev. B_ 101 (2), 024508 (2020) .
  * (68) Belsky, A., Hellenbrandt, M., Karen, V. L. & Luksch, P.  New developments in the inorganic crystal structure database (icsd): accessibility in support of materials research and design.  _Acta Crystallogr.B Struct. Sci._ 58 (3), 364–369 (2002) .
  * (69) Yamase, H., Sakurai, Y., Fujita, M., Wakimoto, S. & Yamada, K.  Fermi surface in la-based cuprate superconductors from compton scattering imaging.  _Nat. Commun._ 12 (1), 1–7 (2021) .
  * (70) An, Y. _et al._ Superconductivity and topological properties of mgb 2-type diborides from first principles.  _Phys. Rev. B_ 104 (13), 134510 (2021) .
  * (71) Inumaru, K., Nishikawa, T., Nakamura, K. & Yamanaka, S.  High-pressure synthesis of superconducting molybdenum nitride $\delta$-mon by in situ nitridation.  _Chem. Mater._ 20 (14), 4756–4761 (2008) .
  * (72) Wang, S. _et al._ The hardest superconducting metal nitride.  _Sci. Rep._ 5 (1), 1–8 (2015) .
  * (73) Choi, H. J., Louie, S. G. & Cohen, M. L.  Prediction of superconducting properties of cab 2 using anisotropic eliashberg theory.  _Phys. Rev. B_ 80 (6), 064503 (2009) .
  * (74) Zbasnik, J. & Toth, L.  Electronic structure of vanadium carbide.  _Phys. Rev. B_ 8 (2), 452 (1973) .
  * (75) Vaitheeswaran, G., Kanchana, V. & Rajagopalan, M.  Structural phase stability and superconductivity of lan.  _Solid State Commun._ 124 (3), 97–102 (2002) .
  * (76) Chu, C., Bucher, E., Cooper, A. & Maita, J.  Electronic transitions in near-equiatomic vanadium-ruthenium alloys.  _Phys. Rev. B_ 4 (2), 320 (1971) .
  * (77) Katsura, Y. _et al._ On the possibility of mgb2-like superconductivity in potassium hexaboride.  _Phys. C: Supercond. Appl_ 470, S633–S634 (2010) .
  * (78) Yan, D. Y. _et al._ Superconductivity in centrosymmetric topological superconductor candidate tac.  _Supercond. Sci. Technol._ 34 (3), 035025 (2021) .
  * (79) Li, Q. & Chen, L.  Superconducting atmospheric structure and pressure-induced novel phases of cobalt mononitride.  _Comput. Mater. Sci._ 174, 109464 (2020) .
  * (80) Hong, D., Zeng, W., Liu, F.-S., Tang, B. & Liu, Q.-J.  Structural, electronic, elastic and mechanical properties of nbc-based compounds: First-principles calculations.  _Physica B Condens. Matter._ 558, 100–108 (2019) .
  * (81) Sernetz, F., Lerf, A. & Schöllhorn, R.  Superconductivity in hydrated ternary tantalum sulfides.  _Mater. Res. Bull._ 9 (12), 1597–1602 (1974) .
  * (82) Ren, Z.-A., Kuroiwa, S., Tomita, Y. & Akimitsu, J.  Structural phase evolution and superconductivity in the non-stoichiometric intermetallic compound niobium diboride.  _Physica C: Superconductivity_ 468 (5), 411–416 (2008) .
  * (83) Kodess, B. N.  Superconductivity of some new w compounds.  _Phys. Status Solidi (a)_ 4 (2), K109–K112 (1971) .
  * (84) Gokhale, A. & Abbaschian, G.  The sc- si (scandium-silicon) system.  _Bull. alloy phase diagr._ 7 (4), 333–336 (1986) .
  * (85) Yan, L. _et al._ Theoretical dissection of superconductivity in two-dimensional honeycomb borophene oxide b2o crystal with a high stability.  _NPJ Comput. Mater._ 6 (1), 1–7 (2020) .
  * (86) Tutuncu, H. M. & Srivastava, G.  Theoretical examination of superconductivity in the cubic antiperovskite cr3gan under pressure.  _J. Appl. Phys._ 114 (5), 053905 (2013) .
  * (87) Ramakrishnan, S., Nigam, A. & Chandra, G.  Resistivity and magnetoresistance studies on superconducting a15 v 3 ga, v 3 au, and v 3 pt compounds.  _Phys. Rev. B_ 34 (9), 6166 (1986) .
  * (88) Chen, Y., Cai, X., Wang, H., Wang, H. & Wang, H.  Novel triadius-like n4 specie of iron nitride compounds under high pressure.  _Sci. Rep._ 8 (1), 1–7 (2018) .
  * (89) Gražulis, S. _et al._ Crystallography open database–an open-access collection of crystal structures.  _J. Appl. Crystallogr._ 42 (4), 726–729 (2009) .
  * (90) Choudhary, K., DeCost, B. & Tavazza, F.  Machine learning with force-field-inspired descriptors for materials: Fast screening and mapping energy landscape.  _Phys. Rev. Mater._ 2 (8), 083801 (2018) .
  * (91) Gurunathan, R., Choudhary, K. & Tavazza, F.  Rapid prediction of phonon structure and properties using an atomistic line graph neural network (alignn).  _Preprint at https://arxiv.org/abs/2207.12510_ (2022) .
  * (92) Bardeen, J., Cooper, L. N. & Schrieffer, J. R.  Theory of superconductivity.  _Phys. Rev._ 108 (5), 1175 (1957) .
  * (93) Kresse, G. & Furthmüller, J.  Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set.  _Phys. Rev. B_ 54 (16), 11169 (1996) .
  * (94) Kresse, G. & Furthmüller, J.  Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set.  _Comput. Mater. Sci._ 6 (1), 15–50 (1996) .
  * (95) Klimeš, J., Bowler, D. R. & Michaelides, A.  Chemical accuracy for the van der waals density functional.  _J. Phys. Condens. Matter_ 22 (2), 022201 (2009) .
  * (96) Giannozzi, P. _et al._ Quantum espresso: a modular and open-source software project for quantum simulations of materials.  _J. Phys. Condens. Matter_ 21 (39), 395502 (2009) .
  * (97) Giannozzi, P. _et al._ Quantum espresso toward the exascale.  _J. Chem. Phys._ 152 (15), 154105 (2020) .
  * (98) Kawamura, M., Gohda, Y. & Tsuneyuki, S.  Improved tetrahedron method for the brillouin-zone integration applicable to response functions.  _Phys. Rev. B_ 89 (9), 094515 (2014) .
  * (99) Poncé, S., Margine, E. R., Verdi, C. & Giustino, F.  Epw: Electron–phonon coupling, transport and superconducting properties using maximally localized wannier functions.  _Comput. Phys. Commun._ 209, 116–133 (2016) .
  * (100) Baroni, S., Giannozzi, P. & Testa, A.  Green’s-function approach to linear response in solids.  _Phys. Rev. Lett._ 58 (18), 1861 (1987) .
  * (101) Gonze, X.  Perturbation expansion of variational principles at arbitrary order.  _Phys. Rev. A_ 52 (2), 1086 (1995) .
  * (102) Garrity, K. F., Bennett, J. W., Rabe, K. M. & Vanderbilt, D.  Pseudopotentials for high-throughput dft calculations.  _Comput. Mater. Sci._ 81, 446–452 (2014) .
  * (103) Perdew, J. P. _et al._ Restoring the density-gradient expansion for exchange in solids and surfaces.  _Phys. Rev. Lett._ 100 (13), 136406 (2008) .
  * (104) Allen, P. B. & Mitrović, B.  Theory of superconducting tc.  _Solid State Phys._ 37, 1–92 (1983) .
  * (105) Margine, E. R. & Giustino, F.  Anisotropic migdal-eliashberg theory using wannier functions.  _Phys. Rev. B_ 87 (2), 024505 (2013) .
  * (106) Marsiglio, F., Schossmann, M. & Carbotte, J.  Iterative analytic continuation of the electron self-energy to the real axis.  _Phys. Rev. B_ 37 (10), 4965 (1988) .
  * (107) Ke, G. _et al._ Lightgbm: A highly efficient gradient boosting decision tree.  _Adv Neural Inf Process Syst_ 30 (2017) .
  * (108) Xie, T. & Grossman, J. C.  Crystal graph convolutional neural networks for an accurate and interpretable prediction of material properties.  _Phys. Rev. Lett._ 120, 145301 (2018) .
  * (109) Schütt, K. T., Sauceda, H. E., Kindermans, P.-J., Tkatchenko, A. & Müller, K.-R.  Schnet – a deep learning architecture for molecules and materials.  _J. Chem. Phys._ 148 (24), 241722 (2018) .
  * (110) Chen, C., Ye, W., Zuo, Y., Zheng, C. & Ong, S. P.  Graph networks as a universal machine learning framework for molecules and crystals.  _Chem. Mater._ 31 (9), 3564–3572 (2019) .
  * (111) Paszke, A. _et al._ Pytorch: An imperative style, high-performance deep learning library.  _Adv Neural Inf Process Syst_ 32, 8026–8037 (2019) .
  * (112) Wang, M. _et al._ Deep graph library: Towards efficient and scalable deep learning on graphs.  _Preprint at https://arxiv.org/abs/1909.01315v2_ (2019) .

##  9 Figure captions

Fig. 1 Schematic showing the steps involved in identifying high-$T_{C}$
superconductors. a) statistical distribution of Debye temperature (K) and b)
statistical distribution of electronic density of states (states/eV/total
number of electrons) at Fermi level from the JARVIS-DFT database, c)
probability that compounds containing a given element have $\theta_{D}>$300 K.
The flow chart shows the application of BCS-inspired screening, density
functional theory calculations and deep-learning training.

Fig. 2 Effect of different k-point and q-point selection on EPC parameters for
VC. a) $\lambda$, b) $\omega_{log}$, c) $T_{C}$. Broadening dependent scaled
d) $\lambda$, e) $\omega_{log}$ and f) $T_{C}$ for the all the systems in the
database. The scaling was done with respect to the corresponding values at
broadening of 0.05 Ry. The legend in panel a corresponds to the number of
k-points and q-points, e.g. 10x10x10 k-points and 10x10x10 q-points (denoted
by black star) is same as 1 k-point per q-point; and 20x20x20 k-point and
2x2x2 q-point correspond to 1000 k-point per q-point (denoted by left purple
triangles).

Fig. 3 Broadening convergence with respect to experimental $T_{C}$ data for 14
superconductors.

Fig. 4 Effect of $\mu^{*}$ on transition temperatures. a) Change in $T_{C}$ as
we increase $\mu^{*}$ for a few representative superconductors, b) change in
$T_{C}$ with respect to $\mu^{*}$ 0.09 for the entire dataset, c) change in
percenatge materials with $T_{c}$ greater equal to 5 K. Here $\mu^{*}$ 0.09 is
shown as 100 %. As we decrease $\mu^{*}$, we add materials and vice-versa.

Fig. 5 Relation between Electron-phonon coupling parameters and EPC function
of some of the potential candidate superconductors. a) $\omega_{log}$ vs
$\lambda$, b) $T_{C}$ vs $\lambda$, c) MoN, d) VC, e) KB6, and f) VTe.

Fig. 6 Atomistic line graph neural network based deep-learning (DL) regression
model performance on 5 % test set for a) Debye temperature and b) DOS.
Classical force-field descriptor (CFID) (c,d,e) and DL (f,g,h) based
regression model performance on 5 % test set for DFT calculated transition
temperature ($T_{C}$), EPC parameter $\omega_{log}$, and EPC parameter
$\lambda$. In (f), we show performance with direct $T_{C}$ prediction (red
color), $T_{C}$ prediction with direct prediction of wlog and lambda and then
using eq. 7 (green color) and $T_{C}$ prediction with Eliashberg function and
then using eq. 5-8 (black color).

Fig. 7 Prediction of Eliashberg function with ALIGNN for the 5 % test set.
ALIGNN can capture peak-positions and heights reasonably well.