1. Introduction

Thermal routing refers to the ability to direct heat energy from a source to a selected drain^[1,2]. By exploiting asymmetric transport mechanisms, it offers significant potential for advanced thermal management scenarios where time-varying heat reallocation is required, such as reconfigurable mode switching between daytime cooling and nighttime overcooling protection in satellites (Figure 1a)^[3,4]. To date, phonon band engineering has been extensively explored as a strategy to realize nonreciprocal heat transfer^[5-8]. However, its intrinsically limited heat-transfer power constrains immediate deployment. In parallel, while thermal nonreciprocity with pronounced flux throughput has been demonstrated with the rising of macroscopic thermal metamaterials^[9-11], they typically rely on mechanically moving media to carry energy, leading to integration challenges and long-term reliability concerns such as mechanical fatigue^[12-15]. Consequently, how to achieve high-power thermal routing in a fully solid-state platform remains an outstanding and unresolved challenge^[16,17].

Display Full Size

Download

Figure 1. Schematic diagram of a solid-state thermal router. (a) Promising application scenario of the programmable thermal routing in satellite. The heat current coming from the source can be directed to the chosen drain by varying the external electric current I; (b) A reciprocal heat flow loop based on internal liquid circulation; (c) A non-reciprocal heat flow loop based on interfacial thermoelectric effect; (d) An enhanced non-reciprocal heat flow loop based on distributive thermoelectric effects.

In recent decades, thermal diode has been under the spotlight of thermal scientists, and the controllable rectification mechanism behind it serves to be a bedrock of thermal routing^[18,19]. One counterintuitive yet influential argument in this field asserts that advection alone cannot make a thermal diode^[20-22], because the time-reversal asymmetry introduced by advection is cancelled by the closed heat-flow loop, as constrained by the mass conservation of macroscopic motion (Figure 1b). However, such argument further provokes a thought: does mass conservation of working media fundamentally prohibit all kinds of the external bias mechanisms from realizing a thermal diode due to heat-flow compensation in the loop? The answer could certainly be “no”, if the binding relationship between mass conservation and heat-energy conservation of the working media could be relaxed. Motivated by this insight, we employ external electric field bias as an analog advection and use solid-state thermoelectric material rather than moving liquid as the working media. By doing so, the initially tightly bounded conservation relationship between external momentum bias and heat flux is loosened due to the existence of entropy-change strategies, enabling thermal rectification and laying a foundation for solid-state thermal routing metadevice^[23-26].

Thermoelectric (TE) effects enable reversible heat–electricity conversion and underpin applications such as cooling^[27,28], energy harvesting^[29,30], and thermal camouflage^[31,32]. While the Thomson effect has received less attention than the Seebeck and Peltier effects due to its weaker magnitude, growing evidence indicates that it is promising to influence heat transport behaviors^[33-35]. In terms of symmetry breaking of heat transfer, the Thomson and Peltier effects are just two sides of the same road whose synergy enables more efficient solid-state thermal routing (Figure 1c,d), yet their combined role remains insufficiently understood.

In this work, we systematically elucidate how the Thomson and Peltier effects break thermal symmetry in distinct but cooperative ways using a thermal diode model. Guided by this analysis, we experimentally realize a thermal metadevice that achieves solid-state thermal routing effect^[36]. Driven by drift-biased carrier–phonon coupling, the compact device delivers a directed heat flux over 300 mW with a near-unity split ratio under vacuum condition. While the Peltier effect governs interfacial heat absorption and release, the Thomson effect plays an auxiliary yet important role in enhancing rectification efficiency, as explained by the evolutionary path of Seebeck coefficient^[37,38]. This work bridges the mutually exclusive high power and solid-state property of thermal metadevice, offering new insights for reconfigurable unidirectional heat transfer and heat information processing.

2. Methods

To elucidate the fundamental physics underlying the routing effect, we begin by analyzing a single-slab model. As illustrated in Figure 2a, the central TE slab is connected to two copper ports, which in turn are attached to external heat reservoirs maintained at constant temperatures T_h and T_c, serving as the heat source and the drain, respectively. Furthermore, two electrodes are applied to the copper ports to apply an external electric field bias to the slab. Within this device, the Peltier effect occurs at the interfaces, while the Thomson and Joule effects contribute to the bulk phenomena. The steady-state governing equation for one-dimensional heat transfer in the slab can be expressed as:

Display Full Size

Download

Figure 2. Demonstration of the thermal rectification effects. (a) The forward case of the single slab thermal rectifier; (b) The asymmetric temperature distribution on the single Bi₂Te₃ slab attributed to Thomson effect; (c) A comparison made between the forward and backward heat flux as a function of electric current I; (d) The rectification ratio of this device. Phases 1 and 2 of the inner Figure represent the parameter domains where the device performance improves or decreases after considering the Thomson effect; (e, f) The S-T diagram of the forward working process in the thermal rectification device with a constant (e) and a temperature-dependent (f) Seebeck coefficient. S-T: Seebeck coefficient–temperature.

(1)$\kappa \frac{d^{2}T}{d{x}^2}-J\tau \frac{dT}{dx}+\frac{J^2}{\sigma }=0$

where κ, J, σ, and τ are thermal conductivity, electric current density, electric conductivity, and Thomson coefficient^[39], respectively. After some algebra, the analytical expression of temperature can be obtained (Section S1).

Here we focus on thermal rectification, quantified by the difference between the forward heat flux q_f and backward heat flux q_b when the source and drain are interchanged. As shown in Figure 2b, in the absence of external bias, heat transport is driven solely by the temperature gradient, producing symmetric temperature profiles in both directions. Once an external electric field bias is applied, however, the profiles become asymmetric due to the odd parity of the Thomson effect with respect to electric current direction: reversing the relative orientation of current and temperature gradient induces either heat absorption or generation within the material. This mechanism accounts for the asymmetry highlighted by the blue curves in Figure 2b. Although the Thomson coefficient is typically small, its influence on heat flux modulation is physically meaningful: as shown in Figure 2c, compared with the case using a temperature-averaged Seebeck coefficient, the forward heat flux increases while the backward flux decreases when the Thomson effect is considered, thus inducing a greater rectification.

To be quantified, the rectification ratio is defined as $\gamma =\frac{q_f-q_b}{max(\left|q_f\right|,\left|q_b\right|)}$, where q_f and q_b are the forward and backward heat fluxes, respectively. Its dependence on electric current is shown in Figure 2d. As current increases, γ first rises, peaks near I = 3 A, and then declines. This behavior reflects the competition between thermoelectric effects (linear scaling) and Joule heating (quadratic scaling). Notably, as long as the constant Seebeck coefficient in the benchmark Peltier-only model is not overestimated, the inclusion of the Thomson effect can enhance the rectification ratio γ over the full current range, as indicated by the Phase 1 region in the inset of Figure 2d. To reveal the underlying physics responsible for the rectification and enhancement effect, we present an intuitive yet in-depth physical picture from the perspectives of the “Seebeck coefficient-Temperature” (S-T) diagram of the working process in the device.

At the microscopic level, a temperature difference ∆T induces an electrochemical potential difference ∆μ_e, giving the Seebeck coefficient $S=-\frac{1}{e}\frac{\Delta {\mu }_e}{\Delta T}$ from the viewpoint of thermopower generation, which couples temperature and electrochemical potential. From the constitutive relation $q=-k\frac{dT}{dx}+JST$, the heat flux q consists of two parts: a purely diffusive term and an electrochemical term, the latter reflecting heat transport by charge carriers. In this sense, the Seebeck coefficient can be viewed as a “specific heat of electricity”^[40,41]. Based on this framework, we first analyze the model in Figure 2a using a temperature-averaged Seebeck coefficient S_av. When considering only the electrochemical term of heat transport, the evolution of the matter state along the current path can be represented in the S-T diagram (Figure 2e). This diagram comprises four processes: (A→B) current flows through the TE slab from hot to cold, transporting heat; (B→C) at the cold end, heat is released as current enters the external electrode, reducing S from S_av to S₀; (C→D) current is driven by the external circuit from the cold end back to the hot end of the rectifier; and finally (D→A) heat is absorbed as current re-enters the TE slab at the hot end, closing the cycle.

According to the cycling process described above, the area enclosed by the path in the graph is exactly the Seebeck voltage generated, which is given by $V_{AB}=-{\int }_{T_c}^{T_h}{S}_B\left(T\right)-S_A\left(T\right)dT$, in response to a temperature difference. Therefore, when multiplied by the electric current, the shaded red region corresponds precisely to the heat released at the cold end during the interfacial process B→C. Similarly, the combined area of the blue and red regions proportionally represents the heat absorbed at the hot end during the process D→A, while the blue area alone corresponds to the net electrical work required to overcome the Seebeck voltage established within the system during one evolution. Since the focus here is on a thermal rectifier, our primary interest lies in the heat transported to the cold end, represented by the red area, which reflects the thermoelectric contribution to thermal transport. Figure 2e shows the S-T evolution path for a device with a constant Seebeck coefficient. In contrast, if the Seebeck coefficient varies with temperature, the corresponding S-T evolution is depicted in Figure 2f. In this case, process A→B can no longer be a horizontal line due to the continuous variation (typically degradation of electrochemical potential) of the Seebeck coefficient along the temperature gradient. As a result, in addition to the original work components W₁ and W₂, an extra contribution W₃ (marked by the pink area) emerges, representing the Thomson effect, which strengthens the heat transported to the drain.

Figure 2e,f represent the forward operation mode of the thermal rectifier; for comparison, the backward operation mode is given in Section S2^[42]. The key distinction between the two cases lies in the direction of the evolution path (clockwise or counter-clockwise), determined by the polarity of the electric current applied to the device. A comparison between them provides an intuitive geometric interpretation: when the temperature dependence of the Seebeck coefficient is taken into account, both heat flux transferred toward the drain end and heat flux pumped from the drain end can be enhanced, leading to the enhanced rectification performance.

3. Results and Discussion

Based on the method above, here we analyze the performance of the solid-state thermal router, where both a theoretical and experimental validation are presented. In general, a thermal router consists of one heat source and multiple drains, enabling directed heat transport from the source to a selected drain. Here, we design a three-port thermal router coupled with three central channels. The ports are fabricated from pure copper. Channels 1 and 3 are made of Bi₂Te₃, while Channel 2 combines copper and Bi₂Te₃. As shown in Figure 3a, the structure is symmetric about the central axis of the heat source. To realize directional heat transport via thermoelectric effects, positive and negative electrodes are attached to the hot ends of Channels 1 and 3, respectively. These two channels are thermally connected but electrically insulated from the source port, achieved using specialized adhesives. This configuration allows current to circulate clockwise or counter-clockwise through the central triangular splitter, thereby routing heat from the source to a chosen drain.

Display Full Size

Download

Figure 3. Performance evaluation of thermal router. (a) Schematic of the device; (b) The split ratio η as a function of the electric current I. The inner figure is the energy flux density on the whole device simulated by COMSOL; (c) Coefficient of performance analysis; (d) The method used to extract the temperature contribution from Thomson effect under the DC excitation. Temperature data are extracted from pairs of symmetric points at the central part of TE channel 1 and channel 3. Then, the current-odd thermoelectric temperature signal can be extracted as T_odd = [T₁(J) - T₃(-J)]/2; (e) Decomposition of the contributions from each effect; (f) The thickness study of the copper sheet coated on channel 2. DC: direct current; TE: thermoelectric; COP: coefficient of performance.

Assuming the source and drains exchange heat with their respective reservoirs at a rate of h₁, and neglecting air convection (i.e., considering a vacuum environment), the mathematical expression of the temperature field distribution can be solved by applying both temperature and heat flux continuity conditions for the coupled equations (Section S3). In the absence of external field bias, the resulting temperature distribution is symmetric with respect to the central axis of the source, characteristic of pure thermal diffusion. By contrast, when a circulating electric current is introduced, a directional heat conveyor belt emerges due to thermoelectric effects. Without loss of generality, consider a counter-clockwise current. As the current flows from channel 1 (Bi₂Te₃) to drain 1 (copper), the electrochemical potential of carriers drops sharply at the junction, releasing heat at a rate of J(S₁ - S₂)*T, where S₁ and S₂ are the Seebeck coefficients of Bi₂Te₃ and copper, respectively, and T is the junction temperature. The heat released is efficiently absorbed by drain 1 due to copper’s high thermal conductivity. Conversely, at the junction between drain 2 and channel 3, carriers move from copper to Bi₂Te₃, where they must absorb heat to traverse from a low to high Seebeck coefficient region. This suppresses the natural heat diffusion toward drain 2. Consequently, heat flow is preferentially routed to drain 1, while drain 2 is effectively isolated.

To quantify routing performance, the split ratio is defined as η = q_D1/(q_D1 + q_D2), where $q_{D1}=q_{D1}^{\text{in}}-q_{D1}^{\text{out}}+q_{\text{joule}}+E_{1}J$ and $q_{D2}=q_{D2}^{\text{in}}-q_{D2}^{\text{out}}+q_{\text{joule}}+E_{2}J$ denote the heat fluxes into drains 1 and 2, respectively. An ideal router yields η = 1 or 0, meaning all heat is directed to a single drain. Detailed derivations are provided in Section S4. As shown in Figure 3b, varying the external current from 0 to 10 A produces a non-monotonic dependence of η, peaking at ~3 A. The inset of Figure 3b illustrates the simulated energy flux distribution at this optimal point: the arrows reveal that nearly all heat is directed to drain 1, with drain 2 almost perfectly isolated. At this condition, the split ratio reaches 0.996, demonstrating excellent routing performance. Besides, since split ratio alone is insufficient to quantify the working performance of an active device like this, an analysis of the coefficient of performance (COP) of electric energy is shown in Figure 3c. It is found that near the optimal working point of thermal routing, a COP of 0.8 can be achieved, which verifies the worth of the energy investment. More details on this can be found in Section S5.

To further reveal the Thomson effect, its contribution to the bulk temperature variation has been isolated from other forms of heating, based on the fact that, among mechanisms directly manifested as bulk temperature changes, only the Thomson effect exhibits an odd-parity response to the direction of electric current (Figure 3d)^[43]. The extracted Thomson-induced temperature contribution increases linearly with electric current, ranging from 1 to 5 A, which agrees well with the governing relation of Thomson heat production rate Jτ∇T. Furthermore, to demonstrate the Thomson contribution to the heat isolation more intuitively, a quantitative analysis of the relative contributions of Peltier, Thomson, and Joule effects is given. As shown in Figure 3e, by considering the Thomson effect, the heat eliminated at the targeted isolation drain can be improved from about 95.4 % to about 97.7 % at the optimal working point, thus proving its auxiliary but unignorable role. Besides, at higher electric current conditions, the improvement brought by Thomson effect can be higher. More details on this part are given in the Section S6.

In addition, it is worth noting that, to prevent the heat released at drain 1 from being further dragged counter-clockwise to drain 2, channel 2 is designed as a composite of both Bi₂Te₃ and copper^[44]. Figure 3f presents the split ratio as a function of both electric current and copper volume fraction, where a trade-off between minimizing parasitic Peltier effects and limiting unwanted thermal conduction has been made (Section S7).

In the previous section, a vacuum environment was assumed to evaluate the ideal thermal split ratio of the device. This assumption is reasonable for aerospace application scenarios where the ambient air is sufficiently rarefied to be neglected, such as in satellite thermal management. To further assess the performance of the thermal router under terrestrial conditions, we conducted experiments at room temperature (~300 K) under atmospheric pressure. In the experiment, a heat reservoir maintained at about 340 K was attached to the source port as the heat supply, while two reservoirs held at about 280 K were attached to drain 1 and 2 as heat sinks. A direct electric current was then applied to flow counter-clockwise through the thermal router (Figure 4a).

Display Full Size

Download

Figure 4. Experimental demonstration of thermal splitting effects. (a) Experimental setups; (b) The thermo diagram of experimental test at I = 5 A; (c) The extracted one-dimensional temperature distribution along channel 1, channel 2, and channel 3 of the thermal router when I = 0 A (upper panel) and I = 5 A (lower panel); (d) The dependence of split ratio on electric current.

The temperature distribution across the device was recorded using an external infrared thermal camera (Figure 4b). To highlight the effect of thermoelectric contributions, we extracted one-dimensional temperature profiles along the three transport channels under two conditions: zero current (upper panel, Figure 4c) and high current (lower panel, Figure 4c). Unlike the idealized case, the experimental profiles at zero current are neither linear nor monotonically convex, reflecting the limited thermal conductivity of Bi₂Te₃ and the influence of natural air convection. This suggests the presence of parasitic heat leakage through the device walls. Therefore, heat exchange with the surroundings must be included in both the theoretical and numerical models for realistic predictions.

To minimize accidental influences, three independent groups of experiments were performed on the device. The raw temperature data recorded by the thermal camera are provided in Section S8, together with the theoretical procedure for quantifying the thermal split ratio from experimental temperature fields. As presented in Figure 4d, when the applied current is relatively small (1-3 A), the measurement results remain highly stable, with nearly negligible variance. However, as the current increases to 4-5 A, significant deviations appear. This discrepancy mainly arises from overwhelming Joule heating at high current and nonlinear heat exchange with the environment at high temperature that combine together to increase the measurement uncertainty. In addition, when compared with the idealized theoretical predictions in Figure 3 from an overall perspective, the experimentally measured split ratio is notably lower. This is understandable because under room-temperature conditions, natural convection of ambient air indiscriminately heats or cools the device, introducing background thermal noise that obscures the thermoelectric contribution to thermal splitting and degrades experimental performance. A more detailed evaluation is provided in Section S9, where it is found that compared with thermal interface resistance, natural convection contributes to most of the performance degradation.

However, despite these environmental effects, the thermal router still achieves a split ratio of approximately 0.75 at 5 A. When these external influences are taken into account, the theoretical, simulated, and experimental results show good agreement, thereby confirming the robustness of the thermal routing functionality. To be noted, while it is confusing that the splitting effect can still emerge even when the temperature distribution appears nearly symmetric with respect to the source port, this counterintuitive behavior originates from the fact that, in thermoelectric transport, the total heat flux is not solely determined by the Fourier conduction term. Therefore, the temperature field alone cannot fully represent the heat-flux distribution, where the relationship between asymmetric temperature distribution and asymmetric heat flux transport can be decoupled to some extent. Finally, while a current of 5 A may seem impractically large for practical applications, the optimal operating current of the device is inherently scale-dependent. Its thermal routing capability is not fundamentally constrained by specific geometric scales or current ranges (Section S10).

4. Conclusion

To summarize, we have demonstrated a solid-state thermal router with both high power and strong routing contrast. The odd-parity changes of electrochemical potential in thermoelectric materials enable the system to break time-reversal symmetry in a mechanically stationary way while effectively directing heat from the source to a selected drain. Through analyzing the either synergistic or competitive roles of hybrid thermal effects, an intuitive yet in-depth physical insight have been provided for the robust operation and optimization of the metadevice^[45-49]. With the solid-state property and a directed heat flux several orders higher than previous solutions, it holds promise for the integrated path-selective heat treatment of advanced electronics, especially for applications where mechanical reliability and operation-mode programmability are required. However, the present work is a proof-of-concept demonstration rather than a product-level reliability test. Therefore, systematic long-term switching tests and package-level durability optimization, including contact resistance and thermal-stress management, are still required before practical deployment, which remain important topics for future device optimization.

Supplementary materials

The supplementary material for this article is available at: Supplementary materials.

Acknowledgments

The authors would like to acknowledge the support from the Research Platform for Energy and Environmental Nanotech, National University of Singapore (Suzhou) Research Institute.

Authors contribution

Ju R: Conceptualization, methodology, investigation, writing-original draft, writing-review & editing.

Wang Y: Experiment, validation.

Xu Z: Investigation.

Sun K: Data curation.

Yang S, Wang D, Qi M, Cao PC: Writing-review & editing.

Chen H, Qiu CW, Li Y: Conceptualization, writing-review & editing, resources, supervision, funding acquisition.

Conflicts of interest

The authors declare no conflicts of interest.

Ethical approval

Not applicable.

Consent to participate

Not applicable.

Consent for publication

Not applicable.

Availability of data and materials

The data supporting the findings can be found in supplementary information published with this paper, and obtained from the corresponding authors.

Funding

This study was supported by the National Natural Science Foundation of China (NNSFC) (Grant Nos. T2550093 and 12475040); the Zhejiang Provincial Natural Science Foundation of China (Grant Nos. LZ24A050002 and LR26A050001); the Ministry of Education, Republic of Singapore (Grant No. A-8002978-00-00); the National Research Foundation Singapore (NRF) under NRFs Medium Sized Centre: Singapore Hybrid-Integrated Next-Generation μ-Electronics (SHINE) Centre funding programme; the Science and Technology Project of Jiangsu Province (Grant No. BZ2022056).

Copyright

© The Author(s) 2026.