1. Introduction

Hydrodynamics phenomena originate from the strong interactions of microscopic (quasi)particles, where the particle number, (quasi)momentum and energy are conserved and a local equilibrium exists during the scattering process^[1-5]. For example, fluid hydrodynamic phenomena are ubiquitous in daily life, such as turbulence, vortices, Poiseuille flow, wave and so on. These universal hydrodynamic behaviors do not depend on the types of quasiparticles and their microscopic interaction^[2,5]. Theoretically, they can be described by the mesoscopic Boltzmann transport equation (BTE)^[6-10] or macroscopic hydrodynamic equations in the continuum limit^[11-13].

However, hydrodynamic behaviors of electrons in crystal solids at room temperature are not easy to be observed because the crystal momentum is usually not conserved during the electron scattering processes^[1,2]. In the 1960s, Gurzhi^[14,15] theoretically proposed that, as a signature of electron hydrodynamics, electrical resistance of a conductor within certain size and temperature ranges could decrease when the temperature increases. The electron hydrodynamic phenomena can be predicted in a certain low temperature window when the momentum-conserving (MC) electron-electron scattering process is much more sufficient than the other momentum-relaxing (MR) scattering processes, including scattering with impurities, phonons and so on. This pioneering work opened the door to electron hydrodynamics in solid materials, which was expected to significantly modify or improve the electric and thermal transport coefficients^[16]. Although electron hydrodynamics have been studied theoretically for over half a century^[17-19], only a little experimental progress has been made in the 20th century^[17,20].

In the past decade, the study of electron hydrodynamics enjoyed a renaissance^[1,2,21,22] due to the discovery of new materials, e.g., graphene^[9,23-25], in which the electron-phonon coupling is weak and the MC scattering process is much stronger. Many electron hydrodynamic phenomena have been measured^[26-32], such as the electron Poiseuille flow^[10], negative nonlocal resistance or vortices/whirlpool^[33-38], Hall viscosity^[39], violation of the Wiedemann-Franz law^[11,40], super-ballistic flow^[41] and so on. The difference in macroscopic behaviors between hydrodynamic, ballistic and diffusive electron transport is widely studied in various materials or geometries, both theoretically^[8,42-45] and experimentally^{[10,32,33,35-37,41]}. Levitov and Falkovich used the linearized electronic Navier-Stokes equation to study the current vortices and viscosity^[44]. Results show that vortices could appear in the hydrodynamic regime while they disappear in the diffusive regime. Similar current vortices, viscosity or whirlpool phenomena have also been predicted in various rectangles^{[12,38,46-48]} or Corbino disk geometries^[4,49-51] by applying voltage in different positions to control current flow. For instance, Shavit et al.^[13] reported that the nonlocal relation between the current and electric field due to momentum-conserving interparticle collisions leads to a total or partial field expulsion from such flows, which results in freely flowing currents in the bulk and a boundary jump in the electric potential at current-injecting electrodes. Li et al. obtained thermoelectric coefficients of the system in the crossover region between charge neutrality and high electron density regime^[52]. The thermal conductance exhibits a sensitive Lorentzian dependence on the electron density. Gall et al. found that local temperature and electric potential are discontinuous at the interfaces with the leads as well as the device resistance in neutral graphene^[53]. They also reported the results of a comprehensive study of the interplay of viscosity, disorder-induced scattering, recombination, energy relaxation, and interface-induced dissipation based on the full consistent hydrodynamic description^[50].

Macroscopic models have been widely used in the above studies and made a success in capturing the qualitative and, in many cases, quantitative features of the electron hydrodynamics. However, these macroscopic approaches rely on the continuum or local equilibrium assumption as well as classical constitutive relations, such as the Navier-Stokes or Ohm-Stokes equations. They cannot access the full momentum-resolved information encoded in the electron distribution function^[7,54], nor can they seamlessly interpolate between the diffusive and hydrodynamic regimes without invoking additional assumptions^[12,44,55]. A mesoscopic kinetic approach, based on the BTE, offers a complementary and more fundamental perspective^[2,8-10,54]. It explicitly describes the evolution of the electron distribution function under electric and magnetic fields, a temperature gradient, MR/MC scattering processes and boundary conditions, providing direct insight into the microscopic origins of macroscopic transport phenomena. Such an approach naturally encompasses all transport regimes within a unified numerical framework, without requiring a priori assumptions about the validity of hydrodynamic or continuous descriptions. In addition, compared to the electric behaviors of hydrodynamic electron transport^[1,2,44,56], its thermal properties have received less attention^{[11,40,52,57,58]}. Actually, apart from the electrical properties, including mobility and conductivity^[17,20], thermal performance is also critical for the thermal management problem^[57,59,60] in semiconductor devices. Electrical and thermal transport affect each other^{[11,40,52,58]}; therefore, it is necessary to simultaneously study electrical and thermal properties of hydrodynamic electron transport.

In this study, the thermal behaviors of hydrodynamic electron transport in a homogeneous Corbino disk geometry under the magnetic field are studied. The electron BTE with a dual-relaxation-time Callaway model^[2,8-10,54] is adopted to systematically study electron transport in a graphene Corbino disk. An implicit discrete ordinate method is developed to solve it numerically across the hydrodynamic and diffusive regimes. The Newton method is used to solve the nonlinear relationships between the Fermi-Dirac distribution, chemical potential and temperature. The competition effects on hydrodynamic transport among magnetic field, electric field and temperature gradient are investigated within different scattering rates. Electric and thermal properties of both diffusive and hydrodynamic electron transport in the spatial domain are compared and discussed. The rest of the paper is organized as follows. In Section 2, the kinetic model and associated numerical methods are introduced. Then a systematical discussion is conducted in Section 3. Finally, a conclusion is made in Section 4.

2. Model Equation and Method

2.1 Electron boltzmann transport equation

Electron transport at steady state can be described by the electron Boltzmann transport equation (eBTE) in the semiclassical limit^[2,6,8-10],

(1)$v\cdot \frac{\partial f}{\partial x}+F\cdot \frac{\partial f}{\partial p}=-\frac{f-f_{mc}^{eq}}{{\tau }_{mc}}-\frac{f-f_{mr}^{eq}}{{\tau }_{mr}}$

which describes the evolution of the electron distribution function f with the spatial position x and momentum $p=ħk$ under the external force F = e(E + v × B) of electric field E = -∇φ and magnetic field B, where e = -1.6 × 10^-19 C is the electron charge, $ħ$ is the reduced Planck constant and φ is the electric potential. Electron group velocity v = ∂ε/∂p depends on the electronic energy band structure ε = ε(p). The actual scattering process is extremely complex, involving interactions between electrons-electrons, electrons-phonons, and electrons-impurities, etc. To simplify the calculation, the classic double relaxation time Callaway model is adopted and the phonon transport or scattering mechanisms are ignored in this paper. The Callaway model has played a significant role in previous studies on the hydrodynamic electron transport, which splits the entire scattering mechanism into two parts. One part is the MC scattering process that satisfies momentum conservation, and the other part is the MR scattering process that does not satisfy momentum conservation, including defect, isotope, disorder, and impurity scattering. Both of these scattering processes satisfy the conservation principles of energy and particle number density, and are expressed by the relaxation time approximation with τ_mc and τ_mr respectively. The MR scattering relaxes the electron distribution to the Fermi-Dirac distribution^[61],

(2)$f_{mr}^{eq}(\epsilon ,\mu ,T)=\frac{1}{\exp \left(\frac{\epsilon -\mu }{k_{B}T}\right)+1}$

where T is the temperature, μ is the chemical potential, and k_B is the Boltzmann constant. The equilibrium state of the MC scattering process is

(3)$f_{mc}^{eq}(\epsilon ,\mu ,u,T)=\frac{1}{\exp \left(\frac{\epsilon -\mu -p\cdot u}{k_{B}T}\right)+1}$

where u is the macroscopic drift velocity.

In this paper, we assume that the variance of the electron number density throughout the entire system is very small, and the spatial distribution of the electric potential is almost determined by the boundary conditions and the Laplace equation ∇²φ ≈ 0. We adopt a relative reference frame based on the conduction band bottom. In this way, the band does not change with the spatial position. The chemical potential μ is determined by the electron number density, while the electric potential φ determines the spatial distributions of the electric field. The combination of these two can yield the electrochemical potential, which will not be extensively discussed in this article.

Macroscopic variables, including particle number density n, electric current j, energy density U and heat flux q, can be obtained by taking the momentum of the distribution function in the framework of the kinetic model,

(4)$n=⟨f⟩,j=e⟨vf⟩,U=⟨\epsilon f⟩,q=⟨v(\epsilon -\mu )f⟩$

where ⟨⟩ represents the integral over the whole momentum space. The local temperature T and chemical potential μ could be updated by assuming a local equivalent equilibrium,

(5)$⟨f⟩=⟨f_{mr}^{eq}(\mu ,T)⟩,⟨\epsilon f⟩=⟨\epsilon f_{mr}^{eq}(\mu ,T)⟩$

Particle number density and energy density are conserved for both MC and MR scattering processes, and momentum is also conserved for the MC scattering process so that

(6)$⟨\frac{f-f_{mr}^{eq}({\mu }_{mr},T_{mr})}{{\tau }_{mr}\left(T\right)}⟩=⟨\epsilon \frac{f-f_{mr}^{eq}({\mu }_{mr},T_{mr})}{{\tau }_{mr}\left(T\right)}⟩=0$

(7)$⟨\frac{f-f_{mc}^{eq}({\mu }_{mc},T_{mc},u)}{{\tau }_{mc}\left(T\right)}⟩=⟨\epsilon \frac{f-f_{mc}^{eq}({\mu }_{mc},T_{mc},u)}{{\tau }_{mc}\left(T\right)}⟩=⟨p\frac{f-f_{mc}^{eq}({\mu }_{mc},T_{mc},u)}{{\tau }_{mc}\left(T\right)}⟩=0$

where μ_mr, T_mr, μ_mc, T_mc are invoked to ensure the conservation principle of scattering processes^[8]. When τ_mc and τ_mr are independent of momentum, T_mc = T_mr = T, μ_mr = μ_mc = μ. The Newton method is used for the above nonlinear equations.

2.2 Implicit discrete ordinate method

To solve the stationary eBTE iteratively, a classical implicit discrete ordinate method is used, where the semi-implicit scheme is used for the scattering terms and the fully-implicit scheme is used for the gradient of the distribution function in both the spatial and momentum spaces,

(8)$v\cdot \frac{\partial f^{m+1}}{\partial x}+F\cdot \frac{\partial f^{m+1}}{\partial p}=-\frac{f^{m+1}-f_{mc}^{eq,m}}{{\tau }_{mc}^m}-\frac{f^{m+1}-f_{mr}^{eq,m}}{{\tau }_{mr}^m}$

where m is the iteration step. The above equation is reformulated into

(9)$\frac{\Delta f^{m+1}}{{\tau }_{mc}^m}+\frac{\Delta f^{m+1}}{{\tau }_{mr}^m}+v\cdot \frac{\partial \Delta f^{m+1}}{\partial x}+F\cdot \frac{\partial \Delta f^{m+1}}{\partial p}=re{s}^m$

where Δf^{m+ 1} = f^{m+ 1} - f^m is the increment of the distribution function between two adjacent iteration steps and the mesoscopic residual res is defined as

(10)$re{s}^m=-\frac{f^m-f_{mc}^{eq,m}}{{\tau }_{mc}^m}-\frac{f^m-f_{mr}^{eq,m}}{{\tau }_{mr}^m}-v\cdot \frac{\partial f^m}{\partial x}-F\cdot \frac{\partial f^m}{\partial p}$

when res→0, the iteration converges and Δf→0. In other words, the formula of the left hand side of Eq. (9) does not influence the final convergent results. Therefore, the inexact Newton method^[62] is used to solve Eq. (9). A first-order upwind scheme is used for the numerical discretization of the left hand side of Eq. (9) and a higher and more accurate numerical discretization method is used for the right.

The whole phase space is discretized into a lot of small pieces, for example, f_i,k, where i and k are the indexes of the discretized spatial position and momentum, respectively. The finite volume method is used to discretize the spatial space

(11)$v\cdot {\left(\frac{\partial f}{\partial x}\right)}_i=\frac{1}{V_i}\sum _{j\in N_i}{S}_{ij}v\cdot n_{ij}{f}_{ij}$

(12)$v\cdot {\left(\frac{\partial \Delta f}{\partial x}\right)}_i=\frac{1}{V_i}\sum _{j\in N_i}{S}_{ij}v\cdot n_{ij}\Delta f_{ij}$

where V_i is the volume of the cell i in the x space, N(i) denotes the sets of face neighbor cell of cell i, ij denotes the interface between cell i and cell j, S_ij is the area of the interface ij, n_ij is the unit normal vector of the interface ij directing from cell i to cell j. For the distribution function at the cell interface,

(13)$v\cdot n_{ij}{f}_{ij}=a_1(f_i+{\sigma }_i(x_{ij}-x_i))+a_2(f_j+{\sigma }_j(x_{ij}-x_j))$

(14)$v\cdot n_{ij}\Delta f_{ij}=a_1\Delta f_i+a_2\Delta f_j$

where a1 = 0.5(n_ij · v + |n_ij · v|), a2 = 0.5(n_ij · v - |n_ij · v|), σ_i is the gradient of the distribution function in the cell i calculated by the van Leer limiter^[63], x_ij is the center of interface ij. A similar strategy is implemented for the numerical discretization of the momentum space,

(15)$F\cdot {\left(\frac{\partial f}{\partial p}\right)}_k=\frac{1}{V_k}\sum _{l\in N_k}{S}_{kl}F\cdot n_{kl}{f}_{kl}$

(16)$F\cdot {\left(\frac{\partial \Delta f}{\partial p}\right)}_k=\frac{1}{V_k}\sum _{l\in N_k}{S}_{kl}F\cdot n_{kl}\Delta f_{kl}$

where V_k is the volume of the cell k in the momentum space, N(k) denotes the sets of face neighbor cell of cells k, kl denotes the interface between cell k and cell l, S_kl is the area of the interface kl, n_kl is the unit normal vector of the interface kl directing from cell k to cell l. For the distribution function at the cell interface,

(17)$F\cdot n_{kl}{f}_{kl}=b_1(f_k+{\sigma }_k(p_{kl}-p_k))+b_2(f_l+{\sigma }_l(p_{kl}-p_l))$

(18)$F\cdot n_{kl}\Delta f_{kl}=b_1\Delta f_k+b_2\Delta f_l$

where b₁ = 0.5(n_kl · F + |n_kl · F|), b₂ = 0.5(n_kl · F - |n_kl · F|), σ_k is the gradient of the distribution function in the cell k and p_kl is the center of interface kl. Under the discretized spatial and momentum spaces, Eq. (9) becomes

where A is a huge coefficient matrix for the discretized phase space. The aforementioned large-scale linear equation system is solved by the lower-upper symmetric-Gauss-Seidel method^[64].

In the hydrodynamic regime, $f\to f_{mc}^{eq}$, hence such approximation can be used

(21)$F\cdot \frac{\partial f}{\partial p}\approx F\cdot \frac{\partial f_{mc}^{eq}}{\partial p}=f_{mc}^{eq}(f_{mc}^{eq}-1)F\cdot \frac{v-u}{k_{B}T}$

Then Eq.(19) becomes

The above approximations can also be used in the diffusive regime with $f_{mc}^{eq}\to f_{mr}^{eq}$ and u→0. The classical implicit discrete ordinate method is not limited by the electron energy band structure or scattering process. Electron transport in various materials can be described by this numerical iterative framework as long as the model equation is valid and the input parameters can be obtained regardless of whether they are obtained from first-principle calculations or empirical formulas.

In this paper, the electron transport and thermal behaviors in single-layer suspended graphene materials are investigated. The basic physical parameters are as follows. An isotropic linear dispersion is assumed, i.e., k = |k|s and $\epsilon =ħ|k|v_F$, where s = (cosα, sinα) is the unit directional vector. The integral over the momentum space in graphene materials can be reformulated into^[8],

(24)$⟨M⟩=\frac{4}{(2\pi \overline{h}{)}^2}\int Mdp=\frac{1}{2\pi }\iint \frac{2\epsilon }{\pi {\left(\overline{h}{v}_F\right)}^2}Md\epsilon d\Omega$

where dΩ is the integral over the solid angle space. The doping concentration is n_D = 10¹² cm^-2. Corresponding Fermi velocity and Fermi energy are v_F = 10⁶ m/s^[8] and $E_F=ħv_F\sqrt{\pi n_D}$^[65], respectively. The discrete details of the momentum space are as follows. The solid angle α ∈ [0, 2π] is discretized equally into N_α = 20 - 40 pieces in the numerical simulations. Electron energy is discretized equally into N_B = 20 - 40 pieces in the range (a, E_F + N_kk_BT₀), where a is the maximum value between 0 and E_F - N_kk_BT₀, T₀ is the background temperature. N_k = 10 in this simulation. M = 201 uniform cells are used for the spatial discretization in each direction, and the van Leer limiter^[63] is used for the spatial gradient of the distribution function. Numerical validations and grid independence tests have been conducted in Supplementary material. Results show that the present scheme and discretizations are acceptable to ensure the numerical simulation accuracy.

2.3 Boundary conditions

The isothermal boundary condition means that all particles hitting the boundary are absorbed^[6], and particles emitting from the boundary to the computational domain follow the Fermi-Dirac equilibrium distribution with boundary temperature T_b and chemical potential μ_b,

(25)$f\left(x_b\right)=f_{mr}^{eq}(x_b,{\mu }_b,T_b),s\cdot n_b>0$

where x_b represents the spatial coordinates of the boundary, n_b refers to the unit normal vector that points from the boundary to the computational domain. The periodic boundary condition means that when an electron leaves one boundary, and meanwhile another electron with the same velocity and energy enters the computational domain from the associated periodic boundary. Each electron distribution deviates from its equilibrium state by the same amount at the corresponding periodic boundaries

(26)$$ f(x_{b1},\epsilon ,s)-f_{mr}^{eq}(T_{b1},{\mu }_{b1},\epsilon )=f(x_{b2},\epsilon ,s)-f_{mr}^{eq}(T_{b2},{\mu }_{b2},\epsilon )$$

where b1 and b2 denote the corresponding periodic boundaries, respectively. The diffusely reflecting boundary condition indicates that the net energy and particle number that penetrate the boundary interface are zero. Besides, electrons reflected from the boundary are isotropic,

(27)${\int }_{s\cdot n_b<0}vfdp={\int }_{s\cdot n_b>0}v{f}_{mr}^{eq}({\mu }^{\prime },T^{\prime })dp$

(28)$$ {\int }_{s\cdot n_b<0}v\epsilon fdp={\int }_{s\cdot n_b>0}v\epsilon f_{mr}^{eq}({\mu }^{\prime },T^{\prime })dp$$

The unknown temperature T’ and chemical potential μ’ at the boundaries can be obtained by solving the above two equations with the Newton method. The specular reflecting boundary condition is

(29)$f(x_b,\epsilon ,\mu ,s)=f(x_b,\epsilon ,\mu ,s^{\prime }),s\cdot n_b>0$

where s’ = s – 2 (s · n_b) n_b.

2.4 Dimensionless analysis

Choosing some reference variables including the Fermi velocity v_F, characteristic length L, background temperature T₀ and Fermi energy E_F, a dimensionless analysis of eBTE is conducted,

(30)$s\cdot \frac{\partial f^{\ast }}{\partial x^{\ast }}+(\frac{-\partial (e\phi {)}^{\ast }}{\partial x^{\ast }}+s\times B^{\ast })\cdot \frac{\partial f^{\ast }}{\partial p^{\ast }}=\frac{f_{mc}^{eq,\ast }-f^{\ast }}{{\tau }_{mc}^{\ast }}+\frac{f_{mr}^{eq,\ast }-f^{\ast }}{{\tau }_{mr}^{\ast }}$

where

(31)$\begin{array}{rlrlrl}{\tau }_{mr}^{\ast }& =\frac{{\tau }_{mr}}{L/v_F},& {\tau }_{mc}^{\ast }& =\frac{{\tau }_{mc}}{L/v_F},& x^{\ast }& =\frac{x}{L}\\ f^{\ast }& =\frac{f}{f_0},& (e\phi {)}^{\ast }& =\frac{\left(e\phi \right)}{E_F},& B^{\ast }& =\frac{B}{E_F/\left|e{v}_{F}L\right|}\\ {\mu }^{\ast }& =\frac{\mu }{E_F},& p^{\ast }& =\frac{p}{E_F/v_F},& T^{\ast }& =\frac{T}{T_0}\\ f_{mr}^{eq,\ast }& =\frac{f_{mr}^{eq}}{f_0},& f_{mc}^{eq,\ast }& =\frac{f_{mc}^{eq}}{f_0},& u^{\ast }& =\frac{u}{v_F}\end{array}$

where f₀ = 1/2. It can be found that the electron transport is mainly determined by these dimensionless parameters ${\tau }_{m{r}^{\prime }}^{\ast }{\tau }_{m{c}^{\prime }}^{\ast }(e\phi {)}^{\ast }$ and B* for given material properties. The relationships between the material size and the length scale of particle transport, such as momentum-conserving and momentum-relaxing mean free paths, are reflected by dimensionless parameters.

3. Results and Discussions

Consider a homogeneous graphene disk^[24,52], where the ratio of the inner to outer radii of the disk is 1 to 5, and the thermal effects under the electric and magnetic fields are studied, as shown in Figure 1a. The outer diameter of the disk is regarded as the system characteristic length. Initially, the temperature and chemical potential inside the domain are equal to the background temperature and the Fermi energy, respectively. The initial electron distribution function satisfies the Fermi-Dirac equilibrium distribution. Temperature and chemical potential in both the inner and outer ends of the disk are fixed at the background temperature T₀ = 300 K and Fermi energy, respectively. Initial particle number density is N₀. An electric potential difference Δφ is applied and a magnetic field B is applied perpendicular to the graphene disk and pointing out of the paper.

Display Full Size

Download

Figure 1. (a) Schematic of a homogeneous Corbino disk geometry driven by electric and magnetic fields; (b, c) Temperature contour and heat flux streamline with a magnetic field B* = 0.8560 perpendicular to the graphene sheet, where (b) ${\tau }_{mr}^{\ast }=0.1$, ${\tau }_{mc}^{\ast }=\infty$ and (c) ${\tau }_{mr}^{\ast }=\infty$, ${\tau }_{mc}^{\ast }=0.1$.

Temperature contour and heat flux streamline in the diffusive and hydrodynamic regimes are plotted in Figure 1b,c, respectively, where (eΔφ)* = (eΔφ)/E_F = 0.0856 and B^* = 0.8560. It can be found that the temperature near the inner boundary is higher than that near the outer boundary, and temperature stays constant tangentially along the radial direction due to symmetry. Heat flows from outside to inside following the direction of the temperature gradient. Under the same electromagnetic field, the temperature rise in the hydrodynamic regime is larger than that in the diffusive regime. More importantly, the heat flux flows roughly along the radial direction in the diffusive regime, but there is a certain degree of deflection in the direction of heat flow in the hydrodynamic regime, namely, the heat flux no longer flows only along the radial direction. There is heat flux in the tangential direction corresponding to curved streamlines.

In order to understand the macroscopic distributions in more details, the spatial distributions of particle number density, temperature and the deflection angles θ of heat flux and electric current along the radial direction from the inside out are plotted in Figure 2, where tan θ = q_t/q_r or j_t/j_r, subscript t and r represent the tangential and radial components of heat flux or electric current, respectively. Please note that if the deflection angle is defined as θ = arctan (q_t/q_r) or arctan (j_t/j_r), it will be difficult to reflect the characteristic that the direction of heat flux is opposite to that of the electric current. Heat flux deflection phenomenon is suppressed by the MR scattering process and promoted by the MC scattering process, as shown in Figure 2a,b. Electrons are constantly moving from the outside to the inside under the electric field forces, which causes electrons to accumulate near the inner circle. As the electrons move inward, they feel the Lorentz forces from the magnetic field and interact with other electrons. The radial temperature gradient is established by the competition between Joule heating (which is stronger in the inner region where the current density is higher) and thermoelectric cooling (or Peltier effects). In the hydrodynamic regime, reduced momentum relaxation leads to lower dissipation and hence a larger temperature rise compared to the diffusive case. The heat flux contains both a dissipative Fourier-like component (opposing the temperature gradient) and a thermoelectric component (driven by the electric current). In the present simulations, the latter dominates in the hydrodynamic regime, resulting in heat flow from the outer (cooler) boundary toward the inner (hotter) boundary. This is a classic thermoelectric response and does not violate Fourier’s law.

Display Full Size

Download

Figure 2. Macroscopic distributions along the radial direction driven by electric potential gradient and magnetic field B* = 0.8560, where normalized distance is In(r/r_in)/In(r_out/r_in), r is the distance from the Corbino disk center. Deflection angles of (a) electric current and (b) heat flux; (c) Normalized particle number density n/N₀; (d) Normalized temperature.

In the diffusive regime, sufficient MR scattering process happens and the frequent energy and momentum exchange among electrons drives the distribution function into a local Fermi-Dirac distribution. On one hand, frequent MR scattering process prevents electrons from constantly moving from the outside in and gathering around the inner circle, so that the difference in particle number density between the inside and the outside becomes smaller. On the other hand, the acceleration tendency of electrons under the action of electric field forces is also destroyed by frequent MR scattering process. Specifically, although the electron is accelerated by the electric field force to obtain the corresponding momentum increment in a mean free path range, this momentum increment is immediately consumed by the MR scattering, which results in a small momentum or energy increment of the electron so that the temperature rise in the whole region is small. When ${\tau }_{mr}^{\ast }\to 0$ and ${\tau }_{mc}^{\ast }\to \infty$, a first-order Chapman-Enskog expansion of the distribution function leads to

(32)$\begin{array}{rl}{f}^{\ast }& =f_{mr}^{eq,\ast }-{\tau }_{mr}^{\ast }s\cdot \frac{\partial f_{mr}^{eq,\ast }}{\partial x^{\ast }}+{\tau }_{mr}^{\ast }\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }}\cdot \frac{\partial f_{mr}^{eq,\ast }}{\partial p^{\ast }}-{\tau }_{mr}^{\ast }(s\times B^{\ast })\cdot \frac{\partial f_{mr}^{eq,\ast }}{\partial p^{\ast }}+o\left({\tau }_{mr}^{\ast ,2}\right)\\ ⟹⟨s{f}^{\ast }⟩& \approx -{\tau }_{mr}^{\ast }⟨ss\cdot (\frac{\partial f_{mr}^{eq,\ast }}{\partial T^{\ast }}\frac{\partial T^{\ast }}{\partial x^{\ast }}+\frac{\partial f_{mr}^{eq,\ast }}{\partial {\mu }^{\ast }}\frac{\partial {\mu }^{\ast }}{\partial x^{\ast }})⟩+{\tau }_{mr}^{\ast }⟨s\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }}\cdot \frac{\partial f_{mr}^{eq,\ast }}{\partial {\epsilon }^{\ast }}\frac{\partial {\epsilon }^{\ast }}{\partial p^{\ast }}⟩\\ & \approx -{\tau }_{mr}^{\ast }⟨ss\cdot (\frac{\partial f_{mr}^{eq,\ast }}{\partial T^{\ast }}\frac{\partial T^{\ast }}{\partial x^{\ast }}+\frac{\partial f_{mr}^{eq,\ast }}{\partial {\mu }^{\ast }}\frac{\partial {\mu }^{\ast }}{\partial x^{\ast }})⟩+{\tau }_{mr}^{\ast }⟨s\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }}\cdot (-\frac{\partial f_{mr}^{eq,\ast }}{\partial {\mu }^{\ast }})s⟩\\ & \approx -{\tau }_{mr}^{\ast }⟨ss\cdot \frac{\partial f_{mr}^{eq,\ast }}{\partial T^{\ast }}⟩\frac{\partial T^{\ast }}{\partial x^{\ast }}+{\tau }_{mr}^{\ast }⟨-ss\frac{\partial f_{mr}^{eq,\ast }}{\partial {\mu }^{\ast }}⟩\cdot (\frac{\partial {\mu }^{\ast }}{\partial x^{\ast }}+\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }})\end{array}$

The first-order magnetic term disappears due to $(s\times B^{\ast })\cdot \frac{\partial f_{mr}^{eq,\ast }}{\partial p^{\ast }}=0$. Note that the gradients of temperature, chemical potential and electric potential are along the radial direction due to symmetry. Hence the deflection of macroscopic flux (Eq. (32)) tangentially along the radial direction in the diffusive regime is much smaller and on the order of $o\left({\tau }_{mr}^{\ast ,2}\right)$ regardless of heat flux or electric current.

In the hydrodynamic regime, sufficient MC scattering process results in a local equilibrium state with a nonzero drift velocity u^*. Frequent scattering process results in a much smaller resistance, which makes it easier for electrons to flow from the outside in so that the particle number density inside increases significantly, as shown in Figure 2c. On the other hand, small resistance indicates small dissipation and a larger temperature rise under the same electric potential difference, as shown in Figure 2d. When ${\tau }_{mr}^{\ast }\to \infty$ and ${\tau }_{mc}^{\ast }\to 0$, the first-order Chapman-Enskog expansion of the distribution function gives

(33)$\begin{array}{rl}{f}^{\ast }& =f_{mc}^{eq,\ast }-{\tau }_{mc}^{\ast }s\cdot \frac{\partial f_{mc}^{eq,\ast }}{\partial x^{\ast }}-{\tau }_{mc}^{\ast }(-\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }}+s\times B^{\ast })\cdot \frac{\partial f_{mc}^{eq,\ast }}{\partial p^{\ast }}+o\left({\tau }_{mc}^{\ast ,2}\right)\\ ⟹⟨s{f}^{\ast }⟩& \approx ⟨s{f}_{mc}^{eq,\ast }⟩-{\tau }_{mc}^{\ast }⟨ss\cdot (\frac{\partial f_{mc}^{eq,\ast }}{\partial T^{\ast }}\frac{\partial T^{\ast }}{\partial x^{\ast }}+\frac{\partial f_{mc}^{eq,\ast }}{\partial {\mu }^{\ast }}\frac{\partial {\mu }^{\ast }}{\partial x^{\ast }}+\frac{\partial f_{mc}^{eq,\ast }}{\partial u^{\ast }}\cdot \frac{\partial u^{\ast }}{\partial x^{\ast }})⟩\\ & -⟨{\tau }_{mc}^{\ast }s(-\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }}+s\times B^{\ast })\cdot \frac{\partial f_{mc}^{eq,\ast }}{\partial {\epsilon }^{\ast }}(\frac{\partial {\epsilon }^{\ast }}{\partial p^{\ast }}-u^{\ast })⟩\\ & \approx ⟨s{f}_{mc}^{eq,\ast }⟩+{\tau }_{mc}^{\ast }⟨-ss\cdot \frac{\partial f_{mc}^{eq,\ast }}{\partial T^{\ast }}⟩\frac{\partial T^{\ast }}{\partial x^{\ast }}+{\tau }_{mc}^{\ast }⟨-ss\frac{\partial f_{mc}^{eq,\ast }}{\partial {\mu }^{\ast }}⟩\cdot (\frac{\partial {\mu }^{\ast }}{\partial x^{\ast }}+\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }})\\ & +{\tau }_{mc}^{\ast }⟨-\frac{\partial f_{mc}^{eq,\ast }}{\partial {\mu }^{\ast }}s(-\frac{\partial \left(e{\phi }^{\ast }\right)}{\partial x^{\ast }}+s\times B^{\ast })⟩\cdot u+{\tau }_{mc}^{\ast }⟨-\frac{\partial f_{mc}^{eq,\ast }}{\partial {\mu }^{\ast }}ss\cdot p^{\ast }⟩\cdot \frac{\partial u^{\ast }}{\partial x^{\ast }}\end{array}$

It can be found that the nonzero drift velocity is significantly affected by the electromagnetic field. The Chapman-Enskog expansion analysis reveals that the flux deflection originates from the emergence of a collective drift velocity in the electron distribution function in the hydrodynamic limit. This electron drift, subject to the Lorentz force, acquires a tangential component in the disk system that cannot be compensated by electric forces in the bulk region. Consequently, the flux deflection is a genuine bulk phenomenon and is largely insensitive to the detailed choice of boundary conditions.

Numerical results under various magnetic field forces and MR/MC scattering processes are shown in Figure 3. It can be found that the deflection angle decreases with the magnetic field force under the same scattering rates. The deflection angle increases when the MC scattering process increases but decreases when the MR scattering process increases. Based on dimensionless analysis of eBTE (Eq. (30)), it can be found that the local flux in the tangential direction of the radius depends on the temperature gradient, electric and magnetic field force and the local electron scattering strength. When MR scattering dominates, there are large momentum relaxation so that the acceleration of electrons by the electric field is less effective. The electron mobility is reduced, and consequently the effect of the magnetic field force is naturally reduced, making it difficult to change the direction of electron motion. Hence, electron motion in the tangential direction tends to zero, and the deflection phenomenon is weak. When the MC scattering process dominates, there is little momentum relaxation and the electric field force drives electrons to move from the outside to the inside. During this motion, it feels a tangential Lorentz force. Different from the Hall effect in cuboid geometry, the Lorentz force in disk geometry is usually not in the radial direction. Therefore, it is not easy to be compensated by the electrostatic field force which is always in the radial direction due to the symmetry of particle number density. Consequently, the macroscopic flux deflection phenomenon appears.

Display Full Size

Download

Figure 3. Macroscopic distributions along the radial direction driven by electric potential gradient and various magnetic field, where normalized distance is In(r/r_in)/In(r_out/r_in), r is the distance from the Corbino disk center. The magnetic field in the first, second, third row is B^* = 0.0856, B^* = 0.428, B^* = 0.856, respectively. Deflection angles of (a, d, g) electric current and (b, e, h) heat flux; (c, f, i) Normalized temperature.

Secondly, the thermal effects in a graphene Corbino disk under the temperature gradient and magnetic field are studied, as shown in Figure 4a. In this case, the electric and chemical potentials at the inner and outer ends of the disk are the same. The dimensionless temperature T^* at the inner and outer boundaries is 1.01 and 0.99, respectively.

Display Full Size

Download

Figure 4. (a) Schematic of a homogeneous Corbino disk geometry driven by temperature gradients and magnetic fields; (b, c) Temperature contour and heat flux streamline with a magnetic field B^* = 0.8560 perpendicular to the graphene sheet, where (b) ${\tau }_{mr}^{\ast }=0.1$, ${\tau }_{mc}^{\ast }=\infty$ and (c) ${\tau }_{mr}^{\ast }=\infty$, ${\tau }_{mc}^{\ast }=0.1$.

Numerical results are shown in Figure 4, Figure 5, and Figure 6. The particle number density changes slightly inside the disk due to small temperature variance, which indicates that the electric field force is much smaller. It can be found that there are larger temperature slips near the inner boundary than at the outside, because the inner radius is smaller so that the ballistic effects inside are more serious. Similarly, there is little heat flux deflection in the diffusive regime due to the sufficient MR scattering process and the heat flux deflection is obvious in the hydrodynamic regime. Opposite to the above results shown in Figure 1b,c driven by electromagnetic fields, the heat flux flows from inside out and electric current flows from outside in under the temperature gradient.

Display Full Size

Download

Figure 5. Macroscopic distributions along the radial direction driven by temperature gradient and magnetic field B* = 0.8560, where normalized distance is In(r/r_in)/In(r_out/r_in), r is the distance from the disk center. Deflection angles of (a) electric current and (b) heat flux; (c) Normalized particle number density n/N₀; (d) Normalized temperature.

Display Full Size

Download

Figure 6. Macroscopic distributions along the radial direction driven by temperature gradient and various magnetic field, where normalized distance is In (r/r_in)/In (r_out/r_in), r is the distance from the Corbino disk center. The magnetic field in the first, second, third row is B* = 0.0856, B* = 0.428, B* = 0.856, respectively. Deflection angles of (a, d, g) electric current and (b, e, h) heat flux; (c, f, i) Normalized temperature.

We elucidate the underlying physical mechanisms in terms of the evolution of the electron distribution function in the momentum space. In the hydrodynamic regime, electron transport presents collective behavior, namely, the electron distribution function is generally deviated in a certain direction with drift velocity u^* in the momentum space. When an electric field is applied, the electron suffers from a radial outside-in electric field force, which causes its distribution function to deviate from the equilibrium state towards the center of the disk. Therefore, the electron will also be driven by a clockwise magnetic field force to move clockwise under the action of a magnetic field, which leads to the electric current flows anticlockwise and the heat flux flows clockwise. Different from the electric field force, the temperature gradient drives the electron distribution function to deviate from equilibrium in a direction away from the center of the disk. Therefore, the electron will be driven by a counterclockwise magnetic field force to move counterclockwise under the action of a magnetic field, which leads to the electric current flowing clockwise and heat flux flowing anticlockwise.

In a word, Figure 2 and Figure 3, Figure 5 and Figure 6 systematically examine the magnetic-field-tuned electron transport in a Corbino disk under electric-field and temperature-gradient drives, revealing the essential distinction between hydrodynamic and diffusive regimes in their magnetic response. The deflection angle θ increases monotonically with B^*, yet its magnitude is governed by the dominant scattering mechanism: enhanced MC scattering (reduced ${\tau }_{mc}^{\ast }$) sustains a collective drift that enables the Lorentz force to generate a persistent transverse momentum, yielding a pronounced deflection as a hallmark of hydrodynamics; enhanced MR scattering (reduced ${\tau }_{mr}^{\ast }$) dissipates momentum within each mean free path, quenching the Lorentz accumulation and restoring diffusive transport with θ→0. The deflection chirality reverses upon switching the drive, inward drift under an electric field produces clockwise deflection, while outward diffusion under a temperature gradient flips the Lorentz force direction. Thermally, the hydrodynamic regime under an electric drive exhibits a substantial temperature rise due to suppressed dissipation, whereas the diffusive regime remains nearly isothermal. Under a temperature-gradient drive, the temperature rise is weak and scattering-insensitive, as the heat flux is boundary-controlled rather than dissipation-dominated. Collectively, the deflection amplitude and temperature rise serve as independent diagnostics for hydrodynamic versus diffusive transport, displaying clear scaling separation under the competition between scattering rates and magnetic field.

4. Conclusion

A kinetic simulation of electron thermal and charge transport in a graphene Corbino disk has been performed, with the magnetic field strength and MC/MR scattering rates being systematically scanned under the framework of the BTE. Three physical distinctions between the hydrodynamic and diffusive regimes are identified. Firstly, flux deflection by the Lorentz force is found to be a hydrodynamic signature, pronounced when electron-electron scattering dominates and quenched when momentum-relaxing scattering prevails. In the diffusive limit, the first-order magnetic correction to the flux vanishes by symmetry, while in the hydrodynamic limit, the finite drift velocity allows the Lorentz force to generate a sustained tangential component. Secondly, significant temperature rise is produced in the hydrodynamic regime under electric-field driving due to reduced energy dissipation, whereas a weaker thermal response is induced under temperature-gradient driving, which is less sensitive to the electron MC/MR scattering rates. Thirdly, the deflection chirality under temperature-gradient driving is observed to be opposite to that under electric-field driving. The current results indicate that thermodynamic physical quantities, especially the heat flow deflection caused by magnetic fields, can also be used as an approach to distinguish hydrodynamic electron flow and classical diffusive transport. However, the spatially resolved, high-resolution measurement^[5,66,67] of heat flux or temperature is much more difficult than that of electric current or potential. Perhaps with the development of future experimental techniques, it will be possible to detect the spatial distribution of heat flux.

Supplementary materials

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

Acknowledgements

The authors acknowledge Beijing PARATERA Tech CO., Ltd. for the HPC resources.

Authors contribution

Zhang C: Supervision, conceptualization, investigation, methodology, formal analysis, funding acquisition, writing-original draft.

Lian M: Investigation, methodology, formal analysis, writing-review & editing.

Liang H, Li X: Formal analysis, writing-review & editing.

Guo Z: Conceptualization, methodology, formal analysis, writing-review & editing.

Lü JT: Supervision, conceptualization, investigation, formal analysis, writing-review & editing.

Conflicts of interest

Jing-Tao Lü is an Editorial Board member of Thermo-X. The other 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 that support the findings of this study are available from the corresponding author upon reasonable request.

Funding

Chuang Zhang acknowledges the support of the National Natural Science Foundation of China (Grant No. 52506078) and the Zhejiang Provincial Natural Science Foundation of China under Grant No. LMS26E060012.

Copyright

© The Author(s) 2025.