Analyzing Hydrothermal Wave Transitions through Rotational Field Application Based on Entropy Production

Abstract

In this study, we evaluated the nonlinear dynamics of convection flow using the thermodynamic variational principle, focusing on scenarios where multiple external forces, such as a thermal gradient and rotational field, are applied to a shallow annular pool. We observed that with the increase in the thermal gradient, the flow changed from an axial flow to a rotational oscillatory flow with the wave amplitudes aligned. Further increasing the temperature difference led to a rotational oscillatory flow characterized by alternating wave generation and annihilation. Our analysis of the flow, considering heat fluxes orthogonal to the thermal gradient, allowed us to describe the flow state as a phase at equilibrium. The state transition of the flow was accompanied by a discontinuous jump in the heat flux, which occurred at the intersection of the entropy production curves. The first transition occurred at a temperature difference $\Delta T=12.4 \mathrm{K} \left(\mathrm{M}\mathrm{a}\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{g}\mathrm{o}\mathrm{n}\mathrm{i} \mathrm{n}\mathrm{u}\mathrm{m}\mathrm{b}\mathrm{e}\mathrm{r},\mathrm{M}\mathrm{a}=1716\right)$ and the second at ΔT = 16.3 K $\left(\mathrm{M}\mathrm{a}=2255\right)$. Analysis based on entropy production could accurately predict the observed transition points.

1. Introduction

Understanding the phenomena of transitions in the spatiotemporal structure of spontaneous convection under an applied external field is crucial in fundamental physics and industry. In fundamental physics, it relates to key areas such as nonlinear dynamics and nonequilibrium thermodynamics related to fluid instabilities [1,2,3,4]. In industry, these transitions are crucial for manufacturing processes using transport phenomena [5,6,7,8,9]. Rayleigh–Bénard convection—a classic example of thermal convection—has been extensively studied both theoretically and experimentally for many years [10,11,12,13,14,15,16,17]. Chandrasekhar [18] stated that the symmetry breaking of Bénard convection is sustained by the competition between the viscosity-induced dissipation of kinetic energy and the internal energy released by buoyancy. He obtained stable conditions for convection generation. Hohenberg and Swift [19] analyzed the complex behavior of Bénard convection under temporally modulated temperature differences using the Lorenz model. They obtained a bifurcation diagram that illustrates the transition process not only from conductive heat transfer to Bénard convection but also from a hexagonal pattern to a roll pattern.

Recently, attempts have been made to theoretically analyze the evolution of nonequilibrium states in systems where various irreversible processes are coupled, using the variational principle of maximum entropy production (MEPP) [20,21,22,23,24]. The MEPP shows that nonequilibrium states evolve into more complex states that generate higher entropy. Here, “complex” refers to the interference between two or more irreversible processes [25]. When two fluids with different viscosities are displaced in a Hele–Shaw cell, the interface becomes destabilized, transitioning from a stable diffusion interface to a finger-like interface due to the Saffman–Taylor instability as the viscosity ratio increases [26]. With a further increase in the viscosity ratio, the flow becomes more complex owing to the splitting of the viscous fingers [27,28]. The flow at the growing interface evolves such that entropy production is maximized [29]. The transition of Bénard convection from a hexagonal to a rolled pattern owing to temperature modulation transitions to a pattern with higher entropy production [30]. The bistability of the hexagonal and rolled patterns arises because, while the rolled pattern exhibits higher heat flux, the hexagonal pattern still produces higher entropy.

In addition to predicting the state transitions, the variational principle can be used to derive the governing equations for complex transport phenomena. Doi et al. derived the time evolution equations for the contact line and droplet shape in the dynamics of an evaporating droplet, driven by the Marangoni effect [31]. They also derived time evolution equations for deposition patterns, such as coffee rings and mountain shapes, during the drying process of the liquid film [32]. Thus, the MEPP has significantly advanced our understanding of various nonequilibrium phenomena. However, it remains unknown whether MEPP can be applied to complex transport phenomena involving multiple external forces. In this study, we quantitatively analyzed the state transitions of thermal convection in a system under a rotational field using the MEPP.

2. Equations of Motion and Numerical Details

Figure 1 shows the schematic of an annular pool with an applied temperature difference. A shallow annular pool with an inner wall of 15 mm $\left(R_{\mathrm{i}}\right)$, an outer wall of 50 mm $\left(R_{\mathrm{o}}\right)$, and a depth of 3 mm$\left(b\right)$ was filled with a silicon melt. The inner wall is fixed at the melting point of silicon, $T_c(=1683 \mathrm{K})$, and the outer wall is set at a higher temperature than the inner wall, $T_h(=T_c+\Delta T)$. The upper liquid surface is a free interface, where Marangoni convection is used. No-slip conditions were imposed on the fluid velocities on all solid walls (inner and outer walls and bottom).

Assuming that the silicon melt behaves as an incompressible Newtonian fluid and the Boussinesq approximation is applicable, the continuity, Navier–Stokes, and energy equations can be expressed as

$${\displaystyle \frac{\partial v_i}{\partial x_i}}=0,$$

(1)

$${\displaystyle \frac{\partial v_i}{\partial t}}+v_i·{\displaystyle \frac{\partial v}{\partial x_i}}=-{\displaystyle \frac{1}{\rho }}{\displaystyle \frac{\partial p}{\partial x_i}}+\nu {\nabla }^2{v}_i,$$

(2)

$${\displaystyle \frac{\partial T}{\partial t}}+v_i·{\displaystyle \frac{\partial T}{\partial x_i}}=\alpha {\nabla }^{2}T,$$

(3)

where $v_i, \rho ,\nu ,$and$\alpha$are the velocity, density, kinematic viscosity, and thermal diffusivity of the fluid, respectively. Additionally, $p$ is the pressure, $T$ is the temperature, and $t$ is time. We assume that if a subscript appears more than once in the same term in the equation, summation over the repeated subscripts is implied. In this thermally driven system, spontaneous convection, resulting from the density difference, and Marangoni convection, resulting from the surface tension difference, coexist. However, the former can be neglected because the density changes owing to temperature differences are small. Therefore, it is unnecessary to consider the buoyancy term in Equation (2) [33]. Initially, there was no convection in the annular pool, and the fluid temperature was set uniformly at 1683 K, which is the melting point of silicon. Convection is induced by the Marangoni effect owing to the presence of a thermal gradient on the free surface, according to the following equation:

$$\eta {\displaystyle \frac{\partial v_x}{\partial z}}={\gamma }_{\mathrm{T}}{\displaystyle \frac{\partial T}{\partial x}},\eta {\displaystyle \frac{\partial v_y}{\partial z}}={\gamma }_{\mathrm{T}}{\displaystyle \frac{\partial T}{\partial y}},$$

(4)

where $\eta$ and ${\gamma }_{\mathrm{T}}$ are the viscosity and the surface tension coefficient of the fluid, respectively. The fluid properties and dimensionless parameters used in the numerical calculations are listed in

Table 1. Physical properties of silicon melt.

Physical Properties	Value
Thermal conductivity (λ)	64 Wm−1K−1
Viscosity (η)	7.0 × 10−4 Pa s
Density (ρ)	2530 kgm−3
Thermal expansion coefficient (β)	1.5 × 10−4 K−1
$\mathrm{Surface} \mathrm{tension} \mathrm{coefficient} ({\gamma }_T$)	−7.0 × 10−5 Nm−1K−1
$\mathrm{Heat} \mathrm{capacity} (C_p$)	1000 Jkg−1K−1
$\mathrm{Melting} \mathrm{temperature} (T_c$)	1683 K

and

Table 2. Numerical parameters of the simulation.

Parameter	Value
Grashof number ($\mathrm{G}\mathrm{r}=g{\rho }^2\beta ∆T{d}^3/{\mu }^2$)	$1.19\times {10}^4$
Prandtl number ($\mathrm{P}\mathrm{r}=\eta C_p/\lambda$)	$1.09\times {10}^{-2}$
Marangoni number ($\mathrm{M}\mathrm{a}=-{\gamma }_T\Delta T\left(R_{\mathrm{o}}-R_{\mathrm{i}}\right)/\eta \alpha ,\alpha =\lambda /\rho C_p$)	$1.38\times {10}^2-2.76\times {10}^3$
Rayleigh number ($\mathrm{R}\mathrm{a}=\mathrm{G}\mathrm{r}\times \mathrm{P}\mathrm{r}$)	$1.19\times {10}^2$

. The relative strength between natural convection due to gravity and Marangoni convection caused by temperature differences was found to be very small $\left({{\mathrm{G}\mathrm{r}}^{1/2}/\mathrm{M}\mathrm{a}}^{2/3}=2.70\times {10}^{-4}\right)$[34]; we therefore neglected the effect of convection due to density changes. The numerical analysis was performed in OpenFOAM 2.4.x using the PISO method, with the discretization of the equations through the finite volume method. We applied the Euler, QUICK, and Gaussian linear schemes to the terms involving the time derivative, convection, and Laplacian in the governing equations, respectively. The grid numbers were 81, 180, and 21 in the radial, $\mathrm{r}$, circumferential, $\mathsf{\theta }$, and vertical directions, $\mathrm{z}$, respectively. We previously reported that the numerical calculations depended on the mesh grid size, and good grid convergence was guaranteed when$\mathrm{r}\times \mathsf{\theta }\times \mathrm{z}$were finer than 60 × 160 × 20, respectively [34,35,36,37]. This grid configuration sufficiently captured hydrothermal waves (HTWs) in fluids with a low Prandtl number (Pr). Thus, a turbulence model was not adopted. We have also validated the simulation code in the absence of external forces, i.e., without rotation. The temperature and velocity profiles showed similar trends to those observed previously [33].

Increasing the temperature difference between the inner and outer walls causes the flow state in a system to transition from an axially symmetric flow (ASF) to rotating oscillatory flows, known as a HTW [38]. The former is a radially symmetric flow, while the latter causes convection in the circumferential direction. The transition point of the flow state in the system is at approximately $\Delta T=7 \mathrm{K}$. Figure 1b shows a typical ASF and HTW displayed in terms of surface temperature fluctuation $\delta T=T-1/2\pi {\int }_0^{2\pi }T\mathrm{d}\theta$. For the ASF, the distribution is almost uniform, showing no significant temperature difference. In contrast, for the HTW, a radial region with alternating hot and cold regions appears and rotates in the circumferential direction.

In a fluid undergoing an irreversible process, the rate of entropy production, $\sigma$, can be expressed via the entropy produced per unit volume (s) and the conserved quantity (${\alpha }_k$) using the following equation [39].

$$\sigma ={\displaystyle \frac{\partial s}{\partial t}}={\displaystyle \frac{\partial s}{\partial {\alpha }_k}}{\displaystyle \frac{\mathrm{d}{\alpha }_k}{\mathrm{d}t}},$$

(5)

$${\alpha }_k=\left\{\begin{array}{c}\epsilon :\mathrm{energy}\\ \rho v :\mathrm{momentum} \\ {\rho }_m :\mathrm{mass}\end{array}\right\}$$

Subscript k represents the conserved quantity. Energy $\epsilon$ is transported by heat conduction, momentum $\rho v$ is transported by fluid movement, and matter ${\rho }_m$ is transported by mass transfer [40,41]. Entropy is produced in the fluid because of the irreversible transport of conserved quantities through the driving force. Thus, the entropy production can be expressed as

$$X_k={\displaystyle \frac{\partial s}{\partial {\alpha }_k}}, J_k={\displaystyle \frac{\mathrm{d}{\alpha }_k}{\mathrm{d}t}}$$

(6)

$$\sigma =X_k{J}_k.$$

(7)

Here, $X_k$ represents the local driving force, and $J_k$ represents the transport of conserved quantities in response to the driving force. In thermal convection, two irreversible processes interact: the transport of energy via heat conduction and the transport of momentum via fluid movement. Because the entropy production of the latter is significantly smaller than that of the former, only the entropy production by the heat flux is considered here [33]. Therefore, subscript k is omitted. Experiments have demonstrated that even with complex flow patterns, such as those of a Bénard cell, the heat flux is linear with respect to the driving force [42,43]. This leads to the following equation:

$$X_i={\displaystyle \frac{\partial }{\partial x_i}}\left({\displaystyle \frac{1}{T}}\right), J_i=-\lambda {\displaystyle \frac{\partial T}{\partial x_i}}$$

(8)

$$\sigma =X_i{J}_i$$

(9)

where $\lambda$ is the thermal conductivity of the fluid. Entropy production is a scalar quantity that represents the amount of heat produced per unit volume due to the thermal gradient in each direction. However, because certain components are highly related to the symmetry with respect to the driving force that keeps the system nonequilibrium, we will denote the entropy production separately for each gradient in this study.

$${\sigma }_x=X_x{J}_x,{\sigma }_y=X_y{J}_y,{\sigma }_z=X_z{J}_z$$

(10)

As the system is far from equilibrium, both $X_i$ and $J_i$ are time-dependent and spatially distributed. To characterize the nonequilibrium state of the system through space- and time-dependent thermodynamic quantities, local thermodynamic quantities at point A = (30 mm, 0 mm, 3 mm) were calculated by time averaging.

$$\left|X_i\right|=1/\Delta t{\int }_{t_o}^{t_1}\sqrt{{\left({\displaystyle \frac{\partial }{\partial x_i}}\left({\displaystyle \frac{1}{T}}\right)\right)}^2}\mathrm{d}t, \left|J_i\right|=1/\Delta t{\int }_{t_o}^{t_1}\sqrt{{\left(\lambda {\displaystyle \frac{\partial T}{\partial x_i}}\right)}^2}\mathrm{d}t,$$

(11)

where $t_o$ = 200 s and $t_1$ = 300 s. Additionally, $\Delta t$ represents the time interval $t_1-t_o$. The entropy production was calculated as the product of $\left|X_i\right|$ and $\left|J_i\right|$, as shown in Equation (9). In addition to the thermal gradient, a counterclockwise rotational field with $\mathsf{\Omega }=0, 0.5, 1,$and $3 \mathrm{r}\mathrm{p}\mathrm{m}$ was applied to the annular pool to assess changes in the flow state and to measure the heat flux and entropy production.

3. Results and Discussion

Figure 2 shows the changes in temperature at point A on the fluid surface for different rotational speeds at $\Delta T=14 \mathrm{K}$. In the absence of a rotating field, the HTW exhibits oscillations with a fluctuation width of approximately 0.7 K and a period of 6.7 s, with slightly modulated amplitudes. For $\mathsf{\Omega }=0.5 \mathrm{r}\mathrm{p}\mathrm{m}$, the HTW is a combination of long-term oscillations, with a period of approximately 75 s, and short-term oscillations, with a period of approximately 2.5 s. The temperature fluctuation of the long-term oscillation was comparable to that without rotation, whereas that of the short-term oscillation was approximately one-third that of the long-term oscillation. At $\mathsf{\Omega }=1.0 \mathrm{r}\mathrm{p}\mathrm{m}$, the long-term oscillation disappeared, and the temperature fluctuation resulted in a constant-amplitude HTW with a period of approximately 1.6 s. For $\mathsf{\Omega }=3.0 \mathrm{r}\mathrm{p}\mathrm{m}$, the temperature fluctuation disappeared, and the ASF appeared. The HTW consists of long-term and short-term oscillations. The Fourier transform of the temperature time series data in the absence of rotation reveals three distinct peaks: $f_{\mathrm{S}}, f_{\mathrm{M}} \mathrm{a}\mathrm{n}\mathrm{d} f_{\mathrm{L}}$. Figure 2b shows that the arithmetic mean of the low and high frequencies, $\left(f_{\mathrm{H}}+f_{\mathrm{L}}\right)/2$, is exactly the same value as the middle frequency, $f_{\mathrm{M}}$. This implies resonant oscillations due to the interactions between oscillators with similar frequencies. In addition, half of the difference between the low and high frequencies, $\left(f_{\mathrm{H}}-f_{\mathrm{L}}\right)/2$, has a value of 0.01667 Hz. By looking closely at the magnified low-frequency region (inset of Figure 2b), a small peak can be observed at the point corresponding to this frequency, which is consistent with the frequency of the temperature oscillation beat. Hence, in the absence of rotation, the HTW can be expressed as a carrier wave of frequency $\left(f_{\mathrm{H}}+f_{\mathrm{L}}\right)/2$ whose amplitude is modulated by an envelope wave of frequency $\left(f_{\mathrm{H}}-f_{\mathrm{L}}\right)/2$. Converting the frequency of the envelope wave to a period corresponds to a wave with a period of 1 s. As the rotational speed of the annular pool is gradually increased, the undulation of the temperature oscillation becomes more pronounced. Moreover, when the period of the envelope wave coincides with the rotational period, only short-wave oscillations emerge. Similar resonant oscillations can be observed in the rotational motion of droplets due to the solutal Marangoni effect [44]. Furthermore, increasing the rotational speed of the annular pool increases the centrifugal force, thereby stabilizing the flow [18].

Figure 3 shows the entropy production ${\sigma }_i$ in each direction, calculated using Equations (10) and (11) as a function of the thermodynamic force $X_i$. As $X_i$ increases, ${\sigma }_i$ increases parabolically. The largest entropy production, ${\sigma }_x$, occurs in the x-direction due to heat flux, as it aligns with the thermal gradient on the inner and outer walls, which serves as the driving force. This is followed by ${\sigma }_y$ and ${\sigma }_z$. If the flow effect is ignored and entropy production is considered to solely result from heat conduction according to Fourier’s law, it can theoretically be expressed as ${\sigma }_i=\lambda T^2{X}_i^2$ [33]. When the theoretical values were plotted using the average temperature $T_{av}=\left(T_c+T_h\right)/2$, the acquired data closely matched the theoretical predictions. Even when the flow state changed from ASF to HTW, the entropy production data exhibited the same curve. This result accurately reflects the increase in entropy production owing to local driving forces but also indicates that changes in state cannot be distinguished by expressing entropy production solely as a function of thermodynamic forces.

Determining the function according to changes in nonequilibrium states is a critical issue. Previous experimental and computational examples have shown that using a function of the global driving force that maintains the nonequilibrium state of the system is effective [45,46]. The emergence of new flow states owing to symmetry breaking in response to the global driving force is often observed in various transport phenomena, such as Bénard convection, Marangoni convection, and the Saffman–Taylor instability [29,30,33,44]. A common feature of this phenomenon is the production of excess entropy in the direction orthogonal to the global driving force. This excess entropy is produced by thermodynamic flux in the direction of symmetry breaking. This thermodynamic flux is referred to as the emergent flux. When the emergent flux is considered as a function of the global driving force, the entropy production similarly becomes a function of the global driving force.

$$J_i=L\prime \left(X_G-\Theta \right)$$

(12)

$$\sigma =L{\left(X_G-\Theta \right)}^2$$

(13)

where $L^{\prime }, L$, and $\Theta$ are phenomenological coefficients. In particular, $\Theta$ is significant as it is related to the free energy required to transition from one state to another and is also related to the degree of coherence of two or more irreversible processes [45,47]. Because the phenomenological coefficients are different for the different states of the system, the states of the system can be described by different lines. Using this relationship between the emergent flux and the global driving force, we can identify the state of the system. The entropy production for each nonequilibrium state can then be obtained using Equation (13). According to the MEPP, the thermodynamic flux changes such that the system transitions to a state with higher entropy production. For the Bénard convection, the global driving force is the thermal gradient provided above and below the vessel [30,42]. As the driving force increased, the thermodynamic flux changed discontinuously with each transition from conduction to convective heat transfer and from hexagonal to roll patterns. Discontinuity represents the transition point of the nonequilibrium state, which is represented by the intersection of the entropy production curves. The MEPP, using the relationship between the emergent flux and global driving force, was used to quantitatively describe complex transitions in thermal convection under a rotating field. Additionally, it was used to predict the transition point.

In this thermally driven system, the thermal gradient provided on the inner and outer walls is the driving force that changes the state of the thermal convection inside. We express the thermal gradient in terms of the dimensions of the thermodynamic force, referring to it as the global driving force $X_G$.

$$X_G={\displaystyle \frac{1}{R_{\mathrm{o}}-R_{\mathrm{i}}}}\left({\displaystyle \frac{1}{T_c}}-{\displaystyle \frac{1}{T_h}}\right)$$

(14)

The heat fluxes in the x- and y-directions at each rotational speed are shown in Figure 4. The lower and upper abscissae represent the global driving force and the corresponding Marangoni number $\mathrm{M}\mathrm{a}$, respectively. The numerical values of $X_G$ are similar to those of $\Delta T$, except for an order of magnitude. $J_x$ increased linearly with respect to the global driving force regardless of the rotational speed. If the flow effects are disregarded and only heat conduction according to Fourier’s law is considered, the heat flux is theoretically expressed as $J_x=\lambda T_c{T}_h\left(R_{\mathrm{o}}-R_{\mathrm{i}}\right)X_G/\left\{r\ln \left(R_{\mathrm{o}}/R_{\mathrm{i}}\right)\right\}$[33]. The heat flux in the x-direction was almost consistent with the theoretical equation and was independent of the rotational speed. As the direction of the global thrust was radial, it agreed with the thermal gradient in the x-direction at the measurement point. Marangoni convection occurred in the x-direction at the free surface; however, the effect of the flow on the heat flux in the x-direction was negligible because the conduction heat transfer was more dominant. The heat flux in the y-direction, where the symmetry was broken by the global driving force, exhibited a nontrivial trend. In the absence of a rotational field, $J_y$ can be described by different straight lines for the ASF and HTW using Equation (12). For ASF, $J_y$ passes through the origin and $\Theta =0$, indicating the complete absence of interference with other irreversible processes. For the HTW, $\Theta$ exhibits a finite value, indicating the degree of interference from other irreversible processes. In the presence of a rotational field, $J_y$ was on the same line as in the ASF, independent of the rotational speed; however, the value of the global driving force required to generate the HTW increased significantly with the increasing rotational speed. This indicates that the centrifugal force generated by the rotational field suppresses HTW generation. As the global driving force further increased, a discontinuous change in the heat flux was observed, and the slope of the straight line became even larger. This implies the existence of a different convection pattern.

The emergent heat flux results suggest the presence of a new HTW flow pattern. Figure 5 illustrates the surface temperature fluctuations for $\Delta T=15 \mathrm{K}$ (top panels) and $\Delta T=17 \mathrm{K}$ (bottom panels) at $\mathsf{\Omega }=1.0 \mathrm{r}\mathrm{p}\mathrm{m}$. For $\Delta T=17 \mathrm{K}$, time progresses from the left to the right panels at an interval of 0.6 s. Waves consisting of eight pairs of high- and low-temperature regions are generated from the inner wall, expanding radially while rotating counterclockwise. When the tip of the growing wave reaches half of the annular pool, the wave detaches from the inner wall and diffuses toward the outer wall. New waves reappear from the inner wall. This results in heat waves propagating in the circumferential direction, with waves being repeatedly generated and annihilated. A normal rotating oscillating flow with a constant width and length that rotates in one direction is represented by HTW$\alpha$. On the other hand, a rotating oscillating flow that repeatedly generates and annihilates waves is represented as HTW$\beta$. In the presence of a rotating field, the heat flux region represented by the straight line with the steepest slope in Figure 4a corresponds to the region where HTW$\beta$ is observed.

The transitions from ASF to HTW$\alpha$ and from HTW$\alpha$ to HTW$\beta$ are discussed using the MEPP. Figure 6 shows the emergent heat flux and the corresponding entropy production at $\mathsf{\Omega }=1.0 \mathrm{r}\mathrm{p}\mathrm{m}$. To obtain the entropy production curve, we first determined the differential coefficient $\mathrm{d}{J}_y/\mathrm{d}{X}_G$ at each point of the heat flux. The heat flux was divided into three regions based on the differential coefficient values. In the regions divided by differential coefficients, the heat flux was fitted using Equation (12), providing three straight lines with different values of L’ and Θ. The entropy production in each region was fitted using Equation (13) with Θ in common, resulting in three entropy curves. The fitting parameters for each region are listed in

Table 3. Fitting parameters using Equations (12) and (13) at $\omega =1.0 {\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$.

	Region	$\mathbf{d}{\mathit{J}}_{\mathit{y}}/\mathbf{d}{\mathit{X}}_{\mathit{G}}$	$\mathit{L}\prime$	$\Theta$	$\mathit{L}$	$\Delta {\mathit{T}}_{\mathit{c}}$	Ma
1. 2. 3.	$0\le T\le 7$ $7<T\le 16$ $16<T\le 18$	$31.8\pm 3.3$ $253.5\pm 56.6$ $1719.2\pm 209.8$	$30.9$ $278.3$ $1593.5$	$0$ $11.1\times {10}^{-5}$ $15.5\times {10}^{-5}$	$5.21\times {10}^{-6}$ $5.61\times {10}^{-4}$ $2.33\times {10}^{-2}$	$12.4$ $16.3$	$1716$ $2255$

.

The entropy production curves of the ASF and HTW$\alpha$ intersect at $\Delta T_c=12.4 \mathrm{K},$ while those of HTW$\alpha$ and HTW$\beta$intersect at $\Delta T_c=16.3 \mathrm{K}$. This is nearly consistent with the point where the emergent heat flux changes discontinuously. Additionally, at $\omega =0.5 {\mathrm{m}\mathrm{i}\mathrm{n}}^{-1}$, the transition point was accurately predicted using the MEPP. The stability of each state may be related to the value of $\Theta$ in Equation (12), particularly that of the first state, which may be highly dependent on the value of $\Theta$ of the second state. If this value is large, the first state will inevitably exist in a wider driving force region, as a higher driving force will be required for the new state to emerge, and more energy will need to be introduced into the system. Predictions based on the MEPP rely on a function relating the thermodynamic flux to the driving force, which would be possible if a one-to-one relationship existed between the functions. Even if there is nonlinearity between the two, the transition point can be predicted because the entropy production can be easily calculated [48]. In general, the most stable state of the system is when the intensive variables are uniformly distributed everywhere. However, when these variables are distributed in the system, forced by external forces applied to the boundaries, the system develops such that the state of the system becomes more uniform, i.e., producing higher entropy. In a thermal convection system with a weak external force, the system tries to homogenize through heat conduction, and the flow in the direction of the external force generates an axisymmetric flow. As the external force increases and the temperature distribution widens, flux emerges in the direction perpendicular to the external force to further magnify the uniformity in the system and improve the heat transfer efficiency. This state corresponds to HTW$\mathsf{\alpha }$, where the rotational oscillatory flow is gently decaying toward the outside, as shown in Figure 5a. As the external force is further increased, the generation and annihilation of the rotational oscillatory flow repeat, and two rotational oscillatory flows are generated inside and outside the system, corresponding to HTW$\mathsf{\beta }$, as shown in Figure 5b. These two rotating oscillatory flows are out of phase, which further promotes the homogenization of the system and increases the efficiency of heat transfer.

Three key observations can be derived from these results. First, the MEPP effectively predicts the transition point of the state, even in a nonequilibrium system with multiple external forces, including the addition of thermal gradients and rotational fields. Second, the relationship between the emergent heat flux and the global driving force leads to the discovery of new flow states. Third, the relationship between the emergent heat flux and the global driving force allows the nonequilibrium states to be specified as if they were a phase of an equilibrium system.

For example, in equilibrium, if the state of water can be categorized into solid, liquid, and vapor phases, each state can be represented by a distinct function, with the chemical potential $\mu$ expressed as a function of the temperature (Figure 7a). The entropy per mole $S_m$ can be obtained from the slope as ${\left(\partial \mu /\partial T\right)}_p=-S_m$. Because $S_m>0$, the chemical potential decreases as the temperature increases. The entropy increases in the order $S_{m.\mathrm{s}\mathrm{o}\mathrm{l}\mathrm{i}\mathrm{d}}<S_{m.\mathrm{l}\mathrm{i}\mathrm{q}\mathrm{u}\mathrm{i}\mathrm{d}}<S_{m.\mathrm{g}\mathrm{a}\mathrm{s}}$, leading to steeper gradients in the same order. For $T_{\mathrm{m}}<T$, the extrapolation of $\mu$ for each phase shows that the solid phase is the most stable because it has the lowest μ. For $T_{\mathrm{m}}<T<T_{\mathrm{b}}$, the liquid phase is stable because it has the lowest μ, while, for $T>T_{\mathrm{b}}$, the gas phase is stable because it has the lowest μ. The phase transition point can be observed at the intersection of the $\mu$ of each phase. These results are logical consequences of the principle of maximum entropy.

In nonequilibrium, if the emergent flux $J_{\perp }$ is expressed as a function of the global driving force $X_G$, each state can be described by a different function, allowing the nonequilibrium state to be defined as a phase, analogous to how the chemical potential defines phases in equilibrium systems. As $X_G$ increases, $J_{\perp }$ increases, which can be accompanied by a discontinuity change at the transition point and an increase in slope. Therefore, the system evolves into a state with higher entropy production. Therefore, to achieve a higher entropy production rate, the system interferes with other irreversible processes and discontinuously changes into a state with higher flux. This leads to a state transition from the ASF to HTW$\alpha$ and further to HTW$\beta$.

4. Conclusions

This study quantitatively analyzed the state transitions of the flow in a thermally convective system under an applied thermal gradient and rotational field using the MEPP. Using a heat flux orthogonal to the applied thermal gradient, the three states of thermal convection, ASF, HTW$\alpha$, and HTW$\beta$, were described by distinct straight lines. The transition points between states are identified by the intersections of the entropy production curves for each state. This suggests that, in nonequilibrium, state transitions are driven by changes in flux, as these changes lead to faster entropy production. The MEPP may have certain universality, because it is found not only in thermal convection phenomena but also in various transport phenomena, such as Marangoni convection and Saffman–Taylor instability. Because of the analogy with equilibrium thermodynamics in terms of phase description and stability, a theorem similar to the equilibrium system may be derived for the nonequilibrium system. The verification of the theorem will lead to the further development of the MEPP.