1. Introduction

Anomalous cooling or heating phenomena are always attractive as exploring the underlying physics can help our understanding of the universe and developing energy-related technologies, such as non-Fourier phonon transport^[1], ballistic-hydrodynamic-diffusion crossover^[2], nonreciprocal or directional thermal radiation^[3-5], to name a few. Among them, Mpemba effect describes a counterintuitive cooling phenomenon that hot water can be cooled and even frozen faster than warm water^[6,7]. In physics, it characterizes an anomalous relaxation process where a system initiated at high temperature relaxes faster than an identical system initiated at low temperature to an even colder state. Such a phenomenon was first reported by Aristotle about 2300 years ago but named after a high school student E. B. Mpemba in the 1960s when he surprisingly observed that hot water freezes faster than cold water^[6]. Since then, Mpemba effect has been frequently debated about its complex mechanism and occurrence conditions, and similar phenomena were reported in broad systems beyond water such as magnetic system^[8,9], colloidal system^[10-12], polymers^[13], granular fluid^[14], spin glasses^[15], quantum system^[16-21], and so on. Analogously, the inverse Mpemba effect can also be observed when a cold system can heat up faster than the identical system with warmer initial temperature to an even hotter state^[8,12,22]. There have been no widely accepted explanations for the Mpemba effect, let alone in different systems^[23-26]. The explanations for the Mpemba effect of water have been attributed to evaporation^[27], convection^[28], supercooling^[29], hydrogen bonding^[30], etc. While Mpemba effect in colloids has been explained by a double well potential that confines the movement of the colloidal particles with clear experimental verification^{[10-12,31,32]}. However, most of these explanations rely on microscopic frameworks, such as Markovian dynamics analyzed through probability distributions^[8-12,24-26] or thermomajorization theory^[33], and a macroscopic counterpart is still lacking.

Both Mpemba effect and its inverse involve three temperatures, T_h, T_w, and T_c to denote the initial hot, warm, and cold states of the system that may vary in different systems but are relatively easy to define for anomalous heating or cooling. But how to define the occurrence of the Mpemba effect has suffered from vague definitions^[33]. In most microscopic analysis involving Markovian dynamics, there exists a specific distance measure D with monotonicity and convexity to quantify the relaxation speed and how close toward the thermal equilibrium state π. In the beginning, $D({\pi }_0^h\to \pi )>D({\pi }_0^w\to \pi )$ must be satisfied but when a crossover $D({\pi }_t^h\to \pi )<D({\pi }_t^w\to \pi )$ happens at a specific time t, the Mpemba effect occurs^[33]. Such definition of occurrence of the Mpemba effect is understandable but how to choose the distance measure D is ambiguous though the total variation distance and the Kullback-Leibler (KL) divergence have been frequently adopted mathematically with ambiguous physical landscape^[8,9]. Therefore, we choose the specific time t as the criterion to define the occurrence of Mpemba effect if $t({\pi }_0^h\to \pi )<t({\pi }_0^w\to \pi )$^[11]. Such definition is unambiguous with clear temperature and time without involving variable system sizes, ingredients, impurities, etc. In the following calculation, we use such definition to analyze the occurrence of Mpemba effect and its inverse.

In this study, we report the Mpemba effect and its inverse and explain it by an inductor-resistor-capacitor (LRC) network that is composed of a typical thermoelectric module (TEM). The evolution of temperature and current can be pure damping or oscillatory, which is further demonstrated to be the key for anomalous heating and cooling, and enables the occurrence of Mpemba effect in p-type TEM and inverse Mpemba effect in n-type TEM. The dependences of Figure-of-Merit (ZT) and inductance (L) are analyzed with different initial temperature difference, which further helps sketch the occurrence domain of the Mpemba effect by taking either minimum t(π₀→π) as the criterion or taking the maximum T_b(t) against reservoir temperature T_r during the temperature oscillation as the more stringent criterion. Our results provide an alternative system to observe and understand the Mpemba effect and its inverse, particularly with a macroscopic perspective.

2. Methods

Considering a TEM in Peltier mode connects a body with temperature T_b and a thermal reservoir at temperature T_r, as shown in the subfigure in Figure 1b, and heat conduction occurs between the TEM and the thermal reservoir due to the initial temperature difference ∆T₀ = T_b0 - T_r. By ignoring the heat dissipation to the environment, we can build a 1-D model to characterize the temperature profile T(x) along the length of the TEM with a constant cross-sectional area according to Domenicali’s equation as^[34]

Display Full Size

Figure 1. Occurrence of Mpemba effect in LRC network consisted of Peltier element connecting a finite body at temperature T_b and a thermal reservoir at temperature T_r. (a) Evolution of temperature T_b with dimensionless time t/τ* for different initial temperatures ΔT₀ of 10 K and 50 K, respectively. Evolution of (b) dimensionless temperature (T_b(t) - T_r)/ΔT₀ and (c) dimensionless current I(t)/I₀ with t/τ*. t_min denotes the first time T_b approaches T_r, which is dependent on the initial temperature ΔT₀. The first minima and first maxima of T_b(t), T_bm1 and T_bm2, occur when the derivative of dT_b(t)/dt vanishes. The subfigures in (b) illustrate the system setup, and in (c) denote the LRC network consisting of a TEM with internal resistance R, thermal conductance k in a closed circuit with an ideal inductance L, and a finite body with thermal capacity C and the thermal reservoir at constant temperature T_r. LRC: inductor-resistor-capacitor; TEM: thermoelectric module.

(1)$\frac{\partial }{\partial x}\left(\kappa \frac{\partial T}{\partial x}\right)=-I^{2}R+IT\frac{\partial S}{\partial x}$

where I is the electric current, and the temperature-dependent parameters κ(T), R(T), and S(T) are thermal conductivity, electrical resistivity, and combined Seebeck coefficient of a typical TEM, respectively. According to Kelvin’s relation, the Peltier heating or cooling power could be characterized by the Peltier coefficient π at a rate of $\dot{Q}=\Pi I=S\left(T\right)TI\left(t\right)$ depending on the direction of electric current. At the same time, heat will be conducted irreversibly between the TEM and the reservoir spontaneously. Then, the total heat flow per unit area along the TEM could be calculated as $q\left(x\right)=\dot{Q}-\kappa \partial T/\partial x$. Combining with thermal boundary conditions T(x = 0) = T_r, T(x = l) = T_b, and assuming the Joule heating power RI² is equally dissipated to both sides of the TEM, resulting in

(2)${\dot{Q}}_b=q\left(l\right)=-S{T}_{b}I+\frac{1}{2}R{I}^2-k(T_b-T_r)$

(3)${\dot{Q}}_r=q\left(0\right)=S{T}_{r}I+\frac{1}{2}R{I}^2+k(T_b-T_r)$

where R is the internal resistance of the TEM and k is the thermal conductance between the TEM and the body and the thermal reservoir. Similarly, the electric profile along the TEM can also be calculated with a combination of the reversible Seebeck effect and the irreversible Ohm’s law, yielding

(4)$L\dot{I}+RI=S(T_b-T_r)$

where L is the TEM electrical inductance, without which Equation (4) has been a standard protocol to describe the heat flow and electric current in a typical TEM^[35,36]. It will be demonstrated later that the introduction of the inductance L is the key in this study, and the influence of either neglecting L or not will be discussed later. With L, the system can be considered as an LRC network, as shown in the subfigure in Figure 1c. Considering the thermal reservoir is infinite with a constant T_r, and the body to be connected through a TEM with the reservoir has a finite temperature T_b and a finite heat capacity C and assuming ${\dot{Q}}_b=C{\dot{T}}_b$, we can obtain a damped harmonic oscillator equation for I(t) as

(5)$LC\ddot{I}+(RC+kL)\dot{I}+(kR+S^2{T}_r)I\approx 0$

after neglecting the Joule heating power RI² and the magnetic energy $LI\dot{I}$ as the summation of these two is rather smaller than S²T_rI, which is discussed in details in Supporting Information. The solution type to Equation (5) is dependent on the corresponding characteristic equation LCx² + (RC + kL)x + kR + S²T_r = 0, especially the sign of Δ = (RC - kL)² - 4LCS²T_r. When Δ ≥ 0, the solution takes the form of I(t) = I₀[exp(-t/τ₁) - exp(-t/τ₂)], where τ₁ and τ₂ are the two damping constants from two distinct solutions of the corresponding characteristic equation as ${\tau }_1=2LC/(RC+kL-\sqrt{\Delta })$ and ${\tau }_2=2LC/(RC+kL+\sqrt{\Delta })$. When ∆ < 0, the solution of interest is I(t) = I₀[exp(-t/τ)sinωt, where τ and ω are the damping and oscillatory constants respectively from the complex solutions of the corresponding characteristic equation. The oscillatory behavior will be demonstrated to enable the overshooting of T_b(t) below T_r, otherwise there is only a damping behavior from the initial T_b(0) to the thermally-balanced T_r eventually. When the TEM is thermally connected with the body and the reservoir, we build the initial conditions for t = 0 as I(0) = 0, T_b(0) - T_r = ∆T₀, and $\dot{I}\left(0\right)=S\Delta T_0/L$. Based on these initial conditions, I₀ will be determined thereafter. For the oscillatory case, τ = 2LC/(RC + kL) and $\omega =\sqrt{\left|\Delta \right|}/2LC$. Note that when L = L* = RC/k, ∆ < 0 will always be fulfilled for any value of S, and corresponding time constant τ* and ω* will be determined thereby as τ* = C/k and ${\omega }^{\ast }=\sqrt{S^2{T}_{r}k/\left(R{C}^2\right)}$, which correspond to the weakest possible damping of I(t) with the maximum of oscillatory frequency ω_max = ω*. These critical values will be specifically discussed later by introducing the standard definition of the dimensionless figure of merit of a TEM at T = T_r as ZT = S²_rT_r/kR = (ω*τ*)^2[35]. The evolution of T_b(t) is therefore as

(6)$T_b\left(t\right)=\{\begin{array}{cc}{T}_r+I_0[(R-L/{\tau }_1)e^{-t/{\tau }_1}-(R-L/{\tau }_2)e^{-t/{\tau }_2}]/S,& \Delta >0\\ T_r+\Delta T_0{e}^{-t/\tau }\cos (\omega t-\delta )/\cos \delta ,& \Delta <0\end{array}$

With tanδ = (R - L/τ)/Lω. In our following calculations, we keep R, k, and C as constant as 0.22 Ω, 0.0318 W/K, and 4.96 J/K, respectively, which are drawn from the experiments^[35].

3. Results and Discussion

Now we achieve the Mpemba effect based on such setup with a p-type TEM in Figure 1. In such a setup, holes will be transferred from the high-temperature body T_b(t) to the low-temperature reservoir T_r, yielding a positive I(t) and S. When ∆ < 0, the body temperature T_b(t) evolutes with both oscillatory and pure damping behaviors. It is seen in Figure 1a that T_b(t) drops first until achieving the minimum T_bm1 even below the reservoir temperature T_r, and then increases to the maximum T_bm2 above T_r, and eventually approaches thermal equilibrium with T_r. Such oscillatory behavior depends on the initial temperature difference ΔT₀. It is seen in Figure 1a that when ΔT₀ = 50 K, T_b(t) drops faster to T_r than that with ΔT₀ = 10 K, which implies a hot object cools faster than the warm one to the thermally balanced state, just as the Mpemba effect states. Such phenomenon can be observed more clearly with the dimensionless temperature (T_b(t) - T_r)/ΔT₀ under the dimensionless time t/τ* in Figure 1b. The shortest time to achieving T_r, defined as t_min, can be characterized by t_min = (π/2 + δ)/ω ≈ π/2ω, will be influenced by the initial temperature difference ΔT₀, temperature-dependent S, ZT, and L. Such t_min can be used to define the criterion of occurrence of Mpemba effect. It is seen that t_min @ΔT₀ = 50 K is smaller than t_min @ΔT₀ = 10 K, which validates the occurrence of the Mpemba effect. The oscillatory behavior can also be observed from the evolution of dimensionless current I(t) in Figure 1c. It is seen that I(t) increases firstly and then drops gradually and even becomes negative, which means the current in the LRC circuit is inversed, and finally I(t) vanishes in thermal equilibrium state.

The Mpemba effect is influenced by both ZT and L, as illustrated in Figure 2 and Figure 3. It is seen in Figure 2 that with the increase of ZT, T_b(t) drops faster with larger derivative of dT_b(t)/dt. Each time when T_b(t) approaches T_r, the dimensionless temperature (T_b(t) - T_r)/ΔT₀ changes its sign, but the excess temperature over T_r become smaller and smaller, just as |T_bm1 - T_r| > |T_bm2 - T_r|. The red dash line in Figure 2a describes a corresponding pure damping relaxation e^-t/τ* with a time constant τ* without oscillation. t_min, a criterion that characterizes the Mpemba effect, is also marked. It is seen that with the increase of ZT, t_min becomes smaller, implying the earlier occurrence of Mpemba effect. Meanwhile, both |T_bm1 - T_r| and |T_bm2 - T_r| become larger, implying the stronger oscillation of T_b(t) in the LRC network. When L = L*, the dimensionless time $t_{min}/{\tau }^{\ast }=\pi /2/\omega \ast \tau \ast =\pi /2/\sqrt{ZT}$, validating the ZT-dependent t_min. It is also seen in the corresponding subfigure in Figure 2a that with the increase of ZT, t_min becomes smaller significantly. Meanwhile t_min is slightly dependent on L. With the increase of L, t_min becomes larger. The evolution of I(t) in Figure 2b shows a similar oscillatory behavior as (T_b(t) - T_r)/ΔT₀. Note that I(t) is highly dependent on ZT and may become negative, implying the direction of current may be inversed. The oscillatory current I(t) is generated by the change of voltage difference at the two ends of the TEM due to the temperature difference ΔT = T_b(t) - T_r, which ultimately induces the voltage difference $L\dot{I}$ through the TEM. Figure 3 shows the L-dependent Mpemba effect. When L ≈ 0, there is no inductance in the LRC network, and thus T_b(t) and I(t) are not oscillatory but just a damping process, which are described by Equation (6). In this case, Δ > 0 will always be fulfilled and there will be no Mpemba effect, as seen in Figure 3a that there is no intersection point between the curves of T_b(t) and T_r. When L = L* and L = 2L*, the oscillatory T_b(t) and I(t) will be observed with t_min marked. It is seen that with the increase of L, t_min becomes larger, verifying the L-dependent Mpemba effect and echoing the subfigure in Figure 2a.

Display Full Size

Figure 2. ZT dependence of Mpemba effect. Evolution of (a) dimensionless temperature (T_b(t) - T_r)/ΔT₀ and (b) dimensionless current I(t)/I₀ with t/τ* under varying ZT between 0.25 to 2.5 in steps of 0.25. The red dash line in (a) denotes a corresponding damping relaxation with a time constant τ* without oscillation. The subfigure in (a) describes the ZT and L dependence of t_min/τ*.

Display Full Size

Figure 3. L dependence of Mpemba effect. Evolution of (a) dimensionless temperature (T_b(t) - T_r)/ΔT₀ and (b) dimensionless current I(t)/I₀ with t/τ* under varying L.

Then one may question about the occurrence domain of Mpemba effect in such LRC system. As shown in Figure 4, we plot the domain of Mpemba effect by varying ZT, L, and ΔT₀. If characterized by t_min, Mpemba effect happens when t_min@ΔT is larger than t_min @(ΔT + ξ), where ξ > 0 is a small temperature increment. It is seen that when ΔT₀ and L tend to be smaller, more cases are pure damping behavior with normal cooling. On the contrary, with the increase of ΔT₀ and L, the domain of Mpemba effect (anomalous cooling) becomes larger, as shown in the enhanced area of blue dots. Fundamentally, it is the oscillatory behavior that enables the occurrence of anomalous Mpemba effect, which will also enable the heat flow from cold temperature to the high temperature temporarily with both theoretically and experimentally demonstrations^[35]. However, one may notice that T_bm1 is below T_r but T_bm2 may be above T_r, as shown in Supporting Information, and question about the definition of Mpenba effect when T_bm2 > T_r. So, we can further define more stringent occurrence criteria of Mpemba effect by both t_min@ΔT > t_min@(ΔT + ε) and T_bm2 < T_r which is more in line with the physical definition of cooling. With such stringent definition, the domain of Mpemba effect is reduced as shown that only a small number of blue dots can turn into red stars in Figure 4c,d. The corresponding T_b(t) and I(t) of black dot, blue dot, and red star can be found in Figure 1, Figure 2, and Figure 3 accordingly. More discussions on the occurrence domain of Mpemba effect with L ≈ 0 and L = L* can be found in Supporting Information.

Display Full Size

Figure 4. Occurrence domain of Mpemba effect. (a-b) Domain of pure damping (black dot) and Mpemba effect (blue dot) characterized by t_min if t_min @ΔT is larger than t_min @(ΔT + ξ); (c-d) Domain of pure damping (black dot) and Mpemba effect characterized by T_bm2 with T_bm2 > T_r (blue dot) and T_bm2 < T_r (red star). All subfigures are with varying ZT between 0.25 to 2.5 in steps of 0.25, and varying ΔT₀ between 1 to 100 in steps of 11. The given electrical inductance L is 1 (a,c) and 2 (b,d), respectively.

With thermoelectric LRC network, heat flow becomes oscillatory, which is the key to achieving Mpemba effect in this study fundamentally speaking. Such processes seem anomalous but do not violate the second law of thermodynamics^[35,37]. It can be understood by considering two successive reversible combined power-refrigeration cycles with four stages in a full-period of oscillatory cycle of T_b(t) alternatively and the TEM acts as a thermoelectric generator (stage 1), a cooler (stage 2), a generator (stage 3), and a heater (stage 4), respectively. In the thermoelectric generator stage 1, the TEM generates electric power by absorbing heat Q₁ from the body and storing the energy 0.5LI(t)² in the inductor during the heat conduction from the body to the reservoir until Tb(t) approaches T_r descendingly. In the thermoelectric cooler stage 2, the TEM is driven by the accumulated inductive energy to continue absorbing heat Q₂ from the body to the reservoir with its temperature from T_r to T_bm1 descendingly until the stored energy is consumed. In these two stages, the reservoir T_r acts as the heat sink. Due to the oscillatory behavior, in the successive thermoelectric generator stage 3, the body T_b(t) acts as the heat sink in reverse and the TEM absorbs heat Q₃ from the reservoir and stores energy 0.5LI(t)² in the inductor again until T_b(t) approaches T_r from T_bm1 to T_r ascendingly, and continues absorbing heat Q₄ from the reservoir with temperature from T_r to T_bm2 ascendingly in thermoelectric heater stage 4. These four stages constitute a full period of T_b(t) with two power-refrigeration cycles, and the transformations of Q₁ from T_b0 to T_r and Q₄ from T_bm2 to T_r are positive, while the transformations of Q₂ and Q₃ for T_bm1 to T_r are negative, respectively. Then, neglecting the energy loss and according to the theorem of equivalence of transformation^[37], we obtain

(7)$\{\begin{array}{c}{\int }_{T_{b0}}^{T_r}C(\frac{1}{T_r}-\frac{1}{T_b})d{T}_b+{\int }_{T_r}^{T_{bm1}}C(\frac{1}{T_r}-\frac{1}{T_b})d{T}_b\approx 0\\ {\int }_{T_{bm1}}^{T_r}C(\frac{1}{T_r}-\frac{1}{T_b})d{T}_b+{\int }_{T_r}^{T_{bm2}}C(\frac{1}{T_r}-\frac{1}{T_b})d{T}_b\approx 0\end{array}$

yielding

(8)$T_r=\frac{T_{b0}-T_{bm1}}{\ln T_{b0}-\ln T_{bm1}}=\frac{T_{bm2}-T_{bm1}}{\ln T_{bm2}-\ln T_{bm1}}$

which implicitly gives the relationship between T_bm1 and T_bm2. From Equation (6), it is seen that larger T_b(0) and lower T_r can give lower T_bm1 and larger T_bm2, implying larger positive transformation of Q₁ and Q₄ to compensate for the larger negative transformation of Q₂ and Q₃, respectively.

Finally, let us consider an n-type TEM in the LRC network with T_r higher than T_b(0) inversely. In such setup, electrons will be transferred from the high-temperature reservoir T_r to the low-temperature body T_b(t), yielding a negative I(t) and S. It is seen in Figure 5 that both T_b(t) and the dimensionless temperature (T_b(t) - T_r)/ΔT₀ exhibit oscillatory behaviors with the dimensionless time t/τ*. When ΔT₀ = -50 K, T_b(t) increases faster to T_r than that for ΔT₀ = -10 K, which implies a colder object can heat up faster than the warm one to the thermally balanced state, just as the inverse Mpemba effect states. Such anomalous heating phenomena can also be characterized by t_min with smaller t_min @ΔT₀ = -50 K than t_min @ΔT₀ = -10 K, and the oscillatory behaviors in the LRC network can be also attributed to it.

Display Full Size

Figure 5. Occurrence of inverse Mpemba effect in LRC network. (a) Evolution of temperature T_b with dimensionless time t/τ* with different initial temperature ΔT₀ of -10 K and -50 K, respectively; (b) Evolution of dimensionless temperature (T_b(t) - T_r)/ΔT₀ with t/τ*. t_min denotes the first time T_b approaches T_r, which is dependent on the initial temperature ΔT₀. The subfigures in (b) illustrate the system setup. LRC: inductor-resistor-capacitor.

4. Conclusion

Through a TEM connecting a body and a reservoir at different temperatures to constitute a LRC network, we demonstrate the pure damping and oscillatory behaviors of both temperature T_b(t) and current I(t), which correspond to the normal and anomalous heating/cooling phenomena. Such a system is built from the macroscopic level rather than the traditional microscopic level, establishing an alternative framework for understanding the Mpemba effect for the first time, complementing existing microscopic studies. Further, we analyze the ZT- and L-dependent occurrence of the Mpemba effect with full consideration of the initial temperature ΔT₀, revealing the occurrence domain of the Mpemba effect by two criteria: t_min and T_bm2 < T_r. Our study not only deepens the understanding of the Mpemba effect and its inverse, but also offers an alternative system to analyze the anomalous heating and cooling phenomena from the macroscopic level. Further experiments may be conducted with very efficient TEM, resistanceless superconducting coils, active gyrator-type circuits, etc.^[35]. With such a macroscopic system, more flexible and dynamic ways for thermal management and energy conversion can be achieved by taking advantage of the oscillatory temperature, current, and heat flow.

Acknowledgments

Deepseek R1 was used solely for language polishing. The authors take full responsibility for the final manuscript.

Authors contribution

Wang Z: Conceptualization, methodology, investigation, writing-original draft, writing review & editing.

He J: Conceptualization, methodology, investigation.

Kim H: Validation, methodology.

Choi: W: Investigation, writing review & editing.

Hu R: Conceptualization, validation, methodology, writing-original draft, writing review & editing.

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 and materials could be obtained from the corresponding author upon request.

Funding

This work is supported by the National Key Research and Development Program of China (Garnt Nos. 2024YFB4104701, 2022YFA1203104), the National Natural Science Foundation of China (Garnt Nos. 52422603, 52511540065, 52521008, W2521166 and 92463311), the Science and Technology Program of Hubei Province (Garnt No. 2023AFA072), the Interdisciplinary Research Program of HUST (Garnt No. 5003120094), the Shenzhen Technology Project (Garnt No. JCYJ20241202123700001), and the Open Research Fund of Suzhou Laboratory (Garnt No. SZLAB-1508-2024-TS016).

Copyright

© The Author(s) 2026.