Accepted Manuscript

Experimental study of dynamic porosity and its effects on simulation of the coffee
beans roasting

Nadia O. Oliveros, J.A. Hernández, F.Z. Sierra-Espinosa, R. Guardián-Tapia, R.

Pliego-Solórzano

PII: S0260-8774(16)30455-1
DOI: 10.1016/j.jfoodeng.2016.12.012
Reference: JFOE 8737

To appear in: Journal of Food Engineering

Received Date: 26 September 2015

Revised Date: 2 December 2016
Accepted Date: 17 December 2016

Please cite this article as: Oliveros, N.O., Hernández, J.A., Sierra-Espinosa, F.Z., Guardián-Tapia, R.,
Pliego-Solórzano, R., Experimental study of dynamic porosity and its effects on simulation of the coffee
beans roasting, Journal of Food Engineering (2017), doi: 10.1016/j.jfoodeng.2016.12.012.

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to
our customers we are providing this early version of the manuscript. The manuscript will undergo
copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please
note that during the production process errors may be discovered which could affect the content, and all
legal disclaimers that apply to the journal pertain.
 ACCEPTED MANUSCRIPT

1 Experimental study of dynamic porosity and its effects on

2 simulation of the coffee beans roasting
3 Nadia O. Oliverosa, J.A. Hernándezb, F.Z. Sierra-Espinosab*, R. Guardián-Tapiab, R.
4 Pliego-Solórzanoa

PT
5
a
6 Postgraduate in Ingeniería y Ciencias Aplicadas de la Universidad Autónoma del Estado

RI
7 de Morelos (UAEM), Av. Universidad No.1001, Col. Chamilpa, Cuernavaca, Morelos,
8 C.P. 62209, México.

SC
b
9 Centro de Investigaciones en Ingeniería y Ciencias Aplicadas (CIICAp), Universidad
10 Autónoma del Estado de Morelos (UAEM), Avenida Universidad No. 1001, Col.
Chamilpa, Cuernavaca, Morelos, C.P. 62209, México.

U
11
12 *Corresponding author: F.Z. Sierra: fse@uaem.mx; Phone number: +52 777 3297084
AN
13
14
M

15 ABSTRACT

16 We investigated coffee beans porosity to demonstrate its effects on predicting heat and
D

17 mass transfer during roasting. Calculated porosity from experimental data from electronic
TE

18 microscopy pore size and distribution for green bean and roasted beans were evaluated.

19 Porosity of 9.8% and 34.2% for green and roasted beans were found. These values defined
EP

20 the limits of a dynamic model, which was developed and introduced as subroutine into a 3-

21 D modeling of coffee roasting. Hot air at 473 K and 503 K roasted single and group beans
C

22 during 200 s for fixed Reynolds number. Modeled and experimental data of temperature
AC

23 from Arabica bean agreed well. It was observed that porosity affected heat and moisture

24 loss rate during roasting. The results can help to understanding the roasting process, or

25 studying processes of food products that require roasting by means of hot air, like grains or

26 seeds.
1
 ACCEPTED MANUSCRIPT

27

28

29 Keywords:
30 Coffee beans roasting, dynamic model coffee porosity, heat transfer, numerical prediction

PT
31
32 1. Introduction

RI
33 Coffee is one of the main traded products in international markets. Obtaining high quality

SC
34 coffee depends on the roasting process. As outlined elsewhere (Baggenstoss et al. 2007;

35 Bottazzi et al. 2012; Scheker et al. 2002; Raemy & Lambelet 1982), during their roasting

U
36 coffee beans develop specific organoleptic properties that determine flavor, aroma and
AN
37 color. Roasting can be divided into two main parts: the drying stage, when bean´s

38 temperature is below 433 K; then roasting, where this temperature increases up to 533 K.
M

39 Therein, when coffee temperature reaches 463 K, pyrolytic reactions start, causing

40 oxidation, reduction, hydrolysis, polymerization, decarboxylation and many other chemical

D

41 changes. It is believed that the bean temperature at which these phenomena occur is the
TE

42 main factor that defines the color, flavor and aroma of the beans, and thus the quality of the

43 coffee (Ortolá et al. 1998; Hernández-Pérez 2002). However, these studies reveal that void
EP

44 porosity coffee beans as property has not been considered into modeling the roasting
C

45 process. This is important because porosity may play a main role during roasting by
AC

46 affecting heat and mass transfer.

47 Roasting involves both complex energy (from the roaster to the bean) and mass transfer

48 (water vapor and volatile compounds from the bean to the environment) leading to change

49 coffee beans in terms of weight, density, moisture, color and flavor (Eggers. & Pietsch

50 2001; Eggers & Pietsch 2001). Endothermic process and reactions (water evaporation)
2
 ACCEPTED MANUSCRIPT

51 occur in the first stage of roasting, whilst the undesirable exothermic pyrolysis of

52 saccharides may occur at the latter roasting stages (Fabbri A. et al. 2011). This work

53 focused on determining coffee beans porosity by electronic microscopy. Frisullo et al.

54 (2012) reported values of coffee bean porosity for different times during roasting. We

PT
55 found that porosity of green coffee beans was 9.8%, while roasted beans were 34.2%

RI
56 porous. It was considered that these values were in good agreement with data from Frisullo

57 et al. (2012). Once implementing a dynamic model of porosity into numerical simulation,

SC
58 results showed that void porosity affected the roasting process by modifying the conduction

59 and sensible heat rates inside the beans. The main reason is that air remains in void pores,

60
U
which has a thermal conductivity many times smaller than coffee bean solid material.
AN
61 Therefore, it is expected that effects like delaying the bean´s moisture loss can be detected.

62 We propose porosity must be part of modeling coffee roasting because to the best of our
M

63 knowledge, porosity has not been yet considered in coffee roasting simulation.
D

64 Consequently, the objective of present paper is to demonstrate that consideration of

TE

65 porosity as a property of coffee beans affects the modeling of heat and mass transfer during

66 coffee beans roasting. The remainder is devoted to emphasize the effects of porosity on
EP

67 prediction of temperature and moisture loss during the coffee beans roasting process.

68 Modeling the quality of coffee during the roasting process is a difficult task. Schwartzberg
C

69 (2002), and Hernández et al. (2007) reported that bean´s temperature was a very important
AC

70 factor on its quality, while Hernández-Pérez (2002), Hernandez et al. (2008a, 2008b)

71 showed that bean´s color, which was a function of roasting temperature, can determine the

72 quality of the product. Therefore, it is extremely important for the industry to obtain

73 numerical tools to determine any of these variables: bean temperature, color, aroma and
3
 ACCEPTED MANUSCRIPT

74 flavor. Alonso-Torres et al. (2013) modeled and validated the heat and mass transfer for an

75 individual coffee bean during roasting using computational fluid dynamics (CFD). The

76 authors used a global heat and mass transfer model proposed by Schwatzberg (2002), which

77 determined temperature distributions in roasted coffee beans. Santos et al. (2015) proposed

PT
78 a real-time on-line monitoring of coffee roasting by near infrared spectroscopy to calculate

RI
79 acidity of coffee. The results were compared against analytical methods with good

80 agreement.

SC
81 Thus, in this investigation results of coffee roasting with CFD in 3D and time-dependent

82 modeling were proposed. Modeling of coffee beans roasting was addressed by considering

83
U
static beans without interacting with each other as a function of time and space, instead of
AN
84 multiple loose grains. A computational model that represents exactly the reality in a coffee

85 roaster is not available at this moment. In the meantime, a 3-D modeling of three fixed
M

86 coffee beans allowed getting experience for future trying to represent the real phenomenon.
D

87 A real modeling increases time for convergence and computer resources. However,
TE

88 compared against a 2D model, the present effort in 3D allowed for accurate predictions of

89 moisture loss. A 2-D model is far from predicting moisture loss; it is also unable of
EP

90 predicting the three-dimensional effects of convective streams around the beans, which

91 develop due to non-symmetry of coffee beans geometry. This topic alone can be discussed
C

92 in full in a separate paper, given its impact in heat transfer analysis between coffee beans
AC

93 and hot air stream.

94 Commercial and well tested computational multi-physics codes are available nowadays,

95 allowing an exceedingly flexible simultaneous numerical solution of the energy, mass and

96 moment equations (Scott & Richardson 1997) and are adopted in the food engineering field

4
 ACCEPTED MANUSCRIPT

97 (Bin & Da-wen 2002; Fabbri 2005). However, despite this huge development, properties

98 and structure of coffee materials like porosity and in-homogeneity remain unattended.

99 Indeed, modeling efforts focused towards basic understanding of transport phenomena

100 inside and around beans. Differential equations of mass, momentum and energy to solve

PT
101 inside outside the beans during roasting are now familiar (Fabbri et al. 2011). In parallel,

RI
102 the use of empirical correlations for mass and energy transfer are still in practice (Bottazi et

103 al. 2012). Therefore, an enhanced model includes conductive and convective terms to link

SC
104 interchange between coffee bean and hot air, to knowing the temperature distribution.

105 Moreover, present model resolves temperature changes among a group of coffee beans

106
U
during roasting instead of one bean alone. Finally, roasting process under effects of
AN
107 dynamic porosity was simulated. By considering porosity in the roasting simulation, heat

108 transfer at bean´s interior varied through the change from pure coffee material, to combined
M

109 coffee material plus air. Thermal conductivity, mass diffusion, specific heat; all these
D

110 properties inside the beans varied because air was present in the pores. Therefore, heat and
TE

111 mass transfer rates varied when introducing porosity. The term of heat conduction in energy

112 equations varied due to properties for air and coffee beans, in a porous media. The pores
EP

113 assumed filled with air, which introduces convection in the space between walls of the

114 pores. Connectivity among air-solid-air assumed diffusion for mass transport and
C

115 convection and conduction for heat transport. Mass transfer rate varied because as beans
AC

116 roast proceeds, mass diffusion of moisture changes as function of temperature and time.

117 The present model represents the phenomenon assuming a dynamic value of porosity along

118 the whole roasting process, based on experimental data of pores, their homogenous

119 distribution and size obtained by electron microscopy in this work and reported in

5
 ACCEPTED MANUSCRIPT

120 literature. Pore size estimated statistically allowed getting a percentage of void pore area

121 inside the beans.

122 Both, experiments and simulation involved injection of hot air to the coffee roaster at air

PT
123 temperatures of 473 and 503 K, during a time interval of over 200 s, for Reynolds number

124 between 10x103 and 40x103. The Reynolds number range in this work was set according to

RI
125 previous studies of coffee bean roasting, and these on time were defined following roaster

126 operating conditions already published in literature. Numerical modeling was validated by

SC
127 comparing against experimental results of interior temperature in Arabica coffee beans,

Pliego-Solórzano (2010). Therefore, a validated model lead to evaluate the influence of a

U
128

dynamic porosity model on moisture release in the bean during roasting.

AN
129

130
M

131 2 Material and methods

132 2.1 Measurement of porosity

D

133 Coffee beans change their aspect and volume during roasting as can be observed in Fig. 1.
TE

134 The volume of coffee beans increases due to porosity increment (Schenker et al., 2000;

135 Frisullo et al., 2012). Schenker et al., (2000) reported increased volume of roasted beans up
EP

136 to 1.8 times compared to green bean volume, for high temperature roasting. In addition, the
C

137 darker the bean the more fragile. This was because the materials become carbonized, which
AC

138 makes the properties different: color, hardness, taste. The roasting process can lead to burnt

139 beans, darker than the right hand side photograph of Fig. 1. Regarding the increment of

140 bean volume, it may be attributed to the presence of void pores, which changes the

141 composition to non-homogeneous solid material. This was confirmed by direct microscopy

142 observation of coffee beans structure. Slices of coffee beans analyzed by electronic
6
 ACCEPTED MANUSCRIPT

143 microscopy allowed measuring porosity. The word slice is used in this work for defining a

144 fragment of coffee, obtained by cracking the beans. The bean fragments were not actually

145 slices, but irregular small pieces of coffee bean material. Therefore, there was no a slice

146 thickness, but small pieces of irregular thickness according to requirements of microscopy

PT
147 observation. The microstructure was not affected by cracking the beans, because the force

RI
148 to produce the cracking was applied externally on the beans, such that irregular fragments

149 of 1-3 mm from the inner regions resulted and remained unaltered. These fragments were

SC
150 mounted on pins, which were mounted on a sampler to receiving a metalizing process

151 through sublimation of gold, Au/Pd (sputter coater). Metalizing the surface of samples

152
U
makes them ready for conduction of electrons, which allows observation of micro-structure.
AN
153 Cracking does not interact with bean´s material at micro-scale. In agreement with the scale

of observation found in the coffee beans, pores ∼20 µm, the scanning electronic microscopy
M

154

155 (SEM-LEO 1450 VP) was chosen for conducting the study of porosity. Microscopy
D

156 requirements in the sample preparation required taking care to not perturbing or deform the
TE

157 objective, their number and their size, which in this case were pores. Special attention

158 focused to the secondary electron detector, which provides a detailed description of pores
EP

159 topography through images that can resolve discontinuities, morphology, and homogeneity

160 of the objects (Moellmer et al. 2010; Thommes 2004).

C

According to Baxes (1994), a shape factor
can be calculated by measuring one pore
AC

161

162 ̅ and its perimeter  in a micrograph, which was identified by a

average surface area 

163 dark region, as:

 

= (1)

164 

7
 ACCEPTED MANUSCRIPT

165 A shape factor
= 1 means the pore represented a perfect circle. The measured surface area

166 of each pore  was thus integrated to obtain a value that represents the measured average

167 ̅ :
pore area, 

̅ =


∑
 

PT
(2)

168

169 Where  was the surface area of the ith pore located in one micrograph from each fragment

RI
170 of bean. To determine a characteristic size of pores,  we used the following expression:

SC

171  =  (3)


U
172 ̅ was normalized by the total surface area investigated in each micro-image, ̅ , in order

AN
173 to obtain void fraction porosity (ASTM 1996; Moellmer et al. 2010; Juang et al. 1986;

174 Pereyra et al. 2017):

M

̅
175 γ = ̅ (4)

D

176 The results presented in Table 1 show the main variables considered in equations (1-4).
TE

177 Two conditions were considered: first, plain green unroasted beans, left hand side of Fig. 1;
EP

178 and second, roasted coffee beans. Microscope scanning was divided into bean´s external

179 surface and interior. Both conditions were studied by cracking the beans to produce sample
C

180 pieces, as shown in Fig. 2a. Bean pieces were mounted on the support shown in Fig. 2b for
AC

181 receiving evaporation and then metalizing as part of the process for subsequent microscopy

182 analysis. A total of N = 650 observed data per each condition roasted/unroasted was

183 investigated to obtain an average estimation of the void fraction porosity.

8
 ACCEPTED MANUSCRIPT

184 Notice the bean´s coffee surface shown in Fig. 3a and 3b, un-roasted and roasted,

185 respectively. The difference of darkness of the external surface of the sectioned grains

186 shown in Figs. 3c and 3d confirm the different condition of roasted darker bean. The scale

187 of these images was 1 mm. However, a glance into micro scale of 20 µm distinguish

PT
188 aspects like the size of flaps and distances between them, The images of Figs. 3e and 3f for

RI
189 un-roasted and roasted grains, respectively, denote that roasted beans have increased

190 distances between pores and flaps in general.

SC
191 Micro images from internal coffee beans confirm previous results. See Fig. 4a for un-

U
192 roasted and 4b for roasted conditions, for a scale of 100 µm, which let observing more and
AN
193 larger black zones like pores in the roasted bean compared against un-roasted. Finally, by

194 increasing the scale to 20 µm, the materials structure look quite different. At this scale,
M

195 pores identified with circles marked in Figs. 4c and 4d, served for quantify a property of

196 coffee named void fraction. Black zones were several times larger in roasted beans,
D

197 confirming that the roasting process modifies coffee material additionally to influencing its
TE

198 chemistry. A summary of porosity results from equations (1) to (4) presented in Table 1

199 compares the porosity of both conditions. Last column indicates that roasted coffee porosity
EP

200 γ = 34.24%, was more than three times compared with un-roasted beans, γ = 9.88%.
C

201 Porosimetric curves were built using electron microscopy by Schenker at al. (2000), and
AC

202 later using X-ray by Frisullo et al. (2012). These studies extended the range from starting

203 the roasting process to roasted dark grains, as shown in Fig. 4d. Frisullo et al. (2012)

204 measured porosity for coffee Arabica, applying high temperature roasting (493 K) to a

205 sample of 100 g green coffee beans. They employed a laboratory scale roaster model BR 74

206 by Probat, fluidized by hot air. Samples of five beans were extracted at 2, 3, 4 and 5 min
9
 ACCEPTED MANUSCRIPT

207 during roasting process for determining the roasting degree. There is consensus on

208 considering that high-temperature roasted coffee has a statistically significantly higher

209 macropore area than low-temperature roasted products, Schenker at al. (2000). Therefore,

210 roasting temperature becomes the most important condition influencing the porosity of final

PT
211 product. Also, of great importance is the medium used in the roaster like direct flame, water

RI
212 steam or hot air. The roaster type is important too, although the heating medium defines

213 what devise is used. In addition, the pores change from origin of grain to grain Schenker at

SC
214 al. (2000). Conditions used by Frisullo et al. (2000) and by present investigation during

215 coffee roasting support the similarity of porosity values found for green and roasted grains:

216
U
grain type Arabica and high temperature for fluidized by hot air roasting. This is why
AN
217 Frisullo et al (2000) data were considered to complement the present results of porosity, in

218 order to simulating dynamic porosity change during the whole roasting process. Therefore,
M

219 the following expression for green through roasted beans porosity,  , is proposed:
D

220   = 0.0002 # + 0.0227 + 9.4573 (5)

TE

221 Where t represents the roasting time. Equation (5) was used to building a user defined
EP

222 function, UDF (see Apendix), which was feed into the main program, ANSYS (2011),

223 which constitutes a dynamic porosity model representing porosity variation as function of
C

224 time and temperature during roasting. We envisaged that porosity vary from one kind of
AC

225 coffee bean to another, from Arabica to Robusta for instance, which is still not addressed. If

226 true, a full predictive model will be even more complex. Furthermore, the chemical changes

227 during roasting may lead to consider products other than air, which represents another

228 challenge for modeling. Present simulation of roasting process assumed pore composition

10
 ACCEPTED MANUSCRIPT

229 as air, based on the scope of this investigation. Void fraction porosity influenced the

230 numerical computation of heat and mass transfer during roasting modeling as described in

231 the following sections.

PT
232 2.2 Experimental procedure for Arabica coffee roasting

233 Valdovinos-Tijerino (2005) and Heyd et al. (2007) carried out investigations on roasting

RI
234 100 g of green beans in a vertical column at 470-533 K, while Pliego-Solórzano (2010)

SC
235 developed vertical a roaster where beans remained suspended, encouraging contact with air

236 and reducing temperatures and roasting times. A fluidized bed was feed by upwards air

U
237 stream in a vertical borosilicate glass cylinder for easy viewing the process. The bed of
AN
238 stainless steel grid of 2 mm size elements hold the coffee beans in the bottom of the

239 column. Hot air supplied by an appliance turbo-ventilator provided adjusted flow and
M

240 temperature, for controlled heat source. Equaled column inner diameter and outer diameter

of appliance nozzle ensured good fluidization of the beans. An 11 cm height column was
D

241

242 according to the amount of beans to roast. Figure 5a contains a schematic of apparatus and
TE

243 main dimensions.

EP

244 A thermal couple type T of 1 mm diameter was inserted 4 mm into one coffee bean, which

245 allowed monitoring its temperature during roasting, as shown in Fig. 5b. Each test of fixed
C

246 airflow and temperature included measuring the air speed at the entrance of the column
AC

247 using a Hygro Thermo Anemometer, by Extech Instruments, Inc., shown in Figs. 5c and

248 5d.

249 Experimental data recorded using a serial data logger through an Agilent 34970A with HP

250 34901A multiplexer, connected to five analogic channels where: channel one recorded
11
 ACCEPTED MANUSCRIPT

251 temperature inside the bean; channels 2 and 3 recorded the column temperature; and

252 channels 4 and 5 recorded input and output air temperature, respectively. A GPIB Bus

253 interface, General Purpose Interface, connected the data logger and computer via a USB

254 port. A commercial software by Benchlink Data Logger 3, version 3.10, ease read out data

PT
255 through in a PC Toshiba Intel Pentium, Dual 1.73 GHz.

RI
256 Arábica coffee green beans were roasted for the conditions given in Table 2. Selected data

results obtained for tests at 473 K and 503 K served to compare and validate the simulation

SC
257

258 results, as described below.

U
259 2.3 Computational modeling
AN
260 We obtained converged solutions for time dependent transport differential equations of

261 mass, momentum and energy using the finite volume approach (Patankar (1980); Versteed
M

262 and Malalasekera (1995)), in a commercial program Fluent (ANSYS, 2011). Turbulent
D

263 flow based on the Reynolds average Navier-Stokes equations, RANS, was resolved, where
TE

264 pressure and velocity were coupled by the SIMPLE algorithm. Equations expressing

265 conservation of mass and momentum for an incompressible Newtonian fluid were:
EP

266 *·, =0 (6)

/,
-. + , · *,1 = −* + 3* # , − ∇5-,
6 7 ´,
68 ´9 + 
C

(7)
/0
267
AC

268 where u represented the velocity vector with velocity components u, v, w, in the coordinate

269 directions, x, y, z, respectively. In equation (7) P, 3 and Si read for pressure, viscosity and a

270 sink term in porous media explained below.

12
 ACCEPTED MANUSCRIPT

271 In RANS approximation the viscous term leads to an extra term, the Reynolds stress tensor,
272 ,
6 7 ´,
6 8 ´, which represents the turbulent viscous dissipation due to the velocity
273 fluctuations, MacComb (1990):

=> =>8
274 ,
6 7 ´,
6 8 ´ = :; <=? 7 + =? @ (8)
8 7

PT
275 where :; represents an eddy viscosity, according to the Boussinesq approximation. We

RI
276 adopted the renormalized group theory model for resolving turbulence; this RNG model

SC
277 had two transport equations: the turbulent kinetic energy (k) and its dissipation (ε):

_ EF
-* · A C = * · [ 3 + *C] + pk − (9)
GH
278

_
UE
-* · A L = * · [ 3 + F *L] +
N
PN pk − PN# P 
AN
(10)
GM O
279

280 where PN was calculated through the following expressions:.

M

281 PN = 1.42 − S (11)

D

T 
S = U.VW
V (12)
X 
282
TE

YZ
=
[\ V
283 (13)
EP

284 and PN# , ]N , ]O were constants while pk represented the production of energy in turbulent
regime due to viscous force. 30 was expressed as:
C

285
AC

O
30 = PE -  (14)
N
286

287 where PE was a semi-empirical constant of value 0.009 (ANSYS, 2011).

13
 ACCEPTED MANUSCRIPT

288 During the roasting of coffee beans, heat was transferred mainly by forced convection from

289 air to the bean surface and by conduction from the surface to the bean interior (Fabbri et al.

290 2011). Numerical simulations assumed that:

291 (a) initial temperature of the bean was uniform (293K);

PT
292 (b) sensible heat due to moisture loss was considered;

RI
293 (c) bean coffee volume as constant;

294 (d) void fraction porosity as function of time, t;

SC
295 (e) coffee beans roasted under forced convection condition.

296 Bases for assumptions (a) to (e) were: a) roasting process does start with beans at ambient

297
U
temperature, getting hotter as time runs, according to real situation. Although forced air at
AN
298 high temperature feeds the roaster, the beans were load initially at ambient temperature; b)

299 refers to the capability of the simulations to taking into account for moisture changes into
M

300 the beans; c) this assumption represents an ideal condition of coffee beans because the real
D

301 situation includes beans volume change as a function of time and temperature. The
TE

302 modeling of coffee roasting is evolving for including the bean´s volume expansion during

303 roasting, leading to real industrial applications, but it is not yet available; d) Dynamic
EP

304 porosity (see equation (5)) was introduced in modeling coffee beans roasting by the first

305 time. According to Table 1, porosity value of 9.8% corresponds to un-roasted beans, and
C

306 34.24% for roasted beans, while intermediate data were taken from literature, see Fig. 4e; e)
AC

307 this assumption follows from previous works, which allows for comparing the results

308 against data from literature.

309 Heat transfer into the beans from hot air allowed phase change latent heat to produce

310 moisture from former humidity contents of green beans set as boundary condition to 12%.

14
 ACCEPTED MANUSCRIPT

311 This model for porous media included coffee bean porosity, by introducing a momentum

312 loss  into equation (7), due to superficial phase velocity in the walls of pores, as part of

313 volumetric flow rate porous region:

^ 
 = − . ,7 + `a ρ|,|,7 1 (15)

PT
_ #
314

315 Where the right hand side included the viscous and inertial loss terms. α represented

RI
316 permeability and C2 was an inertial resistance factor.

In modeling the moisture loss to the airflow, a convection-diffusion equation for the local

SC
317

318 mass fraction of species ci was set as (Rogers and Mayhew (1996)):

U
/γρc EF i
+ * ∙ γρ,e  = −* ∙ f.−γρg, + 1 *e − g, j (16)
/0 F 
319
AN
320 Diffusion was due to concentration and temperature gradient, expressed in the first and

321 third terms on the right hand side of equation (16), respectively. This way, moisture loss
M

322 was affected by porosity γ, which is multiplied by the mass diffusion coefficient Di,m, for
D

323 species i in turbulent regime. Sct was the turbulent Schmidt number set 0.7, and DT,i was
TE

324 the thermal diffusion coefficient.

325 Similarly, the energy equation for a porous media read as:
EP

/
γρk lk + 1 − γρm lm  + ∇ ∙ , ρk lk + n =
/0
326

327 ∇ ∙ 5−Cokk ∇p − ∑rOs ℎ t + τu ∙ ,9 + kv (17)

C

328 where h and Ef were enthalpy and total fluid energy, as well as Es was a total energy for the
AC

329 solid medium, and Cokk was the thermal conductivity for solid-fluid porous media:

330 Cokk = γCk + 1 − γCm (18)

331 where Ck accounted for the fluid phase thermal conductivity (including the turbulent

332 contribution C0 ); and Cm for the solid medium thermal conductivity, which becomes

15
 ACCEPTED MANUSCRIPT

333 affected by the porosity; in equation 17 kv represented a fluid enthalpy source term and

334 t was the diffusion flux of species i (moisture) given as:

E i
335 t = f.−γρg, + F 1 *e − g, 
j (19)
F

PT
336 For forced convective flow around the coffee beans and for conjugate solid-liquid heat

337 exchange, Patankar (1980) gave a general form of the energy equation as:

RI
/
-P + ∇ ∙ , -l + n + ∇ ∙ 5−Cokk ∇p + ∇ ∙ ,9 = 0 (20)
/0
338

SC
339 Where Cokk was an effective thermal conductivity, Cokk = C + C0 , C0 was a turbulent

340 thermal conductivity given by the model for turbulence, and E represented energy through

341 enthalpy h, as:

U
AN
€ ,
l =ℎ− + (21)
 #
342
M

343

344
D

345 2.3.1 Boundary conditions for modeling

TE

346 The computational domain considered the geometry of three zones; one was the fluid

347 represented by hot air that flows through the tube roaster and surrounds the coffee beans. A
EP

348 second zone, defined by coffee beans material, was initially declared as solid porous

349 material with porosity from equation (5), and temperature of 293 K. Finally, the third zone
C

350 represented a region around the beans wall, which was declared as shell. This shell linked
AC

351 the solid porous conductive material of beans to the convective fluid surrounding. Once the

352 computational domain was defined (see Fig. 6a), mass, momentum and energy transport

353 inside the coffee beans were invoked with the model for turbulence in the momentum

354 equations, the RNG-k-ε. The equation for species transport was activated to modeling

16
 ACCEPTED MANUSCRIPT

355 moisture loss as well. Properties and conditions given in tables 3 and 4 feed the program by

356 declaring three materials: coffee bean as solid porous with dynamic porosity as function of

357 time; hot air that roasted the beans, and the humidity that will get out from the beans

358 interior during roasting. Roaster inlet and outlet boundary conditions were set out as given

PT
359 in Table 3, while all walls declared as adiabatic. The patch option was set up to link coffee

RI
360 beans external surface to air stream. A schematic of boundary conditions was given in Fig.

361 6d.

SC
362
363 2.3.2 Model geometry

364
U
Coffee bean´s geometry simulated the shape and size of Arabica according to Fabbri et al.
AN
365 (2011) and Hernandez-Díaz et al. (2008). One coffee bean represented by a semi-elliptical

366 curve with longitudinal diameter of 0.0098 m, equatorial diameter of 0.007 m and depth of
M

367 0.0037 m was built and meshed in the program Gambit. The mesh of one bean was
D

368 reproduced for each of three beans placed inside the roasting column as shown in Fig. 6a,
TE

369 which presents one mid-plane that cuts axially the roasting column with three beans. The

370 mesh of the roaster was a 3-D domain of 781,591 cells. One bean´s computational mesh
EP

371 was represented by 22,508 cells, shown in Fig. 6b. The beans were marked G1, G2 and G3

372 as shown in Fig. 6c for further discussion. A mesh refinement process as shown in Table 5
C

373 demonstrated independence of results, where data of each mesh allowed calculating a grid
AC

374 convergence index, GCI, according to Roache (1998). Results of temperature for one

375 position in the reactor: (x = 0.0, y/y = 0.7, z = 0.0) were used in the calculation of GCI,

376 which decreased as refinement increased, improving mesh independence. For the present

377 work grid 5 was judged appropriate.

17
 ACCEPTED MANUSCRIPT

378
379 3 Results

380 A comparison of predictions against experimental data was made through profiles of

381 temperature and moisture loss within the coffee beans, as function of time as follows. This

PT
382 comparison provided a basis for validation of the numerical model used in subsequent

383 simulations including void fraction porosity (equation (5)). Afterwards, flow dynamics

RI
384 around the coffee beans and temperature distribution in the beans continued the discussion.

SC
385 Finally, Nusselt number as function of Reynolds number of roasting flow illustrated the

386 heat transfer around the coffee beans.

U
387 Figure 7 presents the results of present experimental data and predicted results from present
AN
388 models without porosity, and with dynamic porosity. Results in this figure showed the

389 temperature variation inside the coffee bean for the case of air roasting temperature equal to
M

390 473 K from start up to 400 s. As observed, the trend of all curves was similar: a good
D

391 agreement of the time at which saturation roasting temperature was achieved, which was
TE

392 approximately equal to 250 s. Small differences of heating time of bean interior, which

393 were faster for present model without porosity than dynamic porosity model and
EP

394 experimental data, after 50 s from start. From model with dynamic porosity, we observed

395 that the trend of simulation is closer to the experimental data.

C

396 A calculation of error from predicted results against the experimental data was tabulated in
AC

397 Table 6, in order to give a measurement of the performance for each numerical model. The

398 heating time was in the first column divided into periods of 25 s. For each period, a

399 calculated error completed the table, following columns from model by Fabri et al. (2011),

18
 ACCEPTED MANUSCRIPT

400 to models without and with dynamic porosity. For comparison, a set of results for the

401 numerical predictions of Fabri et al. (2011) was include in the second column.

402 The results indicated good performance of present model without considering porosity

403 compared with Fabri et al (2011), with average error of 2.2 against 5.8. Our predictions are

PT
404 closer to experimental data than results from Fabri et al. (2011) due to differences of

RI
405 boundary conditions and roaster design used by them. In addition, we observe good

406 performance of our modeling using dynamic porosity compared against modeling without

SC
407 porosity. This was because in the first seconds of roasting the real process is like pure solid

408 material. Even porosity is almost null, the dynamic model (equation (5) allows the

409
U
computation getting a fine adjustment as time advances. Thereafter, using a model without
AN
410 porosity worst from 50 to 100 s, while results with dynamic porosity improves. An

411 interpretation of difference from a time of roasting of 50 s and so on until completing the
M

412 process over 200 s is that porosity is not the unique property of coffee beans that changes
D

413 with time, as discussed above.

TE

414 Returning back to the model without porosity, from a roasting time of 100 s its

415 performance remains good until ending the roasting. This is because the temperature
EP

416 gradient between the bean interior and exterior grows drastically, which favors heat transfer

417 rate. Furthermore, the results in Table 6 indicate that modeling with dynamic porosity leads
C

418 to good performance all along the total roasting time of 400 s, except for the period from 50
AC

419 to 120 s. We attribute this disagreement of dynamic model with experimental data to the

420 importance of releasing products other than moisture, which is considered in this study.

421 Instead, predictions without porosity worst compared against the model with dynamic

422 porosity from 50 s of roasting and after until a time of 120 s. As time advances, the

19
 ACCEPTED MANUSCRIPT

423 predicted results indicate that without taking into account the porosity the heating was

424 much faster than taking porosity into account, confirming again the role of thermal

425 conductivity of air pores, and the temperature gradient between bean´s interior and exterior.

426 Porosity represents a delay of heating, because heat conductivity of air in pores was much

PT
427 smaller that coffee bean solid material. This result was entirely attributed to the presence of

RI
428 pores, which represent natural obstacles for heat conduction because thermal conductivity

429 of air contained in the pores is reduced compared against coffee´s material alone (see Table

SC
430 4 for values of thermal conductivity between air and solid coffee bean). The present

431 modeling reproduces the experimental data nicely from a time of roasting of 120 s on until

432
U
saturation as can be observed in Fig. 7. As conclusion, modeling with porosity affects the
AN
433 predictions of heating in coffee beans roasting.

434
M

435 4 Discussion
D

436 4.1 Group of beans temperature

TE

437 Bean temperatures of a group instead of a coffee bean alone were analyzed. First, it is

438 worth to have a glance on the surface condition of the bean compared to its interior. Figure
EP

439 8a shows the temperature profiles plotted against the time, obtained for each of three beans

440 classified as G1, G2 and G3. As observed, there was an average time retarding = 10 s
C

441 between the interior and the exterior of the beans, which represents how long it takes to
AC

442 heating the bean´s interior from its surface. Arrowed insert served to highlight the results in

443 a range from 30 s to 60 s allows to observe the difference of temperature among the beans

444 during this short period. The insert revels that heating inside and outside the beans is

445 characterized by a time delay of 10 s, as mentioned. This varies little from the beginning up

20
 ACCEPTED MANUSCRIPT

446 to reaching the saturation temperature, except for the last fifty seconds where the curves

447 present a strong gradient before the asymptotic behavior was achieved. Although there was

448 not a constant slope, the trend of these curves was very smooth during the first 50 s. In

449 addition, a temperature difference between bean´s interior and surface was of 3 K,

PT
450 attributed to the influence by position with respect to the wall of the cylinder. Furthermore,

RI
451 the difference from one grain to another one was about 1 degree. This changes if we

452 observe G1 and G2 for a time of 70 s for instance: the temperature among them was 1

SC
453 degree, while less than 1 degree existed between G2 and G3. This change was attributed to

454 a slightly uneven distribution of streamlines as will be discussed below. In summary, heat

455
U
transfer for a centrally placed bean G1 gets faster compared to the other two grains over the
AN
456 first 80 s.

457 Figure 8b shows the difference of heating process between beans interior and exterior, for
M

458 constant and dynamic porosity. In this case, two arrowed inserts were used to emphasize
D

459 the difference of heating. It is observed that beans interior reach a temperature of 430 K in
TE

460 20 s, while modeling with dynamic porosity takes longer, 90 s. Therefore, by considering

461 porosity leads to delaying the bean´s heating up to 70 s. However, an insert in the range
EP

462 from 50 s to 100 s lets one to observe that small temperature difference existed among the

463 grains surface and interior, which are due to uneven flow distribution and position to the
C

464 roaster wall, as discussed above.

AC

465 Next we examine the temperature distribution inside three coffee beans (through the mid-

466 plane of each bean) for several time intervals of roasting. A contour time sequence of

467 temperature distribution reveals the characteristics of the heating process. Figure 9a, shows

468 the temperature contours after five seconds from start of roasting process. At this time there

21
 ACCEPTED MANUSCRIPT

469 was no influence of the hot air stream on the interior condition of the beans. Therefore, the

470 time delay of 10 s between bean inner and outer heating deduced from temperature profiles

471 in Fig. 8 is confirmed.

472 After ten seconds from start, the beans present some heating in the surface and little from

PT
473 the extremes towards the center, as observed in Fig. 9b. A relevant effect of hot air takes

RI
474 place after 20 s, which was indicated by a strong temperature gradient of 100 degree

475 between the bean center and the surface as observed in Fig. 9c.

SC
476 In Fig. 9d the temperatures in the center of beans raised 150-160 degree above their initial

477 value, while the surface is at least 200 K above the initial value. The process ends after

478
U
eighty seconds of roasting; at this point the temperature on the center equals the
AN
479 temperature of surface of the beans, as shown in Fig. 9e. This description confirms the

480 results reported by Schwartzberg (2002), where the coffee bean temperature increases
M

481 rapidly due to the speed with which the bean loses water, which will be discussed ahead in
D

482 this work.

TE

483
484 4.2 Velocity contour, Nusselt number and moisture loss
EP

485 The shape of grains, their size and the distance from each other affected the flow of air

486 around the beans. In combined effects, the flow condition in the proximity of the bean
C

487 surface, in turn, affects the heat transfer towards-from the bean and the moisture loss from
AC

488 its interior. Therefore, it was important having a glance over the flow distribution. The

489 results presented in the form of velocity contours, are shown in Fig. 10a. This figure shows

490 an arrow indicating the axial direction of air flow and the distribution of coffee grains, in a

491 mid-plane intersecting the body of the coffee beans and roaster. As observed, hot air

22
 ACCEPTED MANUSCRIPT

492 distributes around beans, with streams of higher velocity located in between one coffee

493 bean and another and with respect to the roaster wall; because a reduction of space through

494 which air must pass through. A slight uneven flow distribution was observed between G2

495 and G3, with faster fluid flow around G2, which is responsible of small temperature

PT
496 difference among these grains as roasting advances. Flow variations around the beans were

RI
497 used for calculating heat convection-conduction rate through the Nusselt number, Nu, a

498 measure of convection heat transfer at the surface. Results of Nu were plotted against the

SC
499 Reynolds number, Re, which characterizes the airflow regime. Both, Nu and Re were

500 calculated following Seghir-Ouali et al. (2006). Re number varied for describing heat

501
U
transfer in presence of the beans as observed in Fig.10b. The plot indicates that heat
AN
502 transfer increases as Re number increases for three positions in the following order: near

503 the roaster wall; near the bean wall; and in between two coffee beans. This was in
M

504 agreement with data published in the literature by Seghir-Ouali et al. (2006) who addressed
D

505 the case of heat transfer in a stationary cylinder as function of axial flow of air. The present
TE

506 results indicate that Nu number was higher for the space from one bean to another, which

507 was observed in Fig. 10a, due to the increment of velocity. It suggests that a roaster
EP

508 roasting pimply fluid velocity was constantly increasing due to the proximity between

509 beans having an influence on quality of coffee roasting beans. Contours of Fig. 10a also
C

510 indicate the formation of vortices and their location, therefore in a rotating cylinder drum
AC

511 roaster vortex are also part of the process. A consequent uneven distribution of flow

512 velocity around the beans has an influence on the roasting process. The results in Fig. 10b

513 indicated, however, that the influence of axial velocity on heat transfer tends to be

514 negligible as Re increases indefinitely.

23
 ACCEPTED MANUSCRIPT

515 Roasting in a porous media is a process thermally governed. As part of this process,

516 moisture loss, treated as diffusion of species has molecular effects. Therefore, the terms that

517 represent convective-conductive heat transfer mechanisms in equations (17) couples with

518 diffusive terms of equation (16), where thermal terms dominate because they are faster than

PT
519 diffusive terms. Although diffusion coefficient increases with increment in porosity, the

RI
520 heat transfer still dominates. Why moisture slows when porosity increases? It is because

521 heat transfer rate varies according to properties involved: As observed in Table 4, thermal

SC
522 conductivity in air is one order smaller than thermal conductivity in coffee solid material.

523 Therefore, by changing porosity we reduce heat transfer rate because it is function of

524
U
thermal conductivity: convective-conductive heat transfer rate in air-solid takes longer
AN
525 compared to pure conductive heat transfer in solid coffee material.

526 Contours of bean moisture contents in the central grain G1 corresponding to different times
M

527 of roasting form a sequence shown in Fig. 10c. The contours show details of resident
D

528 contents of water before it evaporates as moisture, for condition without porosity. We
TE

529 plotted it as a continuous process modeling with and without porosity, in Fig. 10d. The plot

530 described the moisture loss as function of time, which represented the sensitivity of
EP

531 moisture loss with change in porosity. The introduction of this property affects drastically

532 by delaying the bean moisture loss up to 160 s. After this time, the effect of dynamic
C

533 porosity makes the moisture loss reducing to near zero. A possible explanation is based on
AC

534 theory as stated by Cengel (2007): the primary driving mechanism of mass diffusion is the

535 concentration gradient, and mass diffusion due to a concentration gradient is known as the

536 ordinary diffusion. However, diffusion may also be caused by temperature gradients in a

537 medium causing thermal diffusion known as Soret effect, which is negligible, unless the

24
 ACCEPTED MANUSCRIPT

538 gradients are very large. This is the case of bean´s heating conditions observed in Figs. 9a,

539 9b and 9c. In principle, we can´t neglect Soret effect. However, porosity affects diffusion,

540 because when the pores of a porous solid are smaller than the mean free path of the gas

541 molecules, the molecular collisions may be negligible and a free molecule flow may be

PT
542 initiated. This is known as Knudsen diffusion. See that porous of green beans are smaller

RI
543 than 10 µm, Fig. 4a and 4c. Therefore, diffusion in green beans may be a combination of

544 Knudsen diffusion-Soret effects, raising the loss of moisture as observed in the dotted line

SC
545 for modeling without porosity in Fig. 10d. Alternatively, when the size of the gas molecules

546 is comparable to the pore size in solids, adsorbed molecules move along the pore walls in a

547
U
process known as surface diffusion, a retarded process, see the curve line for dynamic
AN
548 porosity in Fig. 10d.

549
M

550 5 Conclusions
D

551 Heat and mass transfer in coffee beans during the roasting process were investigated under
TE

552 the effects of void fraction bean porosity. This important coffee bean property was part of

553 modeling coffee roasting process. Electronic microscopy was used to characterizing
EP

554 unroasted and roasted beans in terms of pore size and distribution. As a result, void fraction

555 porosity was calculated for green unroasted and roasted coffee beans, respectively. The
C

556 results were complemented by data from literature to building a dynamic model of porosity,
AC

557 which was incorporated into a numerical 3-D model for simulating roasting of single and

558 group of 3 coffee beans. The model computes heat and mass transfer rates between the

559 surrounding stream and coffee beans during the roasting process. Roasting temperatures of

560 473 K and 503 K were studied for roasting times up to 400 s. The predicted results were

25
 ACCEPTED MANUSCRIPT

561 validated using experimental data of temperature from Arabica coffee beans roasting in a

562 fluidized bed roaster. We observed that constant porosity modified the temperature

563 distribution in the beans interior during the heating time, compared with a no porosity

564 prediction. In addition, porosity as function of time affected moisture loss by delaying its

PT
565 contents in coffee beans, and decreasing its value to zero along the time of roasting. As

RI
566 conclusion, modeling with porosity affects the predictions of heating and species transport

567 in coffee beans roasting. This procedure can be useful for understanding the roasting

SC
568 process, or studying processes of food products that require roasting by means of hot air,

569 like grains or seeds.

U
AN
570 Acknowledgments

571 The authors thank the Consejo Nacional de Ciencia y Tecnología (CONACyT) for the
M

572 financial support received for the development of this research. First author acknowledges

573 the scholarship provided by CONACyT during her PhD studies.

D

574
TE

575
C EP
AC

26
 ACCEPTED MANUSCRIPT

576 References

577 ASTM D4641 – 12 Standard Practice for Calculation of Pore Size Distributions of

578 Catalysts and Catalyst Carriers from Nitrogen Desorption Isotherms- Committee D32 on

579 Catalysts (1996).

PT
580 ASTM E. 112 Standard test methods for determining average grain size (1996).

RI
581 Alonso-Torres B., Hernández-Pérez J.A., Sierra-Espinosa F., Schenker S., and Yeretzian

582 C., (2013). Modeling and Validation of Heat and Mass Transfer in Individual Coffee Beans

SC
583 during the Coffee Roasting Process Using Computational Fluid Dynamics (CFD). CHIMIA

584 International Journal for Chemistry, 67 (4), 291-294.

585
U
ANSYS, FLUENT User's Guide., ANSYS, Inc., 2011.
AN
586 Baxes G. A. (1994). Digital Image Processing. p. 157, J. Wiley & Sons Ed. USA,

587 Baggenstoss, J., Poisson, L., Luethi, R., Perren, R., Escher, F., (2007). Influence of water
M

588 quench cooling on degassing and aroma stability of roasted coffee. Journal of Agricultural
D

589 and Food Chemistry, 55, 6685–6691.

TE

590 Bin X., Da-wen S., (2002). Applications of computational fluid dynamics (CFD) in the

591 food industry: a review. Computers and Electronics in Agriculture, 34, 5–24.
EP

592 Bottazzi D., Eggers, R., Pietsch, A., (2001). Roasting. In: Clarke, R.J., Vitzthum, O.G.

593 (Eds.), Coffee: Recent Developments. Blackwell Science, London, 90–107.

C

594 Cengel Y.A., (2007). Heat and mass transfer: A practical approach. 3rd Ed. Mc Graw Hill,
AC

595 NY, pp. 773-816.

596 Fabbri A., (2005). Applications of numerical simulation methods in agricultural

597 engineering. In Proceedings of AIIA Meeting, Catania, 273–300.

27
 ACCEPTED MANUSCRIPT

598 Fabbri A., Cevoli C., Alessandrini L., Romani S., (2011). Numerical modeling of heat and

599 mass transfer during coffee roasting process. Journal of Food Engineering, 105, 264-269.

600 Farina S., Milani M., Montorsi L., (2012). A numerical approach for the analysis of the

601 coffee roasting process. Journal of Food Engineering, 112, 243–252.

PT
602 Frisullo P., Barnabà M., Navarini L., Del Nobile M.A., (2012). Coffea arabica beans

RI
603 microstructural changes induced by roasting: An X-ray microtomographic investigation,

604 Journal of Food Engineering, 108, 232-237.

SC
605 Hernández-Díaz W.N., I.I. Ruiz-López, M.A. Salgado-Cervantes, G.C. Rodríguez-Jiménez,

U
606 M.A. García-Alvarado, (2008). Modeling heat and mass transfer during drying of green
AN
607 coffee beans using prolate spheroidal geometry. Journal of Food Engineering, 86, 1-9.

608 Hernández-Pérez J.A., (2002). Étude de la torréfaction: modélisation et détermination de

M

609 degré de torréfaction du café en temps réel. PhD. Thesis. École Nationale Supérieure des

610 Industries Agricoles et Alimentaires, Paris.

D

611 Hernández J.A., Heyd B., Irles C., Valdovinos B. and Trystram G., (2007). Analysis of the
TE

612 heat and mass transfer during coffee batch roasting. Journal of Food Engineering 78 (4),

613 1141-1148.
EP

614 Hernández J.A., Heyd B., and Trystram G., (2008a). Prediction of brightness and surface
C

615 area kinetics during coffee roasting. Journal of Food Engineering, 89 (2), 156-163.
AC

616 Hernández J.A., Heyd B., and Trystram G., (2008b). On-line assessment of brightness and

617 surface kinetics during coffee roasting. Journal of Food Engineering, 87 (3), 314-322.

618 Heyd B., Broyart B., Hernández J.A., Valdovinos-Tijerino B. and Trystram G., (2007).

619 Physical model of heat and mass transfer in a spouted bed coffee roaster. Drying

620 technology, 25 (7-8), 1243-1248.

28
 ACCEPTED MANUSCRIPT

621 Juang C. H., Holtz, R. D. (1986). A probabilistic permeability model and the pore size

622 density function. Int. J. for Num and Analytical Methods in Geomechanics. 10, 543- 553.

623 McComb W.D. (1990). The physics of fluid flow, Oxford Univ. Press.

PT
624 Moellmer J., Celer E.B., Luebke R., Cairns A.J., Staudt R., Eddaoudi M., Thommes M.,

625 (2010). Insights on Adsorption Characterization of Metal-Organic Frameworks: A

RI
626 Benchmark Study on the Novel soc-MOF, Microporous Mesoporous Mater. 129 (3), 345–

SC
627 353.

628 Ortolá M.D., Londoño L., Gutiérrez C.L.,Chiralta., (1998). Influence of roasting

629
U
temperature on physicochemical properties of different coffees. Food Science and
AN
630 Technology International, 4, 59-66.

631 Patankar S.V. (1980). Numerical Heat Transfer and Fluid Flow, Taylor and Francis, New
M

632 York.
D

633 Pliego-Solórzano R., (2010). Heat and mass transfer of coffee roasting by fluidized bed.
TE

634 MSc Thesis, Universidad Autónoma del Estado de Morelos (CIICAp), UAEM, Morelos.

635 Raemy A., Lambelet P.(1982). A calorimetric study of self-heating in coffee and chicory,
EP

636 Journal of Food Technology,17, 451-460.

637 Roache P.J. (1998). Fundamentals of Computational Fluid Dynamics, Hermosa, New
C

638 Mexico.
AC

639 Rogers G. F. C. and Mayhew Y. R. (1996). Engineering Thermodynamics: Work and Heat

640 Transfer. Ed. Longman, Essex.

29
 ACCEPTED MANUSCRIPT

641 Santos J.R., Lopo M., Rangel A.O.S.S. & Lopes J.A., Exploiting near infrared spectroscopy

642 as an analytical tool for on-line monitoring of acidity during coffee roasting, Food Control

643 (2015), doi: 10.1016/j.foodcont.2015.08.007.

644 Scheker S., Heinemann C., Huber M., Pompizzi R., Perren R., Escher F., (2002). Impact or

PT
645 Roasting Conditions on the Formation of Aroma Compounds in Coffee Beans. Food

RI
646 Engineering and Physical Properties, 67, 60-66.

647 Schwartzberg H.G., (2002). Engineering and Food the 21st Century. Chapter Modeling

SC
648 bean heating during batch roasting of coffee beans. CRC Press LLC. Boca Raton, Florida.

649 Scott G., Richardson P., (1997). The application of computational fluid dynamics in the

650
U
food industry. Trends Food Science Technology 8, 119–124.
AN
651 Schenker S., Handschin S., Frey B., Perren R. and Escher F., (2000). Pore Structure of
M

652 Coffee Beans Affected by Roasting Conditions, J Food Science. 65 (3), 452-457

653 Seghir-Ouali, S., Saury, D., Harmand, S., Phillipart, O., Laloy, D., (2006). Convective heat
D

654 transfer inside a rotating cylinder with an axial air flow, Int. J. of Thermal Sc., 45, 1166-
TE

655 1178.

656 Severini, C., Nicoli, M.C., Mastrocola, D., Lerici, C.R., (1992). Influence of heating rate on
EP

657 some physical physico-chemical properties of coffee beans. In Proceedings of 14th

C

658 International Science Colloquium on Coffee. ASIC, San Francisco, 641–648.

AC

659 Thommes M. (2004). Physical Adsorption Characterization of Ordered and Amorphous

660 Mesoporous Materials, In: Nanoporous Materials, Science & Engineering (Lu G.Q. & Zhao

661 X.S., Eds), Chapter 11, Imperial College Press, London.

30
 ACCEPTED MANUSCRIPT

662 Valdovinos-Tijerino B.E., (2005). Étude de la torréfaction du café: Modélisation et

663 développement des outils pour maitriser la qualité du produit en ligne. PhD. Thesis.

664 Universite Pierre and Marie Curie, Paris.

665 Versteeg H. K. and Malalasekera W. (1995). An Introduction to Computational Fluid

PT
666 Dynamics, The Finite Volume Method , Ed. Longman, Essex.

RI
667

668

U SC
AN
M
D
TE
C EP
AC

31
 ACCEPTED MANUSCRIPT

669 Appendix
670 We present the UDF which was introduced into the program Fluent (ANSYS, 2011), in

671 order to simulate dynamic change of porosity γ (t), during roasting of coffee beans:

672 #include"udf.h"

PT
673 #define a1 0.0002
674 #define a2 0.0227
675 #define a3 9.4573
676

RI
677 DEFINE_PROFILE(porosity,thread,index)
678
679 {

SC
680 cell_t c;
681
682 begin_c_loop(c,thread)
683 {
684 real t=RP_Get_Real("flow-time");

U
685 C_PROFILE(c,thread,index)=(a1*t*t)+(a2*t)+a3;
686 }
AN
687 end_c_loop(c,thread)
688 }
M
D
TE
C EP
AC

32
 ACCEPTED MANUSCRIPT

List of Figures
Fig.1 Photograph of un-roasted and roasted conventional coffee beans.

Fig. 2 Samples of coffee beans for surface and interior observation of material structure; a)
Plan view of beans on the support plate; b) Lateral view of beans on the support plate inside
the evaporator machine.

PT
Fig. 3 Aspect of the surface of coffee beans under microscope, scales of 1 mm and 20 µm;
a) Flat side of unroasted bean; b) Flat side of roasted bean; c) Convex side of unroasted

RI
bean; d) Convex side of roasted bean; e) Micrograph of convex surface unroasted bean; f)
Micrograph of convex surface roasted bean.

SC
Fig. 4 Micrographs of interior of coffee beans; a) Un-roasted coffee, scale 100 µm; b)
Roasted coffee, scale 100 µm; c) Un-roasted coffee, scale 20 µm; b) Roasted coffee, scale
20 µm; e) Porosity as function of time for coffee bean roasting.

U
Fig. 5 Schematic of experimental column used for roasting the coffee beans; a) Cylindrical
column made of Pyrex with main dimensions; b) detail of thermal couple insertion into the
AN
coffee bean; c) Overview of source of hot air, column with bean inside and thermal couples
connected with interface and computer; d) Photograph of mounted apparatus.

Fig. 6 3-D computational grid of the roasting column and three coffee beans; a ) Mid-plane
M

view of coffee beans and roasting column; b) Coffee bean´s grid by flat face view; c) Three
coffee beans group, numbered for further discussion; d) Schematic of boundary conditions
used in computational modelling.
D

Fig. 7 Validation of dynamic model of porosity by direct comparison of present numerical

TE

results against experimental data, for air temperature = 473 K.

Fig. 8 Time history of temperature during the column roasting process by profiles in the
center and surface of the coffee beans; a) Without porosity; b) Sensitivity of temperature
EP

with porosity for coffee bean´s interior.

Fig. 9 Sequence of temperature contours inside the coffee beans along the roasting process;
C

a) 5 s; b) 10 s; c) 20 s; d) 50 s; e) 80 s from start.
AC

Fig. 10 Flow dynamics, heat transfer and moisture loss in the coffee beans during roasting;
a) Flow distribution around the coffee grains; b) Nusselt number as function of axial
Reynolds number; c) Time sequence of moisture loss in the central bean; d) Sensitivity of
moisture loss as function of time without porosity and with dynamic function of porosity.
 ACCEPTED MANUSCRIPT

PT
Table 1 Coffee bean porosity before and after roasting.
Grain Sample Scanned Pore average Void Void fraction

RI
condition size área size surface porosity

̅   ̅
 γ

SC
2
mm µm2 µm2 µm2 %

Un-roasted 1.10 4.97x104 167.20 4.51x103 9.88

U
Roasted 1.10 5.53x105 AN 1306.23 1.89x105 34.24
M
D
TE
C EP
AC
 ACCEPTED MANUSCRIPT

Table 2 Conditions for experimental coffee roasting.

Test Temperature Air flow velocity Comments
K m s-1

PT
1-3 473 2.8

RI
4-6 483 2.8 total amount of
green coffee
7-9 493 2.8 beans was fixed
to 5.1 g ± 0.1 g

SC
10-12 503 2.8
ambient
temperature was
13-15 513 2.8
299 ± 1.7 K
16-18 523 2.8

U
AN
M
D
TE
EP
C
AC
 ACCEPTED MANUSCRIPT

Table 3 Boundary conditions for numerical model

Zone Inlet Outlet Wall Material

velocity pressure

PT
m s-1 Pa

Roaster 0.02 0 no slip fluid

RI
Beans static no slip solid shell conduction

Shell no slip solid 0.0001 m thick

U SC
AN
M
D
TE
C EP
AC
 ACCEPTED MANUSCRIPT

Table 4 Properties of air and coffee beans.

PT
Density Thermal Void fraction Humidity Mass- Specific heat Molecular viscosity
conductivity porosity contents diffusivity
kg m-3 W m-1K-1 % % m2 s-1 J kg-1 K-1 kg m-1 s-1

RI
Air
1.225 0.0242 ----- ----- ----- 1006.43 1.7894x10-5

Coffee
1723 0.131 γ (t) 12 2.88x10-9 1256 -----

U SC
AN
M
D
TE
C EP
AC
 ACCEPTED MANUSCRIPT

Table 5 Grid independence test data.

T(x, y/y=0.7, z, t) Maximum Grid

Grid Number Size ratio

PT
volume of temperature GCI
number of cells r
cell (m3) relation, ε

RI
1 362593 0.457 4.37E-10 - - -
2 436790 0.507 3.12E-10 1.40E+00 -0.09861933 7.43E-01

SC
3 461636 0.527 3.40E-10 9.18E-01 -0.03795066 -1.40E+00
4 850288 0.597 2.96E-10 1.15E+00 -0.11725293 2.39E+00
5 931877 0.549 2.05E-10 1.45E+00 0.087431694 5.84E-01

U
AN
M
D
TE
C EP
AC
 ACCEPTED MANUSCRIPT

Table 6 Calculated error from results of different models

against experimental data of bean´s heating.
Roasting time Numerical Numerical
Numerical
period without with dynamic
Fabbri et al.
s porosity porosity

PT
(2011)
0-25 5.99 2.56 0.95
25-50 10.30 3.90 2.01

RI
50-75 9.14 4.12 2.21
75-100 7.15 2.86 1.79
100-125 5.30 1.89 1.38

SC
125-150 3.88 1.09 0.89
150-175 2.80 0.59 0.52
175-200 2.04 0.21 0.20

U
AN
M
D
TE
EP
C
AC
 ACCEPTED MANUSCRIPT

PT
RI
U SC
Fig. 1
AN
M
D
TE
EP
C
AC
 ACCEPTED MANUSCRIPT

PT
RI
U SC
AN
a b
M
D

Fig. 2
TE
EP
C
AC
 ACCEPTED MANUSCRIPT

PT
RI
a b

U SC
AN
M
D

c d
TE
EP
C
AC

e f

Fig. 3
 ACCEPTED MANUSCRIPT

PT
RI
a b

U SC
AN
M
D

c d
TE
EP
C
AC

Fig. 4
 ACCEPTED MANUSCRIPT

PT
RI
SC
cavity size

U
AN
insulation

thermocouple
M

a b
D
TE

coffee
thermocouple
grain
EP

source of hot air

C
AC

c d

Fig. 5
 ACCEPTED MANUSCRIPT

PT
a b

RI
U SC
c
AN
M
D

Pressure outlet
(P=0 Pa)
TE
EP

Adiabatic walls
-1
V=0 ms
C

Coffee bean interior

Ti(x, y, z, t=0 s) = 293 K;
AC

Hi = 11%;
γ = 9.8%

z d

x
y Inlet velocity
-1
Vi = 2.8 ms Fig. 6
Ti = 473 K; 503 K
 ACCEPTED MANUSCRIPT

PT
RI
U SC
AN
M
D
TE

Fig. 7
EP
C
AC
 ACCEPTED MANUSCRIPT

PT
RI
U SC
AN a
M
D
TE
EP

b
C

Fig. 8
AC
 ACCEPTED MANUSCRIPT

PT
K

RI
473

455

SC
437

419
a b

U
401
AN
383

365

347
M

c d
329

311
D

293
TE

e
EP

e
C
AC

Fig. 9
 ACCEPTED MANUSCRIPT

-1
ms
0.04

0.0361

0.0321

PT
0.0281

0.0241

RI
002021

0.0161

SC
0.012

0.008

0.0004

U
a
AN
M
D
TE
EP
C
AC

b
 ACCEPTED MANUSCRIPT

H20

0.13
0.12
0.11

PT
0s 10 s 20 s
0.1
0.09
0.08
0.07

RI
0.06
0.05
0.04
50 s 80 s 100 s

SC
0.03
0.02
0.01
0

U
AN
1200 150 s
c
M
D
TE

Numerical dynamic
EP

porosity
C
AC

Fig. 10
 ACCEPTED MANUSCRIPT

Highlights

• Dynamic porosity model used to simulating coffee beans roasting.

• Roasted and green coffee bean porosity determined by electronic microscopy.

PT
Effect of porosity on roasting modeling change temperature and moisture loss.

• Experimental data of coffee beans roasting validated the dynamic numerical model.

RI
• The results can be useful for understanding roasting of other grains.

U SC
AN
M
D
TE
C EP
AC