Introduction

Rapidly growing population and infrastructure development are behind the rising energy demand, which will further increase from 23 Terawatts in 2030 to 30 Terawatts in 20501. According to the Global Renewable Energy Community (REN21) statistics, nearly 80% of the total energy relies on conventional energy resources and renewable energy sources add only the rest 20%2. An excessive dependence on fossil fuels causes global warming and destructive environmental issues3. There has been a worldwide push to find sustainable and clean alternatives to fossil fuels to counter such problems4. Among the natural renewable energy sources, hydrogen is an ideal sustainable energy source due to its high energy density and environmentally benign5. However, precious and less abundant metal-based catalysts have been used for hydrogen production6, impeding their widespread utilization7. Therefore, exploring novel and metal-free catalysts is a viable route for the mass production of hydrogen8,9. Thermoelectric generators are also excellent alternatives for clean and renewable energy resources, considering the abundance of waste heat accompanied by infrequent maintenance and long device life, as no moving parts are involved in the technology10,11. Although Bi2Te3 has been widely used in thermoelectric generators, the toxicity and scarcity of tellurium restrict their usage12. In addition, bipolar conduction suppresses the figure of merit of Bi2Te3 above 450 K due to its narrow bandgap12. Hence abundant and non-toxic materials with reasonable bandgap would be an optimal choice.

Since the experimental realization of graphene13, tremendous attention has been devoted to other two-dimensional (2D) carbon allotropes due to their peculiar physical properties14, topological states15, massless Dirac cones16,17, and semiconducting behavior18,19. Porous nitrogenated holey graphene20, polyaniline21, phagraphene22, naphyne23, graphtetrayne24, and biphenylene25 are a few examples of experimentally synthesized 2D carbon allotropes.

Among these, graphyne, with varying sp and sp2 hybridized carbon atoms constitutes one of the biggest families of graphene allotropes26. They possess exceptional flexibility, high carrier mobility, a Dirac cone characterized electronic band structure, efficient adsorption of ions and molecular selectivity due to porous structures, and reduced thermal conductivity due to acetylenic bonds with sp state27,28,29,30,31,32. Recently, the bottom-up technique has been employed to synthesize an ultrathin 2D carbon allotrope named holey graphyne33. The nanosheet shows excellent mechanical, thermal, and dynamic stability. Unlike graphene, it is a direct bandgap semiconductor with high carrier mobility (promising for applications in optoelectronics) and possesses sp and sp2 hybridized carbon atoms uniformly distributing the porous architectures (favorable for gas separation, water desalination, energy storage, and catalysis)34. Holey graphyne may also be considered an anchoring material in metal-sulfur batteries, like other materials with similar crystal structures previously researched for this purpose35,36. However, to the best of our knowledge, none of these applications of holey graphyne has been uncovered so far. Herein, we conducted a comprehensive study to explore its potential in thermoelectricity and H2 production. By the presence of unique distinct bonding, κl is anticipated to be lower in holey graphyne compared to the other flat 2D materials from the graphene family, resulting in an enhanced thermoelectric figure of merit. On the other hand, variation in bond charge density accompanied by a highly porous plane, which increases the number of reactive sites, makes it an excellent choice for catalysis.

Computational details

The Vienna Ab-initio Simulation Package (VASP)37,38 is used to perform the density functional theory (DFT) calculations by employing the Perdew − Burke − Ernzerhof and Heyd − Scuseria − Ernzerhof (HSE06) hybrid exchange–correlation functionals39,40. A plane wave cutoff energy of 550 eV and a Γ-centred k-mesh of 9 × 9 × 1 (30 × 30 × 1) is used to sample the first Brillouin zone for self-consistent (non-self-consistent) calculations. The crystal is optimized until Hellmann–Feynman forces drop below 10–4 eV/Å. A vacuum of 15 Å, perpendicular to the sheet, is taken to eliminate the interactions between adjacent layers. The structural formation energy is calculated using Eform. = (Etot. − n.EC)/n formula, where Etot. and EC are the total energies of holey graphyne and an isolated C atom, respectively, and n is the total number of atoms in a unit cell. The phonon dispersion, thermal transport coefficients, and lattice thermal conductivity are calculated using the Hiphive41, BoltzTraP242, and the ShengBTE43 codes, respectively. This sophisticated method has been used to analyze the transport characteristics of numerous materials44,45,46. As inputs to ShengBTE, the 2nd and 3rd-order force constants are computed using a 3 × 3 × 1 supercell. The 3rd-order force constants are calculated by displacing atoms up to the tenth nearest neighbors. A dense q-mesh 30 × 30 × 1 is utilized to obtain well-converged lattice thermal conductivity. The produced lattice thermal conductivity is well converged, as the value at 300 K varies by less than 4% and 1% from the result obtained with displacing atoms up to the ninth nearest neighbors and 20 × 20 × 1 q-mesh, respectively. The change in Gibbs free energy change (ΔGH) is calculated using the relation, ΔGH = ΔEH + ΔEZPE − TΔS, where ΔEH, ΔEzpe, T, and ΔS stand for the hydrogen adsorption energy, change in the zero-point energy, temperature (298.15 K), and change in entropy, respectively47. The value of ΔEZPE – TΔS is equal to 0.24 eV, giving rise to ΔGH = ΔEH + 0.2447. The adsorption energy of hydrogen is obtained as ΔEH = Etotal − Epristine − 1/2EH2, where Etotal, Epristine, and EH2 are the total energies of holey graphyne with adsorbed H atom, pristine holey graphyne, and single H2-molecule in the gas phase, respectively.

Results and discussions

Holey graphyne contains two rings comprised of six and eight carbon atoms, with a big pore, as shown in Fig. 1. Each primitive cell has an optimized in-plane lattice constant of 10.84 Å, in agreement with the previous study33. The C–C bonds have different bond lengths of 1.23 Å (d1; sp hybridized triple bond), 1.41 Å (d2; sp hybridized single bond), 1.40 Å (d3; sp2 hybridized double bond), and 1.46 Å (d4; sp2 hybridized single bond), as shown in Fig. 1b. The formation energy of holey graphyne is − 8.46 eV/atom, comparable to graphene (− 9.22 eV/atom) and h-BN (− 8.80 eV/atom) monolayers48. It is worthwhile to point out that holey graphyne has been experimentally synthesized using the bottom-up approach33. Besides that, ab-initio molecular dynamics simulations have determined that holey graphyne is thermally stable even at higher temperatures of 900 K.33 Fig. 2 shows the calculated direct bandgaps of 0.50 eV (PBE) and 1.0 eV (HSE06), which agrees with the previously calculated value and is also close to the experimental value of 1.10 eV.33.

Figure 1
figure 1

The side (a) and top (b) views of the optimized crystal structure of holey graphyne unit cell. di (i = 1–4) represents the bond length. Dashed lines enclose the unit cell.

Full size image
Figure 2
figure 2

Calculated electronic band structures of holey graphyne. The orange and red curves correspond to valence and conduction bands, respectively. The blue dashed line represents the Fermi level.

Full size image

The calculated phonon dispersion of holey graphyne has a total of 72 phonon modes without imaginary frequencies, assuring its dynamic stability (see Fig. 3). The out-of-plane flexural acoustic (ZA) phonons have the lowest frequency among the acoustic modes, followed by the in-plane transverse acoustic (TA) and longitudinal acoustic (LA) phonon modes. The first optical phonon mode has a noticeably small frequency of 2.55 THz, favoring the high phonon–phonon scattering owing to the coupling between ZA and the lowest optical phonons favors high phonon–phonon scattering49,50. Also, the less dispersive nature of optical phonons results in their small phonon group velocities. These distinct characteristics collectively contribute to low κl, compared to other carbon allotropes.

Figure 3
figure 3

Calculated phonon dispersion of holey graphyne.

Full size image

Figure 4 reveals the electronic transport coefficients with respect to varying carrier concentrations and temperatures. The electronic part of electrical (σ/τ) and thermal (κe/τ) conductivities increases and the Seebeck coefficient (S) declines linearly with carrier concentration51. The room temperature |S| (300 µV/K at 1 × 1012 cm−2) is approximately 6 times greater than that of graphene (~ 50 µV/K) at the same carrier concentration52. The p-type (n-type) |S| turns out to be 114 µV/K (110 µV/K) even at the highest considered carrier concentration (2 × 1013 cm−2) and at 600 K. The p-type S2σ/τ rises with growing carrier concentration until it approaches 5 × 1011 W/mK2s (9 × 1011 W/mK2s) at 300 K (600 K) and then drops as the doping concentration elevates. Such a remarkably high S and S2σ/τ suggest that investigating holey graphyne for thermoelectrics is worthwhile. It is further notable that identical dispersions of the valence band maxima and the conduction band minima result in similar variance in p-type and n-type electronic transport coefficients.

Figure 4
figure 4

Calculated p-type and n-type electronic transport coefficients as a function of carrier concentration (ρ) at 300 K and 600 K.

Full size image

Further, we have looked into how lattice vibrations affect heat conduction in holey graphyne (see Fig. 5a). The calculated room temperature κl turns out to be 29.3 W/mK, lower than graphene (3000 W/mK)53, γ-graphyne (76.4 W/mK)54, C3N (128 W/mK)53, and C2N (82.22 W/mK)55. The calculated κl at 300 K differs by 0.01% from the results obtained using 20 × 20 × 1 q-mesh (i.e., 29.4 W/mK), leaving no relevant effect on our conclusion. The atom displacement to the ninth nearest neighbors yields κl of 28.05 W/mK, which deviates by 4% from the value obtained by displacing atoms up to the tenth nearest neighbors. The ultralow κl of holey graphyne, which is appealing from the thermoelectric perspective, is attributed to its bonds being less stiff than those of the above-mentioned flat materials. The C–C bond stiffness is determined by the spring constant (K) of holey graphyne, which is calculated as the trace of the harmonic force constant tensor between the closest adjacent atoms. It is written as \(K={\Phi }_{CC}^{xx}+ {\Phi }_{CC}^{yy}+{\Phi }_{CC}^{zz}\), where \({\Phi }_{MX}^{\alpha \alpha }\) is the second derivate of energy with respect to displacement of atoms along the Cartesian axis α. The bond stiffness of d1, d2, d3, and d4 are 87 eV/Å2, 36 eV/Å2, 44 eV/Å2, and 31 eV/Å2, which are substantially smaller than that of graphene (10,105 eV/Å2). Such a bond feature inhibits heat transfer via lattice vibrations, leading to low κl values. Figure 5a shows that as temperature increases κl decreases due to the pronounced phonon–phonon scattering and follows the relation κl ∝ 1/T.

Figure 5
figure 5

(a) Calculated and fitted (κl ∝ 1/T) κl as a function of temperature, (b) cumulative lattice thermal conductivity (κcl) as a function of phonon mean free path (MFP), (c) phonon scattering rates, and (d) phonon group velocities at 300 K as a function of phonon frequency. The black arrows in (c and d) represent frequency 2.55 THz of first optical phonon mode at Γ point.

Full size image

Nanostructuring of materials, where materials are composed of nanometer size grains and further with nanoscale internal structures56, can reduce κl without affecting σ. To better comprehend the scope of phonon engineering, the cumulative lattice thermal conductivity (κcl) as a function of the phonon mean free path (MFP) is investigated and presented in Fig. 5b. As MFP decreases, the phonons scattering increases reducing the heat transfer. The contribution of phonons with different MFPs to lattice thermal conductivity is studied by calculating κcl, thereby deducing phonons most relevant to thermal conduction43. At 300 K (600 K) 75% of κl is by phonons having MFP 335 nm (110 nm), which implies reducing κl through nanostructuring is a viable strategy for holey graphyne. In other words, a sample size of 335 nm (110 nm) could help reduce the inherited value of κl at 300 K (600 K) by one-fourth. κcl increases (decreases) with rising MFP (temperature) and shows a plateau above 20,092 nm (7924 nm) at 300 K (600 K). A lower plateau at 600 K is due to the stronger phonon scattering at elevated temperatures10. Given that phonon transports are predominantly dependent on phonon scattering rates and phonon group velocities, our calculated results for acoustic and optical phonon modes are illustrated in Fig. 5c,d. The highest scattering rate for acoustic phonons at room temperature is 2.74 ps−1, comparatively higher than C3N monolayer (2 ps−1)53. It is due to the coupling of acoustic and optical phonon modes, which leads to increased scattering rates of acoustic phonons, leading to a considerable reduction of κl50. The phonon group velocities are calculated to validate the above analysis and presented in Fig. 5d. The highest phonon group velocity of acoustic phonon mode at room temperature is 13.9 km/s, much lower than that of graphene (~ 22 km/s)57 and nitrogenated holey graphene (18.48 km/s)58. Such low group velocity is a consequence of flat phonon modes (see Fig. 3)50. Thus, smaller κl of holey graphyne arises from low phonon group velocity and high scattering rates.

The figure of merit (ZT) depends on the relaxation time (τ) and is vital for showcasing the potential of a material for its thermoelectric applications. In this work, we have adapted the value for τ from the deformation potential theory used in the previous study33. The relation τT = 300*τ300/T is engaged to yield its value at 600 K. The values of τ for holes (electrons) turn out to be 3.27 ps (1.16 ps) and 1.64 ps (0.58 ps) at 300 K and 600 K, respectively. The calculated ZT in Fig. 6 follows an upward trend similar to S2σ/τ. The p-type ZT is higher than the n-type and is counter to the trend seen in electronic transport coefficients (see Fig. 4). It is a result of the fact that holes and electrons have different relaxation times. The relation ZT = S2T/(κe/σ + κl/σ), where κl/σ is influenced by τ, explains that larger τ values result in greater ZT. The ZT secures peak values of 3.36 (1.50) and 1.82 (0.71) at 600 K (300 K) for p-type and n-type dopings, respectively. In the previous study33, ab-initio molecular dynamics simulations determined that holey graphyne is thermally stable even at higher temperatures of 900 K, as established based on the ab initio molecular dynamics simulations. The room temperature p-type ZT of holey graphyne is higher than that of holey graphene (1.13)59 and significantly higher than pristine graphene (0.55 × 10–3)60 and γ-graphyne (0.48)54.

Figure 6
figure 6

Calculated figure of merit (ZT) of holey graphyne.

Full size image

Since active and cost-effective metal-free catalysts, as opposed to costly Pt-based catalysts, for HER are necessary to attain sustainable energy systems. We have also uncovered the response of holey graphyne to HER. As the standard descriptor for HER activity is the Gibbs free energy (ΔGH), which is derived from the hydrogen adsorption energy. In the first instance, the H atom is adsorbed on all the possible 9 adsorption sites (s1–s9 in Fig. 7a). The H atom finds the s2 site the most favorable. The corresponding ΔGH of s2 site is 0.20 eV, which is much smaller than that of graphene (1.41 eV)61, phosphorene (1.25 eV)62, and C3N4 (0.58 eV)63 and comparable to biphenylene (0.29 eV)61, see Fig. 7b. As a matter of fact, ΔGH = 0 stands out as an optimal value for HER. However, a value of |ΔGH|< 0.2 eV signifies the better catalytic performance of materials for HER activity64. We further engaged strain engineering to evaluate its impact on ΔGH. It is observed that 2% compressive strain improves the catalytic performance of holey graphyne by reducing ΔGH to 0.03 eV. Moreover, the-ab initio molecular dynamics simulation assures the thermal stability of hydrogen adsorbed holey graphyne concerning the low energy fluctuations alongside the intact H and the sheet structure after a time of 8 ps (see Fig. 7c). To address the photocatalytic performance of holey graphyne, the positions of valence band maxima (VBM) and conduction band minima (CBM) relative to the vacuum level along with H+/H2 reduction and (O2/H2O) oxidation potentials for water splitting are presented in Fig. 7d. It is evident that CBM is more positive than the H+/H2 potential, suggesting that holey graphyne is a potential material for photocatalytic hydrogen production. In contrast, the VBM is higher than the O2/H2O oxidation potential; hence, the holey graphene is unsuitable for oxidizing H2O to O2. However, an appropriate band engineering to shift the VBM downwards (e.g., doping, application of an external bias, heterojunction, etc.) may enable the water splitting into H2 and O2.

Figure 7
figure 7

(a) The unit cell of holey graphyne with the possible adsorption sites (s1–s9) of hydrogen. (b) Calculated Gibbs free energy (ΔGH) of hydrogen adsorption on holey graphyne. The corresponding values on graphene [Ref. 61], C3N4 [Ref. 63], and biphenylene [Ref. 61] are also included for comparison. (c) The fluctuation of total energy during ab-initio molecular dynamics simulations at 300 K alongside the initial and the final hydrogen-adsorbed holey graphene structures. The brown and green spheres represent C and H atoms, respectively. (d) Calculated band edge positions of pristine and strained (2% compressive strain) holey graphene relative to the vacuum level.

Full size image

Conclusions

Using the first-principles computations, we investigated the thermal transport properties and HER catalytic activity of holey graphyne, a semiconducting material with a direct bandgap of 1.0 eV. The closely packed phonons display positive frequencies across the entire Brillouin zone, assuring the dynamic stability of holey graphyne. The room temperature |S| is as high as 300 μV/K at 1 × 1012 cm−2 (nearly sixfold than graphene). Flexural acoustic phonons couple with the lowest optical phonon mode to enhance phonon–phonon scattering and eventually decrease κl. The estimated room temperature κl of 29.3 W/mK is ultralow compared to that of graphene (3000 W/mK), attributed to its low bond stiffness. With the aid of nanostructure engineering, κl further reduces by a factor of four for a crystal size of around 335 nm (110 nm) at 300 K (600 K). The room temperature p-type ZT achieves a maximum value of 1.50, largely greater than graphene (0.55 × 10–3), demonstrating its superior thermoelectric performance. Lastly, holey graphyne having low overpotential and more positive CBM than the H+/H2 potential, is also capable of catalyzing HER.