 Research article
 Open Access
 Published:
Optimizationbased investigations of a twophase thermofluidic oscillator for lowgrade heat conversion
BMC Chemical Engineering volume 1, Article number: 20 (2019)
Abstract
Background
The noninertivefeedback thermofluidic engine (NIFTE) is a twophase thermofluidic oscillator capable of utilizing heat supplied at a steady temperature to induce persistent thermalfluid oscillations. The NIFTE is appealing in its simplicity and ability to operate across small temperature differences, reported as low as 30 ^{∘}C in early prototypes. But it is also expected that the NIFTE will exhibit low efficiencies relative to conventional heat recovery technologies that target highergrade heat conversion. Mathematical modeling can help assess the full potential of the NIFTE technology.
Results
Our analysis is based on a nonlinear model of the NIFTE, which we extend to encompass irreversible thermal losses. Both models predict that a NIFTE may exhibit multiple cyclic steady states (CSS) for certain design configurations, either stable or unstable, a behavior that had never been hypothesized. A parametric analysis of the main design parameters of the NIFTE is then performed for both models. The results confirm that failure to include the irreversible thermal losses in the NIFTE model can grossly overpredict its performance, especially over extended parameter domains. Lastly, we use the NIFTE model with irreversible thermal losses to assess the optimization potential of this technology by conducting a multiobjective optimization. Our results reveal that most of the optimization potential is achievable via targeted modifications of three design parameters only. The Pareto frontier between exergetic efficiency and power output is also found to be highly sensitive to these optimized parameters.
Conclusions
The NIFTE is of practical relevance to a range of applications, including the development of solardriven pumps to support smallholder irrigation in the developing world. Given its low capital cost, potential improvements greater than 50% in efficiency or power output are significant for the uptake of this technology. Conventional heat recovery technologies are known to have higher efficiencies than those reported in this work, but they also have more complex designs and operations, higher capital costs, and may not even be feasible for the temperature differences considered herein. Future work should focus on confirming this modelbased assessment via dedicated experimental campaigns and on investigating design modifications to mitigate irreversible thermal losses.
Background
Hightemperature heat engines are ubiquitous in the power generation industry, stationary distributed power (co)generation, the transport and many other sectors. This demand for highgrade heat has led to ever increasing rates of fossilfuel consumption and resource depletion, while at the same time concerns exist in relation to the emissions associated with the combustion of these fuels. There is therefore a pressing need for exploiting alternative heat sources.
Broadly speaking, heat available at a relatively low temperature and for which there is no direct use is classified as “lowgrade”. Such heat is abundant in the process and chemical industries, where it is often “wasted” by release into the environment without further utilization. Large amounts are also generated and available in solar collectors and geothermal wells, among other applications. Unfortunately, all technologies for lowgrade heat conversion to power are associated with low thermal efficiencies [1, 2], a restriction imposed by the second law of thermodynamics [3].
Thermofluidic oscillators are a class of engines for converting heat into useful work via the creation of oscillatory fluid motion. The noninertivefeedback thermofluidic engine (NIFTE) is a particular type of twophase thermofluidic oscillator [4–6]. The description of the configuration and stepbystep operation of the NIFTE are given as part of Fig. 1. The working fluid exists in both the liquid and vapor phases, and changes phase periodically during the operation (evaporation and condensation). It can be tailored to suit the available external heat source and sink temperatures. A NIFTE is simple to construct and has few moving parts in comparison with conventional engines, making it both a reliable and low capitalcost technology, that is capable of converting lowgrade heat.
It has been established experimentally on early prototypes [7] that a NIFTE is capable of operating with a temperature difference as small as 30 ^{∘}C between the heat source and sink, and more recent versions of the technology may be able to operate with temperature differences even smaller than this. An early prototype was reported as achieving thermal efficiencies up to 1% and exergetic efficiencies—the thermal efficiency as a fraction of the Carnot efficiency—between 1–2%, when operated with a heat source at 65–90 ^{∘}C, a heat sink at 4–12 ^{∘}C, and npentane as the working fluid [6, 8].
Such efficiencies can be considered low in comparison to those of other heat conversion technologies that are typically designed and intended for highergrade heat sources, such as organic Rankine cycle (ORC) systems, Stirling and thermoacoustic engines [3]. In general, thermal efficiencies in the approximate range of 10–20% are expected from ORC systems, also depending on the system size, the lower end achievable at heat source temperatures of about 100–150 ^{∘}C and the upper end at up to 400–500 ^{∘}C [9]; while Stirling and thermoacoustic engines have much lower efficiencies at those temperatures, but can reach thermal efficiencies between 20–30% at heatsource temperatures of 600–800 ^{∘}C [3]. In one example [10], a smallscale regenerative ORC power system comprising a scroll expander and with R123 as the working fluid was reported as having a thermal efficiency upwards of 7% at a heatsource temperature of 165–183 ^{∘}C and an ambient air temperature of 22.5 ^{∘}C cooling the condenser. In a travelingwave thermoacoustic engine with air as the working gas [11], a thermal efficiency of 14.5% was reported under 25 bar, with a heatsource temperature of 650 ^{∘}C and a cooling temperature of 12 ^{∘}C. Importantly, it can be seen that these conventional power generation systems are neither technically suitable nor economically feasible when the heatsource temperature drops below 100 ^{∘}C. Furthermore, the construction and operation of these systems are invariably more complicated than those for a thermofluidic oscillator, which has no or few moving parts or mechanical components and dynamic seals.
An alternative type of twophase thermofluidic oscillator, known as the UpTHERM engine, has also been recently investigated for heat conversion applications [12, 13]. A labscale UpTHERM engine had a predicted exergetic efficiency of about 0.8%, when operating with a heat source and sink at 450 ^{∘}C and 10 ^{∘}C, respectively, and npentane as the working fluid [14]. It is also worth noting that an optimized pilotscale UpTHERM engine with the same working fluid was predicted to yield promising exergetic efficiencies between 10–40% for a heat source temperature in the range 200–400 ^{∘}C [13], which is close to the efficiency of its ORC counterparts and could even be higher with optimally selected working fluids [15].
Nearly all of the mathematical models proposed so far to analyze the dynamical behavior of the NIFTE feature linear dynamics. Smith [4–6] first developed a linear model, based upon earlier work on thermoacoustic engines [16] and electrical analogies for thermal and fluid processes [17, 18]. In this early work, the NIFTE components were described by linear equations in analogy with electrical components—namely, resistors to account for fluid drag and thermal resistance, and capacitors to represent hydrostatic pressure and vapor compressibility—and the heatexchange rate between the hot and cold elements and the working fluid was assumed to be either linear or constant. A complete engine model then connects these electrical components into a circuit according to the actual thermalfluid technology layout. Inductors were later introduced in the electrical analogy to account for fluid inertia [7, 19], which allowed improvements to the accuracy of the predictions, and an effort to consider performance over wide parameter ranges was also made [20].
Notwithstanding their ability to provide quantitative predictions of the NIFTE performance via numerical simulation at marginal stability, these linear models cannot describe the sustained and robust periodic oscillations of the working fluid as observed in practice, which in essence is a nonlinear phenomenon. Subject to a steady temperature difference between the two heatexchangers, the engine reaches a cyclic steady state (CSS) after a sufficiently long settling time, where the amplitude and frequency of the oscillations remain constant. When variations in the heatexchanger temperature difference or other external disturbances occur, which is inevitable in practice, the engine will eventually again return to a CSS after some time. This is the signature of an asymptotically stable limit cycle in the dynamics, which may only occur in nonlinear dynamic systems.
Inspired by earlier work on exploiting the nonlinearity of core physical elements to predict the efficiency and robustness of thermodynamic engines [21–24], Markides and coworkers [25] were the first to propose a nonlinear inertive model of the NIFTE, referred to as the NTP model. The key nonlinearity in this model corresponds to the physical saturation in heat transfer rate between the hot and cold elements and the working fluid, whereas the rest of the engine components are described using the same linear electrical analogy used in earlier work.
Failure to account for various irreversible thermal losses in the early NIFTE models can lead to an overestimation of the engine’s efficiency and work output by an order of magnitude or more. Smith [6] identified the two predominant exergy losses to be: (i) shuttle losses, associated with the irreversible twoway heat transfer between the working fluid and the cylinder walls; and (ii) entrance and retrograde condensation losses, caused by the alternating parasitic condensation and evaporation of the working fluid. Solanki and Markides [8, 26] were the first to account explicitly for these losses in a linear NIFTE model, which led to much more reliable efficiency predictions. However, a similar extension of the NTP model to encompass exergy losses has not yet been reported in the literature.
This paper describes a detailed, modelbased assessment of the performance of the NIFTE technology. We base our assessment on the existing NTP model, and propose a model extension, referred to as NTPTL, that encompasses the irreversible thermal losses. Our methodology entails a modelbased analysis of key design parameters—via the resistors, capacitors, inductors and other gain values in the electrical analogy—in terms of exergetic efficiency and average power output under CSS operation, that leads to new insight into the operation and limits of the NIFTE. We also perform a multiobjective optimization based on the NTPTL model to explore the optimal tradeoff between efficiency and power output, and gain further insight into NIFTE’s optimization potential.
Mathematical modeling
The two main thermalfluid transport processes that govern the operation of a NIFTE device are: (i) pressure differences, driving workingfluid flows through fluid components; and (ii) temperature differences, responsible for heat flows in thermal components. We provide a summary of the models developed to describe these processes in the fluid and thermal domains first, followed by a complete statement of both the NTP and NTPTL device models, and a description of the performance indicators considered for the assessment. A nomenclature of the notation used for mathematical modeling can be found in Table 1. For further details about the model development, refer to the original paper [25].
Fluid domain
The main components ➁, ➅, ➆, ➇, ➈ of the engine (cf. labels in Fig. 1) are modeled based on electrical analogies. The drag experienced by the working fluid in the load ➁ (due to viscosity) and in the feedback connection ➅ (due to small diameter and/or a valve) are represented by resistors, R_{i}. In particular, the model assumes that the fluid flow in these elements is akin to viscous laminar flow in a smooth pipe. The inertive effect experienced by the oscillating liquid flow in the load ➁, feedback connection ➅, power cylinder ➆, and displacer cylinder ➇ are modeled by inductors, L_{i}. The hydrostatic pressures at the bottom of the power cylinder ➆ and displacer cylinder ➇, and the adiabatic compression/expansion in the vapor region ➈ are represented by capacitors, C_{i}. The processes of vapor compression and expansion undergone by the working fluid in the vapor region ➈ are assumed to be isentropic (i.e. both adiabatic and reversible).
These R,L,C parameters may be expressed in terms of the physical and geometric properties of the NIFTE, namely density, ρ and dynamic viscosity, μ of the working and pumped fluids (indexed with _{wf},_{pf}, respectively), and diameter, d, length, l and crosssectional area, A of the tubes and valves (indexed with _{l},_{f},_{p},_{d},_{ad} in reference to the components ➁, ➅, ➆, ➇, ➈, respectively):
with V^{∘} and P^{∘}, the equilibrium (timeaveraged) volume and pressure in the gas phase; γ, the adiabatic index; and g, the gravitational acceleration. All of the physical properties, and therefore the R,L,C parameters themselves, are assumed to be temperature independent. But for a particular fluid, the R,L,C parameter values can be adjusted by varying the geometric properties (length l, diameter d, and/or crosssectional area A) of the corresponding elements, which provides the basis for the parametric analysis and optimization in the Results section below.
The electrical analogy assumes linearized firstorder governing equation between the volumetric flowrate of the working fluid, U_{i} through a given component i and the pressure difference, ΔP_{i} across that component:
with Z_{i} the impedance of component i. The impedance for each of the main components of the engine is expressed in terms of resistance, inductance and capacitance in the Laplace domain, using the principle of superposition to describe multiple effects occurring within the same component:
These R,L,C elements, as arranged in a complete NIFTE circuit, may be visualized in Fig. 2, for which further details will be provided below.
Thermal domain
The thermal domain, which comprises the hot and cold heatexchangers ➃, ➄ within the displacer cylinder ➇, converts heat into pressurization and fluid flow. We adopt the usual assumptions that: (i) the spatial temperature profile along the vertical walls of the heatexchangers is externally imposed as a boundary condition to the device and does not vary with time; and (ii) the convective heat transfer associated with phase change dominates over the forced convection taking place away from an active region near the vaporliquid interface of the working fluid within which phasechange heat transfer occurs. Under these conditions, the heat transferred to the working fluid, \(\dot {X}_{\text {th}}\) is given by:
with \(\dot {S}_{\text {th}}\), the associated rate of change of entropy; T^{∘}, the equilibrium (timeaveraged) temperature; h_{th}, the (constant) phasechange convective heat transfer coefficient; A_{th}, the surface area over which the heat transfer takes place; T_{ad}, the temperature of the working fluid; and T_{th}, the heatexchanger wall temperature experienced by the working fluid at the vaporliquid interface ➂, y.
The connection between Eq. 3 and the fluid domain proceeds by: (i) relating the entropy rate of change, \(\dot {S}_{\text {th}}\) to an equivalent volumetric flowrate of vapor (due to phase change), U_{th} via
with ρvap∘, the equilibrium (timeaveraged) density of the vapor, and Δs_{vap}, the specific entropy of vaporization; and (ii) relating temperatures to pressures as
with \(\left (\frac {\mathrm {d}T}{\mathrm {d}P}\right)_{\text {sat}}\), the rate of change in saturation temperature of the working fluid with respect to the saturation pressure. Using the same electrical analogy as in Eq. 2, we obtain:
where the thermal impedance Z_{th} consists of the thermal resistance between the working fluid and the walls of the heatexchangers,
The heatexchanger wall temperature T_{th}, and the pressure P_{th} likewise, vary with the vertical position y of the vaporliquid interface ➂, and can be related to the pressure variable P_{d}=ρ_{wf} g y. In particular, the function f in Eq. 4 may taken as either linear [6, 7] or nonlinear [25] depending on modeling assumptions.
In addition to the thermal resistance R_{th}, two important thermal loss mechanisms in the displacer cylinder ➇ are the shuttle loss and the entrance and retrograde condensation loss [6]. The former is associated with irreversible alternating conductive heat transfer between the working fluid and the cylinder walls, which may be modeled by a heat diffusion equation with either an isothermal or adiabatic boundary condition. The latter is caused by the alternating parasitic condensation and evaporation of the working fluid via convective heat transfer. Following the derivations in Refs. [8, 26], a total impedance representing these thermal losses can be solved in the Laplace domain, to arrive at the expressions:
where the thermal loss parameters A_{tl},B_{tl} and R_{tl,0} are given by
with α_{gl},κ_{gl},l_{gl} and Σ_{gl}, the thermal diffusivity, thermal conductivity, length and surface area of the glass component, respectively. We reiterate that the model assumes the physical properties to be temperature independent.
In order to facilitate expressing the dynamic system in the time domain, we consider an RCcircuit approximation for the nonlinear expression of Z_{tl} herein, in the form:
where n denotes the number of parallel RC components in series (see RC circuit depicted in red on Fig. 2). The resistance and capacitance values R_{tl,i} and C_{tl,i} may be estimated to match Z_{tl} over a suitable frequency range.
Complete NIFTE models
A graphical depiction of the R,L,C elements within each component of both the fluid and thermal domains is shown in Fig. 2. The original NTP model [25] comprises the components shown in black only. Addition of the components in red describing irreversible thermal losses, gives rise to the extended NTPTL model.
A mathematical formulation of the NTP model (Eqs. 7a–7e) from its electrical circuit analogue is obtained by applying Kirchhoff’s voltage and current laws. A minimum of five independent thermodynamic variables is required to fully describe the system, here taken as: the pressure in the adiabatic vapor region, \(\hat {P}_{\text {ad}}\); the hydrostatic pressures in the displacer cylinder, \(\hat {P}_{\mathrm {d}}\) and in the power cylinders, \(\hat {P}_{\mathrm {p}}\); and the volumetric flowrates in the feedback tube, \(\hat {U}_{\mathrm {f}}\) and in the power cylinder \(\hat {U}_{\mathrm {p}}\). The resulting system of ordinary differential equations (ODEs) is furthermore nondimensionalized (denoted with \(\hat {\,}\)), by using P^{∘},U^{∘} and τ as reference values for the pressures, flowrates and time, respectively.
The derivation of the NTPTL model (Eqs. 7a’–7f) proceeds similarly to that of the NTP model above, yet with n extra thermodynamic variables \(\hat {P}_{\text {tl},i}\) representing the pressures associated with each thermal loss capacitance C_{tl,i}. Eq. 7a’ describing the dynamics of pressure \(\hat {P}_{\text {ad}}\) also comprises an extra term compared with Eq. 7a, to account for the thermal losses.
Observe that both NTP and NTPTL models consist of linear ODEs in the state variables, with the exception of the term \(f(\hat {P}_{\mathrm {d}})\) in Eqs. 7a and 7a’ describing the transfer function between \(\hat {P}_{\mathrm {d}}\) and \(\hat {P}_{\text {th}}\). Early models of the NIFTE [4–7, 19] have considered linear (static) transfer functions of the form \(\hat {P}_{\text {th}}=K_{\text {th}}\,\hat {P}_{\mathrm {d}}\), where K_{th} is a proportional gain. An alternative approach that also preserves linearity uses the (dynamic) transfer function \(\hat {P}_{\text {th}}\,s=K_{\text {th}}\,\hat {P}_{\mathrm {d}}\) [8, 27]. But because they are fully linear, none of these models can capture the robust periodic oscillations corresponding to a (stable) limit cycle nonetheless.
To account for the fact that the temperature on the heatexchanger walls T_{th} cannot increase or decrease indefinitely with the liquid level ➂ in the heatexchanger block, one can introduce a nonlinearity in the form of a saturation function as the liquid level in the displacer cylinder moves away from the equilibrium position (which lies halfway between HHX ➃ and CHX ➄). Markides and coworkers [25] proposed the saturation function
where the parameter \(\hat {P}^{\text {max}}_{\text {th}}\) sets the maximum (saturation) temperature amplitude on the heatexchanger walls; and the nonlinear gain K_{th} determines the maximum rate of heat transfer and is a function of the physical distance between the hot and cold heat elements. More specifically, K_{th} is inversely proportional to the distance between the HHX and the CHX; and directly proportional to the temperature difference between the heat source and sink. The resulting NIFTE models are comprised of nonlinear ODEs, which have been shown to feature limit cycles [25].
Key performance indicators
An important indicator of the performance of a NIFTE is the average useful power output from the load, \(\dot {W}_{\mathrm {l}}\) given by
where P_{l} and U_{l} denote the pressure and volumetric flowrate in the load, respectively; and the average value is calculated over a steadystate cycle of the engine.
A second key performance criterion is the exergetic (or, second law) efficiency of the engine, η_{ex} defined as
where P_{th} and U_{th} are the input thermal pressure and volumetric flowrate associated with the phase change, respectively; and \(\dot {X}_{\text {th}}\) the exergy power input into the device.
For given temperatures of the heat source and sink, recall that the Carnot efficiency, η_{ca} is constant, and therefore η_{ex} is proportional to the thermal efficiency of the engine, η_{th}=η_{ex} η_{ca}. The extra pressure and volumetric flowrate states in Eqs. 9–10 are readily determined from the application of Kirchoff’s circuit laws in the electrical circuit analogue (Fig. 2). Further information on the definition and interpretation of these performance indicators may be found in Reference [7].
Methods
Cyclic steadystate detection and analysis
The search for a CSS of the NIFTE can be automated by solving the following constrained dynamic optimization problem,
where T denotes a cycle period; the ODEs righthand side, f corresponds to either the NTP model (Eqs. 7a–7e) with the state vector x comprised of \(\left (\hat {P}_{\text {ad}}, \hat {P}_{\mathrm {d}}, \hat {P}_{\mathrm {p}}, \hat {U}_{\mathrm {f}}, \hat {U}_{\mathrm {p}}\right)\), or the NTPTL model (Eqs. 7a’–7f) with x comprised of \(\left (\hat {P}_{\text {ad}}, \hat {P}_{\mathrm {d}}, \hat {P}_{\mathrm {p}}, \hat {P}_{\text {tl}}, \hat {U}_{\mathrm {f}}, \hat {U}_{\mathrm {p}}\right)\); and p is the vector of model parameters (cf. Table 2).
Any solution \(\overline {\mathbf {x}}(\cdot,\mathbf {p})\) to the problem (11) with T>0 corresponds to a CSS of the NIFTE. However, this optimization problem is overparameterized, in the sense that there exists infinitely many initial conditions x_{0} describing a given CSS with T>0. This issue is easily circumvented – and without loss of generality – by forcing one of the initial states to zero; e.g., \(\hat {P}_{\text {ad}}(0)=0\).
The optimization problem (11) is nonconvex and may therefore present multiple (local) optima. These optima may correspond to different CSS or equilibrium points of the engine, but they could also be “parasitic” local optima corresponding to multiple periods, e.g. 2T,3T,…, of the same CSS. Herein, we use a simple singleshooting approach as implemented in our inhouse optimization toolkit CRONOS^{Footnote 1} to compute such local optima, and we apply a multistart search heuristic based on Sobol’ sampling to detect the presence of multiple local optima. CRONOS provides an interface to the nonlinear programming (NLP) solver IPOPT [28] coupled with the ODE integrator CVODES in SUNDIALS [29] for computing the objective/constraint functions and their gradients. A local optimization typically converges within a few seconds of CPU time.
The stability of the limit cycle corresponding to a solution \(\overline {\mathbf {x}}(\cdot,\mathbf {p})\) can be characterized by analyzing the socalled monodromy matrix Φ_{p}(T,0), which is the fundamental matrix solution to the variational equation [30, 31]
In practice, the monodromy matrix can be conveniently computed as \(\Phi _{\mathbf {p}}(T,0) = \frac {\partial \overline {\mathbf {x}}(T,\mathbf {p})}{\partial \mathbf {x}_{0}}\), e.g. using the sensitivity analysis capabilities of the integrator CVODES.
One (at least) of the eigenvalues of the monodromy matrix shall always be equal to 1. The limit cycle \(\overline {\mathbf {x}}(\cdot,\mathbf {p})\) is asymptotically stable (attractor) if all of the remaining eigenvalues lie within the open unit disk,
If any of these eigenvalues has a modulus greater than 1, the limit cycle is otherwise unstable (repeller).
Bifurcation occurs when a small change made to the values of a system’s parameters p causes a sudden qualitative change in its asymptotic behavior [31]. Of particular interest in the NIFTE are socalled Hopf bifurcations, whereat the equilibrium position x=0 (origin) loses stability and a periodic solution arises. Such bifurcations can be detected upon repeating the solution of the optimization problem (11) at a large number of parameter values. One can then use dedicated codes to confirm the presence of a bifurcation and refine its location, such as MatCont [32, 33] which implements a numerical continuation approach in MATLAB. The presence of multiple Hopf bifurcations for a particular parameter can be an indication for the coexistence of multiple limit cycles for certain engine configurations.
Multiobjective optimization
Finding the design parameters that maximize both the exergetic efficiency and the average power output of a NIFTE leads to a multiobjective optimization (MOO) problem in the form of
with \(\dot {W}_{\mathrm {l}}\) and η_{ex} as defined in Eqs. 9–10; and using the same notation as in the problem (11). Notice that the initial conditions x_{0} and cycle time T are decision variables alongside the design variables p in order to enforce optimal performance under CSS operation. The solution set to the problem (14) gives rise to a socalled Pareto optimal frontier [34], whereby each point is nondominated in the sense that no better feasible solution exist in terms of the two objectives simultaneously: betterment of the exergetic efficiency compared to a point on the Pareto frontier results in worsening of the average power output, and vice versa.
Herein, we approximate the Pareto frontier of the problem (14) via scalarization, and we proceed by solving the following singleobjective optimization problems for a collection of scalars α∈[0,1]:
Clearly, the problems with α=0 and α=1 correspond to the two extremes of maximizing the average power output and the exergetic efficiency, respectively. It is noteworthy that such a simple scalarization approach fails to describe nonconvex portions of a Pareto frontier in general [34], but our results later on will confirm that this situation does not occur for the Pareto frontiers of interest herein.
In order to tackle the singleobjective optimization problem (15), we furthermore apply a twolevel decomposition: for a given design vector p, the inner problem determines the values of x_{0} and T that yield the CSS with best performance, e.g. in case multiple CSS coexist, or detects that no (nontrivial) CSS is possible; whereas the outer optimization problem operates on the design vector p only, and solves the inner optimization problem at each iteration in order to calculate the performance under CSS operation for a particular instance of p. The inner optimization problems are solved in the same way as the CSS detection problem (11) earlier, i.e. using a combination of local search and multistart search heuristic. As for the outer problem, we implement a simple NelderMead simplex algorithm [35], a derivativefree approach with slow convergence but tolerant to numerical noise, with a termination tolerance on the minimum simplex diameter (set to 10^{−3} in logscale here). In order to rule out convergence to local (suboptimal) solutions and initialize the NelderMead simplex algorithm, the search space is probed using lowdiscrepancy sampling based on Sobol’ sequences [36] (with 10,000 points here), and the NIFTE performance under CSS operation is computed at each sampling point. Algorithm 1 provides a summary of the multiobjective optimization procedure. The parameter range is set to ± 3 orders of magnitude around their nominal values (see below), and the set of weights considered in the scalarized cost of the problem (15) are \(\left \{0,\frac {1}{2000},\frac {1}{1000},\frac {1}{500},\frac {1}{200},\frac {1}{100},\frac {1}{50},\frac {1}{20},\frac {1}{10},\frac {1}{5},1\right \}\). Each run of the NelderMead simplex algorithm in order to determine one point on the Pareto frontier takes up to a few hundred iterations to converge.
Results and discussion
Case study description
The assessment presented in this section is based on a working prototype of the NIFTE. A temperature difference of 80 ^{∘}C is applied across the heatexchanger, with the cold source at ambient temperature and a corresponding Carnot efficiency of about 20%. The pumped and working fluids in this prototype are water and npentane, respectively. The load corresponds to a delivery pipe of 1 m in length and with a diameter of 0.01 m that would in practice be connecting a supply tank of water to the destination tank. Further details concerning this prototype may be found in Refs. [6, 7, 19].
Table 2 reports the nominal values for the electrical analogues R_{i},L_{i} and C_{i} in the NTP model (Eqs. 7a–7e) as well as the parameters \(K_{\text {th}}, \hat {P}_{\text {th}}^{\text {max}}\) in the heat transfer function (Eq. 8) and the reference values P^{∘},U^{∘},τ. These were calculated based on the prototype’s physical properties, and it is noteworthy that the NTP model was tested against experimental data with these parameter values (referred to as SetI in Ref. [25]). In spite of its simplicity compared to complexity of the actual device, the NTP model is capable of earlystage engineering predictions of the NIFTE’s operation and performance.
The same parameter values are used in the NTPTL model (Eqs. 7a’–7f). Values for the extra parameters n, R_{tl,i},C_{tl,i} in the approximate thermalloss impedance (Eq. 6) were determined by fitting the thermalloss impedance with isothermal boundary condition (Eq. 5a), using standard MATLAB tools. In order to obtain accurate gain and phase predictions in the frequency range of interest, 0.1– 100 rad s^{−1} for the NIFTE prototype at hand, a minimum of n=4 RCcircuits in series is necessary here (Fig. 3). The fitted values of R_{tl,i} and C_{tl,i} for i=1…4 are reported in Table 2. Notice that the values of the parameters A_{tl},B_{tl},R_{tl,0} in Eq. 5 used for the calibration were also calculated based on the prototype’s physical properties, reported in Table 2 for completeness. Interestingly, the use of RCcircuit in series to approximate the thermalloss impedance with isothermal boundary condition also yields an excellent agreement for its counterpart with adiabatic boundary condition (Eq. 5b) in the frequency range 0.2– 100 rad s^{−1}. Such close agreement suggests that both types of boundary conditions in Eq. 5 to express the thermalloss impedance should lead to minimal differences in quantitative predictions by the NTPTL model. Both types of boundary conditions are indeed relevant as they are representative of different applications of the NIFTE technology.
Cyclicsteadystate multiplicity and stability
This subsection investigates the asymptotic behavior of a NIFTE based on the NTP and NTPTL models. It has already been established [25], based on the NTP model, that a given engine may or may not exhibit a limit cycle. A key finding herein is that a given engine could in fact exhibit multiple CSS, which could be either stable or unstable, an observation that had never been reported before.
To illustrate this behavior, we first consider a NIFTE with all its parameters at their nominal values (Table 2), except the feedback inductance L_{f}=4.74×10^{5} kg m^{−4} s^{−1} (10fold smaller than the nominal value). When using the NTP model and solving the optimization problem (11) as described earlier, three possible asymptotic behaviors are identified for this NIFTE configuration. The first one is a trivial equilibrium point at the origin, which is unstable. The other two are limit cycles with periods T=1.30 s (CSS1) and T=0.88 s (CSS2). The spectra of the monodromy matrices (Eq. 12) at CSS1 and CSS2 are reported in Table 3. Since all the eigenvalues of the monodromy matrix at CSS1 lie within the unit disk (Eq. 13), this CSS describes a stable limit cycle. By contrast, the monodromy matrix at CSS2 has two of its eigenvalues outside the unit disk (shown in red), and so this is an unstable limit cycle.
The pressurevolume phase diagrams in Fig. 4 are obtained by simulating the NIFTE over an extended time horizon, from different initial states. They represent the pressure P_{l} against the volume V_{l} of the working fluid, given by \(\dot {V}_{\mathrm {l}}=U_{\mathrm {l}}\), and so the perimeter of the cycle limit is proportional to the work output, W_{l} over one complete cycle:
It is clear that CSS1 acts as an attractor for the trajectories, and the origin as a repeller (Fig. 4a). Furthermore, trajectories initialized near CSS2 eventually escape away from that limit cycle and are again attracted by CSS1, thereby showing two clusters of trajectories (Fig. 4b). All of these observations are in agreement with the numerical stability results in Table 3.
The comparison of performance indicators in the bottom half of Table 3 reveals that a 70%greater average useful power \(\dot {W}_{\mathrm {l}}\) (Eq. 9) is generated by operating the NIFTE at the (unstable) CSS2 compared with the (stable) CSS1. Because the engine is attracted to CSS1 when initialized near the origin (Fig. 4a), reaching (and then remaining around) an alternative, potentially better limit cycle such as CSS2 would require a tailored control strategy. This could be done, for example, by manipulating the valve in the feedback connection ➅ (Fig. 1), a dynamic control that has never been applied in practice however. Notice also that CSS2 has a 30%shorter period T and an 8fold greater exergetic efficiency η_{ex} (Eq. 10) than CSS1. The latter is not only the result of the greater useful power output at CSS2, but also a 5fold lower exergy power input \(\dot {X}_{\text {th}}\) (Eq. 10) into the engine from the heat source. Coincidentally, the combination of a shorter period and a larger power output at CSS2 results in about the same work output during a complete cycle of CSS1 (W_{l}=1.03 J) or CSS2 (W_{l}=1.07 J). This can seen on Fig. 4b, where the limit cycles for CSS1 and CSS2 appear to have identical perimeters.
Unlike the NTP model, a simulation of the same NIFTE with the NTPTL model in problem (11) predicts a unique, stable limit cycle (CSS1’). This CSS has a slightly longer period than CSS1, while both the exergetic efficiency and power output are reduced by about 3fold (Table 3). This dramatic reduction is of course attributed to the extra thermal losses in NTPTL, confirming that their inclusion is paramount to the accuracy of the model predictions. The perimeter of the limit cycle CSS1’ on the pressurevolume phase diagrams in Fig. 4c is accordingly much shorter in comparison with that of CSS1 on Fig. 4a.
In the following subsections, a parametric analysis of the NIFTE is carried out by varying (combinations of) the design parameters in the feedback, displacer and power cylinders and heatexchanger, and analyzing the performance of the corresponding CSS using both the NTP and NTPTL models. A multiobjective optimization is also conducted with the NTPTL model in order to investigate the optimal tradeoffs between exergetic efficiency and power input when varying three or more parameters.
Parametric analysis of (R _{f},L _{f})
Previous studies [6, 7] have shown that the performance of NIFTE is highly sensitive to the feedback resistance R_{f} and the feedback inductance L_{f}. These parameters can be varied in a physical device by modifying the geometry of the feedback tube (e.g., l_{f},d_{f} or A_{f}; see Eq. 1) or the position of a feedback valve thereof. This subsection presents a joint parametric analysis of the feedback parameters R_{f} and L_{f} using both models, while all other parameters are held constant at their nominal values (Table 2). We use the superscript ^{+} to denote the ratio between a parameter and its nominal value.
NTP model
Both the exergetic performance η_{ex} and useful power output \(\dot {W}_{\mathrm {l}}\) predicted by the NTP model under CSS operation are shown in Fig. 5, for a wide parameter range spanning several orders of magnitude around the nominal parameter values. By repeating the solution of the optimization problem (11) over this range, we could identify up to two limit cycles, denoted by CSS1 and CSS2 in agreement with previous notation. The stability of each of these CSS is further determined by the eigenvalues of their monodromy matrices (Eq. 12), and we report the separation between stable and unstable CSS on each plot in Fig. 5. Notice that multiple CSS only coexist in the subdomain \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right)\in [0.001,10]\times [0.044,0.145]\), with either one stable and one unstable CSS, or even two stable CSS.
As far as the effect of R_{f} is concerned, a reduction of the feedback resistance with all other parameters fixed leads to a larger useful power output and a larger exergetic efficiency, for both CSS1 (Fig. 5a) and CSS2 (Fig. 5b). This behavior is consistent with previous parametric analyses conducted with the NTP model [25]. The parameter L_{f} is also found to have a high impact on the exergetic performance η_{ex} and useful power output \(\dot {W}_{\mathrm {l}}\), but analyzing these variations is more intricate. In the case of CSS1, an overall trend is for the useful power output to increase as L_{f} is lowered, whereas the exergetic efficiency is maximal around the nominal L_{f} value \(\left (L_{\mathrm {f}}^{+}=1\right)\). These different trends, alongside the wavy pattern on Fig. 5a, is attributed to a shift in phase between oscillations in the load and displacer cylinders, as a result of the variation in gasliquid interface and hydraulic pressure due to the phase change. Notice also that nearly all of CSS1 points are stable, apart from those corresponding to the smallest L_{f} values that turn into unstable before vanishing. In the case of CSS2 by contrast, a trend is for the exergetic efficiency to decrease as L_{f} is lowered, while the useful power output passes through a maximum that is shifting to lower L_{f} values as R_{f} itself is lowered.
The results in Fig. 5 thus reveal a tradeoff between η_{ex} and \(\dot {W}_{\mathrm {l}}\) in terms of maximizing NIFTE performance. The maximal efficiency η_{ex} is predicted to be as high as 77% at \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right) = (0.001,0.145)\), yet with a very small power output of 0.1 W; whereas the maximum power output is 5.4 W at \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right) = (0.001,0.001)\), but with a lower exergetic efficiency of about 36%. These values should be considered relative to the nominal NIFTE configuration \(\left (R_{\mathrm {f}}^{+}=L_{\mathrm {f}}^{+}=1\right)\), which has a predicted exergetic efficiency and average power output of about 2.6% and 0.86 W, respectively. Whether the primary target is to maximize the efficiency or to produce more power should normally drive the design of a NIFTE in practice. Here, a good compromise could aim for a design near the separatrix between stable and unstable CSS2 in Fig. 5b; for instance at \(L_{\mathrm {f}}^{+} = 0.03\) and \(R_{\mathrm {f}}^{+}<0.01\), where the predicted power output and exergetic efficiency are greater than 5.2 W and 45%, respectively. Because these CSS are at the verge of instability, a control strategy [37] might become necessary in order to stabilize the NIFTE nonetheless.
Although the NTP model is capable of accurate predictions of several key operational characteristics of the NIFTE, such as its oscillation frequency, we shall see in the remainder of the paper that such improvements in efficiency and power output are unfortunately an artifact of this model as it fails to describe irreversible thermal losses.
NTPTL model
We conduct the same parametric analysis of NIFTE, now accounting for thermal losses with the NTPTL model. The predicted exergetic efficiency η_{ex} and average power output \(\dot {W}_{\mathrm {l}}\) under CSS operation are presented in Fig. 6, for a wide parameter range spanning several orders of magnitude around the nominal parameter values. A first notable difference with the NTPTL model is that a unique limit cycle is detected via the repeated solution of the optimization problem (11) for \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right) \in [\!0.001,10]^{2}\), which we refer to as CSS1’ in agreement with Table 3. This CSS1’ is furthermore stable over the explored parameter range.
The variations in exergetic efficiency and power output of CSS1’ do not display the same wavy pattern that could be observed for CSS1 in Fig. 5a either, although CSS1’ presents several minima and maxima across the parameter range. The maximal efficiency η_{ex} is predicted to be around 1.1% at \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right)=(3.16,0.001)\), while a maximal average power output \(\dot {W}_{\mathrm {l}}\) of 0.31 W is predicted at \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right)=(1.38,0.316)\). Thee values are now to be compared with the predicted exergetic efficiency of 0.69% and predicted average power output of 0.24 W for the nominal NIFTE configuration \(\left (R_{\mathrm {f}}^{+}=L_{\mathrm {f}}^{+}=1\right)\). A general trend on Fig. 6 is for both the efficiency and power output to increase as L_{f} is reduced from its nominal value. In particular, these results suggest that enhanced performance of the NIFTE prototype by 30–50% could be achieved via a significant reduction in L_{f}, while maintaining R_{f} close to nominal.
These predictions by the NTPTL model are a major departure from those of the NTP model. The comparison between Figs. 5a and 6a, where the feasible ranges of (R_{f},L_{f}) largely overlap, shows that accounting for thermal losses may cause a 10 to 20fold drop in the predicted exergetic efficiency and average power output, especially for the most promising NIFTE configurations determined by the NTP model. Despite the fact that the NTP and NTPTL model predictions are close for certain operational characteristics of the NIFTE, such as its oscillation frequency, these extreme differences in efficiency and power output predictions shed doubt on whether the NTP model is capable of accurate, quantitative predictions of the NIFTE performance over wide parameter ranges. By contrast, the level of performance predicted by the NTPTL model is in good agreement with direct experimental evidence [6, 7], where the exergetic efficiency of the NIFTE prototype was measured between 0.4–1.6%. It was also observed experimentally that the exergetic efficiency could improve by increasing the feedback resistance R_{f}, a behavior that is not described by the NTP model.
Parametric analysis of (C _{p},L _{p})
Apart from the feedback cylinder, the design of the power cylinder where the working fluid is in direct contact with the pumped fluid in the load (Fig. 1), is also known to impact the NIFTE performance greatly. This subsection carries out a joint parametric analysis of the capacitance C_{p} and inductance L_{p} based on the NTPTL model only, with all the other parameters held constant at their nominal values (Table 2). Both C_{p} and L_{p} could be varied in a physical device by modifying either the length l_{p} or the crosssection area A_{p} of the power cylinder (Eq. 1).
The predicted exergetic efficiency η_{ex} and average power output \(\dot {W}_{\mathrm {l}}\) under CSS operation are presented in Fig. 7, for \(\left (L_{\mathrm {p}}^{+},C_{\mathrm {p}}^{+}\right) \in [0.01,100]^{2}\). By repeating the solution of the optimization problem (11) over this parameter range, we could identify three distinct limit cycles, denoted by CSS1”, CSS2” and CSS3”. Notice that CSS1” is the unique limit cycle in a wide domain around the nominal NIFTE design \(\left (L_{\mathrm {p}}^{+} = C_{\mathrm {p}}^{+} = 1\right)\). The domains of existence of CSS2” and CSS3” overlap with CSS1” for values of \(C_{\mathrm {p}}^{+}\leq 0.025\) and \(L_{\mathrm {p}}^{+}\geq 17\), respectively, yet these domains are much narrower.
By far the most promising limit cycle is CSS1” in Fig. 7a, both in terms of exergetic efficiency and power output. The entire range of CSS1” is furthermore stable. The maximal efficiency η_{ex} is predicted to be slightly below 1% at \(\left (C_{\mathrm {p}}^{+},L_{\mathrm {p}}^{+}\right) = (24.5, 0.105)\), with efficiencies above 0.9% achieved in a wide subdomain of CSS1” with \(C_{\mathrm {p}}^{+}\geq 3\) and \(L_{\mathrm {p}}^{+}\leq 1\). A maximal average power output \(\dot {W}_{\mathrm {l}}\) of 0.31 W is predicted at \(\left (C_{\mathrm {p}}^{+},L_{\mathrm {p}}^{+}\right) = (4.90,0.01)\), with improved power outputs again observed in a wide band with \(3\leq C_{\mathrm {p}}^{+}\leq 10\) and \(L_{\mathrm {p}}^{+}\) around or below nominal. Similar to the sensitivity analysis for (R_{f},L_{f}) earlier, these results suggest that the performance of the NIFTE prototype could be enhanced by 30–50% by increasing C_{p} and reducing L_{p} from their nominal values in the NIFTE prototype at hand.
Parametric analysis of (C _{d},L _{d})
A third sensitivity analysis is conducted using the NTPTL model in terms of the capacitance C_{d} and inductance L_{d} of the displacer cylinder, with all the other parameters held constant at their nominal values (Table 2). Similar to the parameter pair (C_{p},L_{p}) earlier, both C_{d} and L_{d} could be varied in a physical device by modifying either the length l_{d} or the crosssection area A_{d} of the displacer cylinder (Eq. 1).
The predicted exergetic efficiency η_{ex} and average power output \(\dot {W}_{\mathrm {l}}\) under CSS operation are shown in Fig. 8, for \(\left (L_{\mathrm {d}}^{+},C_{\mathrm {d}}^{+}\right) \in [0.01,100]^{2}\). By repeating the solution of the optimization problem (11) over this parameter range, two distinct limit cycles emerged, denoted by CSS1”’ and CSS2”’. Like previously with (C_{p},L_{p}), we find that CSS1”’ is the unique limit cycle in a wide subdomain around the nominal NIFTE design \(\left (L_{\mathrm {d}}^{+} = C_{\mathrm {d}}^{+} = 1\right)\). The domain of existence of CSS2”’ overlaps with that of CSS1”’ over a narrow band with \(C_{\mathrm {d}}^{+} \in [0.046,0.063]\). A majority of the limit cycles are furthermore stable, under either CSS1”’ or CSS2”’.
The most promising limit cycle is CSS1”’ in Fig. 8a, with predicted exergetic efficiencies and power outputs three to fourtimes higher than CSS2”’ in Fig. 8b. A maximal efficiency η_{ex} of about 0.93% is predicted at \(\left (C_{\mathrm {d}}^{+}, L_{\mathrm {d}}^{+}\right) = (0.229, 0.010)\); and the maximum power output \(\dot {W}_{\mathrm {l}}\) is predicted to be 0.32 W at \(\left (C_{\mathrm {d}}^{+}, L_{\mathrm {d}}^{+}\right) = (0.055, 95.5)\). The gain in power output is about identical to the gains predicted by adjusting the feedback parameters (R_{f},L_{f}) or the displacer parameters (C_{p},L_{p}); whereas the gain in efficiency, in the order of 35% compared to the nominal NIFTE, is lower than the gains predicted by adjusting by (R_{f},L_{f}) or (C_{p},L_{p}). A general trend on Fig. 8a is for both the efficiency and power output to increase as C_{d} is reduced to within the range \(C_{\mathrm {d}}^{+}\in [0.15,0.3]\), while the parameter L_{d} exhibits a relatively lower sensitivity in the NIFTE prototype.
Parametric analysis of K _{th}
To complete the assessment, we conduct a final parametric analysis of the nonlinear gain K_{th} using the NTPTL model, with all the other parameters held constant at their nominal values (Table 2). The gain represents the maximum rate of heat transfer across the heatexchanger, which may be varied in a physical device by adjusting either the distance between the HHX and CHX elements (Fig. 1), or the temperature difference thereof.
The topleft plot in Fig. 9 shows the predicted exergetic efficiency η_{ex} versus average power output \(\dot {W}_{\mathrm {l}}\) under CSS operation, for \(K_{\text {th}}^{+} \in [0.09, 10]\). The embedded table reports the values of \(K_{\text {th}}^{+}\) at four selected points (red circles), in addition to the nominal point (green triangle). The general trend is for η_{ex} to increase and \(\dot {W}_{\mathrm {l}}\) to decreases as K_{th} is reduced, i.e. as the HHX and CHX elements are drawn further apart or the temperature difference is reduced. The sensitivity of \(K_{\text {th}}^{+}\) appears to be very low around the nominal value \(\left (K_{\text {th}}^{+}=1\right)\), where increasing \(K_{\text {th}}^{+}\) only leads to marginal variations in average power output and exegetic efficiency. The engine performance becomes notably more sensitive to values of \(K_{\text {th}}^{+}<0.5\), but a CSS ceases to exist below \(K_{\text {th}}^{+} \approx 0.09\) due to the HHX and CHX elements no longer allowing for a sufficient heat transfer rate. The exergetic efficiency there is predicted to increase by over 10% compared to the nominal efficiency, but the NIFTE may no longer be practical because of the dramatic reduction in power output that ensues.
A physical interpretation of this parametric analysis is that the performance of the NIFTE would be seldom affected if either the temperature difference or the distance between heatexchangers were to vary by a factor of 2 or even larger. In particular, these results corroborate experimental evidence and suggest that a NIFTE might indeed remain feasible for temperature differences lower than 30 ^{∘}C, without reducing the power output or exergetic efficiency significantly.
Multiobjective optimization by varying multiple parameters
Results of the parametric analyses conducted in the previous subsections revealed that adjusting the pairs (R_{f},L_{f}),(C_{p},L_{p}) and (C_{d},L_{d}) of design parameters in the different NIFTE cylinders could enhance the work output by up to 30% or the exergetic efficiency by over 50% in the NIFTE prototype. For instance, the topright plot in Fig. 9 shows the Pareto frontier of maximal exergetic efficiency versus maximal power output for the pair (R_{f},L_{f}). This Pareto frontier was constructed using the same parametric sensitivity data as in Fig. 6. Values of \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+}\right)\) corresponding to the labeled points A–E on the frontier are also reported in the embedded table for comparison. Moving along the Pareto frontier from A to E in this case entails increasing the value of \(R_{\mathrm {f}}^{+}\) by 2fold, while simultaneously decreasing the value of \(L_{\mathrm {f}}^{+}\) toward zero. Such Pareto frontiers are helpful to visualize the tradeoffs between objectives, identify an acceptable compromise based on engineering priorities, and assess the sensitivity of Paretooptimal decisions toward the design parameters.
Next, we apply multiobjective optimization to characterize the Pareto frontier for combinations of three or more design parameters, by following the scalarization and bilevel decomposition approach described earlier. The Pareto frontier obtained by adding the design parameter C_{p} to the pair (R_{f},L_{f}) is displayed on the bottomleft plot in Fig. 9. Also shown on this plot are the results of the Sobol’ sampling used to crosscheck the Pareto points determined by the NelderMead simplex algorithm.
A first observation is that the performance of the NIFTE is clearly enhanced by the addition of an extra design parameter to the optimization problem, especially in terms of the maximal power output which improves by a further 15%. Values of the triplet \(\left (R_{\mathrm {f}}^{+},L_{\mathrm {f}}^{+},C_{\mathrm {p}}^{+}\right)\) at points A–E along the frontier (embedded table in bottomleft plot) suggest a radically different design strategy compared with the (R_{f},L_{f}) case on the topright plot. Here, reaching the Pareto frontier entails decreasing the value of \(R_{\mathrm {f}}^{+}\) toward zero and moving along the frontier from A to E by reducing \(L_{\mathrm {f}}^{+}\) from 0.34 to 0.28 and \(C_{\mathrm {p}}^{+}\) from 4.7 to 3.9 simultaneously. Notice, in particular, the high sensitivity to changes in L_{f} and C_{p}, whereby the full frontier from A to E is described by a mere 20% variation in their values. A practical implication of this high sensitivity is that mathematical models capable of highly accurate predictions will likely be paramount to support the detailed design of NIFTE, also calling for dedicated experimental campaigns to further validate these models.
It is important to note, at this point, that the selection of alternative triplets as decision variables in the multiobjective optimization problem, including (R_{f},L_{f},L_{p}) and (R_{f},L_{f},C_{d}), produces similar improvements of the NIFTE performance. This, in turn, motivates the consideration of larger parameter subsets in the multiobjective optimization. The bottomright plot in Fig. 9, for instance, shows the Pareto frontier corresponding to the quadruplet (R_{f},L_{f},C_{p},C_{d}). Comparison with the bottomleft plot shows marginal improvement in terms of exergetic efficiency or power output after adding the parameter C_{d}. The consequence of C_{d} being redundant with the other parameters (R_{f},L_{f},C_{p}) also manifests through the lack of consistent trends in the design parameter values at points A–E along the frontier (embedded table in bottomright plot). The consideration of alternative quadruplets to optimize does not produce further improvement for this NIFTE prototype. And even the optimization of (R_{f},L_{f},C_{p},L_{p},C_{d},L_{d},K_{th}) all together leads to either a maximal efficiency close to 1.2% or a maximal power output of 0.35 W. Overall, these results establish that nearoptimal performance of the NIFTE can be achieved through targeted modifications of three parameters only, such as (R_{f},L_{f},C_{p}). In other words, optimizing the length and/or diameter of the feedback and power cylinders would be sufficient in practice for NIFTE to reach its full potential.
Conclusions
The NIFTE is a simple and robust technology for lowgrade heat conversion to useful power. Earlystage prototypes were found to exhibit relatively low thermal and exergetic efficiencies, but these twophase thermofluidic engines are appealing compared to other lowgrade heat recovery technologies since they can operate with temperature differences as small as (and possibly lower than) 30 ^{∘}C between heat source and sink, and may be technoeconomically viable in applications where conventional heatrecovery technologies are not feasible. The main objective of this paper was, therefore, to conduct a modelbased assessment of the NIFTE technology under cyclic steadystate operation, with a view to better understanding its operation and quantifying its true potential.
A nonlinear dynamic model of the NIFTE describing cyclic steady states as stable limit cycles was used as the basis for our investigations, and we extended it to account for irreversible thermal losses in nominally adiabatic components of the device. For the first time, we could predict the coexistence of multiple limit cycles for certain NIFTE configurations, which may be either stable or unstable.
A parametric analysis of the most sensitive design parameters of the NIFTE device (feedback, displacer and power cylinders, and heatexchangers) was conducted, followed by a systematic multiobjective optimization to characterize the Pareto frontier of maximal exergetic efficiency versus maximal power output. Our results establish that failure to account for irreversible thermal losses can grossly overpredict the exergetic efficiency and power output that may be delivered by a NIFTE. The assessment conducted with the extended model predicts improvement greater than 50% in either the work output or the exergetic efficiency via targeted design modifications to an earlystage NIFTE prototype. The optimization study furthermore revealed that nearoptimal performance of the NIFTE can be realized by adjusting combinations of no more than three parameters, such as the dimensions of the feedback tube and power cylinder.
These improvements are particularly relevant for the development of NIFTE pumps to support smallholder irrigation in the developing world [38]. Though for the NIFTE technology to be adopted for lowgrade heat recovery at larger scale, minimization of the irreversible thermal losses via further design modifications would be worth investigating, given that these losses entail a dramatic reduction in the exergetic efficiency and work output. Other recommended followups to our paper include the execution of dedicated experimental campaigns to further calibrate the NIFTE models, and the implementation of optimized NIFTE designs to verify the modelbased predictions. In doing so, one should keep in mind the simplicity of the electrical analogy models used relative to the complexity of actual NIFTE devices. As more experimental data of NIFTE become available, it would be interesting to attempt the development of hybrid semiparametric models that enable more reliable predictions by combining a physicsbased models with datadriven components [39].
Availability of data and materials
Data sharing not applicable to this article as no datasets were generated or analyzed during the current study.
Notes
 1.
Available from: https://github.com/omegaicl/cronos
Abbreviations
 CHX:

Cold heatexchanger
 CSS:

Cyclic steady states
 HHX:

Hot heatexchanger
 NIFTE:

Noninertivefeedback thermofluidic engine
 NLP:

Nonlinear programming
 NTP:

Nonlinear temperature profile
 ORC:

Organic Rankine cycle
 ODE:

ordinary differential equation
 TL:

Thermal loss
References
 1
Chen H, Goswami DY, K SE. A review of thermodynamic cycles and working fluids for the conversion of low grade heat. Renew Sustain Energy Rev. 2010; 14:3059–67.
 2
Ammar Y, Joyce S, Norman R, Wang Y, Roskilly AP. Low grade thermal energy sources and uses from the process industry in the UK. Appl Energy. 2012; 89(1):3–20.
 3
Markides CN. Lowconcentration solarpower systems based on organic Rankine cycles for distributedscale applications: Overview and further developments. Front Energy Res. 2015; 3:47.
 4
Smith TCB. Power dense thermofluidic oscillators for high load applications. In: Proceedings of the 2nd International Energy Conversion Engineering Conference. Providence (RI): 2004. p. 1–15. https://arc.aiaa.org/doi/book/10.2514/MIECEC04.
 5
Smith TCB. Asymmetric heat transfer in vapour cycle liquidpistion engines. In: Proceedings of the 12th International Stirling Engine Conference and Technology Exhibition. Durham: 2005. p. 302–14. https://www.worldcat.org/title/proceedingsofthe12thinternationalstirlingengineconferenceandtechnologyexhibitiondurham79september2005/oclc/64554313.
 6
Smith TCB. Thermally driven oscillations in dynamic applications. PhD thesis, University of Cambridge, UK. 2006. https://www.repository.cam.ac.uk/handle/1810/284059.
 7
Markides CN, Smith TCB. A dynamic model for the efficiency optimization of an oscillatory low grade heat engine. Energy. 2011; 36:6967–80.
 8
Solanki R, Mathie R, Galindo A, Markides CN. Modelling of a twophase thermofluidic oscillator for lowgrade heat utilisation: Accounting for irreversible thermal losses. Appl Energy. 2013; 106:337–54.
 9
Acha S, Mariaud A, Shah N, Markides CN. Optimal design and operation of distributed lowcarbon energy technologies in commercial buildings. Energy. 2018; 142:578–91.
 10
Peterson RB, Wang H, Herron T. Performance of a smallscale regenerative Rankine power cycle employing a scroll expander. Proc IME Part A J Power Energy. 2008; 222(3):271–82.
 11
Wang K, Sun D, Zhang J, Xu Y, Zou J, Wu K, Qiu L, Huang Z. Operating characteristics and performance improvements of a 500 W travelingwave thermoacoustic electric generator. Appl Energy. 2015; 160:853–62.
 12
Taleb AI, Timmer MAG, ElShazly MY, Samoilov A, Kirillov VA, Markides CN. A singlereciprocatingpiston twophase thermofluidic primemover. Energy. 2016; 104:250–65.
 13
Kirmse CJW, Oyewunmi OA, Haslam AJ, Markides CN. Comparison of a novel organicfluid thermofluidic heat converter and an organic Rankine cycle heat engine. Energies. 2016; 9(7):479.
 14
Kirmse CJW, Oyewunmi OA, Taleb AI, Haslam AJ, Markides CN. A twophase singlereciprocatingpiston heat conversion engine: Nonlinear dynamic modelling. Appl Energy. 2017; 186:359–75.
 15
Oyewunmi OA, Kirmse CJW, Haslam AJ, Müller EA, Markides CN. Workingfluid selection and performance investigation of a twophase singlereciprocatingpiston heatconversion engine. Appl Energy. 2017; 186:376–95.
 16
Backhaus S, Swift GW. A thermoacousticStirling heat engine: detailed study. J Acoust Soc Am. 2000; 107:3148–66.
 17
Ceperley PH. A pistonless Stirling enginethe traveling wave heat engine. J Acoust Soc Am. 1979; 66:1508–13.
 18
Huang BJ, Chuang MD. System design of orifice pulsetube refrigerator using linear flow network analysis. Cryogenics. 1996; 36:889–902.
 19
Solanki R, Galindo A, Markides CN. Dynamic modelling of a twophase thermofluidic oscillator for efficient low grade heat utilization: Effect of fluid inertia. Appl Energy. 2012; 89(1):156–63.
 20
Palanisamy K, Taleb AI, Markides CN. Optimizing the noninertivefeedback thermofluidic engine for the conversion of lowgrade heat to pumping work. Heat Trans Eng. 2015; 36(March 2017).
 21
Daw CS, Kennel MB, Finney CEA, Connolly FT. Observing and modeling nonlinear dynamics in an internal combustion engine. Phys Rev E. 1998; 57:2811–9.
 22
Chen L, Sun F, Wu C. Influence of heat transfer law on the performance of a Carnot engine. Appl Therm Eng. 1997; 17:277–82.
 23
Hoffmann KH, Burzler JM, S S. Endoreversible thermodynamics. J NonEquilib Thermodyn. 1997; 22:311–55.
 24
Sun F, Chen W, Chen L, Wu C. Optimal perfomance of an endoreversible Carnot heat pump. Energy Convers Mgmt. 1997; 38:1439–43.
 25
Markides CN, Osuolale A, Solanki R, Stan GBV. Nonlinear heat transfer processes in a twophase thermofluidic oscillator. Appl Energy. 2013; 104:958–77.
 26
Markides CN, Solanki R, Galindo A. Working fluid selection for a twophase thermofluidic oscillator: Effect of thermodynamic properties. Appl Energy. 2014; 124:167–85.
 27
Solanki R, Galindo A, Markides CN. The role of heat exchange on the behaviour of an oscillatory twophase lowgrade heat engine. Appl Thermal Eng. 2013; 53(2):177–87.
 28
Wächter A, Biegler LT. On the implementation of a primaldual interior point filter line search algorithm for largescale nonlinear programming. Math Program. 2006; 106(1):25–57.
 29
Hindmarsh AC, Brown PN, Grant KE, Lee SL, Serban R, Shumaker DE, Woodward CS. SUNDIALS: Suite of nonlinear and differential/algebraic equation solvers. ACM Trans Math Softw. 2005; 31(3):363–96.
 30
Grass D, Caulkins JP, Feichtinger G, Tragler G, Behrens DA. Optimal Control of Nonlinear Processes: With Applications in Drugs, Corruption, and Terror. Berlin: Springer; 2008.
 31
Teschl G. Ordinary Differential Equations and Dynamical Systems. Graduate studies in mathematics. Providence (RI): American Mathematical Society; 2012.
 32
Dhooge A, Govaerts W, Kuznetsov YA. MatCont: A MATLAB package for numerical bifurcation analysis of ODEs. ACM Trans Math Softw. 2003; 29:141–64.
 33
Govaerts W, Kuznetsov YA, Dhooge A. Numerical continuation of bifurcations of limit cycles in MATLAB. SIAM J Sci Comput. 2005; 27:231–52.
 34
Miettinen K. Nonlinear Multiobjective Optimization. Dordrecht: Kluwer Academic Publishers; 1999.
 35
Conn AR, Scheinberg K, Vicente LN. Introduction to DerivativeFree Optimization. Philadelphia: MOSSIAM Series on Optimization; 2009.
 36
Sobol’ IM, Asotsky D, Kreinin A, Kucherenko S. Construction and comparison of highdimensional Sobol’ generators. Wilmott Mag. 2012; Nov:64–79.
 37
Pyragas K. Continuous control of chaos by selfcontrolling feedback. Phys Lett A. 1992; 170(6):421–8.
 38
Thermofluidics Ltd. http://www.thermofluidics.co.uk. Accessed: 20190408.
 39
Thompson ML, Kramer MA. Modeling chemical processes using prior knowledge and neural networks. AIChE J. 1994; 40(8):1328–40.
Acknowledgements
YW is grateful to the Department of Chemical Engineering at Imperial College London for a doctoral scholarship.
Funding
Not applicable.
Author information
Affiliations
Contributions
YW (yukun.wang11@imperial.ac.uk), CM (c.markides@imperial.ac.uk) and BC (b.chachuat@imperial.ac.uk) all made substantial contributions to the ideas that led to this manuscript and the discussion of the results. YW was in charge of implementing the mathematical models and generating the computational results. YW and BC prepared the initial draft of the paper. All authors read and approved the submitted manuscript.
Corresponding author
Correspondence to Benoît Chachuat.
Ethics declarations
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Competing interests
The authors declare that they have no competing interests.
Additional information
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Received
Accepted
Published
DOI
Keywords
 Heat engine
 Lowgrade heat
 Thermofluidic oscillator
 Cyclic steady state
 Asymptotic stability
 Parametric analysis
 Exergetic efficiency