Abstract

Solar-driven thermochemical energy storage systems are proven to be promising energy carriers (solar fuels) to utilize solar energy by using reactive solid-state pellets. However, the production of solar fuel requires a quasi-steady-state process temperature, which represents the main challenge due to the transient nature of solar power. In this work, an adaptive model predictive controller (MPC) is presented to regulate the temperature inside a tubular solar reactor to produce solid-state solar fuel for long-term thermal storage systems. The solar reactor system consists of a vertical tube heated circumferentially over a segment of its length by concentrated solar power, and the reactive pellets (MgMn2O4) are fed from the top end and flow downwards through the heated tube. A countercurrent flowing gas supplied from the lower end interacts with flowing pellets to reduce it thermochemically at a temperature range of 1000—1500 °C. A low-order physical model was developed to simulate the dynamics of the solar reactor including the reaction kinetics, and the proposed model was validated numerically by using a 7-kW electric furnace. The numerical model then was utilized to design the MPC controller, where the control system consists of an MPC code linked to an adaptive system identification code that updates system parameters online to ensure system robustness against external disturbances (sudden change in the flow inside the reactor), model mismatches, and uncertainty. The MPC controller parameters are tuned to enhance the system performance with minimum steady-state error and overshoot. The controller is tested to track different temperature ranges between 500 °C and 1400 °C with different particles/gas mass flowrates and ramping temperature profiles. Results show that the MPC controller successfully regulated the reactor temperature within ± 1 °C of its setpoint and maintained robust performance with minimum input effort when subjected to sudden changes in the amount of flowing media and the presence of chemical reaction.

1 Introduction

The immense potential of solar energy not only brings life to our planet but also enables us to generate power and fuels. Solar energy is utilized in a variety of energy conversion applications such as water heaters, photovoltaics, desalination units, and others [1,2]. It is possible to store solar energy in the form of “solar fuel.” For example, photo-electrochemical (PEC) is a solar energy conversion method to produce solar fuel via CO2 reduction. Earth-abundant materials are an option to develop efficient photo-electrodes in CO2 reduction with PEC systems [3]. Recently, the performance of PEC was enhanced with electrode improvements. These include nano-coned silicon electrodes, versus planar silicon, or nanoporous TiO2 coatings on CuBi2O4 photocathodes [4], or coating nanoporous CuBi2O4 photocathode with TiO2 [5]. However, the main drawback of the PEC is the limitation to the presence of sunlight; otherwise, no electricity will be generated [6].

Storing concentrated solar energy into a long-duration solar fuel in a solid-state form via thermochemical processes is a promising technology. Endothermic heat for the reaction is supplied by concentrated solar energy to initiate the chemical conversion of reactants inside a solar reactor. Heat captured via concentrated solar energy is transferred to the reactants either directly through an aperture of the reactor or indirectly via a windowless reactor. There are various thermochemical energy storage systems and the use of each type is limited by the temperature range of the charging and discharging process, volumetric energy, and the number of reusable cycles (lifetime) [7].

The cost and thermal stability for long-duration storage play a critical role in the selection process of thermochemical energy storage candidates. For example, BaO2/BaO has typically 20 useable lifetime cycles [8]. Recent research shows that the life cycle of BaO2/BaO has increased to 200 by adding 33 wt% MgO under the temperature range of 600–750 °C [9]. Co3O4 and Co3O4/Al2O3 both have the same lifetime of 100 cycles and react at temperature ranges of 840–940 °C [10] and 700–1000 °C, respectively [11]. At higher temperatures, magnesium manganese oxide (MgMn2O4) reacts between 1000 and 1500 °C with a lifetime of 30 cycles and was proven to show unique thermal characteristics and stability when stored at room temperature [12], and the cyclic lifetime of thermochemical storage media is predicated on good temperature control.

One of the essential parameters in an efficient chemical conversion of reactants in a solar reactor is the “residence time.” Residence time is defined as the duration from the time when the reactants enter the reactor to when they leave, this is a key metric in reduction and oxidation processes during the energy charging and discharging steps. Residence time can be preserved by adjusting the flowrate of the reactants. This assists in maintaining an efficient chemical conversion process inside a solar reactor.

Another important parameter in achieving an efficient chemical conversion process inside a solar reactor is maintaining a stable supply of endothermic heat in the reactor. Given the fact that solar energy is intermittent by nature, it is a challenge to keep the reactor temperature constant. There are two possible ways to adjust incoming direct normal irradiance (DNI) intercepted by the solar reactor. One promising method is via using a mechanically variable aperture [13]. Another approach is by controlling the gas flowrate entering the reactor [14]. A linear proportional integral (PI) controller can be used to regulate solar reactor temperature [15], while an adaptive Model Predictive Control (MPC) provides a more comprehensive adjustment based on the system dynamics [16].

In our previous work [17], we analyzed the thermal characteristics of a horizontal tubular solar reactor with a fixed bed to develop a low-order physical model which was used to design a PI control system to regulate temperature [15]. In this work, we extend our previous study by covering the dynamics of a vertical tubular reactor with a moving bed. This work presents an adaptive model predictive controller (AMPC) to regulate the temperature of a tubular solar reactor. The control system is developed based on a system model where the analysis incorporates a one-dimensional heat transfer model to accurately predict the transient temperature response of a high-temperature plug flow solar reactor. The system includes a vertically oriented tubular reactor surrounded by heating elements in one segment and MgMn2O4 pellets fed into the reactor from the top. The sensible heat of the downward-flowing pellets is recuperated via a countercurrent flowing air such that the pellets/gas enter and exit the reactor tube at approximately room temperature. An in-house code was developed to optimize the performance of the controller by including an updating algorithm for the system identification process. The low-order physical model was validated by an experimental test of different amount of particle mass flowrate and heating levels. The performance of the proposed adaptive MPC was simulated to track different temperature setpoints and different operation conditions including the chemical reaction at high temperatures.

2 Methodology

This research includes numerical and experimental work on the production of a solid-state solar fuel MgMn2O4 in a solar reactor featuring countercurrent flowing gas and solid particles forming a complex heat transfer interference. A simplified one-dimensional (1D) heat transfer model is developed to obtain a temperature distribution profile and to understand system dynamics during the process of solar fuel production. The ratio between mass flowrates for gas and particles needs to be regulated to maintain an efficient heat recuperation process and practically using (cpsm˙s=cpgm˙g) would be used to calculate the amount of gas flowrate for the corresponding amount of particle mass flowrate. To observe the system dynamics and temperature distribution, the phenomenon was monitored by heating a vertical tube housing flow of MgMn2O4 particles and a countercurrent flowing gas.

2.1 Experimental Methodology.

The reactor consists of a 152.4 cm in length alumina tube centered vertically through a 30.5-cm heating zone, and more geometry and specifications are listed in Table 1. [19] The reactor tube is heated circumferentially via adjustable level heat flux by a 7-kW electric furnace. A funnel filled with MgMn2O4 particles is mounted at the top of the reactor tube, and the tube receives 3–3.66 mm MgMn2O4 particles from the funnel. The flowrate of the particles was controlled by a gas-pulsation valve mechanism at the bottom end of the reactor, and the discharged particles were collected inside a tank beneath the setup.

Table 1

List of thermal properties and dimensions

MaterialThermal conductivity (W/mK)Specific heat (J/kg K)Mass density (kg/m3)Dimensions (mm)
Reactor tube (alumina)kw=1.495×1013Tw5+4.123×1010Tw44.47×107Tw3+2.69×104Tw20.11Tw+34.96[19]Cpw=2.696×1016Tw51.69×1012Tw4+3.933×109Tw34.267×106Tw2+2.3×103Tw+0.691[19]3900Do = 57.15, Di = 50.8
Insulation0.315Thickness = 152.4
Gas (averaged)0.087312060.88
Particles (MgMn2O4)ks=5.27×1010Ts31.89*107Ts2+4.72×104Ts+0.3889 [19]870(1 − ɛ)ρs = 2003Size = 3
MaterialThermal conductivity (W/mK)Specific heat (J/kg K)Mass density (kg/m3)Dimensions (mm)
Reactor tube (alumina)kw=1.495×1013Tw5+4.123×1010Tw44.47×107Tw3+2.69×104Tw20.11Tw+34.96[19]Cpw=2.696×1016Tw51.69×1012Tw4+3.933×109Tw34.267×106Tw2+2.3×103Tw+0.691[19]3900Do = 57.15, Di = 50.8
Insulation0.315Thickness = 152.4
Gas (averaged)0.087312060.88
Particles (MgMn2O4)ks=5.27×1010Ts31.89*107Ts2+4.72×104Ts+0.3889 [19]870(1 − ɛ)ρs = 2003Size = 3

As aforementioned, a countercurrent flowing gas was supplied to the reactor from the bottom end. The gas flowrate is controlled by a digital Alicat gas flow controller. The reactor tube is divided into three sections: (1) upper, (2) middle which is the heating zone, and (3) lower section. The upper and middle sections are insulated with fibrous insulation, whereas the lower section is exposed to the ambient. Six B-type thermocouples are installed along the reactor tube length to record wall temperature during the heating and cooling processes as depicted in Fig. 1. Thermocouples T1, T2, and T3 record the temperatures at the lower section, T4 records the temperature in the heating zone, whereas T5 and T6 record the temperatures in the upper section. A GraphTec digital data logger is used to record the temperature every second. The reactor tube sections and the locations of thermocouples are illustrated in Fig. 1.

Fig. 1
Schematic diagram of the experimental setup illustrates thermocouples locations and reactor sections
Fig. 1
Schematic diagram of the experimental setup illustrates thermocouples locations and reactor sections
Close modal

2.2 Numerical Methodology.

The assembly of the reactor tube shown in Fig. 1 was analyzed numerically by modeling the heat transferred via convection, conduction, and radiation between the reactor tube and interpenetrating gas/solid phases. The tube model is subjected to a constant heat flux over a segment of its length centered on the midpoint of the tube length as depicted in Fig. 2. The reactor tube and the flowing media are discretized into finite control volumes, and the model was developed based on the assumption that the superficial velocity remains constant along the flow path.

Fig. 2
Schematic diagram of the heat transfer model of reactor tube and solar fuel particles
Fig. 2
Schematic diagram of the heat transfer model of reactor tube and solar fuel particles
Close modal

Formulation of the temperature profiles was derived from the energy conservation of the reactor tube in addition to gas/solid phases as follows:

  • Energy balance of the reactor wall
    ρwcpwTwt=x(kwTwx)+hgwaw(TgTw)+hswaw(TsTw)+hraw(TsTw)+haw(TTw)+awoqs
    (1)
  • Energy balance of the gas
    ερgcpgTgt+ρgcpgugTgx=x(kgefTgx)+hgsags(TsTg)+hgwaw(TwTg)
    (2)
  • Energy balance of the particles
    (1ε)ρscpsTst+ρscpsvsTsx=x(ksefTsx)+hgsags(TgTs)+hswaw(TwTs)+hraw(TwTs)+R
    (3)
    where s, g, and w stand for solar fuel, gas, and tube wall, respectively. ug and vs are gas and particles flow velocities, respectively, ρ is the mass density, cp is the specific heat capacity, and ɛ is the porosity. The boundary conditions for the particles and gas inlet/outlet temperatures are as follows:
    Atx=0,Tg=(Tg)inandTsx=0
    (4)
    Atx=L,Ts=(Ts)inandTgx=0
    (5)
    The differential equations for the heat transfer model described in Eqs. (1)(3) are discretized into 100 cylindrical-shaped control volumes as illustrated in Fig. 2. The diffusion term was discretized using second-order central differencing, third-order upwind scheme (QUICK) for convection, and Crank–Nicolson method in time.

The reaction of the (Mg–Mn–oxide) particles occurs when the temperature of the solid particles reaches above 1000 °C [12], and the energy balance for the particles will include a nonzero value of the heat sink (R). This represents the total chemical energy stored during the reduction reaction and is defined as [18]
R=A1f1(ΔH,us)+A2f2(ΔH,dαdt)
(6)
where A1 and A2 are constants, ΔH is the change in enthalpy of reaction, us is the velocity of the particles, and dαdt the conversion rate which is defined as [18]
dαdt=f(Po2,Ts,α)
(7)
where PO2 is the oxygen partial pressure, Ts is the temperature of solid particles, and α is the extent of reaction. Following the assumption of constant pressure, a simplified form of Eq. (15) can be obtained by fitting the term dαdt with the corresponding temperature of solid (Ts) within the range of 1000 °C < Ts < 1500 °C. Also, the reaction rate is proportional to the residence time of particles which is governed by the particle mass flowrate along the heating zone. Hence, the simplified form of Eq. (15) can be written as
R=f(m˙s,Ts)
(8)
where f(m˙s,Ts) is a simplified 2D fitted function which can be obtained by testing the system experimentally with different particle mass flowrates m˙s for different steady-state temperatures Ts.

3 Results and Discussions

In order to observe the reactor dynamics, an experiment was conducted through a multi-stage process of heating-cooling and various gas/particle flowrates. To monitor the influence of flowing media on heat transfer and the temperature distribution along the flow path, the mass flowrate of the gas was adjusted to be less than the required amount that preserves the heat recuperation process. The entire experimental run lasted for 14 h, and the exact experimental conditions were simulated. The validity of the numerical model was evaluated by comparing temperature distribution at various periods of time throughout the experimental run.

3.1 Experimental Results.

Temperature distribution along the reactor tube is monitored by testing the system in three stages. In the first stage, the reactor tube is filled with MgMn2O4 pellets, and the system is heated from room temperature to 1000 °C in 3 h with a 5.5 °C/min heating rate and no gas/particle flow. The second stage starts once the temperature at the center of the heated zone reaches 1000 °C, and gas/particle flowrates are initiated with 1 g/s particle mass flowrate and an average gas flowrate of 38 standard litre per minute (SPLM). At the same time, the heating process continues to ramp the temperature up to 1400 °C, after which the furnace holds a steady temperature for 90 min. In the third stage, the gas/particle mass flowrates were stopped and a cooling process of the reactor tube starts until the entire system cools down to room temperature.

Figure 3(a) illustrates the temperature profiles for three sections of the reactor tube during the entire experimental run. During the first stage (0–1000 °C), the system, a packed bed with no gas flow, was heated uniformly. As a result, the temperature distribution along the reactor tube shows a steady increase. It should be noticed that the profile of T3 (in the lower section) increases faster since it is located at the closest position to the heating zone. Once the gas/particle mass flowrates are initiated in the second stage (the shaded period in Fig. 3(a)), the growth of the temperature profile along the upper section (T4 and T5) was suppressed since flowing pellets removed heat from this section toward the particle flow path. In contrast, a significant temperature increase is noticeable in the lower section (T2 and T3), which implies that more heat was transferred to this section via the mass flow of particles. This is reasonable since the average tested gas flowrate was 38 SPLM while the actually required gas flowrate for an efficient heat recuperation process is approximately 54 SPLM. During the third stage, the system was cooled down to room temperature. However, the temperature reading in Fig. 3(a) starts from 200 °C because of the sensitivity range of Type-B thermocouples. As the system cools down naturally, the initial thermal energy stored in the system from the end of stage 2 is redistributed longitudinally across the tube wall. The fork-shaped profiles of upper T3 and lower T2 sections indicate that the hotter section (lower) losses heat towards the cooler section (upper) to reach a steady-state as the process proceeds. Nevertheless, the temperatures along the upper section (T5 and T6) tend to be higher near the steady-state since the upper section is insulated.

Fig. 3
Experimental test: (a) temperature profiles along the reactor tube during three stages, (b) corresponding particle mass flowrate, and (c) corresponding gas flowrate
Fig. 3
Experimental test: (a) temperature profiles along the reactor tube during three stages, (b) corresponding particle mass flowrate, and (c) corresponding gas flowrate
Close modal

3.2 Numerical Results.

The accuracy of the numerical model described in Eqs. (1)(3) was validated with the experimental results by simulating exact boundary conditions and the heating process of the experimental run. The parameters, values, correlations, and thermal properties used in the numerical simulation are listed in Tables 1 and 2. The system was discretized into 100 cylindrical-shaped control volumes, and the experimental run was simulated in three stages. During the first stage, the system initially was at room temperature and heated up to 1000 °C. Thereafter, the second stage started with the final temperature of the first stage as the initial system temperature. Similarly, the initial temperatures of the third stage are the final temperatures of the second stage. For the simulation, stability proposes, the time-step was set to 0.1 ms.

Table 2

List of parameters and correlations

ParameterSymbolValue/CorrelationRef.
Porosityɛ0.375+0.34dpDi[20]
Solid-gas heat transfer coefficienthsg(2+1.1Re12.Pr13)kgdp[21]
Solid-wall heat transfer coefficienthsw10.085+12πPeLkgdp
Wall-gas heat transfer coefficienthwg(0.023Re0.8.Pr0.4)kgDiDerived from Ref. [22]
Radiation heat transfer coefficienthrσ(Ts2+Tw2)(Ts+Tw)1eb+1ew1
Emissivity of the bulkeb0.5(1 + ep)
Emissivity of the tube wallew0.6–0.8[23]
Emissivity of the particlesep0.7
ParameterSymbolValue/CorrelationRef.
Porosityɛ0.375+0.34dpDi[20]
Solid-gas heat transfer coefficienthsg(2+1.1Re12.Pr13)kgdp[21]
Solid-wall heat transfer coefficienthsw10.085+12πPeLkgdp
Wall-gas heat transfer coefficienthwg(0.023Re0.8.Pr0.4)kgDiDerived from Ref. [22]
Radiation heat transfer coefficienthrσ(Ts2+Tw2)(Ts+Tw)1eb+1ew1
Emissivity of the bulkeb0.5(1 + ep)
Emissivity of the tube wallew0.6–0.8[23]
Emissivity of the particlesep0.7

Figure 4 illustrates the experimental and simulation profiles along the reactor tube sections. In general, simulated temperatures show a similar trend when compared to experimental results. In the first stage, the temperature profiles show good agreement with the experimental results. However, during the second stage, the temperature at the upper section (T5) decreases as expected, whereas a slight change is observed in the temperature at the top upper section T6. On the other hand, temperatures of the lower sections T2 and T3 increase simultaneously, although T3 shows deviation with respect to measured data. In the third stage, it is seen that the profiles of the upper and lower section temperatures T2 and T5 simulate the aforementioned phenomenon of energy redistribution during the cooling process. Moreover, it is noticeable that the temperatures in the upper section are higher near the steady-state since the section is insulated.

Fig. 4
Simulation and experimental results of temperature profiles along the reactor sections
Fig. 4
Simulation and experimental results of temperature profiles along the reactor sections
Close modal

To obtain better insights into the comparison between numerical and experimental results, temperature distribution along the reactor tube was traced with time for various temperature levels for both cooling and heating processes. Figure 5 illustrates a comparison between the numerical and experimental results during the heating process for five temperature heating levels, namely 300, 600, 1000, 1200, and 1400 °C. It is clear that the temperature distribution is uniform during heating to 1000 °C (first stage). Nevertheless, when gas/particle flow starts as the system is heated from 1000 °C to 1400 °C (stage 2); temperature distributions indicate an appreciable cooling at the upper section and significant heating at the lower section since the flowing mass of particles transfers heat toward the lower section.

Fig. 5
Numerical and experimental results of temperature distribution along the reactor tube during the heating process
Fig. 5
Numerical and experimental results of temperature distribution along the reactor tube during the heating process
Close modal

The mismatch in numerical results is in the range of 0.7–9.7% at the gas inlet and 1.7–5% at the particles inlet during the heating process. Similar to the heating process, Fig. 6 depicts the time tracing of stage 3 during the cooling down process. Both numerical and experimental results illustrate the energy redistribution across the tube length as the system cools down. However, in the upper section, the temperature tends to be higher as the system reaches a steady-state. This is attributed to the presence of insulation in the upper section as discussed earlier. Two observations can be drawn from the simulation of the heating and cooling processes shown in Figs. 5 and 6. First, the temperature in the upper section (insulated) shows good agreement with the experimental result. Second, at the lower section, it is noticeable that the simulated temperature at the location of T3 deviates from the experimental measurement during the process of heating and cooling, especially when the system temperature is between 1400 °C and 600 °C. The lower section is exposed to the ambient environment, and the heat loss is calculated by using a combined heat transfer coefficient of hloss = 30 W/m2.°C. Since the temperature deviation occurs between 1400 °C and 600 °C, it can be mainly attributed to the heat loss by radiation. The good agreement at the insulated section supports this finding. Although hloss accounts for radiation, a more rigorous heat loss model should be implemented in future work to accurately address this issue.

Fig. 6
Numerical and experimental results of temperature distribution along the reactor tube during cooling
Fig. 6
Numerical and experimental results of temperature distribution along the reactor tube during cooling
Close modal

4 Design of the Model Predictive Controller Controller

The main concept in the development of an MPC control strategy is to utilize system models to predict future outputs over a prediction horizon, Np, by taking into account different control actions which change over a control horizon, Nu. Also, the known past input is a useful parameter that contributes to response optimization by obtaining a set of future control signals that preserve the process output very close to the desired reference trajectory. Among other control methods, the MPC control strategy has many advantages since it utilizes a system model to optimize the controller's performance and handles challenging dynamics. In this work, a Generalized Predictive Controller (GPC) is implemented with the linearized form of the solar reactor. It is given as an autoregressive with extra input (ARX) model which is described as follows [24]:
A(q1)=B(q1)u(k1)+e(k)
(9)
where u(k) and y(k) are the control/input (heat flux) and output (the average temperature at the center of heating zone in °C), respectively. q−1 is the backward shift operator, and e(k) is a white noise (which is negligible in this work) representing the system disturbances and the polynomials which are given as
A(q1)=1+a1q1+a2q2++anaqna
(10)
B(q1)=b0+b1q1+b2q2++bnbqnb
(11)
a1,…, ana and b0, …, bnb are the model parameters with polynomial orders of 3 (na = nb = 3), which to be determined. For linear systems a1,…, ana and b0, …, bnb are fixed constants; however, in the case of the solar reactor with the thermochemical process, the system dynamics is nonlinear, and hence, using fixed model parameters will reduce the accuracy of the control model. Therefore, it is convenient to implement an online system identification process to address the system’s nonlinearity (as shown in Fig. 7) during the reactor operation via a recursive least square (RLS) method with a forgetting factor as follows [13]:
y^=φT(k)θ^(k1)
(12)
where φT represents the column vector of past input and output values which is defined as:
φT(k)=[y(k1),y(k2),y(k3),u(k),u(k1),u(k2),u(k3)]
(13)
where θ^ represents the estimation of model parameters which is defined as
θ^(k)=[a1,a2,a3,bo,b1,b2,b3]
(14)
and calculated by using
θ^(k)=θ^(k1)+P(k1)φ(k)λ(k)+φT(k)P(k1)φ(k)(y(k)y^(k))
(15)
where P(k) represents the covariance matrix of estimation errors which is defined as
P(k)=1λ(k)[P(k1)P(k1)φ(k)φT(k)P(k1)λ(k)+φT(k)P(k1)φ(k)]
(16)
where λ(k) represents the forgetting factor which is calculated from
λ(k)=max[(1y(k)y^(k)1+(y(k)y^(k))2),λmin]
(17)
Fig. 7
Feedback control loop of the adaptive model predictive controller
Fig. 7
Feedback control loop of the adaptive model predictive controller
Close modal

The system identification process was performed by simulating the system response to a pseudorandom generated input power signal that changes over a random period of time for 20 h as shown in Fig. 8. Using fixed model parameters results in a maximum error of 177.8 C in the estimated temperature, whereas the error was eliminated when the model parameters were calculated online with every time-step. Figure 9 illustrates the identification progress over time-steps, where the model parameters initially are zero and then start to converge after 100 time-steps (iteration) where the final values are (−2.4813, 1.9871, −0.50525) for a1,…, a3, and (0.0040861, −0.002985, −2.1876 × 10-5, −0.00093613) for b0,…, b3, The number of iteration can be significantly reduced by choosing a pre-calculated initial guess for b0,…, b3 and a1,…, a3.

Fig. 8
System identification process: (a) simulated temperature output in comparison to the estimated one based on 8-steps ahead prediction and (b) a pseudorandom generated power input
Fig. 8
System identification process: (a) simulated temperature output in comparison to the estimated one based on 8-steps ahead prediction and (b) a pseudorandom generated power input
Close modal
Fig. 9
Estimated parameters of the ARX model: (a) A(q−1) and (b) B(q−1)
Fig. 9
Estimated parameters of the ARX model: (a) A(q−1) and (b) B(q−1)
Close modal
The estimated ARX model was then used to optimize the control action of the system by implementing the GPC algorithm, which minimizes the cost function J defined as [13]
J(Np,Nu)=j=1Np[y^(k+j|k)w(k+j)]2+j=1Nuω(j)[Δu(k+j1|k)]2
(18)
where y^(k+j|k) is ahead output prediction at time interval k, and w(k + j) is the future reference trajectory, which is defined as the desired setpoint. ω(j) is the weighing parameter, which is to be adjusted properly to follow the reference trajectory. The control action is the change in the input power, and the desired output is the average temperature in the reactor's heating zone in °C. To avoid excessive heating rate, the increment in the input power is limited to 100 W, and the total input power is u(tk) = u(tk−1) + Δu(tk), where Δu ≤ |100| W, and the time-step (k) is taken to be 60 s. The MPC controller needs to be tuned in a manner similar to tuning the proportional–integral–derivative (PID) controller. The parameters to be tuned are prediction horizon (Np), control horizon (Nu), and weighting parameter ω. To ensure optimal control performance (minimum control effort, steady-state error, and overshoot), the cost function (J) in Eq. (18) needs to be minimized. The predictive horizon is usually chosen to be greater than the control horizon [25]. The numerical tuning procedure showed that choosing Np = Nu + 3 results in the lowest control effort, where the cost function (J) can be minimized by tuning Nu and ω accordingly.

Figure 10 illustrates the change in the cost function with different values of control horizon (Nu) and fixed value of the weighting sequence (ω = 0.2). It is obvious that the short control horizon results in oscillatory behavior, and increasing the control horizon to Nu = 20 results in a better optimization process. However, increasing the control horizon further to Nu = 30 would not enhance the optimization process, and hence, the value of Nu = 20 will be considered to be optimal for further simulations. The value of the weighting sequence (ω) was tuned in a similar way by fixing the value of the control horizon to Nu = 20 and simulating the performance of the MPC controller with different values of (ω) as shown in Fig. 11. It can be noticed that a small value of the weighting sequence (ω = 0.001) increases the overshot and settling time and the control action results in an input increment greater than (100 W), which resulted in straight lines in the input profile. On the other hand, increasing (ω) to 1.2 results in a longer settling time and greater overshot. Among different values of (ω), the optimal (ω) is 0.2 which obtains optimal performance as shown in Fig. 11(a).

Fig. 10
Normalized cost function (J) with different values of control horizon (Nu) and fixed value of ω = 0.2
Fig. 10
Normalized cost function (J) with different values of control horizon (Nu) and fixed value of ω = 0.2
Close modal
Fig. 11
Tuning of MPC controller with different values of Ω and fixed value of control horizon (Nu = 20): (a) response and (b) input power
Fig. 11
Tuning of MPC controller with different values of Ω and fixed value of control horizon (Nu = 20): (a) response and (b) input power
Close modal

The performance of the control system was tested with different operation conditions by tracking a trajectory of different setpoints namely 1000, 1400, and 500 °C with heating and cooling rates of 5 °C/min. For all runs, the controller system was switched from an open loop to a closed loop after 3 min. First, the MPC controller was tested to track the prescribed setpoints with a fixed bed. The system was heated from room temperature to 1000 °C and maintained steady for 2 h. Thereafter, the system was heated up to 1400 °C and preserved steady for 3 h, and then cooled down to 500 °C for the rest period of the simulation as shown in Fig. 11(a). The controller tracked the setpoints closely with a maximum error of ±1 °C and a maximum overshoot of 2%.

Figure 13 illustrates the performance of the MPC controller including flowing media and the chemical reaction, where the controller was tested to track the prescribed setpoints with a 1 g/s particle mass flowrate (and approximately 50 L/min gas flowrate). Since the (Mg–Mn–oxide) particles start to react in a temperature range of 1000–1500 °C [12]; therefore, the flow starts within the period of 5–11 h where the reactor temperature lies within the reaction temperature range. The MPC controller tracked the ramping setpoint closely and maintain a steady temperature during the reaction temperature.

The system identification algorithm updates the model parameters momentarily which captures the system nonlinearity and results in better performance. Comparing the input power in Figs. 12(b) and 13(b), the input power in Fig. 13(b) is greater within the flowing period, which is attributed to the additional energy required to overcome the effect of flowing media and chemical reactions. Also, the MPC controller was tested further by simulating variable particle/gas mass flowrate for the prescribed temperature setpoints as shown in Fig. 14. The controller was capable to track the setpoint closely and maintain a robust response with very small steady-state error and reasonable overshoot during the reactor operation

Fig. 12
Performance of MPC controller for fixed-bed system: (a) response and (b) input power
Fig. 12
Performance of MPC controller for fixed-bed system: (a) response and (b) input power
Close modal
Fig. 13
Performance of MPC controller with chemical reaction and 1 g/s particle mass flowrate: (a) response and (b) input power
Fig. 13
Performance of MPC controller with chemical reaction and 1 g/s particle mass flowrate: (a) response and (b) input power
Close modal
Fig. 14
Performance of MPC controller with chemical reaction and different particles mass flowrate: (a) response and (b) input power and flow
Fig. 14
Performance of MPC controller with chemical reaction and different particles mass flowrate: (a) response and (b) input power and flow
Close modal

Finally, the performance of the MPC controller was compared to experimental results obtained by the PID controller in 7-kW electric furnace as shown in Fig. 15. The system was heated up from 1000 to 1400 °C within 90 min and maintained steady with a particle mass flowrate of 1.25 g/s for one hour. The particle mass flowrate was then reduced to 0.75 g/s (and approximately 36 L/m gas flowrate) to observe the response of the control system. Unlike the PID controller, the MPC maintained a stable temperature with the sudden change in particle mass flowrate. Moreover, the overshoot in input power is less when compared to the PID controller, which implies the superiority of MPC in utilizing the future prediction of the input control to perform the optimal control decision.

Fig. 15
Performance of MPC controller compared to industrial PID controller: (a) response and (b) input power and flow
Fig. 15
Performance of MPC controller compared to industrial PID controller: (a) response and (b) input power and flow
Close modal

5 Conclusions

An adaptive MPC was designed to regulate the temperature inside a tubular solar reactor that is used to produce solar fuel. The production process includes downward-flowing MgMn2O4 particles passing through a vertical tubular reactor. The reactor tube is heated circumferentially at a limited segment of its length such that the MgMn2O4 particles are subjected to 1000–1500 °C for a residence time of less than 10 min as it flows in order to achieve a chemical reaction. The sensible heat of flowing particles is recuperated via counter-current flowing gas such that the gas/particles enter and exit the reactor tube at approximately room temperature. The design of the control system was based on a simplified system model where the governing equations for the heat transfer model were obtained by using energy balance within a control volume considering conduction, convection, and radiation heat transfer. The numerical model was experimentally validated using a reactor prototype made of a 152.4-cm alumina tube heated by a 7-kW electric furnace. The 1D model presented in this study was utilized as a simplified physical model to design a MPC system. The performance of the MPC controller was tested by tracking prescribed setpoints including ramping up from room temperature to 1000 °C, 1500 °C, and cooling down to 500 °C. Also, the controller was tested in with different operating conditions including the presence of chemical reactions and various particle/gas mass flowrates. The control model captured the system nonlinearity and was able to track the desired setpoints closely with maximum steady-state error of ±1 °C and maximum overshoot of 2%. In the actual solar reactor, the system nonlinearity is greater due to solar concentration uncertainties and more modification needs to be considered to ensure the efficacy of the control system.

Acknowledgment

This material is based upon work supported by the U.S. Department of Energy's Office of Energy Efficiency and Renewable Energy (EERE) under the Solar Energy Technology Office (SETO) Award Number DE-EE0008992.

Conflict of Interest

There are no conflicts of interest.

Data Availability Statement

The authors attest that all data for this study are included in the paper.

Nomenclature

ags =

gas particles area per volume (m−1)

aw =

particle-wall area per volume (m−1)

cp =

specific heat (J.kg−1.K−1)

eb =

emissivity of the bulk

ep =

emissivity of the SoFuel particles

ew =

emissivity of the tube wall

hgw =

gas-to-wall convective heat transfer coefficient (W.m−2.K−1)

hgs =

solid-to-gas convective heat transfer coefficient (W.m−2.K−1)

hsw =

solid-to-wall convective heat transfer coefficient (W.m−2.K−1)

hr =

radiation heat transfer coefficient (W.m−2.K−1)

h =

convection transfer coefficient (W.m−2.K−1)

kgef =

effective thermal conductivity of gas (W.m−1.K−1)

ksef =

effective thermal conductivity of pellets (W.m−1.K−1)

tr =

particles residence time (s)

Di =

tube inner diameter (m)

Do =

tube outside diameter (m)

PeL =

Péclet number

Qin =

input power at heating zone (W)

Re =

Reynolds number

Tg =

gas temperature (K)

Ts =

solar fuel temperature (K)

Tw =

reactor wall temperature (K)

T =

ambient temperature (K)

q =

heat flux (W.m−2)

q−1 =

backward shift operator

ug,vs =

gas and particles velocities (m.s−1)

kp,τ1 =

coefficients of the PI controller

na,nb =

polynomial orders for A and B

Np =

prediction horizon

Nu =

control horizon

Pr =

Prandtl number

Greek Symbols

ɛ =

Porosity

ρ =

mass density (kg.m−3)

References

1.
Azad Gilani
,
H.
, and
Hoseinzadeh
,
S.
,
2021
, “
Techno-economic Study of Compound Parabolic Collector in Solar Water Heating System in the Northern Hemisphere
,”
Appl. Therm. Eng.
,
190
, p.
116756
.
2.
Nyarko
,
F. K. A.
,
Takyi
,
G.
, and
Amalu
,
E. H.
,
2020
, “
Robust Crystalline Silicon Photovoltaic Module (c-Si PVM) for the Tropical Climate: Future Facing the Technology
,”
Sci. African
,
8
, p.
e00359
.
3.
Chen
,
Y.
,
Liu
,
Y.
,
Wang
,
F.
,
Guan
,
X.
, and
Guo
,
L.
,
2021
, “
Toward Practical Photoelectrochemical Water Splitting and CO2 Reduction Using Earth-Abundant Materials
,”
J. Energy Chem.
,
61
, pp.
469
488
.
4.
Pang
,
H.
,
Yang
,
G.
,
Li
,
P.
,
Huang
,
H.
,
Ichihara
,
F.
,
Takei
,
T.
, and
Ye
,
J.
,
2020
, “
Wafer-Scale Si Nanoconed Arrays Induced Syngas in the Photoelectrochemical CO2 Reduction
,”
Catal. Today
,
339
, pp.
321
327
.
5.
Wang
,
Y.
,
Wang
,
H.
, and
He
,
T.
,
2021
, “
Study on Nanoporous CuBi2O4 Photocathode Coated With TiO2 Overlayer for Photoelectrochemical CO2 Reduction
,”
Chemosphere
,
264
, p.
128508
.
6.
Han
,
N.
, and
Ho
,
J. C.
,
2014
, “3—One-Dimensional Nanomaterials for Energy Applications,”
Nanocrystalline Materials (Second Edition)
,
S-C
Tjong
, ed.,
Elsevier
,
Oxford
, pp.
75
120
7.
Li
,
W.
,
Klemeš
,
J. J.
,
Wang
,
Q.
, and
Zeng
,
M.
,
2022
, “
Salt Hydrate–Based Gas-Solid Thermochemical Energy Storage: Current Progress, Challenges, and Perspectives
,”
Renew. Sustain. Energy Rev.
,
154
, p.
111846
.
8.
Fahim
,
M. A.
, and
Ford
,
J. D.
,
1983
, “
Energy Storage Using the BaO2BaO Reaction Cycle
,”
Chem. Eng. J
,
27
(
1
), pp.
21
28
.
9.
Bayon
,
A.
,
Bader
,
R.
,
Jafarian
,
M.
,
Fedunik-Hofman
,
L.
,
Sun
,
Y.
,
Hinkley
,
J.
,
Miller
,
S.
, and
Lipiński
,
W.
,
2018
, “
Techno-economic Assessment of Solid–Gas Thermochemical Energy Storage Systems for Solar Thermal Power Applications
,”
Energy
,
149
, pp.
473
484
.
10.
Hutchings
,
K. N.
,
Wilson
,
M.
,
Larsen
,
P. A.
, and
Cutler
,
R. A.
,
2006
, “
Kinetic and Thermodynamic Considerations for Oxygen Absorption/Desorption Using Cobalt Oxide
,”
Solid State Ionics
,
177
(
1–2
), pp.
45
51
.
11.
Karagiannakis
,
G.
,
Pagkoura
,
C.
,
Halevas
,
E.
,
Baltzopoulou
,
P.
, and
Konstandopoulos
,
A. G.
,
2016
, “
Cobalt/Cobaltous Oxide Based Honeycombs for Thermochemical Heat Storage in Future Concentrated Solar Power Installations: Multi-cyclic Assessment and Semi-quantitative Heat Effects Estimations
,”
Sol. Energy
,
133
, pp.
394
407
.
12.
Randhir
,
K.
,
King
,
K.
,
Rhodes
,
N.
,
Li
,
L.
,
Hahn
,
D.
,
Mei
,
R.
,
AuYeung
,
N.
, and
Klausner
,
J.
,
2019
, “
Magnesium-Manganese Oxides for High Temperature Thermochemical Energy Storage
,”
J. Energy Storage
,
21
, pp.
599
610
.
13.
Abuseada
,
M.
, and
Ozalp
,
N.
,
2020
, “
Experimental and Numerical Study on a Novel Energy Efficient Variable Aperture Mechanism for a Solar Receiver
,”
Sol. Energy
,
197
, pp.
396
410
.
14.
Petrasch
,
J.
,
Osch
,
P.
, and
Steinfeld
,
A.
,
2009
, “
Dynamics and Control of Solar Thermochemical Reactors
,”
Chem. Eng. J.
,
145
(
3
), pp.
362
370
.
15.
AlSahlani
,
A.
,
Randhir
,
K.
,
Ozalp
,
N.
, and
Klausner
,
J.
,
2022
, “
A Forward Feedback Control Scheme for a Solar Thermochemical Moving Bed Counter-current Flow Reactor
,”
ASME J. Sol. Energy Eng.
,
144
(
3
), p.
031004
.
16.
Abedini Najafabadi
,
H.
, and
Ozalp
,
N.
,
2018
, “
Aperture Size Adjustment Using Model Based Adaptive Control Strategy to Regulate Temperature in a Solar Receiver
,”
Sol. Energy
,
159
, pp.
20
36
.
17.
Al Sahlani
,
A.
,
Randhir
,
K.
,
Ozalp
,
N.
, and
Klausner
,
J.
,
2022
, “
A Simplified Numerical Approach to Characterize the Thermal Response of a Moving Bed Solar Reactor
,”
ASME J. Therm. Sci. Eng. Appl.
,
14
(
8
), p.
081010
.
18.
Randhir
,
K.
,
King
,
K.
,
Rhodes
,
N.
,
Li
,
L.
,
Hahn
,
D.
,
Mei
,
R.
,
AuYeung
,
N.
, and
Klausner
,
J.
,
2022
, “
A Continuum Model for Heat and Mass Transfer in Moving-Bed Reactors for Thermochemical Energy Storage
,”
Appl. Energy
,
313
, p.
118842
.
19.
Hayes
,
M.
,
Masoomi
,
F.
,
Schimmels
,
P.
,
Randhir
,
K.
,
Klausner
,
J.
, and
Petrasch
,
J.
,
2021
, “
Ultra-high Temperature Thermal Conductivity Measurements of a Reactive Magnesium Manganese Oxide Porous Bed Using a Transient Hot Wire Method
,”
ASME J. Heat Transfer-Trans. ASME
,
143
(
10
), p.
104502
.
20.
Wen
,
D.
, and
Ding
,
Y.
,
2006
, “
Heat Transfer of Gas Flow Through a Packed Bed
,”
Chem. Eng. Sci.
,
61
(
11
), pp.
3532
3542
.
21.
Huang
,
W.
,
Korba
,
D.
,
Randhir
,
K.
,
Petrasch
,
J.
,
Klausner
,
J.
,
AuYeung
,
N.
, and
Li
,
L.
,
2022
, “
Thermochemical Reduction Modeling in a High-Temperature Moving-Bed Reactor for Energy Storage: 1D Model
,”
Appl. Energy
,
306
, p.
118009
.
22.
Stenberg
,
V.
,
Sköldberg
,
V.
,
Öhrby
,
L.
, and
Rydén
,
M.
,
2019
, “
Evaluation of Bed-to-Tube Surface Heat Transfer Coefficient for a Horizontal Tube in Bubbling Fluidized Bed at High Temperature
,”
Powder Technol.
,
352
, pp.
488
500
.
23.
Touloukian
,
Y. S.
, and
Buyco
,
E. H.
,
1970
,
Specific Heat: Nonmetallic Solids
,
IFI/Plenum
,
New York
.
24.
Abuseada
,
M.
, and
Ozalp
,
N.
,
2020
, “
Experimental and Numerical Study on Heat Transfer Driven Dynamics and Control of Transient Variations in a Solar Receiver
,”
Sol. Energy
,
211
, pp.
700
711
.
25.
Camacho
,
E. F.
, and
Alba
,
C. B.
,
2013
,
Model Predictive Control
,
Springer
,
London
.