Hierarchical chemical bonding and multi-valley band edge-induced high performance in layered Bi₆Ag₂O₆Se₄: a theoretical study

Abstract

Exploring novel intrinsically low lattice thermal conductivity materials has become an effective strategy to obtain high thermoelectric performance. Based on first-principles calculations and Boltzmann transport theory, this study investigates the electronic band structure and thermoelectric properties of the layered oxychalcogenide Bi₆Ag₂O₆Se₄. As an indirect bandgap semiconductor, Bi₆Ag₂O₆Se₄ possesses a sandwich structure composed of positively charged [Bi₂O₂]²⁺ oxide layers and negatively charged Se²⁻ and [Ag₂Se₂]²⁻ layers stacked along the c-axis. This layered architecture imparts a distinct anisotropic crystal structure and transport properties. The hierarchical chemical bonding, weak interlayer interaction, lone pair electrons of Bi and rattler-like behavior of Ag result in strong anharmonicity and intrinsically low lattice thermal conductivity. The Ag-4d and Se-4p hybridization at the valence band maximum induces multi-valley band edges and gives rise to excellent electrical transport properties, especially along the in-plane direction. Benefiting from the synergistic effect of low thermal conductivity and multi-valley band edges, the maximum ZT of p-type doping increases from ∼0.24 at 300 K to ∼2.19 at 900 K with the optimized carrier concentration of 2.94 × 10²⁰ cm⁻³ and 4.92 × 10²⁰ cm⁻³, respectively. This study identifies Bi₆Ag₂O₆Se₄ as a promising high performance thermoelectric candidate. Moreover, understanding the relation between the anisotropic crystal structure and thermoelectric properties enables further optimization of thermoelectric performance.

---

1. Introduction

Against the backdrop of the increasing global energy demand and escalating environmental pollution challenges, the development of high-efficiency clean energy conversion technologies has emerged as a pivotal endeavor for sustainable development. Thermoelectric materials, as functional materials capable of direct thermoelectric conversion, exhibit substantial application potential in addressing energy scarcity and environmental challenges.^1–4 The energy conversion efficiency of thermoelectric materials is typically characterized by the dimensionless figure of merit ZT, defined as ZT = S²σT/κ_tot, where S, σ, T, and κ_tot refer to the Seebeck coefficient, electrical conductivity, operating temperature in Kelvin, and total thermal conductivity, respectively. Superior thermoelectric conversion efficiency is inherently linked to materials exhibiting a high dimensionless figure of merit ZT over a wide temperature range. Considering the inherent coupling among these parameters, achieving simultaneous enhancement of the Seebeck coefficient and electrical conductivity while reducing the total thermal conductivity remains a pivotal challenge in elevating the ZT value.^5,6 The attainment of this objective requires that the materials satisfy two critical prerequisites: first, enhancing the power factor (PF = S²σ) by improving the Seebeck coefficient and electrical conductivity; second, effectively suppressing the total thermal conductivity to sustain a substantial temperature gradient within the material.

In recent years, a series of strategies have been developed to improve the thermoelectric performance and have achieved significant progress.^7–11 These include enhancing the power factor through regulation of carrier concentration,^12,13 electronic band structure engineering,¹⁴ construction of heterostructures or superlattices,^15,16 and nanoscale structural modulation of electrical transport;^17–19 as well as strengthening phonon scattering to reduce lattice thermal conductivity (κ_lat) by introducing multi-scale defects,^20,21 constructing nanoscale precipitates for scattering,^22,23 designing porous structures,^24,25 regulating lattice anharmonicity,^26–30 and developing novel materials with intrinsically low κ_lat.

Layered oxychalcogenides possess a unique structural characteristic in which strong intralayer covalent bonds coexist with weak interlayer interactions, providing an optimal platform for balancing electrical transport and thermal conduction characteristics. Continuous intralayer electronic states ensure high carrier mobility, whereas weak interlayer coupling suppresses the lattice thermal conductivity through enhanced phonon scattering, thus making them emerge as potential candidates for thermoelectric performance regulation.^31–34 Within this system, as a typical representative, BiCuOCh (Ch = Se and S) has been extensively studied. The power factor of n-type BiCuOSe can be enhanced by 54% and that of n-type BiCuOS by 74% through tensile strain.³⁵ Layered LaCuSeO, featuring diverse chemical bonds, exhibits multi-valley behavior in the band structure near the Fermi level. The hybridization of La-p and O-p orbitals within layers enhances electrical transport performance, and weak interlayer coupling effectively suppresses lattice thermal conductivity, enabling the optimal ZT to reach 1.46 for n-type LaCuOSe.³⁶ Bi₆Cu₂Se₄O₆, analogous to BiCuSeO, was modified by doping transition metal elements (Ti, Zr, and Ce) to increase the electron carrier concentration, achieving a peak ZT value of ∼0.16 at 873 K.⁸ The NbSe₂ introduced resonant level effectively enhanced the peak ZT values to ∼0.40 at 873 K.³⁷ However, although Bi₆Ag₂O₆Se₄ exhibits an interlayer heterostructure analogous to that of Bi₆Cu₂Se₄O₆, its intralayer bonding characteristics and carrier transport mechanisms remain insufficiently investigated. Ag's larger atomic radius and more flexible electronic configuration endow the material with regulatory space for carrier transport and lattice thermal conductivity not available in Bi₆Cu₂Se₄O₆. Thus, a systematic investigation of Bi₆Ag₂O₆Se₄ is essential for advancing high-performance layered oxychalcogenide thermoelectric materials, as it provides a viable strategy for balancing high thermoelectric efficiency and operational reliability.

This study systematically investigates the thermoelectric transport properties of the layered oxychalcogenide Bi₆Ag₂O₆Se₄ by combining density functional theory (DFT) and the Boltzmann transport equation (BTE). Regarding electrical transport, the intralayer covalent bond network forms highly delocalized electronic states, constructing low-resistance channels for carrier transport. The valence band maximum is derived from the hybridization of Ag-4d and Se-4p orbitals, exhibiting significant anisotropic dispersion and multi-valley characteristics, which are essential for excellent electrical transport properties. For thermal transport, weak interlayer coupling, hierarchical chemical bonding strength, lone pair electrons of the Bi atoms, and three-phonon scattering effectively lower the lattice thermal conductivity. By combining phonon transport and various carrier scattering mechanisms, the anisotropic thermoelectric transport properties of n-type and p-type Bi₆Ag₂O₆Se₄ were evaluated. The results show that p-type doping exhibits higher ZT values along both the in-plane and out-of-plane directions, which originate from the multi-valley band edge below the Fermi level. This work verifies low lattice thermal conductivity and multi-valley band edges of Bi₆Ag₂O₆Se₄ and identifies it as a promising layered thermoelectric material.

2. Computational methodology

First-principles calculations were performed based on Density Functional Theory (DFT) via the Vienna Ab initio Simulation Package (VASP),³⁸ using the Generalized Gradient Approximation (GGA) and the Projector Augmented Wave (PAW) method.^39,40 The cutoff energy was set to 450 eV, with a Brillouin zone sampling grid of 13 × 13 × 2. Both lattice parameters and atomic coordinates were fully relaxed. The electronic band structure was calculated using the Tran–Blaha modified Becke–Johnson (mBJ) functional.⁴¹ Interatomic force constants were obtained via the Moment Tensor Potentials (MTP) method, as implemented in the Machine Learning Interatomic Potentials (MLIP) software.⁴² The training sets were generated from ab initio molecular dynamics (AIMD) simulations within the NVT ensemble at various temperatures (50–800 K),⁴³ as displayed in Fig. S1. Phonon dispersion and density of states were computed using the Phonopy package.⁴⁴ Lattice thermal conductivity was calculated based on the phonon Boltzmann transport equation using the ShengBTE software⁴⁵ with a converged q-grid of 13 × 13 × 4. Carrier relaxation times and transport properties were calculated using the AMSET package⁴⁶ with three scattering mechanisms considered: ionized impurity (IMP) scattering, acoustic deformation potential (ADP) scattering, and polar optical phonon (POP) scattering. This software is widely applied to both the two-dimensional and bulk structural thermoelectric materials.^47,48 Detailed calculation methods are provided in the SI.

3. Results and discussion

3.1 Crystal structure and chemical bonding

As depicted in Fig. 1(a), the layered oxychalcogenide Bi₆Ag₂O₆Se₄ crystallizes in the tetragonal P4/nmm space group, which is constructed by the alternating stacking of positively charged [Bi₂O₂]²⁺ oxide layers and negatively charged Se²⁻ and [Ag₂Se₂]²⁻layers along the c-axis, thereby forming a typical sandwich structure through electrostatic attraction. The optimized lattice constants are a = b = 3.91 Å and c = 21.51 Å. In the [Ag₂Se₂]²⁻ chalcogenide layers, Ag atoms are typically located at the center of tetrahedra and are coordinated with Se atoms. The tetrahedra of Ag atoms are connected through shared edges. However, the Bi atom is located at the top of the square formed by four O atoms and adjacent square pyramids are connected by sharing a common edge. The Se²⁻ layer inserts between two [Bi₂O₂]²⁺ layers and forms Bi–Se bonds with adjacent Bi atoms. This multi-layered crystal structure usually possesses anisotropic transport properties, i.e. high carrier mobility along the intralayer direction and lower thermal conductivity along the interlayer direction.

Fig. 1 (a) Schematic diagram of the Bi₆Ag₂O₆Se₄ crystal structure. (b) Calculated electronic band structure of Bi₆Ag₂O₆Se₄ along high-symmetry directions in the first Brillouin zone. (c) Atomic projected density of states (PDOS) of Bi₆Ag₂O₆Se₄. (d) The projected crystal orbital Hamilton population (pCOHP) analysis. The bonding length and −ICOHP values (within parentheses) are also listed.

The band gap is a critical parameter in the evaluation of the transport properties in semiconductors. The Tran–Blaha modified Becke–Johnson (MBJ) local potential was adopted to obtain the band gap,^49,50 which can obtain a band gap comparable to that obtained using expensive hybrid functions.⁵¹ Clearly, the valence band maximum (VBM) is concentrated along the Γ–M direction, whereas the conduction band minimum (CBM) is located at the Γ point, illustrating the indirect wide band gap (∼1.19 eV) in Bi₆Ag₂O₆Se₄, as displayed in Fig. 1(b). There are two valleys that are located in X and R directions within the energy ∼0.2 eV lower than the VBM. However, for the conduction band, the energy difference ∼0.02 eV can be found between the CBM and the nearest valley located at the Z point. Such a small energy difference at the band edge illustrates the multi-valley band structure in both the valence and conduction bands, suggesting the potential high degeneracy and larger Seebeck coefficient in Bi₆Ag₂O₆Se₄. Fig. 1(c) shows the orbital projected density of states (PDOS). The valence band is mainly contributed by the strong hybridization of Ag-4d and Se1-4p orbitals. A sudden increase of PDOS at the valence band edge indicates the larger effective mass of holes, as displayed in Table S1, which is a key factor for a large Seebeck coefficient. While the conduction band mainly originates from Bi-6p orbitals, the steady increase of PDOS illustrates the strong dispersion of the conduction band (Fig. 1(b)), low effective mass and potential high carrier mobility of electrons (Table S1).

The projected crystal orbital Hamiltonian population (pCOHP) analysis was conducted to understand the bonding characteristics, as displayed in Fig. 1(d). A positive −pCOHP value indicates a bonding state, while a negative −pCOHP value corresponds to an antibonding state. Clearly, Bi–Se1, Bi–Se2 and Bi–O exhibit antibonding states around the VBM, originating from the partially oxidized Bi³⁺ ions with lone pair 6s² electrons. However, for Ag–Se, a wide range of antibonding states can be observed, which results from the hybridization of Ag-4d and Se1-4p orbitals. These occupied antibonding states are evidenced to enhance lattice anharmonicity and lower the lattice thermal conductivity,^52–54 which is essential for high performance thermoelectric materials. Moreover, the antibonding states in the VBM can lead to a large defect tolerance, which is beneficial for defect modulation in thermoelectric materials. The integrated crystal orbital Hamiltonian population (ICOHP) quantifies bond energy strength, with a larger absolute −ICOHP value signifying stronger bonding. As shown in Fig. 1(d), the intralayer Bi–O bonds exhibit strong bonding characteristics, alongside enhanced atomic orbital hybridization, which endows the material with excellent structural rigidity and thermal stability. In contrast, the interlayer Bi–Se2 bonds display relatively weak bonding, illustrating the weak interaction between [Bi₂O₂]²⁺ and adjacent [Ag₂Se₂]²⁻ and Se²⁻ layers. These diverse bonding lengths and strengths suggest the existence of hierarchical bonding in Bi₆Ag₂O₆Se₄, which is an essential factor for the low thermal conductivity.

3.2 Phonon transport properties

Both the lone pair electrons from Bi³⁺ and hierarchical chemical bonds along different directions suggest the damped phonon transport in Bi₆Ag₂O₆Se₄, which is fundamental in high performance thermoelectric materials. Fig. 2(a) reveals the calculated phonon spectrum of Bi₆Ag₂O₆Se₄ along high symmetry lines. The absence of imaginary frequencies suggests the dynamic stability of Bi₆Ag₂O₆Se₄. The spectrum consists of 3 acoustic branches and 51 optical branches, due to the relatively larger unit cell with 18 atoms. A relatively low frequency of ∼0.2 THz can be observed in acoustic phonons at the Z point. The low-lying optical branches observed at Γ, M, and A points overlap with acoustic phonon branches and further decrease the acoustic phonon frequency and group velocity, which contributes to the low thermal conductivity. The eigenvectors of the lowest optical branches at Γ, M, and A points suggest contrasting amplitudes of atoms, as displayed in Fig. S2, originating from the different interatomic force constants and coexistence of weak and strong chemical bonds in the atomic pairs. Moreover, clearly rattler-like behavior of Ag and Bi atoms can also be found in the atomic trajectory of Bi₆Ag₂O₆Se₄ at 600 K, as displayed in Fig. S3. An energy gap in the range of 5–7 THz between the lower and upper optical branches can be found, which might originate from the coexistence of the strong and weak bonding strength in the intralayer and interlayer regions (Fig. 1(d)). Such a frequency gap can restrict the phase space for normal three-phonon processes due to energy and momentum conservation constraints.^55–57 However, the overlap between the acoustic branches and low-lying optical branches will reduce the lattice thermal conductivity by strengthening the phonon scattering intensity, as discussed later.

Fig. 2 Calculated (a) phonon dispersion, (b) atomic resolved phonon density of states (PhDOS), (c) anisotropic potential well of Bi, Ag, O, and Se atoms in the Bi₆Ag₂O₆Se₄ compound as a function of the atomic displacements, and (d) temperature-dependent atom displacement parameter (ADP) of Bi₆Ag₂O₆Se₄.

The contributions of each atom to the phonon vibrational modes in Bi₆Ag₂O₆Se₄ were evaluated through the projected phonon density of states (PhDOS), as presented in Fig. 2(b). The low frequency range is mainly contributed from the heavy Bi and Ag atoms. Several relatively sharp peaks can be observed for Bi and Ag atoms in the low frequency range, illustrating the potential local vibration of those two atoms, as presented in Fig. S2 and S3. With a relatively light mass, Ag atoms, in the low-to-mid frequency range, actively participate in phonon vibrational modes due to their flexible coordination environment within the [Ag₂Se₂]²⁻ layers. The high-frequency region is dominantly contributed by the lighter O atoms. The eigenvectors of the two sharp peaks are calculated and are displayed in Fig. S2. The significantly larger amplitudes of eigenvectors of Ag and Bi atoms indicate the presence of local vibrations of Ag and Bi atoms, which can act as the local oscillators and introduce strong anharmonicity.

Fig. 2(c) exhibits the potential wells of Bi, Ag, O, and Se along different directions. The deep potential well can be found for Bi and O atoms along both the in-plane and out-of-plane directions, consistent with the strong chemical bonding in Fig. 1(d). A relatively shallow potential is observed for Ag and coordinated Se1 atoms along the in-plane and out-of-plane directions, respectively, suggesting the weak bonding between Ag and Se1 atoms. Notably, Se2 atoms exhibit the shallowest potential well along the in-plane direction, due to the weakly bonded Se2 atoms in the single layer structure of Se2. The thermally induced anisotropic atomic displacement parameters (ADPs) are also calculated, as displayed in Fig. 2(d). Compared to other atoms, the relatively shallow potential well of Ag atoms along both the in-plane and out-of-plane directions reflects the weak confinement by surrounding atoms and gives rise to a significant ADP over the whole temperature range. Se2 exhibits stronger vibrational modes than Se1 in the [Ag₂Se₂]²⁻ layers, which originates from the weak interaction with the lone pair electrons of the Bi atom along the out-of-plane direction. Meanwhile, Ag atoms also present larger ADP with elevated temperature, as evidenced by the steeper slope. At 300 K, the ADP of Ag atoms is 1.5 times that of Bi atoms, consistent with the flexible coordination and weak bonding environment of Ag in the [Ag₂Se₂]²⁻ layers, as displayed in Fig. 1(d). In contrast, the Bi atom exhibits a deeper potential well in its potential energy curve, which not only restricts the amplitude of its thermal vibrations but also gives rise to strong lattice anharmonicity in its vibrational characteristics.

The aforementioned discussions all suggest the potential of low thermal conductivity of Bi₆Ag₂O₆Se₄; the lattice thermal conductivity was calculated based on the Boltzmann transport equation with the three phonon scattering mechanism being considered. The introduction of fourth-order phonon scattering might lead to further reduction in the lattice thermal conductivity, but the extent of the reduction varies for each material.^58,59 The convergence of lattice thermal conductivity with respect to the q-grid was tested, as shown in Fig. S4. Fig. 3(a) displays the direction-dependent lattice thermal conductivity as a function of temperature. A significant anisotropy of lattice thermal conductivity can be observed, ∼1.91 W m⁻¹ K⁻¹ and ∼0.26 W m⁻¹ K⁻¹ along the a and b axes (in-plane IP) and c axis (out-of-plane, OP) directions at 300 K, which originates from the very weak interaction across the stacked multi-layers. It can be observed that κ_lat gradually decreases with increasing temperature, i.e. from ∼1.91 (∼0.26) W m⁻¹ K⁻¹ at 300 K to ∼0.66 (∼0.09) W m⁻¹ K⁻¹ at 900 K along the IP (OP) direction, resulting from the enhanced phonon scattering with the elevated temperature. The obtained lattice thermal conductivity is comparable to that reported for isostructural compounds such as Bi₆Cu₂O₆Se₄ (ref. 60) and La₂Bi₄Cu₂O₆Se₄.⁹ The grain size is a fundamental factor affecting the lattice thermal conductivity of polycrystalline materials.⁶¹ Fig. 3(b) illustrates the cumulative lattice thermal conductivity as function of the phonon mean free path (MFP) at 300 K. Obviously, the lattice thermal conductivity increases with elevated MFP and achieves its saturation after the MFP becomes larger than ∼75.6 nm along both the in-plane and out-of-plane directions. Thus, the fine grain size less than ∼2.61 nm can effectively suppress the lattice thermal conductivity to 50% along the in-plane direction. However, for out-of-plane directions, the grain size less than ∼19.75 nm is required to reduce the lattice thermal conductivity to 50%. These results suggest that the grain size modulation is a potential effective strategy to further decrease the lattice thermal conductivity in Bi₆Ag₂O₆Se₄. Fig. 3(c) shows the thermal conductivity spectrum and cumulative lattice thermal conductivity as a function of phonon frequency at 300 K. The lattice thermal conductivity is dominantly contributed by the acoustic phonon branches for both the in-plane and out-of-plane directions, which is evidenced by the slope of cumulative lattice thermal conductivity. A sharp decline in the lattice thermal conductivity spectrum is found along the in-plane direction at the acoustic phonon cutoff frequency, illustrating the effective phonon suppression due to the coupling between acoustic and optical phonons.

Fig. 3 Thermal conductivity of Bi₆Ag₂O₆Se₄. (a) Temperature-dependent lattice thermal conductivity. (b) Cumulative κ_lat as a function of the mean free path (MPF) at 300 K. (c) Thermal conductivity spectra dκ_lat/dω (dashed line) and frequency dependence of cumulative κ_lat at 300 K.

The group velocity as a function of frequency at 300 K is displayed in Fig. 4(a). The weak interaction between atoms usually leads to slow group velocity in the low frequency range. Relatively smaller group velocities can be observed in the low frequency range, which originates from the strong coupling between the low-lying optical branches and acoustic branches. However, for the higher frequency range, the increase of group velocity can vanish due to the minor contribution of optical branches as shown in Fig. 3(c). The mode-dependent Grüneisen parameters at 300 K show that the value in the low-frequency region is as high as 20 (Fig. 4(b)), suggesting the strong anharmonicity of phonon vibration, which is attributed to the coexistence of weak and strong interactions and lone pair electrons of Bi atoms. Thus, the small group velocity and strong anharmonicity contributes to the low lattice thermal conductivity in Bi₆Ag₂O₆Se₄.

Fig. 4 Mode-dependent (a) group velocity and (b) Grüneisen parameter (γ) of Bi₆Ag₂O₆Se₄. (c) Three-phonon scattering phase space for Bi₆Ag₂O₆Se₄, including total, combination and splitting processes. (d) Three-phonon scattering rates for Bi₆Ag₂O₆Se₄, including total, splitting and combination processes.

To gain further understanding of the mechanism of the low lattice thermal conductivity, the thermal transport behavior was evaluated in this work. Fig. 4(c) displays the weighted phase space (WP3) of Bi₆Ag₂O₆Se₄ within the three-phonon scattering mechanism. The weighted phase space is an indicator of the available scattering channels among all phonon modes, and the larger weighted phase space illustrates more available scattering channels and stronger phonon scattering.⁶² Strengthening phonon scattering requires relatively larger phase space. A higher weighted phase space observed in this work is comparable to that in Bi₆Cu₂O₆Se₄ (ref. 60) and La₂Bi₄Cu₂O₆Se₄,⁹ suggesting the enhanced phonon scattering in Bi₆Ag₂O₆Se₄. Moreover, the larger WP3 is contributed by the combination process especially in the acoustic phonon mode. However, for the higher frequency range, the WP3 predominately originates from the splitting process. Such strong phonon scattering results in relatively higher scattering rates, as displayed in Fig. S5. Clearly, the 3ph scattering process plays a dominant role compared to the boundary and isotope scattering processes, indicating the necessity of considering three-phonon scattering. Moreover, 3ph scattering rates increase with increasing phonon frequency in the low-frequency range, illustrating strong coupling between acoustic and low-lying optical phonon branches. Fig. 4(d) presents the phonon frequency-dependent 3ph scattering rates that include total, splitting and combination processes at 300 K. Clearly, the combination process plays an important role in the low (<2 THz) and middle frequency range (8–11 THz), illustrating that the strong anharmonicity induced by lone pair electrons of the Bi atom acts as the essential mechanism in the combination process. However, the scattering rates of the splitting process increase with increasing phonon frequency, especially in the low frequency range. Thus, the low group velocity induced by the coexistence of weak and strong interactions and significant 3ph scattering brought by lone pair electrons contribute to the low lattice thermal conductivity in Bi₆Ag₂O₆Se₄, which is required by high-performance thermoelectric materials.

3.3 Thermoelectric transport properties

Considering the anisotropic crystal structure of Bi₆Ag₂O₆Se₄, the thermoelectric transport properties along the IP and OP directions were symmetrically investigated. The carrier scattering rate was evaluated through the inclusion of ADP, POP and IMP scattering mechanisms, as implemented in the AMSET package. The evaluated electrical transport properties, including electrical conductivity (σ), Seebeck coefficient (S), and power factor (PF), for p- and n-type doping along IP and OP directions are presented in Fig. 5 and S6, respectively. As shown in Fig. 5(a), the electrical conductivity increases significantly with the rise in carrier concentration for both n-type and p-type Bi₆Ag₂O₆Se₄ along the IP direction. The σ for n-type doping is obviously higher than that of the p-type doping counterpart, originating from the large mobility induced by the smaller effective mass of electrons than holes (Table S1), as displayed in Fig. S7 and 1(b). Moreover, the decrease of σ with increasing temperature can be observed in both the n- and p-type doping for a given carrier concentration, which is attributed to enhanced carrier scattering with elevated temperature. As observed, the σ along the OP is significantly lower than in the IP, which is attributed to strong interlayer scattering-induced low carrier mobility, as displayed in Fig. S6(a). Different from the dependence of σ on carrier concentration and temperature, the absolute Seebeck coefficient (|S|) increases with increasing temperature at given carrier concentration, while it decreases with rising carrier concentration for a constant, as displayed in Fig. 5(b). Moreover, the |S| for the p-type doping is higher than that for the n-type doping under the same conditions, which resultes from the larger effective mass of hole carriers and smaller energy difference between the VBM and sub-valleys. The contrary dependence relation of σ and S on carrier concentration leads to an optimized power factor (PF = S²σ), as presented in Fig. 5(c). The PF increases with increasing carrier concentration and achieves significantly exceptional values in the p-type doping system, suggesting the advantages of multi-valley contribution to Seebeck coefficient in electrical transport optimization. Furthermore, the PF along the IP direction significantly outperforms that along the OP direction for both the n- and p-type doping, as presented in Fig. S6(c).

Fig. 5 Evaluated thermoelectric properties of n-type (solid symbols) and p-type (hollow symbols) Bi₆Ag₂O₆Se₄ in the in-plane directions at 300, 500, 700, and 900 K. (a) Electrical conductivity (σ), (b) Seebeck coefficient (S), (c) power factor (PF), (d) total thermal conductivity (κ_tot), (e) electronic thermal conductivity (κ_ele), and (f) figure of merit (ZT).

The total thermal conductivity (κ_tot), composed of the electronic thermal conductivity (κ_ele) and the lattice thermal conductivity (κ_lat), is displayed in Fig. 5(d). The κ_tot remains a constant in the low carrier concentration range, while exhibiting a sudden increase in the higher carrier concentration range, which results from the continuous increase of κ_ele (Fig. 5(e)). The κ_ele is obtained using the Wiedemann–Franz law, κ_ele = LσT, where L is the Lorentz number, obtained from the Seebeck coefficient based on the single-parabolic model. Therefore, the κ_ele presents the same trend as σ with respect to carrier concentration and temperature. Obviously, a relatively low κ_tot can be observed along the OP direction, as displayed in Fig. S6d, which is attributed to the smaller carrier mobility (Fig. S7) induced lower electronic thermal conductivity (Fig. S6(e)). With exceptional electrical transport properties and lower κ_tot, an optimized figure of merit (ZT) can be anticipated in Bi₆Ag₂O₆Se₄, especially along the IP direction. Fig. 5(f) shows the temperature dependent ZT as a function of carrier concentration for both the n- and p-type doping. ZT increases first and then decreases with increasing carrier concentration and achieves peak values at optimized carrier concentration. For p-type Bi₆Ag₂O₆Se₄, the peak ZT value increases from ∼0.24 at 300 K to ∼2.19 at 900 K with the optimized carrier concentration of 2.94 × 10²⁰ cm⁻³ and 4.92 × 10²⁰ cm⁻³, respectively. Such a high ZT value is significantly higher than ∼0.57 (900 K) in n-type doping along the IP direction, ∼0.59 (900 K) in p-type doping, and ∼0.48 (900 K) in n-type doping along the OP direction, illustrating the predominant contribution of multi-valley bands and high carrier mobility to ZT values, especially in intrinsically low lattice thermal conductivity materials.

4. Conclusion

In summary, evaluating the performance of materials with intrinsically low lattice thermal conductivity provides valuable insights for exploring novel high performance thermoelectric materials. Based on first-principles calculations and Boltzmann transport theory, this study systematically investigates the bonding structure and thermoelectric transport properties of the layered oxychalcogenide Bi₆Ag₂O₆Se₄. The material exhibits significant anisotropy from strong intralayer Bi–O covalent bonds and weak interlayer Bi–Se2 bonds, laying a unique foundation for directional regulation of electron and thermal transport. The detailed analysis of phonon properties reveals that the lone pair electrons of Bi³⁺ and weak bond interactions give rise to strong anharmonicity and efficient 3ph scattering. Combined with the strong coupling between acoustic and optical branches, these factors collectively lead to low lattice thermal conductivity. To evaluate the electron transport characteristics, the carrier scattering rate was evaluated with the ADP, POP, and IMP scattering mechanisms included. The results show that the high carrier mobility improves the electrical conductivity, and the multi-valley band edge further enhances the Seebeck coefficient. The significant improvements in the electron transport properties and excellent low lattice thermal conductivity result in the optimal p-type doped system achieving a maximum ZT of 2.19 at 900 K (4.92 × 10²⁰ cm⁻³) along the in-plane direction. The high ZT in the in-plane direction highlights the significant advantages of multi-valley band edge and high carrier mobility to ZT values, which contributes to optimizing the performance of layered thermoelectric materials. This study not only demonstrate the potentially high performance of Bi₆Ag₂O₆Se₄ along the in-plane direction but also provide a guideline to optimize the performance of layered thermoelectric materials with intrinsically low lattice thermal conductivity.

Conflicts of interest

The authors declare no conflicts of interest.

Data availability

Data will be made available from the corresponding authors upon reasonable request.

Supplementary information: computational details; effective mass; AIMD silulations; phonon eigenvectors of phonon modes; atomic trajectory; convergence test of lattice thermal conductivity; phonon scattering rate; thermoelectric preoperties along out-of-plane direction; carrier mobility. See DOI: https://doi.org/10.1039/d5ta09379a.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (No. 12204156), the Natural Science Foundation of Henan Province (No. 252300423012), and the China Postdoctoral Science Foundation (No. 2023TQ0315 and 2023M743224).

References

1. S. Liu, S. Bai, Y. Wen, J. Lou, Y. Jiang, Y. Zhu, D. Liu, Y. Li, H. Shi, S. Liu, L. Wang, J. Zheng, Z. Zhao, Y. Qin, Z. Liu, X. Gao, B. Qin, C. Chang, C. Chang and L.-D. Zhao, Science, 2025, 387, 202–208.
2. S. Bai, X. Zhang and L.-D. Zhao, Acc. Chem. Res., 2023, 56, 3065–3075.
3. X. Li, X. Meng, Z. Miao, B. Liang and Q. Zou, Appl. Therm. Eng., 2024, 245, 122889.
4. J. Wang, Q. Gao, Y. Gao, Y. Luo, J. Guo, Q. Sun and E. Liang, Appl. Phys. Lett., 2021, 118, 222105.
5. L. D. Zhao, V. P. Dravid and M. G. Kanatzidis, Energy Environ. Sci., 2014, 7, 251–268.
6. D. Zhao, A. Würger and X. Crispin, J. Energy Chem., 2021, 61, 88–103.
7. H. An, D. Wang and W. He, Appl. Phys. Lett., 2025, 127, 143903.
8. J. Zheng, D. Wang and L.-D. Zhao, Scripta Mater., 2021, 202, 114010.
9. P. Ai, S. Tang, S. Bai, D. Wan, W. Guo, P. Zhang, T. Zheng, H. Wang and T. Yan, J. Energy Chem., 2025, 107, 376–385.
10. K. Zhao, D. Wang, T. Hong, J. Zhu, S. Wang, S. Cheng, X. Gao, C. Shan and L.-D. Zhao, Adv. Sci., 2025, 12, e06963.
11. B. Qin, D. Wang, W. He, Y. Zhang, H. Wu, S. J. Pennycook and L.-D. Zhao, J. Am. Chem. Soc., 2019, 141, 1141–1149.
12. Z. Wei, C. Wang, J. Zhang, J. Yang, Z. Li, Q. Zhang, P. Luo, W. Zhang, E. Liu and J. Luo, ACS Appl. Mater. Interfaces, 2020, 12, 20653–20663.
13. J. Zhu, D. Wang, T. Hong, K. Zhao, B. Qin, Y. Qin, S. Wang, X. Gao, S. Cheng, C. Shan and L.-D. Zhao, Appl. Phys. Lett., 2025, 126, 023903.
14. X. Zhang, Z. Wang, B. Zou, M. K. Brod, J. Zhu, T. Jia, G. Tang, G. J. Snyder and Y. Zhang, Chem. Mater., 2021, 33, 9624–9637.
15. D. Wang, J. Ding, Y. Ma, C. Xu, Z. Li, X. Zhang, Y. Zhao, Y. Zhao, Y. Di, L. Liu, X. Dai, Y. Zou, B. Kim, F. Zhang, Z. Liu, I. McCulloch, M. Lee, C. Chang, X. Yang, D. Wang, D. Zhang, L.-D. Zhao, C.-a. Di and D. Zhu, Nature, 2024, 632, 528–535.
16. Z. Xing, Q. Tian, J. Wei, H. Wu, G. Qin and Z. Qin, Appl. Phys. Lett., 2025, 127, 132402.
17. L. Xu, X. Wang, Y. Wang, Z. Gao, X. Ding and Y. Xiao, Energy Environ. Sci., 2024, 17, 2018–2027.
18. J. Lei, H. Wuliji, Q. Ren, X. Hao, H. Dong, H. Chen, T.-R. Wei, J. Zhang, P. Qiu, K. Zhao and X. Shi, Energy Environ. Sci., 2024, 17, 1416–1425.
19. X. Li, Y. Lou, K. Jin, L. Fu, P. Xu, Z. Shi, T. Feng and B. Xu, Angew. Chem., Int. Ed., 2022, 61, e202212885.
20. J. Zhu, Q. Ren, C. Chen, C. Wang, M. Shu, M. He, C. Zhang, M. D. Le, S. Torri, C.-W. Wang, J. Wang, Z. Cheng, L. Li, G. Wang, Y. Jiang, M. Wu, Z. Qu, X. Tong, Y. Chen, Q. Zhang and J. Ma, Nat. Commun., 2024, 15, 2618.
21. F. L. Madarasz and P. G. Klemens, Int. J. Thermophys., 1987, 8, 257–262.
22. Y. Zhu, Z. Han, B. Han, F. Jiang, X. Wu, C.-G. Han, Y. Deng and W. Liu, ACS Appl. Mater. Interfaces, 2021, 13, 56164–56170.
23. S. Yamasaka, Y. Nakamura, T. Ueda, S. Takeuchi and A. Sakai, Sci. Rep., 2015, 5, 14490.
24. G. Romano and A. M. Kolpak, Sci. Rep., 2017, 7, 44379.
25. K. Zhao, Y. Jiao, Q. Sun, R. He, A. Sanson, Z. Zhong, E. Liang and Q. Gao, Phys. Rev. B, 2024, 110, 064322.
26. B. Qin, D. Wang and L.-D. Zhao, InfoMat, 2021, 3, 755–789.
27. Y. Xia and M. K. Y. Chan, Appl. Phys. Lett., 2018, 113, 193902.
28. M. Yan, H. Geng, P. Jiang and X. Bao, J. Energy Chem., 2020, 45, 83–90.
29. Z. Liu, J. Mao, T.-H. Liu, G. Chen and Z. Ren, MRS Bull., 2018, 43, 181–186.
30. Z. Y. Tian, Q. Y. Zhang, Y. W. Xiao, G. A. Gamage, F. Tian, S. Yue, V. G. Hadjiev, J. Bao, Z. Ren, E. Liang and J. Zhao, Phys. Rev. B, 2022, 105, 174306.
31. S. D. N. Luu and P. Vaqueiro, J. Materiomics, 2016, 2, 131–140.
32. P.-P. Shang, J. Dong, J. Pei, F.-H. Sun, Y. Pan, H. Tang, B.-P. Zhang, L.-D. Zhao and J.-F. Li, Research, 2019, 2019, 9253132.
33. P. Jood, M. Ohta, O. I. Lebedev and D. Berthebaud, Chem. Mater., 2015, 27, 7719–7728.
34. M.-L. Liu, L.-B. Wu, F.-Q. Huang, L.-D. Chen and I. W. Chen, J. Appl. Phys., 2007, 102, 116108.
35. C. Zou, C. Lei, D. Zou and Y. Liu, Materials, 2020, 13, 1755.
36. N. Wang, M. Li, H. Xiao, Z. Gao, Z. Liu, X. Zu, S. Li and L. Qiao, npj Comput. Mater., 2021, 7, 18.
37. J. Zheng, T. Hong, D. Wang, B. Qin, X. Gao and L.-D. Zhao, Acta Mater., 2022, 232, 117930.
38. J. Hafner, J. Comput. Chem., 2008, 29, 2044–2078.
39. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865–3868.
40. P. E. Blöchl, Phys. Rev. B:Condens. Matter Mater. Phys., 1994, 50, 17953–17979.
41. J. A. Camargo-Martínez and R. Baquero, Phys. Rev. B:Condens. Matter Mater. Phys., 2012, 86, 195106.
42. I. S. Novikov, K. Gubaev, E. V. Podryabinkin and A. V. Shapeev, Mach. Learn.: Sci. Technol., 2021, 2, 025002.
43. M. E. Tuckerman, P. J. Ungar, T. von Rosenvinge and M. L. Klein, J. Phys. Chem., 1996, 100, 12878–12887.
44. A. Togo and I. Tanaka, Scripta Mater., 2015, 108, 1–5.
45. W. Li, J. Carrete, N. A. Katcho and N. Mingo, Comput. Phys. Commun., 2014, 185, 1747–1758.
46. A. M. Ganose, J. Park, A. Faghaninia, R. Woods-Robinson, K. A. Persson and A. Jain, Nat. Commun., 2021, 12, 2222.
47. T. Yue, Y. Zhao, J. Ni, S. Meng and Z. Dai, npj Comput. Mater., 2023, 9, 17.
48. D. Wang, K. Zhao, T. Hong, J. Zhu, H. Shi, B. Qin, Y. Qin, G. Wang, X. Gao, S. Cheng, C. Shan and L.-D. Zhao, Mater. Today Phys., 2025, 51, 101654.
49. A. D. Becke and E. R. Johnson, J. Chem. Phys., 2006, 124, 221101.
50. H. Jiang, J. Chem. Phys., 2013, 138, 134115.
51. J. Heyd and G. E. Scuseria, J. Chem. Phys., 2004, 121, 1187–1192.
52. J. He, Y. Xia, W. Lin, K. Pal, Y. Zhu, M. G. Kanatzidis and C. Wolverton, Adv. Funct. Mater., 2022, 32, 2108532.
53. J. Yuan, Y. Chen and B. Liao, J. Am. Chem. Soc., 2023, 145, 18506–18515.
54. A. Das, K. Pal, P. Acharyya, S. Das, K. Maji and K. Biswas, J. Am. Chem. Soc., 2023, 145, 1349–1358.
55. S. Mukhopadhyay, L. Lindsay and D. S. Parker, Phys. Rev. B, 2016, 93, 224301.
56. J. Liu, T. A. Strobel, H. Zhang, D. Abernathy, C. Li and J. Hong, Mater. Today Phys., 2021, 21, 100566.
57. L. Shu, Y. Xia, B. Li, L. Peng, H. Shao, Z. Wang, Y. Cen, H. Zhu and H. Zhang, npj Comput. Mater., 2024, 10, 2.
58. A. R. Natarajan, L. Ponvijayakanthan, M. K. Gupta, R. Mittal, D. J. Singh and V. Kanchana, Phys. Rev. Mater., 2023, 7, 025405.
59. T. Feng, L. Lindsay and X. Ruan, Phys. Rev. B, 2017, 96, 161201.
60. S. Bai, H. Shi, Y. Wen, Y. Hu, J. Zheng, Y. Qin, L. Su, S. Liu, D. Liu, T. Gao, T. Hong, X. Gao, F. Zhu, B. Qin and L.-D. Zhao, Phys. Rev. X, 2025, 15, 031033.
61. C. Zhao, Q. Deng, W. Yuan, X. An, W. Su, Z. He, Y. Xie, Z. Zhao and R. Ang, J. Mater. Chem. A, 2024, 12, 9066–9074.
62. P. Ai, S. Tang, S. Bai, D. Wan, W. Guo, P. Zhang and S. Dong, Chem. Eng. J., 2024, 499, 156458.

---