Ultra-thin Carbon Biphenylene Network as an Anisotropic Thermoelectric Material with High Temperature Stability Under Mechanical Strain

Gözde Özbal Sargın Faculty of Engineering and Natural Sciences, Sabanci University, Istanbul, 34956 Turkey Salih Demirci Department of Physics, Kırıkkale University, Kırıkkale, 71450 Turkey Kai Gong Department of Civil and Environmental Engineering, Rice University, Texas, 77005 USA V. Ongun Özçelik [email protected] Faculty of Engineering and Natural Sciences, Sabanci University, Istanbul, 34956 Turkey Materials Science and Nano Engineering Program, Sabanci University, Istanbul, 34956 Turkey

Abstract

Carbon biphenylene network (C-BPN), which is an ultra-thin material consisting of carbon atoms arranged in square-hexagonal-octagonal (4-6-8) periodic rings, has intriguing properties for nano-scale device design due to its unique crystal structure. Here, using the Landauer formalism in combination with first-principles calculations, we show that C-BPN is a highly stable thermoelectric material at elevated temperatures under mechanical strain, where its thermoelectric efficiency can be anisotropically engineered. Transport calculations reveal that C-BPN’s transmission spectrum has significant degrees of directional anisotropy and it undergoes a metal-insulator transition under strain, which leads to an increase in its Seebeck coefficient. C-BPN’s lattice thermal conductance can be selectively tuned up to 35 $\%$ bidirectionally at room temperature by strain engineering. Enhancement in its power factor and the suppression of its lattice thermal conductance improves the p-type figure of merit up to 0.31 and 0.76 at 300 and 1000 K, respectively. Our findings reveal that C-BPN has high potency to be used in thermoelectric nano devices with selective anisotropic properties at elevated temperatures.

Designing thermoelectric (TE) materials that can efficiently convert energy between heat and electricity is one of the bottlenecks of solid-state energy harvesting problem where the ideal material should have high electrical conductivity and high Seebeck coefficient ( $S$ ). In this field, nanomaterials and two-dimensional (2D) films are gaining interest due to their light weight, ease of production, flexibility, and variety of applications. Among these, 2D allotropes of carbon such as graphyne, penta-graphene, T-graphene, Net-graphene, and $\Psi$ -graphene have been under spotlight of researchers since the successful realization of graphene.Zhang et al. (2015); Liu et al. (2012); Li et al. (2017); Özçelik and Ciraci (2013); Enyashin and Ivanovskii (2011) These pristine carbon-based materials and their functionalized forms have diverse periodic arrangements as well as symmetries in their crystals which lead to unique physical properties such as ultra-high electron thermal conductivity, negative differential resistance, enhanced charge–discharge rate, superconductivity, high storage capacities for Na and Li-ion batteries. Liu et al. (2012); Tong et al. (2022); Bravo et al. (2018); Zhang et al. (2024); Qiao et al. (2022); Yadav and Dhilip Kumar (2021)

Carbon biphenylene network (C-BPN), a two-dimensional (2D) material consisting of carbon atoms arranged in 4-6-8 periodic rings, was previously predicted based on first-principles calculations Wang et al. (2013a); Denis (2014), and recently synthesized in 2021 by Fan et.alFan et al. (2021). This material stands out due to its unique anisotropic $Pmmm$ symmetry (space group: No. 47) and its synthesis was followed by numerous discoveries on its vibrational, mechanical, magnetic and transport properties Demirci et al. (2022a); Mashhadzadeh et al. (2022); Veeravenkata and Jain (2021a); Son et al. (2022); Ying et al. (2022). Type-II Dirac fermionic and magnonic states of C-BPN’s band structure offer novel topological features Liu et al. (2021); Son et al. (2022); Zhang and Yao (2023). Phonon analysis has shown that anisotropy in the mechanical properties of C-BPN may arise from inconsistent group velocities and different intensities of hybridization between two different types of carbon atoms in its structure.Wang et al. (2022) Based on classical molecular dynamics (MD) simulations, C-BPN has been predicted to have high melting point (4024 K) and Young’s modulus (1019.4 GPa) which are comparable to the reported values for graphene.Pereira et al. (2022) Altering the conductivity and mechanical properties of C-BPN has been achieved by strain engineering and doping.Yang et al. (2023); Ren et al. (2023) Recently, Zhou et al. showed that, buckled analogues of biphenylene can be promising candidates for nanomechanics and ferroelectric storage applications.Zhou et al. (2024) In addition to carbon, other group IV elements and group IV-IV, III-V and IV-VI compounds were also successfully integrated in the biphenylene geometry.Demirci et al. (2022a, b); Gorkan et al. (2023); Zhou (2023); Guo et al. (2023); Abdullahi and Ersan (2024); Nemati-Kande et al. (2024)

Low-dimensional materials have an advantage on exhibiting high thermoelcetric performance due to the quantum confinement effect. Dresselhaus and Hicks have shown that confinement length can be used to optimize the electronic part of the figure of merit ( $ZT$ , a measures of the TE performance of materials) and reduce the lattice thermal conductance ( $\kappa_{ph}$ ) value, Hicks and Dresselhaus (1993a, b) which signals that C-BPN has a potential to be used as a thermoelectric material since low dimensional materials benefit from quantum confinement effect. However, pristine C-BPN also has a disadvantage due to its metallic character for utilization in TE applications, since an improvement in the TE figure of merit requires both high electrical conductivity and high $S$ . Therefore, nano-engineering strategies such as cutting nanoribbons, constructing bilayers, introducing defects and alloying that can open the band gap of C-BPN are needed to increase the $S$ and thus enhance the TE performance of C-BPN.Xie et al. (2023); Lv et al. (2023); Chowdhury et al. (2022)

Refer to caption — Figure 1: (a) Top and side views of the C-BPN monolayer. Its unit cell is indicated by the dashed rectangle with lattice constants $a$ =4.52Å and $b$ =3.76Å . The structure has different C-C bond distances with $d_{1}$ =1.46Å , $d_{2}$ =1.41Å , $d_{3}$ =1.44Å and bond angles with $\alpha$ =125°, $\beta$ =110°, $\gamma$ =145°, $\phi$ =125°, and $\theta$ =90°. (b) The phonon dispersion relation of C-BPN under zero strain. The inset shows the fluctuation of C-BPN’s the total energy calculated by AIMD simulations at room temperature. (c) The evolution in the percentage of III-fold coordinated carbon atoms (left axis) and potential energy (right axis) in both unstrained and strained C-BPN as the temperature increases from 1 to 5000 K over 10 ns. The inset figures show the C-BPN structures at selected temperatures. (d) The fluctuations in the total energy (in eV) per unit cell, and (e) the percentage of III-fold coordinated carbon atoms for both unstrained and strained C-BPN during 10 ns of equilibration at 400 K. The inset figures in (e) show the side views of the unstrained and strained C-BPN structures.

Here, we show that monolayer C-BPN maintains its structural stability at elevated temperatures under mechanical stress and its TE efficiency can be selectively engineered biaxially (in the x-and y-directions) by strain engineering. We conduct ballistic transport calculations of stable C-BPN structures using the Landauer formalism combined with ab-initio methods based on density functional theory (DFT).Sup Band gap opening is accomplished at 10 $\%$ -biaxial tensile strain and we show that direction-dependent geometry manifests itself in both the electronic and the phonon transmission spectra. Improvement in the $S$ of C-BPN was achieved due to band gap opening and the emergence of nearly flat bands in its electronic structure under strain. Our results also show that strain engineering at elevated temperatures can enhance the electronic $TE$ coefficient of C-BPN, as well as allowing the selective engineering of its lattice thermal conductance anisotropicly such that, it can be tuned up to 35 $\%$ along the x- and y-directions at room temperature.

Structural stability of C-BPN:

The structural properties of pristine C-BPN and its phonon diagram optimized by DFT are shown in Fig. 1. The phonon diagram indicates stability since all vibrations have positive frequencies. Room temperature stability of C-BPN is tested by calculating the fluctuations of total energy using ab-initio MD simulations, which show that the amplitude of oscillation in total energy (0.08eV) is much lower than the cohesive energy of the material (7.50eV).Demirci et al. (2022a) High-temperature stability tests using a reactive force field MD was used to further justify the material’s structural stability at elevated temperatures and under external mechanical stress. Fig. 1(c) illustrates the evolution of the percentage of III-fold coordinated carbon atoms in both unstrained and strained C-BPN as they are heated from 1 to 5000 K at a rate of 0.5 K/ns, employing a force field MD approach. The unstrained and strained C-BPN exhibit 100% III-fold coordination below a temperature of approximately 540 and 640 K, respectively, indicating stability. The lower threshold for the strained structure suggests that its thermal stability is slightly lower than the unstrained case. Beyond these thresholds, a discernible emergence of non-III-fold coordination in carbon atoms is observed, suggesting the occurrence of C-C bond breaking and formation. This phenomenon is visually depicted in the inset of Fig. 1(c), which show the breaking of C squares and the formation of holes, particularly evident in the case of strained C-BPN.

Following a period of decrease, the percentage of III-fold coordinated carbon atoms begins to rise with increasing temperature. Concurrently, there is a heightened formation of hexagonal carbon structures, as evidenced by the angular distribution evolution.Sup The transformation from C squares and pentagons to C hexagons during this phase is associated with a reduction in potential energy, indicating the energetic favorability of this transformation. To further assess the stability of unstrained and strained C-BPN, we subjected them to a 10 ns heating process at 400K using force field MD. In Fig. 1(d), the fluctuations in the total energy per unit cell for both C-BPN structures oscillate around an average value with an amplitude of 0.03eV per unit cell. This amplitude is significantly smaller than the cohesive energy of the structure (7.50eV), suggesting that these oscillations are insufficient to break apart the monolayer. Coordination analysis in Fig. 1(e) reveals that both unstrained and strained C-BPN maintains 100 % III-fold coordination through the heating process indicating that both of them are stable at 400K. The inset figures in Fig. 1(e) show the straight structure of strained C-BPN, whereas unstrained C-BPN exhibits surface wrinkles.

Anisotropic evolution of the electronic properties:

Having justified the high temperature stability of C-BPN, we next focus on the evolution of its electronic structure under strain. C-BPN shows an intriguing electronic band structure which possesses a type-II Dirac cone along the $\Gamma-X$ (x-direction) path as shown in Fig. 2(a). Here, type-II Dirac point is originated from two crossed bands in which one of the crossing bands is nearly flat and the other one is highly dispersive. Previous studies have shown that, this unique band structure type is present in topological states, type-II Dirac magnons and superconductors.Liu et al. (2021); Zhang and Yao (2023); Liu et al. (2023) As the strain value reaches 10%, the material transforms from a metal to a semiconductor with an indirect band gap of 0.53eV based on DFT calculations using HSE06Heyd et al. (2003); Krukau et al. (2006) functional (0.11eV based on PBEPerdew et al. (1996) calculations) Sup by opening the tilted Dirac cone along the $\Gamma-X$ path.

The effective masses of charge carriers are calculated using the inverse of the local curvature of the electronic band energy by $m^{*}=\hbar\left(\frac{\partial E^{2}(\mathbf{k})}{\partial k^{2}}\right)$ . To calculate the derivative of the band band energy with respect to the k, the finite difference method was used instead of simple parabolic fit in order to obtain a more accurate $m^{*}$ . The effective mass of holes along $\Gamma$ -X (x-direction) and $\Gamma$ -Y (y-direction) are highly anisotropic. The nearly flat valence band along the x-direction results in $m^{*}=8.02m_{0}$ , while the linear-like valence band (VB) along $\Gamma$ -Y gives $m^{*}=0.03m_{0}$ , where $m_{0}$ is the free electron mass. In the lowest conduction band (CB), the parabolic portion has $m^{*}=0.194m_{0}$ along the x-direction and the linear-like portion induces $m^{*}=0.03m_{0}$ along the y-direction. Fig. 2(b) shows that nearly flat portion of the VB gives rise to sharp peak in total density of states (DOS) at the VB edge. This means that more conducting channels will be allowed to participate in the p-type electronic transmission. The main contribution to the total DOS comes from $p_{z}$ orbitals of the C₂ type carbon atoms, while the remaining contribution originated from $p_{z}$ orbitals of C₁ type atoms at the VB edge. Contrarily, at the CB edge, total DOS is smaller than it is at the VB edge as is evident from having more dispersive band shape compared to the VB. Unlike the VB, contribution of the C₁ atom is larger than that of C₂ at the CB edge. In addition, the site projected electronic band structure of strained C-BPN, illustrated in Fig. 2(c), confirms that the valence band maximum and the conduction band minimum are dominated by C₂ and C₁ atoms, respectively.

Phonons and lattice thermal conductivity:

Phonon dispersion relations, which are important for dynamical stability and lattice thermal conductivity, are presented in Fig. 3(a) for pristine and strained C-BPNs. There are no imaginary frequencies in the phonon spectra of pristine and biaxially strained C-BPN. The maximum phonon frequency of unstrained C-BPN reaches to 49.5 THz which is close to that of graphene.Veeravenkata and Jain (2021b) As the biaxial tensile strain increases, the quadratic nature of the ZA mode around the $\Gamma$ point turns into linear along the x- and y-direction, where the group velocity ( $\upsilon_{g,s}(\boldsymbol{q})=\partial\omega_{\alpha}(\boldsymbol{q})/\partial\boldsymbol{q}$ ) of the linearized ZA mode is still smaller than that of the TA and LA modes. At zero strain, ZA phonon mode follows $q^{2}$ dependency up to 3.7THz and 5.8THz for $\Gamma$ -X and $\Gamma$ -Y directions, respectively, which are the highest frequencies of zone boundaries along these directions, as presented in Fig. 3(b). This leads the phonon transmission function of C-BPN along the $\Gamma$ -Y direction to be larger than it is along the $\Gamma$ -X direction at low energies. The group velocities of LA and TA modes along the $\Gamma$ -X ( $\Gamma$ -Y) direction are 18.505km/s (19.801 km/s) and 11.380km/s (11.372km/s), respectively, as $q\rightarrow 0$ . However at the BZ boundary, phonon anisotropy becomes more prominent.Wang et al. (2022) Therefore, it is expected that $\tau_{\mathrm{ph}}(\omega)$ and $\kappa_{\mathrm{ph}}$ of C-BPN to take higher values along the $\Gamma$ -Y direction than the $\Gamma$ -X direction. As presented in Fig. 3(a), phonon frequencies shift down when the interaction between the atoms are weakened due to external strain. The highest phonon frequency reduces to 38.5 THz when the applied strain reaches to 10 $\%$ .

According to the site projected phonon DOS diagram, C₁ and C₂ atoms contribute less to the phonon densities for the strained system than the unstrained case at lower frequencies, which is due to the linearization of the ZA mode under strain as shown in Fig. 3 (c). However, phonon DOS of the strained system becomes dominant at moderate frequencies since the phonon bands become flattened with increasing strain. Phonon DOS of the strained C-BPN goes to zero at low frequencies as also seen from the opening of the phonon band gap between optical phonon modes. Enhancement in the DOS of strained C-BPN contributes to phonon transmission spectrum while decrease in group velocities of LA and TA modes have detrimental impact on phonon transmission spectrum. Since the $\tau_{\mathrm{ph}}(\omega)$ depends on the product $\delta(\omega-\omega_{\alpha}(q))\nu_{\alpha}(q)|_{\parallel}$ Karamitaheri et al. (2013); Datta (2005), competition between the phonon DOS and the group velocity determines the trend of the phonon transmission spectrum.

We also observe anisotropy in the phonon transmission spectrum( $\tau_{\mathrm{ph}}(\omega)$ ) of C-BPN as shown in Fig. 4(a). Phonon transmission decreases with external strain and it gets squeezed into a narrower energy range in line with the phonon dispersion. As strain increases, fewer conducting channels contribute to phonon transport. Due to the higher frequencies of acoustic phonons along the y-direction, $\tau_{\mathrm{ph}}(\omega)$ reaches higher values along the y-direction than the x-direction for all applied strain values. $\kappa_{ph}$ is calculated by integrating the product of phonon window function $\rho(\omega,T)$ and $\tau_{\mathrm{ph}}(\omega)$ over all phonon frequencies. Phonon window function, $\rho(\omega,T)=\hbar\omega(\partial f_{BE}(\omega,T)/\partial T)$ , determines the contribution of phonon modes to the $\kappa_{\mathrm{ph}}$ .

The dependence of $\rho(\omega,T)$ to frequency and temperature can be seen in Fig. 4(c). As temperature increases, the width $\rho(\omega,T)$ enlarges and approaches to 1 even for low frequencies. At low temperatures $\rho(\omega,T)$ filters out the high frequency optical modes, and acoustic modes dominate the $\kappa_{\mathrm{ph}}$ . However for low strain values at 1000 K, $\rho(\omega,T)$ converges to 98 $\%$ of its maximum value which will be achieved with the contribution of optical modes above 1000 K. When there is no strain in the y-direction, $\kappa_{ph}$ has its maximum value while its minimum value appears when there is 10 $\%$ strain in the x-direction. With regard to TE efficiency, low lattice thermal conductance is favorable. Since at high temperatures $\rho(\omega,T)$ approaches to 1, $\kappa_{\mathrm{ph}}$ starts to saturate. The $\kappa_{\mathrm{ph}}$ values are calculated as 1.76 and 2.38 $nW/(Knm^{-1})$ along the x- and y-directions under zero strain at room temperature. Under biaxial tensile strain, $\kappa_{\mathrm{ph}}$ decreases to two-third of its unstrained value along the x-direction. The influence of the biaxial tensile strain on $\kappa_{\mathrm{ph}}$ along the y-direction is relatively weak as compared to the x-direction. The $\kappa_{\mathrm{ph}}$ values drop to 1.14 and 1.90 along the x- and y-directions under 10 $\%$ strain at 300 K. The $\kappa_{\mathrm{ph}}$ values of C-BPN are comparable to those reported for zigzag and armchair graphene nanoribbons (1.64 and 1.15 $nW/K$ ). Xie et al. (2023)

Strain dependent electronic transmission:

Strain-dependent electronic transmission spectra, $\tau_{\mathrm{e}}(E)$ , of C-BPNs under biaxial tensile strain ranging from 0 $\%$ to 10 $\%$ were calculated using their electronic band diagrams (obtained by HSE functionals) where the $\tau_{\mathrm{e}}(E)$ values were normalized according to the width of their transport channels (see Figure S4).Sup . Pristine C-BPN has direction-dependent electron transmission character where the transmission is higher in the y-direction in comparison to the x-direction. Alcón et al. (2022) It should be noted that, due to the existence of available transport channels around its Fermi level, the metallic C-BPN under zero strain has the highest transmission value in comparison to C-BPNs under strain.

When the applied strain reaches 10 $\%$ , $\tau_{\mathrm{e}}(E)$ drops to zero both in the x and y directions around the Fermi level as shown in Fig. 4 (d). In this case, the n-type $\tau_{\mathrm{e}}(E)$ is higher than the p-type $\tau_{\mathrm{e}}(E)$ in the x-direction, while the p-type $\tau_{\mathrm{e}}(E)$ is significantly higher than the n-type $\tau_{\mathrm{e}}(E)$ along the y-direction. As shown in Fig. 4 (d), the p-type $\tau_{\mathrm{e}}(E)$ steeply increases at the VB edge along the y-direction, which is due to the existence of the nearly flat band along the $\Gamma-X$ direction in the electronic band diagram of C-BPN under 10 $\%$ biaxial strain (Fig. 2 (a)). The parabolic CB creates a square root energy dependence of the n-type $\tau_{\mathrm{e}}(E)$ along the y-direction, whereas the linear portion of the CB creates a linear dependence along x-direction. Therefore, the transmission spectra also demonstrate that there is a strong anisotropy in the electronic transport properties of strained C-BPNs. The energy isosurface plots presented in Fig. 4 (e) also show the directional contribution of the conduction and valence bands to $\tau_{\mathrm{e}}(E)$ . When the Fermi level crosses the topmost VB, cigar-shaped isosurfaces elongate in the x-direction. As we go down into deeper energy states, isosurfaces get wider along the y-direction and the cigar-shaped isosurfaces disappear, which indicate that p-type $\tau_{\mathrm{e}}(E)$ along y-direction is higher than p-type $\tau_{\mathrm{e}}(E)$ along x-direction. The three-dimensional visualizations of VBM and CBM presented in Fig. 4 (f) also shows the strong anisotropy discussed above.

Table 1: Reported

ZT

and

\kappa_{\mathrm{ph}}

values of 2D carbon-based materials at room temperature.

2D-Material	ZT(p/n)	$\kappa_{\mathrm{ph}}$
		( $Wm^{-1}K^{-1}$ )
Graphene Reshak et al. (2014)	0.08	14.1
	^a0.56	-
Graphene Anno et al. (2017)	0.55 $\times$ 10^-3	3367.6 $\pm$ 485.0
	^b2.65	-
Graphene Wang et al. (2013b)	-/0.0094	-
Graphyne Wang et al. (2013b)	-/0.157	-
$\gamma$ -Graphyne Jiang et al. (2017)	0.49/0.59(x)	76.4
	0.48/0.56(y)	-
^cPhagrapheneGhosh et al. (2022)	0.58	0.58
CBN Lv et al. (2023)	1.80/-(armchair)	33.19
Graphdiyne Tan et al. (2015)	0.99/3.26	22.3
$\gamma$ -Graphyne Tan et al. (2015)	0.91/1.81	82.3
Penta-graphene Chen et al. (2017)	0.0531/0.019	197.85
^d0.481/0.148	116.73
Graphenylene Zhang et al. (2023)	-/0.48(x)	378
	-/0.39(y)	-
^eC-BPN	0.26 $\times$ 10^-2/0.24 $\times$ 10^-3(x)	1.76
	0.47 $\times$ 10^-2/0.21 $\times$ 10 ${}^{-4}(y)$	2.38
^fC-BPN	0.06/0.04(x)	1.14
	0.31/0.04(y)	1.90

^aH₂S adsorption ^bDefect-ZT/ZT₀ ^cElectronic $ZT$ ^dStrained ^eOur work-unstrained ^fOur work-strained

Anisotropic thermoelectric transport:

Having seen anisotropy in the electronic transmission properties of C-BPNs, we next focus on the anistropic TE transport properties of C-BPN and investigate the effects of strain on the $ZT$ . The electronic TE properties can be derived using the $L_{n}$ integrals, as detailed in the Methods section. In Fig. 5 (a), the electrical conductance, $G_{e}$ , as a function of the chemical potential, $\mu$ , for C-BPN under 10 $\%$ strain is presented in terms of the quantum conductance unit ( $G_{0}=2e^{2}/h$ ) calculated at the HSE level. In the limit $T\rightarrow 0$ , $G_{e}$ approaches to $\tau_{\mathrm{e}}(E)$ while it softens with increasing temperature. A strong anisotropy in the p-type $G_{e}$ along the x-and the y-direction is clearly seen near the Fermi level. Unlike the p-type $G_{e}$ , n-type $G_{e}$ exhibits almost the same trend up to 1eV. Above 1eV, the anisotropic behavior appears with the contribution of the second CB to the transport. In contrast to the high $G_{e}$ , $S$ takes its lowest values in the metallic C-BPN, resulting from the bipolar conduction (shown in Figure S6Sup ). The increase in the $S$ of C-BPN under 10 $\%$ -strain is due to the band gap opening (up to 0.53eV), which is sufficient to avoid bipolar conduction effect up to 600 K. Above 600 K, $S$ starts to decrease due to the bipolar effect. While the increase in the power factor ( $PF$ ) and the $ZT$ accelerate with temperature, it slows down due to the decrease in the the $S$ . Improvement of p-type $S$ along the y-direction is not only due to the band gap opening but also results from the step-like increase in the $\tau_{\mathrm{e}}(E)$ at the the VB edge which is evident from the Mott’s relationCutler and Mott (1969); Paulsson and Datta (2003) given by $S(T,\mu)=-(\pi^{2}k_{B}^{2}T/3|e|)$ [ $\partial$ ln $\tau(E)$ / $\partial E$ ]. Accordingly, as the slope of $\tau_{\mathrm{e}}(E)$ increases, $S$ should also increase. Combining $G_{e}$ and $S$ gives the power factor, which plays a key role in determining TE performance. The p-type $PF$ along the y-direction is significantly higher than it is in the x-direction. The electronic thermal conductance follows the same trend with $G_{e}$ .

The $ZT_{max}$ values as a function of temperature for C-BPN at 10 $\%$ strain calculated at the HSE06 level of theory is presented in Fig. 5(e). (Results at the PBE level of theory are available in the Supporting Information.) The fact that the p-type $ZT_{max}$ along the y-direction is much higher than it is along the x-direction is a another result of the anisotropy of C-BPN, which supports the above discussions. Table 1 summarizes the maximum $ZT$ as well as $\kappa_{\mathrm{ph}}$ values of known 2D carbon-based materials and our work at for 300K. Here $\kappa_{\mathrm{ph}}$ indicates lattice thermal conductance, normalized conductance or conductivity depending on the transport regime. Accordingly, the TE performance of C-BPN is larger than similar carbon-based monolayer semiconductors, such as graphyne, which is a direct semiconductor.Wang et al. (2013b) The p-type $ZT_{max}$ along the y-direction is 0.31 at 300K which is much higher than TE efficiency of graphene as seen from Fig. 5(e). Moreover, $ZT$ of strained C-BPN is comparable to those of $\gamma$ -graphyne, phagraphene, graphenylene and penta-graphene with induced strain as presented in Table 1. The calculated p-type $ZT_{max}$ along the y-direction reaches 0.76 at 1000 K. The p- and n-type $ZT_{max}$ along the x-direction achieves values around 0.15 at 1000 K, while n-type $ZT_{max}$ along the y-direction shows significant increase above 600 K and reaches 0.26. In order to provide a broader perspective into the TE efficiency of C-BPN for both directions under 10 $\%$ -strain, $ZT$ is plotted as a function of temperature and chemical potential in Fig. 5(f). One can see that the difference between the p- and n-type $ZT$ in the y-direction is significantly large, while the n- and p-type efficiencies in the x-direction are close to each other.

In conclusion, using the Landuer formalisim in combination with DFT calculations, we reveal that C-BPN is a highly stable material even at elevated temperatures and external mechanical strain where its TE properties can be selectively engineered bidirectionally using strain engineering. The material undergoes a metal-insulator transition at elevated strain values where the band gap opening improves its $S$ and power factor leading to improved TE performance under strain. In addition, we show that flat band portion of the valence band gives rise to an increase in p-type electronic TE coefficients of C-BPN. Dissimilar effective masses of holes along the x- and y-directions cause a strong anisotropy in the DOS, electronic transmission spectra, and hence electronic TE coefficients. Direction-dependent anisotropy also exists in the phonon dispersion of C-BPN which induces distinct $\kappa_{ph}$ values for x- and y-direction. Lattice thermal conductance $\kappa_{ph}$ is suppressed by 35 $\%$ along x-direction while it is reduced by 20 $\%$ along y-direction at 300 K via strain engineering compared to its unstrained values. The combined effect of the enhancement in its power factor and the suppression of its $\kappa_{ph}$ increases the p-type $ZT$ along y-direction up to 0.31 and 0.76 at 300 and 1000K, respectively. Our findings reveal that ultra-thin carbon biphenylene network is stable and has high potential to be used as a TE material even at elevated temperatures and external mechanical stress.

Methods:

Density functional theory calculations were carried out using the Vienna ab-initio simulation package (VASP).Hafner and Kresse (1997); Kresse and Furthmüller (1996) The projector augmented wave (PAW) method was employed to define the the pseudopotential of Carbon atom.Blöchl (1994); Kresse and Joubert (1999) The exchange-correlation potential has been treated by generalized gradient approximation (GGA) using Perdew–Burke–Ernzerhof (PBE) functionals.Perdew et al. (1996) Energy cutoff value for plane wave basis set was taken as 520eV. Brillouin zone (BZ) is sampled by a $\Gamma$ -centered 12 × 14 × 1 $k$ -meshes.Monkhorst and Pack (1976) Simulation cell was adopted with a vacuum spacing of 20Å along the z direction to eliminate the artificial interactions between adjacent monolayers. During the structural relaxation, the convergence criteria for energy difference between electronic steps was set to $10^{-5}$ eV and the force on each atom was set to $10^{-2}$ eV/Å. Since GGA PBE formalism underestimates the electronic band gaps, electronic band structures were also obtained by applying Heyd–Scuseria–Ernzerhof (HSE06) functionals method.Heyd et al. (2003); Krukau et al. (2006) Dynamical stability analysis were performed by density functional perturbation theory (DFPT) Baroni et al. (2001) as implemented in the PHONOPY package.Togo and Tanaka (2015); Togo (2023); Togo et al. (2023) In phonon dispersion calculations, supercells were taken sufficiently large 4 × 5 × 1. Electronic and phonon band structures are plotted with SUMO toolkit.Ganose et al. (2018)

The thermal stability of the systems was initially tested by ab-initio MD simulations using a microcanonical ensemble method (NVT) Nosé (1984); Hoover (1985) with constant temperatures (300 K), with time steps of 1 fs and a total simulation time of 3 ps. The stability of both unstrained and strained C-BPN structures was further evaluated using force field MD simulations, facilitated by the QuantumATK NanoLab software package. Schneider et al. (2017); Smidstrup et al. (2019) These simulations employed a three-body Tersoff potential, Tersoff (1988a, b) which has been widely used for modeling carbon-based materials such as graphene and carbon nanotubes. This potential explicitly considers the bond order contingent upon the local atomic environment and includes a short-range interaction term delineated by the universal Ziegler-Biersack-Littmark (ZBL) potential. Ziegler and Biersack (1985) The parameters for the force field were taken from previous the literature.Bellido and Seminario (2012) Each structure, consisting of a 20x20 supercell with 2400 atoms, was subjected to an MD run spanning 107 steps, implemented under periodic boundary conditions within an NVT ensemble and controlled by the Nose–Hoover thermostat.Martyna et al. (1992) A timestep of 1fs was employed, culminating in a total simulation time of 10ns for each scenario.

Electronic transmission, $\tau$ (E), is calculated by counting the number of contributing channels to the spectrum at a given energy. This method was shown to be successful in calculating the ballistic TE properties of the 2D transition metal dichalcogenides, group III-VI and group-VA monolayers.Özbal et al. (2019); Çınar et al. (2021); Sevinçli (2017) In order to compute $\tau_{\mathrm{e}}(E)$ accurately, using a dense k-point sampling is a necessary condition. In our electrical transmission calculations, the number of sampling points was set $101\times 101\times 1$ for both PBE and HSE06 based calculations. The electrical conductance $G_{e}=e^{2}L_{0}$ , the Seebeck coefficient (thermopower) $S=(L_{1}/L_{0})/eT$ , the power factor $PF=S^{2}G_{e}$ and the electrons’ contribution to the total thermal conductance $\kappa_{el}=(L_{2}-L_{1}^{2}/L_{0})/T$ , were expressed using the $L_{n}$ integrals within the Landauer formalismSivan and Imry (1986); Esfarjani et al. (2006) using

L_{n}(\mu,T)=-\frac{2}{h}\int\tau(E)(E-\mu)^{n}\left(-\frac{\partial f_{FD}}{\partial E}\right)dE,

(1)

where $\partial f_{FD}(E)/\partial E$ is the derivative of the Fermi-Dirac distribution function. The partial derivative $\partial f_{FD}(E)/\partial E$ is defined as the Fermi window function, which describes the energy window over conducting channels that contribute to the electronic transport. Similarly, the Landauer formalism is used to obtain lattice thermal conductance, which is expressed as

\kappa_{ph}=\frac{1}{2\pi}\int\tau_{ph}(\omega)\hbar\omega\left(\frac{\partial f_{BE}(\omega,T)}{\partial T}\right)d\omega.\cite[cite]{\@@bibref{Authors Phrase1YearPhrase2}{datta1997electronic,rego1998quantized,sevinccli2019green}{\@@citephrase{(}}{\@@citephrase{)}}}

(2)

Here, $f_{BE}(\omega,T)$ is the Bose-Einstein distribution function and $\tau_{ph}(\omega)$ is the phonon transmission spectrum per width, which is obtained using a $101\times 101\times 1$ q-point sampling. Thermoelectric figure of merit, $ZT=S^{2}GT/(\kappa_{el}+\kappa_{ph})$ , is calculated using the electronic and lattice thermal transport parameters.

Acknowledgments:

The authors thank Hâldun Sevinçli for valuable discussions. VOO acknowledges funding from the Scientific and Technological Research Council of Turkey (TUBITAK) under project no 123F264. GÖS acknowledges the support from the TUBITAK project, 123C159. Computing resources used in this work were provided by the National Center for High Performance Computing of Turkey (UHeM) under grant numbers 1013432022 and 1019342024.

References

Zhang et al. (2015) S. Zhang, J. Zhou, Q. Wang, X. Chen, Y. Kawazoe, and P. Jena, Proc. Natl Acad. Sci. USA 112, 2372 (2015).
Liu et al. (2012) Y. Liu, G. Wang, Q. Huang, L. Guo, and X. Chen, Phys. Rev. Lett. 108, 225505 (2012).
Li et al. (2017) X. Li, Q. Wang, and P. Jena, J. Phys. Chem. Lett. 8, 3234 (2017).
Özçelik and Ciraci (2013) V. O. Özçelik and S. Ciraci, J. Phys. Chem. C 117, 2175 (2013).
Enyashin and Ivanovskii (2011) A. N. Enyashin and A. L. Ivanovskii, Phys. Status Solidi B 248, 1879 (2011).
Tong et al. (2022) Z. Tong, A. Pecchia, C. Yam, T. Dumitrică, and T. Frauenheim, Adv. Energy Mater. 12, 2200657 (2022).
Bravo et al. (2018) S. Bravo, J. Correa, L. Chico, and M. Pacheco, Sci. Rep. 8, 11070 (2018).
Zhang et al. (2024) R. Zhang, Y. Hou, X. Guo, X. Li, W. Li, X. Tao, and Y. Huang, Phys. Chem. Chem. Phys. 26, 455 (2024).
Qiao et al. (2022) S.-X. Qiao, C.-H. Sui, L. Yang, Y.-P. Li, Y.-X. Sun, N.-X. Zhang, J.-Q. Bai, N. Jiao, and H.-Y. Lu, Phys. Chem. Chem. Phys. 24, 25767 (2022).
Yadav and Dhilip Kumar (2021) N. Yadav and T. J. Dhilip Kumar, Sustainable Energy Fuels 5, 4060 (2021), URL http://dx.doi.org/10.1039/D1SE00657F.
Wang et al. (2013a) X.-Q. Wang, H.-D. Li, and J.-T. Wang, Phys. Chem. Chem. Phys. 15, 2024 (2013a).
Denis (2014) P. A. Denis, J. Phys. Chem. C 118, 24976 (2014).
Fan et al. (2021) Q. Fan, L. Yan, M. W. Tripp, O. Krejčí, S. Dimosthenous, S. R. Kachel, M. Chen, A. S. Foster, U. Koert, P. Liljeroth, et al., Science 372, 852 (2021).
Demirci et al. (2022a) S. Demirci, Ş. Çallıoğlu, T. Gorkan, E. Aktürk, and S. Ciraci, Phys. Rev. B 105, 035408 (2022a).
Mashhadzadeh et al. (2022) A. H. Mashhadzadeh, M. Z. Dehaghani, F. Molaie, S. Fooladapanjeh, O. Farzadian, and C. Spitas, Comput. Mater. Sci. 214, 111761 (2022).
Veeravenkata and Jain (2021a) H. P. Veeravenkata and A. Jain, Carbon 183, 893 (2021a).
Son et al. (2022) Y.-W. Son, H. Jin, and S. Kim, Nano Lett. 22, 3112 (2022).
Ying et al. (2022) P. Ying, T. Liang, Y. Du, J. Zhang, X. Zeng, and Z. Zhong, Int. J. Heat Mass Transf. 183, 122060 (2022).
Liu et al. (2021) P.-F. Liu, J. Li, C. Zhang, X.-H. Tu, J. Zhang, P. Zhang, B.-T. Wang, and D. J. Singh, Phys. Rev. B 104, 235422 (2021).
Zhang and Yao (2023) M.-H. Zhang and D.-X. Yao, Phys. Rev. B 108, 144407 (2023).
Wang et al. (2022) K. Wang, K. Ren, D. Zhang, Y. Cheng, and G. Zhang, Appl. Phys. Lett. 121 (2022).
Pereira et al. (2022) M. Pereira, W. da Cunha, R. de Sousa, G. A. Nze, D. Galvão, and L. Ribeiro, Nanoscale 14, 3200 (2022).
Yang et al. (2023) G. Yang, Y. Hu, Z. Qiu, B.-L. Li, P. Zhou, D. Li, and G. Zhang, Appl. Phys. Lett. 122 (2023).
Ren et al. (2023) X. Ren, K. Wang, Y. Yu, D. Zhang, G. Zhang, and Y. Cheng, Comput. Mater. Sci. 222, 112119 (2023).
Zhou et al. (2024) M. Zhou, J. Li, J. Li, X. Kong, W.-J. Gong, F. M. Peeters, and L. Li, ACS Appl. Nano Mater. 7, 7906–7915 (2024), eprint https://doi.org/10.1021/acsanm.4c00454, URL https://doi.org/10.1021/acsanm.4c00454.
Demirci et al. (2022b) S. Demirci, T. Gorkan, Ş. Çallıoğlu, V. O. Özçelik, J. V. Barth, E. Aktürk, and S. Ciraci, J. Phys. Chem. C 126, 15491 (2022b).
Gorkan et al. (2023) T. Gorkan, S. Demirci, J. V. Barth, E. Aktürk, and S. Ciraci, J. Phys. Chem. C 127, 8770 (2023).
Zhou (2023) Y. Zhou, Comput. Mater. Sci. 226, 112206 (2023).
Guo et al. (2023) A. Guo, F. Cao, W. Ju, Z. Wang, H. Wang, G.-L. Li, and G. Liu, AIP Advances 13 (2023).
Abdullahi and Ersan (2024) Y. Z. Abdullahi and F. Ersan, Comput. Mater. Sci. 242, 113103 (2024), ISSN 0927-0256, URL https://www.sciencedirect.com/science/article/pii/S0927025624003240.
Nemati-Kande et al. (2024) E. Nemati-Kande, S. Faramarzi, S. Yavari, and M. Shafiee, ACS Appl. Nano Mater. 7, 12431–12444 (2024), eprint https://doi.org/10.1021/acsanm.4c00699, URL https://doi.org/10.1021/acsanm.4c00699.
Hicks and Dresselhaus (1993a) L. D. Hicks and M. S. Dresselhaus, Phys. Rev. B 47, 16631 (1993a).
Hicks and Dresselhaus (1993b) L. D. Hicks and M. S. Dresselhaus, Phys. Rev. B 47, 12727 (1993b).
Xie et al. (2023) Z.-X. Xie, X.-K. Chen, X. Yu, Y.-X. Deng, Y. Zhang, W.-X. Zhou, and P.-Z. Jia, Comput. Mater. Sci. 220, 112041 (2023).
Lv et al. (2023) F. Lv, H. Liang, and Y. Duan, Phys. Rev. B 107, 045422 (2023).
Chowdhury et al. (2022) S. Chowdhury, S. Ghosal, D. Mondal, and D. Jana, J. Phys. Chem. Solids 170, 110909 (2022).
(37) See Supporting Information.
Liu et al. (2023) G.-H. Liu, L. Yang, S.-X. Qiao, N. Jiao, Y.-J. Chen, M.-Y. Ni, M.-M. Zheng, H.-Y. Lu, and P. Zhang, Phys. Chem. Chem. Phys. 25, 2875 (2023).
Heyd et al. (2003) J. Heyd, G. E. Scuseria, and M. Ernzerhof, J. Chem. Phys. 118, 8207 (2003).
Krukau et al. (2006) A. V. Krukau, O. A. Vydrov, A. F. Izmaylov, and G. E. Scuseria, J. Chem. Phys. 125 (2006).
Perdew et al. (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996).
Veeravenkata and Jain (2021b) H. P. Veeravenkata and A. Jain, Carbon 183, 893 (2021b), ISSN 0008-6223, URL https://www.sciencedirect.com/science/article/pii/S0008622321007661.
Karamitaheri et al. (2013) H. Karamitaheri, N. Neophytou, and H. Kosina, J. Appl. Phys. 113 (2013).
Datta (2005) S. Datta, Quantum transport: atom to transistor (Cambridge University Press, 2005).
Alcón et al. (2022) I. Alcón, G. Calogero, N. Papior, A. Antidormi, K. Song, A. W. Cummings, M. Brandbyge, and S. Roche, J. Am. Chem. Soc. 144, 8278 (2022).
Reshak et al. (2014) A. Reshak, S. A. Khan, and S. Auluck, J. Mater. Chem. C 2, 2346 (2014).
Anno et al. (2017) Y. Anno, Y. Imakita, K. Takei, S. Akita, and T. Arie, 2D Mater. 4, 025019 (2017).
Wang et al. (2013b) X.-M. Wang, D.-C. Mo, and S.-S. Lu, J. Chem. Phys. 138 (2013b).
Jiang et al. (2017) P. Jiang, H. Liu, L. Cheng, D. Fan, J. Zhang, J. Wei, J. Liang, and J. Shi, Carbon 113, 108 (2017).
Ghosh et al. (2022) M. Ghosh, S. Ghosal, and D. Jana, Adv. Theory Simul. 5, 2200221 (2022).
Tan et al. (2015) X. Tan, H. Shao, T. Hu, G. Liu, J. Jiang, and H. Jiang, Phys. Chem. Chem. Phys. 17, 22872 (2015).
Chen et al. (2017) C.-P. Chen, C. Liu, L.-L. Liu, L.-S. Zhao, and X.-C. Wang, Mater. Res. Express 4, 105031 (2017).
Zhang et al. (2023) C. Zhang, C. Hou, Y. Lu, L. Zhao, H. Wu, H. Song, J. Rong, L. Yu, and X. Yu, Phys. Chem. Chem. Phys. 25, 31301 (2023).
Cutler and Mott (1969) M. Cutler and N. F. Mott, Phys. Rev. 181, 1336 (1969).
Paulsson and Datta (2003) M. Paulsson and S. Datta, Phys. Rev. B 67, 241403 (2003).
Hafner and Kresse (1997) J. Hafner and G. Kresse, in Properties of Complex Inorganic Solids (Springer, 1997), pp. 69–82.
Kresse and Furthmüller (1996) G. Kresse and J. Furthmüller, Phys. Rev. B 54, 11169 (1996).
Blöchl (1994) P. E. Blöchl, Phys. Rev. B 50, 17953 (1994).
Kresse and Joubert (1999) G. Kresse and D. Joubert, Phys. Rev. B 59, 1758 (1999).
Monkhorst and Pack (1976) H. J. Monkhorst and J. D. Pack, Phys. Rev. B 13, 5188 (1976).
Baroni et al. (2001) S. Baroni, S. De Gironcoli, A. Dal Corso, and P. Giannozzi, Rev. Mod. Phys. 73, 515 (2001).
Togo and Tanaka (2015) A. Togo and I. Tanaka, Scr. Mater. 108, 1 (2015).
Togo (2023) A. Togo, J. Phys. Soc. Jpn. 92, 012001 (2023).
Togo et al. (2023) A. Togo, L. Chaput, T. Tadano, and I. Tanaka, J. Phys.: Condens. Matter (2023).
Ganose et al. (2018) A. M. Ganose, A. J. Jackson, and D. O. Scanlon, J. Open Source Softw. 3, 717 (2018).
Nosé (1984) S. Nosé, J. Chem. Phys. 81, 511 (1984).
Hoover (1985) W. G. Hoover, Phys. Rev. A 31, 1695 (1985).
Schneider et al. (2017) J. Schneider, J. Hamaekers, S. T. Chill, S. Smidstrup, J. Bulin, R. Thesen, A. Blom, and K. Stokbro, Modelling Simul. Mater. Sci. Eng. 25, 085007 (2017), URL https://dx.doi.org/10.1088/1361-651X/aa8ff0.
Smidstrup et al. (2019) S. Smidstrup, T. Markussen, P. Vancraeyveld, J. Wellendorff, J. Schneider, T. Gunst, B. Verstichel, D. Stradi, P. A. Khomyakov, U. G. Vej-Hansen, et al., J. Phys.: Condens. Matter 32, 015901 (2019).
Tersoff (1988a) J. Tersoff, Phys. Rev. B 37, 6991 (1988a).
Tersoff (1988b) J. Tersoff, Phys. Rev. Lett. 61, 2879 (1988b).
Ziegler and Biersack (1985) J. F. Ziegler and J. P. Biersack, in Treatise on heavy-ion science: volume 6: astrophysics, chemistry, and condensed matter (Springer, 1985), pp. 93–129.
Bellido and Seminario (2012) E. P. Bellido and J. M. Seminario, J. Phys. Chem. C 116, 4044 (2012).
Martyna et al. (1992) G. J. Martyna, M. L. Klein, and M. Tuckerman, J. Chem. Phys. 97, 2635 (1992).
Özbal et al. (2019) G. Özbal, R. T. Senger, C. Sevik, and H. Sevinçli, Phys. Rev. B 100, 085415 (2019), URL https://link.aps.org/doi/10.1103/PhysRevB.100.085415.
Çınar et al. (2021) M. N. Çınar, G. Ö. Sargın, K. Sevim, B. Özdamar, G. Kurt, and H. Sevinçli, Phys. Rev. B 103, 165422 (2021), URL https://link.aps.org/doi/10.1103/PhysRevB.103.165422.
Sevinçli (2017) H. Sevinçli, Nano Lett. 17, 2589 (2017), pMID: 28318269, eprint https://doi.org/10.1021/acs.nanolett.7b00366, URL https://doi.org/10.1021/acs.nanolett.7b00366.
Sivan and Imry (1986) U. Sivan and Y. Imry, Phys. Rev. B 33, 551 (1986).
Esfarjani et al. (2006) K. Esfarjani, M. Zebarjadi, and Y. Kawazoe, Phys. Rev. B 73, 085406 (2006).
Datta (1997) S. Datta, Electronic transport in mesoscopic systems (Cambridge University Press, 1997).
Rego and Kirczenow (1998) L. G. Rego and G. Kirczenow, Phys. Rev. Lett. 81, 232 (1998).
Sevinçli et al. (2019) H. Sevinçli, S. Roche, G. Cuniberti, M. Brandbyge, R. Gutiérrez, and L. M. Sandonas, J. Phys.: Condens. Matter 31, 273003 (2019).