Global Stability of the Curzon-Ahlborn Engine with a Working Substance That Satisfies the van der Waals Equation of State

In this paper, we show an analysis of the global stability of a Curzon–Ahlborn engine considering that the working substance of the engine satisfies the Van der Waals equation of state, which is more general than the ideal gas case. We use the Lyapunov stability theory for the case where the engine operates at a maximum power output. We analyze the steady state of the intermediate temperatures as well as the asymptotic behavior of the performance of the engine. Additionally, we study the relationship between the inherent time delay by analyzing the dynamic properties of the system and the stability of the steady state. We present illustrative graphs of the obtained results. Finally, we include a brief discussion of the obtained results and appropriate conclusions.


Introduction
Since the Curzon-Ahlborn paper, published in 1975 [1], the so-called finite-time thermodynamics (FTT) started its development as a branch of classical equilibrium thermodynamics. Diverse ways and methods to qualify the performance of the Curzon and Ahlborn engines and other heat engines have been obtained [2][3][4][5][6][7]. One of the main objectives of FTT is the analysis of the optimal operation of thermal engines to obtain new operating limits of efficiency and power output in thermal devices. That is, there was no thermal equilibrium between the working fluid and the thermal reservoirs in the isothermal branches of the cycle. Curzon and Ahlborn (CA) showed that such a heat engine provides nonzero power, positive entropy production, and a more realistic efficiency than the Carnot engine. In Figure 1, we show the representation of a CA engine, in which T 1 and T 2 are the temperature of the heat reservoirs with T 1 > T 2 , J 1 and J 2 are the heat fluxes of the engine through thermal conductors, both with the same conductance (α) and the same heat capacity (C), and x and y are the hot and cold temperatures of the Carnot cycle isothermal branches, respectively. The efficiency of the CA engine ( Figure 1) at maximum power output is given by 1 √ , (with ⁄ and ). This efficiency is a better approximation for real heat engines than the one proposed by Carnot. Most studies in FTT of thermal systems have focused on their steady-state energy properties, which are important from a design The efficiency of the CA engine ( Figure 1) at maximum power output is given by η CA = 1 − √ τ, (with τ = T 2 /T 1 and T 1 > T 2 ). This efficiency is a better approximation for real heat engines than the one proposed by Carnot. Most studies in FTT of thermal systems have focused on their steady-state energy properties, which are important from a design viewpoint. FTT studies of the dynamic properties [3,4,7] of the system have also been conducted, such as the response to a noisy disturbance or the stability of the steady state of the system. Recently, studies have proposed to analyze the stability and robust dynamics of thermal engine models. In real systems, when modifying the value of a variable, the change in the system is not appreciated immediately, but a certain amount of time must lapse to observe the effect of said change in the dynamic response of the system. These systems are known as time dynamic systems. The robustness of a system can manifest itself in one of two different ways: the system returns to its current dynamic attractor after the disturbance, or it moves to a different attractor that maintains the function of the system.
Recently, the distinct types of stability problems in dynamical systems have been studied through stability theory, which plays a significant role in systems and engineering studies. In 2001, Santillan et al. [8] conducted the first work on the local stability analysis of a CA endoreversible engine working at maximum power output. These authors showed that a clever design generates stability between dynamic robustness and optimal thermodynamic properties. Páez-Hernández et al. [9,10] studied the effects of robustness and time delay on the performance and local stability properties of a non-endoreversible Curzon-Ahlborn engine. Páez-Hernández et al. [11] studied dynamic properties in an endoreversible Curzon-Ahlborn engine using a Van der Waals gas as a working substance. Páez-Hernández et al. [12] studied the effects of time delays in endoreversible and nonendoreversible thermal engines working at different regimes. Chimal-Eguía J. C. et al. [13] analyzed the stability of a Curzon-Ahlborn engine with the Stefan-Boltzmann heat transfer law. Guzman-Vargas et al. [14] performed local stability analysis of the CA engine under different heat transfer laws. Barranco-Jiménez M. A. et al. [15] analyzed the local stability of a thermo-economic model of a Novikov-Curzon-Ahlborn heat engine. Other authors have studied the stability under different operating regimes, such as the regime of maximum efficiency and maximum ecological function [16][17][18][19][20][21][22][23][24]. It is important to note that the finite time thermodynamics continues to address more systems of current interest, as shown in [25][26][27][28][29][30].
Part of the stability theory is the equilibrium points that are generally analyzed by Lyapunov's stability theorem. An equilibrium point is said to be stable if all solutions that start in the vicinity of the equilibrium point remain in the vicinity of the equilibrium point; otherwise, the equilibrium point is unstable. If a stable point fulfills that all the solutions that start in the vicinity of the stable point remain in the vicinity of the equilibrium point, and also tends towards equilibrium as time approaches infinity, then the equilibrium point is said to be asymptotically stable.
Different models of thermal engines can be studied [31][32][33][34] using the Lyapunov stability theory for the analysis of the dynamic properties of the system. In these studies, it is assumed that the working substance of the engine is an ideal gas. However, in real engines, this substance is a non-ideal gas. In this paper, we analyze the influence of assuming a non-ideal gas as the working substance in the engine, particularly a Van der Waals gas. We note that the efficiency of the CA engine, at maximum power output, can be written as the limit at zero of a series of powers of a factor involving the compression ratio, as shown by Gutkowicz-Krusin et al. [35]. In the case of a Van der Waals gas, this efficiency has the same form, and only differs in that the noted factor includes the available volume of one mole of particles; Ladino-Luna D [36]. We use the Lyapunov stability theory and show the global stability analysis of a Curzon-Ahlborn engine with a working substance that obeys the Van der Waals equation of state when the engine operates at maximum power output. We analyze the steady state of the intermediate temperatures, as well as the delay times of the machine. The work is organized as follows. In Section 2, we present the steadystate characteristics of the Curzon and Ahlborn engine with a Van der Waals gas (CAEVW). The stability analysis and effects of the time delay of the CAEVW are presented in Section 3. The global stability of the CAEVW is presented in Section 4. For each important result, we show illustrative graphs of the obtained results. Finally, in Section 5, important conclusions and remarks are presented.

Steady-State Characteristics of the CAEVW
The heat engine, shown in Figure 1, operates between temperatures T 1 and T 2 , with T 1 > T 2 . We can consider the steady-state temperatures as x and y, where Considering a Newton heat transfer law, where the heat exchange is carried out through the thermal conductors, both with conductance α, the irreversible fluxes from T 1 to x and from y to T 2 are J 1 and J 2 , respectively, where J 1 and J 2 represent the steady-state heat fluxes of the engine, and The following inequality follows from the Clausius theorem and the fact that the inner Carnot-like engine works in irreversible cycles: Considering the endoreversibility hypothesis, the above equation results in The system's steady-state power output and efficiency can be defined as and By combining Equations (1)-(3), and (5), we can solve for x and y in terms of T 1 , T 2 and η as where τ = T 2 /T 1 .
In our analysis, we use a working substance that obeys the Van der Waals equation of state, which, unlike the ideal gas equation PV = nRT, considers a correction in the volume of the molecules, which are treated as small hard spheres, leading to a minimum volume b (for 1 mol). That is, in the ideal gas equation V is replaced by V − b. Additionally, the Van der Waals equation considers the attraction between molecules, which tends to decrease the volume, considering a small superficial layer of gas; its attraction for the internal molecules should be proportional to the square of the density.
Thus, the efficiency of the Curzon-Ahlborn heat engine, operating at maximum power output and considering a Van der Waals gas working substance following [36], is given by With where γ is the ratio of the heat capacities, V max and V min are the volumes supplied by the gas in one cycle, and b is a gas-dependent constant. Nevertheless, the b parameter in the Van der Waals equation is so small compared with the spanned volumes of the engine, namely V max and V min . Therefore, we have , and V max V min = r c is the compression ratio of the engine, which has a particular value depending on the chosen cycle: Diesel, Otto, etc. In the case of the Diesel cycle, 12 ≤ r c ≤ 15; for the Otto cycle, 2 ≤ r c ≤ 15 [37]. Considering r c = 12 and substituting in Equation (9), we obtain λ vw = 1.006, for the simplicity calculus we use λ vw = 1, now substituting this value in (8), we obtain the efficiency η VW = η C /2. This substitution of lambda is justified by observing the behavior of the plots in Figure 2. It is even shown an interval of τ, namely τ [0.322, 0.674], in which real plants usually work [38,39]. As we can see, the behavior of η VW and η C /2 is practically the same in that interval; so there is no loss of generality in the present study if we suppose λ vW = 1, and we obtain where η C represents the Carnot efficiency. Substituting Equations (6) and (7) into Equation (10), to express the steady state to the CAEVW and the results and tends to decrease the volume, considering a small superficial layer of gas; its attraction for the internal molecules should be proportional to the square of the density. Thus, the efficiency of the Curzon-Ahlborn heat engine, operating at maximum power output and considering a Van der Waals gas working substance following [36], is given by where is the ratio of the heat capacities, and are the volumes supplied by the gas in one cycle, and b is a gas-dependent constant. Nevertheless, the b parameter in the Van der Waals equation is so small compared with the spanned volumes of the engine, namely and . Therefore, we have ≅ , and is the compression ratio of the engine, which has a particular value depending on the chosen cycle: Diesel, Otto, etc. In the case of the Diesel cycle, 12 15; for the Otto cycle, 2 15 [37]. Considering 12 and substituting in Equation (9), we obtain 1.006, for the simplicity calculus we use 1, now substituting this value in (8), we obtain the efficiency 2 ⁄ . This substitution of lambda is justified by observing the behavior of the plots in Figure 2. It is even shown an interval of , namely 0.322,0.674 , in which real plants usually work [38,39]. As we can see, the behavior of and 2 ⁄ is practically the same in that interval; so there is no loss of generality in the present study if we suppose 1, and we obtain where represents the Carnot efficiency. Substituting Equations (6) and (7) into Equation (10), to express the steady state to the CAEVW and the results For the power output of the CAEVW in the steady state, by combining Equations (1), (5) and (6) into Equation (10), we obtain To perform the stability analysis, we write T 1 , T 2 , and P in terms of x and y, that is

Dynamic Model of the CA Engine
As was noted, Figure 1 is the representation of a CA engine, in which T 1 and T 2 are heat reservoirs with T 1 > T 2 ; J 1 and J 2 are the heat fluxes of the engine through thermal conductors, both with the same conductance (α) and the same heat capacity (C); x and y are the hot and cold temperatures of the Carnot cycle isothermal branches, respectively. From the capacity definition, C = dQ dT , we obtain the dynamics equation, dx dt = 1 C dQ dt . Now, considering the rate change of the internal temperatures x and y, and the thermal conductance α, these temperatures can vary with time, according to the following differential equations: with J 1 the heat flow from x to the working substance and J 2 the flow from the Carnot engine to y. Both equations cancel when x, y, J 1 , and J 2 take their steady-state values. We write J 1 and J 2 in terms of x, y and the output power of the engine, Assuming that the power output of the CA engine is related to the temperatures x and y, as in the case of the power output for steady-state that is related to the temperatures x and y, then we can write the power output of the CA engine out of steady-state in terms of x and y of the form Now, by substituting Equations (19)-(21) into Equations (17) and (18), we obtain the following expressions in terms of the temperatures x, y, T 1 , and T 2 ,

Local Stability Analysis
To perform a local stability analysis of the CAVEW dynamic model, we transform Equations (22) and (23) as the following set of ordinary differential equations: dy dt = g(x, y).
The system's steady state is couple (x, y) that simultaneously satisfies f (x, y) = 0 and g(x, y) = 0. With a little algebra, it is easy to show that the Curzon-Ahlborn engine has Entropy 2022, 24, 1655 6 of 11 a unique steady state, given by Equations (11) and (12). Based on Strogatz [40], the local stability in the steady state is determined by the eigenvalues of the Jacobian matrix: where Following Páez-Hernández et al. [9] step by step, both eigenvalues have negative real parts, which implies that the system's steady state is a stable node, and allows us to define relaxation times t i = −1/λ i , i = 1, 2. The time evolution of a given perturbation from the steady state is generally determined by both relaxation times. However, the perturbation long-term behavior is dominated by the longest relaxation time. In Figure 3, the relaxation times t 1 and t 2 are plotted against τ. Notice that t 1 > t 2 for all values of τ, that t 2 is a monotone-increasing function of τ, and that t 1 has a maximum at τ = 0 and a minimum at τ = 1.
The system's steady state is couple ̅ , that simultaneously satisfies ̅ , 0 and ̅ , 0. With a little algebra, it is easy to show that the Curzon-Ahlborn engine has a unique steady state, given by Equations (11) and (12). Based on Strogatz [40], the local stability in the steady state is determined by the eigenvalues of the Jacobian matrix: where ̅ , Following Páez-Hernández et al. [9] step by step, both eigenvalues have negative real parts, which implies that the system's steady state is a stable node, and allows us to define relaxation times 1/ , i = 1, 2. The time evolution of a given perturbation from the steady state is generally determined by both relaxation times. However, the perturbation long-term behavior is dominated by the longest relaxation time. In Figure 3, the relaxation times and are plotted against . Notice that for all values of , that is a monotone-increasing function of , and that has a maximum at 0 and a minimum at 1.

Dynamic Effects of the Time Delay
Consider a small perturbation near the steady state. It is determined by the linear combination based on exp , i = 1, 2, …, being solution of the characteristic equation, McDonald [41],

Dynamic Effects of the Time Delay
Consider a small perturbation near the steady state. It is determined by the linear combination based on exp(λ i t), i = 1, 2, . . . , being λ i solution of the characteristic equation, McDonald [41], Therefore, we can consider the delay differential system equations where π/2 represents the variable with time delay. Because the characteristic equation has an infinite number of solutions, the stability analysis of a dynamical system with time delay can be complicated. Furthermore, as it is Entropy 2022, 24, 1655 7 of 11 known, a common effect of time delays is to destabilize previously stable steady states, as it induces sustained oscillations. In Páez-Hernández et al. [12], it is shown that no time delay can destabilize the CAEVW. Figure 4 shows that surfaces ρ 1 and ρ 2 , which result from solving Equation (33), do not intersect at any point.
Therefore, we can consider the delay differential system equations where /2 represents the variable with time delay. Because the characteristic equation has an infinite number of solutions, the stability analysis of a dynamical system with time delay can be complicated. Furthermore, as it is known, a common effect of time delays is to destabilize previously stable steady states, as it induces sustained oscillations. In Páez-Hernández et al. [12], it is shown that no time delay can destabilize the CAEVW.

exp
(33) Figure 4 shows that surfaces and , which result from solving Equation (33), do not intersect at any point.

Lyapunov Function Description
In this section, we apply the Lyapunov stability theory to analyze and check the global asymptotic stability for the CA heat engine operating at maximum output power by constructing the Lyapunov function. This function, denoted by , , must satisfy the following conditions: (1) , must be positive definite in a region around the steady state; (2) , must be negative definite, that is, , 0.

Lyapunov Function for CA Engine Model
To build the matrix , we consider the right-hand side of the differential equations given by Equations (22) and (23) and after evaluating the steady-state temperatures in terms of , i.e., Equations (11) and (12).
Considering that the previous matrix is negative definite, the Lyapunov function is written as

Lyapunov Function Description
In this section, we apply the Lyapunov stability theory to analyze and check the global asymptotic stability for the CA heat engine operating at maximum output power by constructing the Lyapunov function. This function, denoted by V(x, y), must satisfy the following conditions: (1) V(x, y) must be positive definite in a region around the steady state; (2) . V(x, y) must be negative definite, that is, .

Lyapunov Function for CA Engine Model
To build the matrix J(x), we consider the right-hand side of the differential equations given by Equations (22) and (23) and after evaluating the steady-state temperatures in terms of τ, i.e., Equations (11) and (12).
Considering that the previous matrix is negative definite, the Lyapunov function is written as It can be verified in the previous equation that in the stationary state V(x, y) = 0, if we consider x and y.
On the other hand, the function fulfills .
x + ∂V ∂y ∂V ∂y Figure 5 shows the behavior of the level curves of the Lyapunov function built for our dynamic system defined by Equations (22) and (23). The behavior indicates that, as the constant k of the level curve decreases, the surface also decreases towards the steadystate values (x, y). Our result is consistent with the work reported by Guzmán-Vargas L et al. [33].

3
It can be verified in the previous equation that in the stationary state ̅ , 0, if we consider ̅ and .
On the other hand, the function fulfills (38) Figure 5 shows the behavior of the level curves of the Lyapunov function built for our dynamic system defined by Equations (22) and (23). The behavior indicates that, as the constant k of the level curve decreases, the surface also decreases towards the steadystate values ̅ , . Our result is consistent with the work reported by Guzmán-Vargas L et al. [33]. In Figure 6, the surface of the Lyapunov function is shown. Comparably to these authors [33], we conduct the numerical study of the decay of the disturbances around the stationary state through Equations (17) and (18), for different initial conditions. In Figure 6, the surface of the Lyapunov function is shown. Comparably to these authors [33], we conduct the numerical study of the decay of the disturbances around the stationary state through Equations (17) and (18) (38) Figure 5 shows the behavior of the level curves of the Lyapunov function built for our dynamic system defined by Equations (22) and (23). The behavior indicates that, as the constant k of the level curve decreases, the surface also decreases towards the steadystate values ̅ , . Our result is consistent with the work reported by Guzmán-Vargas L et al. [33]. In Figure 6, the surface of the Lyapunov function is shown. Comparably to these authors [33], we conduct the numerical study of the decay of the disturbances around the stationary state through Equations (17) and (18), for different initial conditions.

Results
In the numerical analysis shown in Figure 7, the initial state is (x 0 , y 0 ), and we consider the following cases: (1) x 0 > x, y 0 > y; (2) x 0 > x, y 0 < y; (3) x 0 < x, y 0 > y; and (4) x 0 < x, y 0 < y. For all cases analyzed, it is observed that they tend to the value of the steady state of the temperature.

Results
In the numerical analysis shown in Figure 7, the initial state is , , and we consider the following cases: (1) ̅ , ; (2) ̅ , ; (3) ̅ , ; and (4) ̅ , . For all cases analyzed, it is observed that they tend to the value of the steady state of the temperature.

Conclusion and Remarks
In this work, we performed a global stability analysis of a CA engine using a Van der Waals gas as a working substance and by applying the Lyapunov stability theory. For the global analysis, we first studied the dynamics of the CA engine to shape the Lyapunov function, and then we analyzed the global properties. Our results show that the stability of the engine is maintained around the stable state even when variations in parameters of the Lyapunov function are considered, i.e., the surface also decreases towards the steadystate values ̅ , , as reported by other authors. In addition, we studied the local analysis of a CA engine using a Van der Waals gas as a working substance and have considered the engine's inherent time delays. Time delays are present in multiple systems subject to dynamic regulation. In the endoreversible and non-endoreversible Curzon-Ahlborn engine, the inherent time delays are incapable of destabilizing the steady state [9,12]. Thus, they seem not to play a role in the trade-off between energetic and dynamic properties. This does not have to be true for all energy-converting systems, however. They should be taken into account in further studies concerning the relationship between the dynamic and the energetic properties of energy-converting systems.