Solar cooking potential in Switzerland: nodal modelling and 1 optimization

7 Solar cooking is one possible solution to reduce the domination of fossil fuel in the domestic sector and to beneﬁt from renewable energy. This study assesses the solar cooking potential in Switzerland. A nodal model, based on energy balance equations of a box-type solar cooker is implemented in Matlab. Model parameters that cannot be determined experimentally or analytically are evaluated through an optimization procedure based on a Genetic Algorithm (GA). The model is able to predict the temperature of the cooking vessel with an average relative error around 5%. Based on its reliability, the model is simulated over a year for diﬀerent locations in Switzerland in order to determine the solar cooking potential. It is characterized by a metric that represents the number of days in a year the oven could be used to cook potatoes for two persons. It is found that the cooking times of potatoes can be well predicted by an Arrhenius law with an activation energy of 74.14 (cid:2) kJmol (cid:3) . The potato cooking criterion is based on the Arrhenius equation and determines if the pot simulated temperature proﬁle of a particular day allows to cook potatoes. The North-East of Switzerland is the least favourable area for solar cooking with theoretically around 155 cooking days per year. Around 240 days are estimated to be suitable for cooking in the cantons of Valais and Grisons, which represents a signiﬁcant potential for solar cooking in Switzerland.


16
Solar energy is abundant on the Earth's surface but, its low density and the important 17 mass contained in the atmosphere do not allow to naturally reach temperatures high enough 18 to cook food. An early description of a solar cooker or solar oven appears in the work of the Each of the nodes has two physical properties: a temperature T and a capacitance C (Associates, 2000). The capacitance is defined as the product between the mass m i and the specific heat capacity c pi of the corresponding node i : (1) The thermal resistance R ij characterizes the heat exchange between the nodes i and j. It is also common to refer to the conductance G ij , defined as the inverse of the thermal resistance : The governing equations of the nodal model are the heat balance equations applied for each node. The heat transfers from all the connecting nodes j must be taken into account. For node i, it is written as : The time derivative of the temperature is discretized with Euler schemes. To avoid problems with the stability of the solution, an implicit Euler method is preferred : By fixing the initial conditions T 0 i , Equation 4 can be solved to determine the temperature at 102 each time step.

103
The solar oven ULOG is divided into the seven following nodes:

104
• Node 1: Cooking vessel 105 • Node 2: All surfaces corresponding to the inner wall 106 • Node 3: All surfaces corresponding to the outer wall 107 Figure 3: Schematic representation of the ULOG oven and the corresponding nodal network, conductances, energy inputs (red arrows) and energy loss (blue arrows) Table 1: ULOG oven dimensions. The subscript number corresponds to the identified numbers as referred to in Fig.3.  The surface area of the other components can be found in Table 1. The red arrows represent 114 the energy inputs of the system. These inputs are modelled as solar gains G through the glass 115 layers. A fixed fraction of the solar gains, depending on the solar oven geometry, is absorbed by 116 the inner wall and the rest by the cooking vessel. A value of 20% has been calculated for this 117 solar cooker, based on experimental data collected at the laboratory. In order to validate the nodal model, an experiment is set up on the laboratory roof to 125 measure the temperature inside the oven when exposed to solar irradiation (see Fig. 4). The 126 solar oven with the cooking vessel inside is placed to face the South, away from any shadow. The experiments detailed above are performed during transient phenomena, such as the heating up of the oven. But it is also interesting to perform experiments at stationary conditions to evaluate, for instance, the U -value of the oven's walls and window. U is the overall heat transfer coefficient. This coefficient is related to the heatQ going through the element with surface area A by the relation 5 :Q = AU ∆T, where ∆T is the difference between the room temperature and the air temperature inside the 140 oven.

141
Such a test is performed in the laboratory by heating up the air inside the solar oven with 142 two electric resistances that deliver 20 W of heat each (see Fig. 4

152
G 46 and G 24 represent the heat losses due to the thermal bridges of the oven. They are difficult to measure but they can be deduced from the U -value experiment. The heat balance at steady-state is written as : The heat supplied by the heatersQ heater equals the heat losses through the walls, the window 153 and the thermal bridges of the window frame and the oven. ψ f rame and ψ oven are the linear 154 thermal transmittances of the thermal bridges and L f rame and L oven are the length of the thermal 155 bridges. ψ f rame L f rame and ψ oven L oven are respectively equal to G 46 and G 24 and since the other 156 terms of Equation 31 are known, the thermal bridges can be determined. The thermal losses due 157 to the thermal bridges are equally distributed among G 24 and G 46 . The infiltration rate is measured by artificially increasing the concentration of carbon dioxide in the oven and monitoring its decrease over time. Both the infiltration rate and the concentration of carbon dioxide of the room are assumed constant. In that case, the concentration of CO 2 in the oven follows a diffusion equation whose solution corresponds to an exponential decay. By fitting the measurements with a function such as ae −λt , the air leakage mass flowṁ leak can be found through the parameter λ which represents the infiltration rate. The air leakage mass flow is simply given byṁ where ρ air is the air density and V oven is the air volume inside the oven. The conductance G 47 is therefore evaluated as the quantity with c p,air being the specific heat capacity of air. Depending on the physical heat transfer mode they represent, conductances are evaluated in different manners. Conductive conductance is defined as where k is the thermal conductivity, A i the cross-sectional area through which heat flows and L 162 is the length between the two nodes.

163
For convective heat transfers, the following formula is used where h is the convective heat transfer coefficient and A i is the area of the node i in contact 164 with the fluid.

165
The heat transfer coefficient is evaluated through the Nusselt number N u with the relation where L is a characteristic length and k is the conductivity of the fluid, air in the present case. According to Churchill and Chu (1975) with the Prandtl number P r set to 0.7 and the Rayleigh number Ra.
The characteristic length is this time the ratio between the area A and the perimeter P of the horizontal surface. The total heat transfer coefficient between the cooking vessel and the inner air, h 17 is evaluated as a weighted average between the lid and the lateral surfaces contributions, as  Equation 12 holds for the inclined glazing with the difference that the gravity g has to be replaced by its component parallel to the inclined glazing. This is given by the quantity g sin π 6 . From Saha et al. (2007), the inclination and the space between the two layers of the glazing are too high to neglect the convection mode. The Nusselt number characterizing the convection between the two glass layers is given by the following formula from Incropera et al. (2013): It is a function of the aspect ratio, the height H and width L of the enclosure.

177
The formulation for the radiative exchanges is slightly different than the conduction and convection modes. The heat fluxq is indeed dependent on the temperatures to the power four. In order to have analytical solutions available, the surfaces are assumed to be gray diffuse. It means that the radiative properties of the surface such as emissivity, absorptivity, reflectance are independent of the wavelength and the direction of the irradiation (Modest, 2003). For such surfaces, net radiative heat fluxq i from surface i among N surfaces is defined as: E b,i is the black body emissive power of the surface i, evaluated as σ i T 4 i . σ is the Boltzmann's constant equal to 5.67 · 10 −8 W m 2 K 4 and i is the emissivity of the node i. F ij is the view factor between nodes i and j and only depends on the geometry. When N = 2, one equation for eacḣ q i is derived from Equation 17. This may be solved forq 1 by eliminatingq 2 8 18.

179
For simple geometries such as two parallel flat plates, the view factors are equal to 1 and Equation 18 reduces to a simple formQ with the conductance given by Equation 20 is valid for the radiative exchange between the two glass layers. beams while g dif f is defined as the value of g dir at 60 • of incidence.

200
The Zenith angle of the Sun is assumed to be such that the solar beams are always perpendicular to the glass window. Its position in the plan perpendicular to the glazing is then determined by the incident angle γ. g dir is written as a polynomial function of γ: The reference γ = 0 • is determined in the normal direction of the glass window. Experiments  Table 2.

203
The solar irradiation is divided into a direct component i dir and a diffuse component i dif f . The irradiation coming from reflection from the ground is treated with an albedo coefficient A bd of 0.2. Solar gains are therefore the sum of three components: the direct irradiation, the diffuse where A * glass and A glass are respectively the apparent and the real surface area of the glass Sun angle. In that case the apparent surface of the reflector A * ref l is equal to A glass cos π 6 , 221 knowing the glazing and the reflector have the same surface area.

222
The solving method is summarized in Fig. 5.  Single objective optimization is conducted. The objective function to minimize uses the least squares method, which is the square of the difference between the measured and the simulated cooking vessel temperature, T 1,mes and T 1 respectively. In order to increase the validity of the optimization, data of all the available measurements are given to the objective function. This yields the following expression: where the index n refers to the time step.

228
To assess the robustness of the optimization procedure, the same problem is optimized with The only constraint is that the capacitances of the two glass layers C 5 and C 6 must be equal, The cooking process of the potatoes is assumed to follow an Arrhenius type behaviour: The exponential term depends on the activation energy φ of the reaction, the universal gas constant R and the temperature T . The constant B is a characteristic of the chemical reaction. It is common to use this theory in lifetime estimation of items under specific constraints (Schüler et al., 2000) and in cooking processes (Petrou et al., 2002) which both involve chemical processes. In the first case, the rate constant r is replaced by the lifetime and by the cooking time t c in the second case. With potatoes, the main reaction is the "cooking" of the starch, or gelatinization, whose activation energy must be determined. Equation 24 is thus rewritten as: It must be noted that in a 1 T against ln t c plot, Equation 25 becomes a straight line whose slope is  To ensure the potatoes are ready, a cooking criterion CC is defined as: This quantity is similar for the four cooking tests, with variations less than 5%, and is a good indicator of the doneness of the potato. Thus for a day to be cookable, it must fulfill the following requirement: with t 0 and t end being respectively the starting time and ending time of the cooking period. conditions and is determined using Eq. 28.
where T ps is maximum absorber plate temperature, T as is ambient air temperature (at stag-284 nation) and H s is the insolation on a horizontal surface at the stagnation time (in W m −2 ).

285
F 2 and P s are obtained from measurements when the solar oven is loaded with water. F 2 is 286 calculated using Eq. 29. where F is heat exchange efficiency factor, η 0 is optical efficiency, C R is heat capacity ratio,

288
M is the mass of water in kg, C is the heat capacity of water J kg −1 K −1 , A is aperture area, τ 289 is time interval in seconds, T w1 is initial temperature of water, T w2 is final temperature of water,

290
T a is average ambient air temperature and H is the average solar radiation on horizontal surfaces 291 in W m −2 .

292
And P s using Eq. 30: where ∆T is a difference of temperature of 50 • C and ∆t the time needed to obtained it in 294 seconds.

296
Results of two temperature measurements on the lab roof are presented in Fig. 6 and 7.   the walls area A wall which is taken as the mean between the inner and the outer surface areas.

326
In a similar way that for G 23 , an overall heat transfer coefficient for the window U glass of

329
G 46 and G 24 represent the heat losses due to the thermal bridges of the oven. They are difficult to measure but they can be deduced from the U -value experiment. The heat balance at steady-state is written aṡ Q heater = U wall A wall ∆T + U glass A glass ∆T + (ψ f rame L f rame + ψ oven L oven ) ∆T The heat supplied by the heatersQ heater equals the heat losses through the walls, the window and the thermal bridges of the window frame and the oven. ψ f rame and ψ oven are the linear thermal transmittances of the thermal bridges and L f rame and L oven are the length of the thermal bridges.Q heater can be written as function of the global heat transfer coefficient of the oven to get U tot A tot = U wall A wall + U glass A glass + ψ f rame L f rame + ψ oven L oven (32) By recognizing that ψ f rame L f rame and ψ oven L oven are respectively equal to G 46 and G 24 , since the other terms of Equation 32 are known, the following result is available As a first assumption, the thermal losses due to the thermal bridges are equally distributed 330 among G 24 and G 46 . Thus, they are both equal to 0.015 W K . These heat losses due the thermal 331 bridges represent slightly more than 5% of the total heat losses through the envelope.

332
The infiltration rate is found to be equal to 0.5 vol h . Knowing the infiltration rate λ, the air leakage mass flow is simply given bẏ m leak = λρ air V oven 3600 = 4.63 · 10 −6 kg s where ρ air is the air density and V oven is the air volume inside the oven. The air properties are given in the appendices and the oven volume is equal to 0.031 [m 3 ]. The conductance G 47 is therefore evaluated as the quantity with c p,air being the specific heat capacity of air. It can be noted here that the heat losses due 333 to air leakage are assumed to be very small and as such are do not play an important role in the 334 heat balance.

335
It is assumed that the heat transfer between the node 1 and 2 is largely dominated by the conduction through the bottom of the cooking pot compared to the radiative heat transfer. The conduction through the bottom of the cooking vessel A cv,bot is evaluated with the following formula The thermal resistance due to the interface between the pot and the inner wall is taken into       Figure 14: Schematic representation of the nodal network, conductances and energy inputs (arrows) the two curves are quite different. This is obvious in Fig. 13 where the simulated curve heats up 377 too quickly.

378
Based on the results from the optimization, the parameter values given in Table 5 are used 379 in the model.  Table 6, the error is almost divided by 384 two for the first measurement. The others are also slightly smaller than in the scenario with the 385 five measurements given to the objective function. It is interesting to note that the optimized 386 values of the parameters are similar to the previous scenario except C 2 whose optimized value 387 decreases to 2119.6 J K .

388
All the simulated temperatures of the oven nodes are presented in Fig. 17. While the main 389 focus was on T 1 , it is also interesting to analyse the other results of the model simulation. As   transfer coefficients for the confined air present also typical values for such convection processes.

404
Tests have been made with data set without water load in the cooking vessel. Due to the high 405 volatility of the solar gains and the relatively large time step compared to these variations, the 406 fit tends to be more difficult. An example of such results is shown in Fig. 19. The relative error 407 is also quite low and the red curve follows the general trend of the blue curve. The simulated 408 temperature however oscillates too much and reaches too high temperatures. The heat losses 409 due the thermal bridges represent slightly more than 5% of the total heat losses through the      The model seems to be in good accordance with the experimental data in particular as 459 compared to the maximum temperature obtained. This is comparable to results from Zafar with previous studies (Guidara et al., 2017). The focus has mainly been on the results with a 462 water load, which represents the principal situation in which the oven will be practically used.

463
This model with high capacitance is optimized more easily than sensitive models. It tends to Tests have been performed with smaller time steps in order to capture these sharp variations.

470
However, the data set and the computing time become very large and the fitting is quite difficult.

471
It can also be noted in Figs The very good fit of the potato cooking times in Fig. 20 demonstrates that the cooking 479 process is well described by an Arrhenius law. The cooking times of potatoes double when the 480 cooking temperature decreases by 10 • C. In Shiotsubo (1984), an activation energy for the 481 starch gelatinization of around 99 ± 25 kJ mol is presented. The value found in the present study 482 is slightly lower but still in good agreement with the literature. The difference could be explained 483 by the fact that the whole potatoes cooking process is characterized in this study while only the 484 starch gelatinization reaction is considered in (Shiotsubo, 1984). This chemical reaction may be 485 dominant but it is not the only chemical process involved in cooking.

486
The results of the simulation over a year show the good potential for solar oven utilization however necessitate real irradiation data for the whole Swiss territory which is quasi impossible 510 to obtain. The interpolation is useful to get an approximation of such data but its results must 511 be interpreted with care.

513
This study has presented an innovative nodal model for the Swiss ULOG solar oven, which 514 predicts the cooking vessel temperature's profile along the day with solar irradiation and ambient 515 temperature as inputs.

516
The model parameters, depending on the physical properties of the oven, are determined an-517 alytically when correlations are available, experimentally or through an optimization procedure.
518 Overall, the model is able to predict the pot temperature with average relative errors around 519 5%. In an upcoming study, the model will be improved by using an adaptive infiltration rate 520 which will be based on the formulations proposed by Malik (1978); Han et al. (2015) and by 521 incorporating new measured data that will be collected for the solar cooker.

522
Based on these reliable results and thanks to the new metric developed, the model is used to 523 assess the solar cooking potential over the Swiss territory. The North-East of Switzerland is the 524 least favourable region for solar cooking, with theoretically around 155 cooking days per year.

525
On the other hand, almost 240 cooking days are reached in the Valais and the Grisons. Similar 526 results are obtained even for locations at high altitude. This is an encouraging conclusion, which 527 supports the initial assumption that solar cooking is suitable for remote places such as high 528 altitude areas. The study highlights a significant potential for solar cooking in Switzerland.

529
Nowadays, solar cooking is still anecdotal in Switzerland and the potential is far from being 530 fully exploited. Obviously, solar cooking will not completely replace standard cooking, due to its 531 long cooking times, the good solar exposure needed and habits. However solar cooking interest is 532 likely to increase and the corresponding awareness might contribute to reducing the high energy 533 consumption per household in Switzerland.