Bread baking oven design optimisation

Commercial bread-baking industry has been showing increasing interest in computational fluid dynamics (CFD) to model the complex mass and heat transfer processes involved in bread baking, with the aim of reducing cooking times without compromising bread quality.

In forced convection ovens (also called “direct fired” ovens), a common oven typology, hot air is injected in the baking chamber by nozzles located on the ceiling and the floor of the chamber. Bread loaves are placed on a tray located at equal distance H from the ceiling and the floor, as shown schematically in Fig. 1.

Figure 1: Three-zone direct fired oven: a) Overview of the oven b) Simplified schematic showing the oven longitudinal section and the mechanism for distributing air through the nozzles on the ceiling and on the floor

Air temperature uniformity in the baking chamber is a key factor to ensure uniform heat transfer and then product quality. Following Khatir et al. (2011ab), temperature uniformity in an oven baking chamber can be defined by the parameter T defined in Eq. (1):

\sigma_{T}=\sqrt{\dfrac{\int_{V} (T_{i}-T_{zone})^{2}dV}{\int_{V}dV}} \quad  \quad (1)

where T_{zone} is the temperature of the air flowing through the nozzles fitted on the baking chamber walls, V is the baking domain and T_{i} is the air temperature at point i of the baking chamber. By the definition in Eq. (1), the optimisation of the oven design consists in finding the set of design parameters that minimises the root mean square temperature variation \sigma_{T}.


In mathematical terms, the oven design optimisation implies finding a global minimum of the function \sigma_{T}. HyGP was used to generate a global metamodel of \sigma_{T} : due to oven simmetry, \sigma_{T} was modelled only on a portion of the baking chamber. Three parameters defining oven and nozzles characteristics were considered as independent variables, as done in Khatir et al. (2011b):
– the nozzle jet diameter (D),
– the dimensionless nozzle-to-surface distance (H/D)
– the nozzle jet velocity (u_{noz}).

The CFD model of the oven baking chamber and the input variables are described in Fig. 2.

Figure 2: CFD baking chamber model with design variables: nozzle jet diameter D, jet velocity u_noz, and distance H between tray surface and chamber ceiling (nozzle to nozzle distance S=200 mm)$

CFD simulations provided \sigma_{T} values on a Latin hypercube DoE made of 30 points, used as building data set. As by definition \sigma_{T} is positive in the whole domain, such knowledge was exploited to improve HyGP symbolic regression through the use of the penalisation introduced in Armani 2014, Section 5.4, Chapter 5. The sign of \sigma_{T} was checked on an additional data set C made of 120 points generated using a full factorial DoE in the region [5; 20] [2; 10] [8; 40]. For the HyGP experiment, 12 independent evolutions were performed, each using a population of 300 individuals for a maximum number of 200 generations. As functional primitives addition, subtraction, multiplication, power, shift, square, cube, sine, cosine, hyperbolic sine, hyperbolic cosine, hyperbolic tangent, exponential were selected.

The best metamodel found by HyGP is shown in Eq. (2):

\sigma_{T}(D, \tfrac{H}{D}, u_{noz})= 6.54966\times10^{-6}\, {\mathrm{D}}^3\, \mathrm{u_{noz}} + 4.26492\times10^{-6}\, {\mathrm{D}}^2\, {\mathrm{u_{noz}}}^2 - 1.84683\times10^{-4}\, \mathrm{D}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^2\, \mathrm{u_{noz}} + 2.69398\times10^{-5}\, \mathrm{D}\, \left(\mathrm{\tfrac{H}{D}}\right)\, {\mathrm{u_{noz}}}^2 + 4.97275\times10^{-2}\, \mathrm{D}\, \left(\mathrm{\tfrac{H}{D}}\right) - 6.163\times10^{-1}\, \mathrm{D} + 5.77163\times10^{-6}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^5\, {\mathrm{u_{noz}}}^2  + 2.13294\times10^{-8}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^4\, {\mathrm{u_{noz}}}^4 - 1.68813\times10^{-4}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^4\, {\mathrm{u_{noz}}}^2  - 1.07612\times10^{-2}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^4 - 1.76832\times10^{-7}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^3\, {\mathrm{u_{noz}}}^4  + 9.4134\times10^{-4}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^3\, {\mathrm{u_{noz}}}^2 + 2.83469\times10^{-2}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^3\, \mathrm{u_{noz}} + 5.52956\times10^{-4}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^2\, {\mathrm{u_{noz}}}^2 - 2.96942\times10^{-1}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^2\, \mathrm{u_{noz}}  + 5.81804\times10^{-1}\, {\left(\mathrm{\tfrac{H}{D}}\right)}^2 + 5.57406\times10^{-1}\, \left(\mathrm{\tfrac{H}{D}}\right)\, \mathrm{u_{noz}} + 1.234\, \left(\mathrm{\tfrac{H}{D}}\right) + 7.574\times10^{-7}\, {\mathrm{u_{noz}}}^4 + 5.78694\times10^{-4}\, {\mathrm{u_{noz}}}^3 - 4.88598\times10^{-2}\, {\mathrm{u_{noz}}}^2 + 5.79625\times10^{-1}\, \mathrm{u_{noz}} + 2.86242 \quad  \quad (2)

Although the generated metamodel shown above is quite accurate on the building data set (RMSE = 4.62702e-02, R2 = 0.999296), its response was found to be positive only inside the hypercube [5; 20] [2:75; 8] [8; 39:7]. This may be explained considering that the used HyGP aggregated fitness gathers conflicting objectives: RMSE minimisation on building data set may hinder the minimisation of the distance of the metamodel responses from the feasible region. An attempt to visualise the behaviour of the metamodel at the boundaries of the hypercube [5; 20] [2:75; 8] [8; 39:7] is shown in Fig. 3a.

Results of the optimisation

In the second stage of the analysis, a GA was used to find the minimum of the \sigma_{T} metamodel as defined in Eq. (2) in the hypercube [5; 20] [2:75; 8] [8; 39:7]. In Table 1 the minimum found is compared with the minimum of another \sigma_{T} metamodel generated by MLSM using the same input data. The values returned by the metamodels were validated with the responses returned by CFD simulations in the minima locations (Khatir et al. 2011b), reported in fourth column.

Table 1: Minima found using metamodels generated by HyGP and MLSM
technique   {D;H/D; u_{noz}}   \sigma_{T}   \sigma_{T} from CFD   rel. error
{ [mm], [ ], [m/s] }   [K]  [K]   %
MLS + GA      {20.00, 6.82, 38.12}   1.16    1.22   4.9%
HyGP + GA    {13.80, 7.31, 39.70}    1.07   1.14   6.1%

The estimates provided by HyGP and MLSM metamodels are affected by relative errors of the same order of magnitude (see fifth column in Table 1). The small size of the building data set may have been the main cause of the poor accuracy of the metamodels.


Despite its lower accuracy, HyGP metamodel helped discover an optimal design characterised by a value of \sigma_{T} (1.14) lower than the one corresponding to the optimum found on MLSM metamodels (1.22). Moreover, HyGP metamodel explicit expression has the intrinsic advantage that it can be easily used to perform Monte Carlo-based sensitivity analysis.

The \sigma_{T} metamodel defined in Eq. (2) is shown in Fig. 3b, 3c, 3d setting each variable to a constant value: the dot in the figures represents the minimum found by GA.

a) \sigma_T on design space boundaries
b) \sigma_T for D=13.80 mm
c) \sigma_T for H/D=7.31
d) \sigma_T for u_noz=39.70 m/s

Figure 3: Behaviour of \sigma_{T} metamodel (Eq. (2)) on the boundaries of the hypercube [5; 20] [2:75; 8] [8; 39:7] and around the minimum



  • Z. Khatir, J. Paton, H. Thompson, N. Kapur, V.V. Toropov, M. Lawes, and D. Kirk. Computational fluid dynamics (cfd) investigation of air flow and temperature distribution in a small scale bread-baking oven. Applied Energy, In Press, Corrected Proof, 2011a
  • Z. Khatir, H. Thompson, N. Kapur, V.V. Toropov, J. Paton, and M. Lawes. The application of computational fluid dynamics (cfd) and oven design optimisation in the british breadbaking industry. In Proceedings of the Eighth International Conference on CFD in Oil & Gas, Metallurgical and Process Industries (SINTEF/NTNU), Trondheim, Norway, 2011b.
  • U. Armani, Development of a hybrid genetic programming technique for computationally expensive optimisation problems. PhD thesis, School of Civil Engineering, University of Leeds, UK, 2014