Modelling particle random walk in a confined environment for inclusion in fluidised bed applications

Powder Technology 221 (2012) 155–163

Contents lists available at SciVerse ScienceDirect

Powder Technology journal homepage: www.elsevier.com/locate/powtec

Modelling particle random walk in a confined environment for inclusion in fluidised bed applications Mike Vanderroost a, Frederik Ronsse a, Koen Dewettinck b, Jan G. Pieters a,⁎ a b

Department of Biosystems Engineering, Ghent University, Coupure links 653, 9000 Ghent, Belgium Laboratory of Food Technology and Engineering, Ghent University, Coupure links 653, 9000 Ghent, Belgium

a r t i c l e

i n f o

Article history: Received 21 September 2011 Received in revised form 28 November 2011 Accepted 25 December 2011 Available online 29 December 2011 Keywords: Random walk Particulate flow Hexagonal close packed lattice

a b s t r a c t In many particulate systems, e.g., a particulate flow or a fluidised bed, particles are subjected to consecutive particle–particle collisions. From a stochastic point of view, these particle–particle collisions give rise to a random particle motion that in most processes is superposed onto a forced particle motion caused by welldescribed external forces that initialise and/or maintain the overall motion of the particulate system (pressure and concentration gradients, gravity, electromagnetic and mechanical forces). Modelling each individual particle–particle collision requires lots of computational time and power. For process control purposes however, models describing many-particle systems such as fluidised beds need to be able to describe the random particle motion in a relatively short amount of time. Therefore, a new stochastic and discrete method was developed that is able to simulate random walk in a confined and relatively dense environment without considering the underlying physics of the phenomena (e.g., particle–particle collisions in a particulate flow) that cause the random walk. The method was successfully tested for relatively dense particulate systems and will be incorporated in a process control model that is currently being developed for fluidised bed coating processes in which particle–particle collisions occur frequently and for which the description of the motion of individual particles is required. © 2011 Elsevier B.V. All rights reserved.

1. Introduction The fluidised bed coating process is widely used in food and pharmaceutical industry for coating of particles. For this type of coating process, however, relatively little research has been performed in developing a quality model, i.e., a model able to predict the final product quality and the process efficiency. In a fluidised bed coating process, three different phases (droplets, particles and air) can be distinguished that continuously interact with each other [1–3]. In the development of a process control model, a detailed description is required of the three phases. A spray sub-model has already been developed that is able to predict the air and droplet temperatures and the air humidity [4]. In this paper, a next step in the development of such a model is presented: a random walk sub-model, able to describe the random walk of particles in a fluidised bed, caused by consecutive particle–particle collisions. This random walk sub-model will be further developed to be built in a particle motion sub-model that also considers forced particle motion caused by air flow. Random walk is the process of an object (a particle, a data packet, a drunk person, etc.) moving from one position to another in some sequential random order. The walk starts at an initial position, and at

⁎ Corresponding author. E-mail address: [email protected] (J.G. Pieters). 0032-5910/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.powtec.2011.12.056

each step the object moves in a randomly chosen direction to another position. Random walk has been widely discussed in literature in a broad spectrum of scientific domains like economics [5,6], biophysics [7], electrical engineering [8], statistical physics [9]. Random walk has also been applied in the domain of chemical engineering. Vitrac and Hayert [10], for example, used random walk to study the diffusion of small and medium-sized molecules in polymeric materials. Caulkin et al. [11] incorporated random walk in their model to study packed bed reactors with cylindrical particles. Nakamura et al. [12] considered random walk to describe particle motion caused by turbulent gas flow and particle–particle collision in a tumbling fluidised bed coater. From a mathematical point of view, a discrete random walk is a Markov chain which describes the trajectory of a particle taking successive random steps. A discrete random walk in an infinite (unbounded) three-dimensional space is characterised by the fact that the displacement Δr of the particle with respect to its initial position is related to the number of steps Nstep [13]: qffiffiffiffiffiffiffiffiffiffi Δr e Nstep

ð1Þ

A frequently used discrete random walk model is that of a random walk on a regular lattice, where at each step a particle jumps from one lattice node to a neighbouring lattice node according to some

156

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163

probability distribution [14]. In a simple random walk, the probabilities of jumping to any of the neighbouring positions are the same. This is typically applied in complex network theory [14,15] and in modelling the shapes of polymer molecules [16]. In dense particulate systems, such as fluidised beds [17,18], random walk is the result of many successive particle–particle collisions. Describing and modelling each individual particle–particle collision [19–21] by detecting colliding particles in the system, is very elaborative and complex. It requires a discrete element method (DEM) to calculate the trajectories of each particle by solving individual equations of motion [22,23]. After each particle–particle collision, the position of each individual particle needs to be recalculated and stored into the computer's memory. The denser or the larger the system, the more particles need to be considered and the more computational power is required. In cases where random motion of particles needs to be simulated without actual interest in the mechanical properties of and interaction between particles, such as process control models for the fluidised bed coating process, it is desirable to have a less complex modelling approach that is able to handle dense particle systems and that can also overcome and significantly reduce the complexity of describing random walk in dense particulate systems. Therefore, in this paper, a new discrete modelling approach was developed that is able to describe random walk in confined dense environments without having to consider the phenomena that cause the random walk and with reduced computational effort. The model presented in this study was developed to be incorporated in a model that is currently being developed for fluid bed processes in which particle–particle collisions occur frequently and for which the description of the motion of individual particles is required. As such, it is important to emphasise that the model in this work should not be considered a competitor for DEM, which has been applied successfully in many fields during the last decades, but as a less complex and quicker alternative. 2. Model Description The new random walk model statistically describes random walk of a particle in a three-dimensional confined and relatively dense environment. The random walk model can be considered as a new statistical approach of a Lattice Boltzmann Method (LBM) defined on a hexagonal close packed lattice. The statistical approach enables the description of random walk without having to consider individual particle–particle collision events which in dense particle system cause the random walk of the particles. The confined environment – from now on referred to as CE – was characterised by its finite volume, fixed boundaries and by the randomly changing positions of the particles within this environment. The random walk was restricted to a framework formed by the hexagonal close packed lattice that was defined within the boundaries of the confined environment. The lattice is thus used to discretise three-dimensional space, as is also the case in Lattice Boltzmann Methods [24]. As already stated, the new random walk model does not consider the underlying physics of the phenomena that cause the random walk, e.g., the particle–particle collisions in a fluidised bed.

[2] Each node in the confined environment can be occupied by at most one particle. [3] The lattice is only defined within the boundaries of the confined environment, meaning that nodes laying outside these boundaries cannot be jumped to. [4] The particles are divided into ‘real’ particles and ‘phantom’ particles. The latter are not considered in simulations of the random walk, but do affect the motion of the real particles as will be explained further on. 2.1.1. Lattice When spheres are regularly packed, a lattice arises by considering the centres of the spheres, the so-called nodes of the lattice. Since the new random walk model should also be able to model random walk in relatively dense particulate systems, the lattice was to be chosen in a way that the lattice packing density, ρlp, defined as the volume fraction of particles when each node in the lattice is occupied with a particle, was sufficiently high. When equally sized spheres are poured in a large container, this will result in a random spherical packing with ρrsp ≈ 0.64. Therefore, a lattice was to be chosen that fulfilled ρlp ≥ ρrsp. It is well known that hexagonal close packing (HCP) is – together with the cubic close packing (CCP) – the most efficient way to stack spheres in a three-dimensional space (Fig. 1). The lattice packing corresponding with the HCPL reaches an average voidage of ρlp ¼ pπffiffiffiffi ≅0:74 [25], 18 which makes it a suitable candidate for being introduced into the random walk model that has to be able to handle dense particulate systems with particle volume fraction of 0.6 or more. In the HCPL, two types of layers (A and B in Fig. 1) are alternately stacked upon each other. Each node in a HCPL has twelve closest neighbours: six in the layer containing the node and three in the adjacent layers above and below, respectively. To denote the position of the nodes in the lattice, a three-dimensional orthogonal lattice coordinate system was defined in the confined environment (Fig. 2) with coordinates (i, j, k). The first layer – by convention an A-layer with k = 0 (grey spheres in Fig. 2) – was defined in such a way that the nodes with the same j-value, are connected by a horizontal zigzag line. In the same layer, nodes having the same i-value are connected by a straight line which is perpendicular to the horizontal zigzag line. The second layer – by definition a B-layer with k = 1 (white spheres in Fig. 2) – is translated with respect to the first layer in a way that the sphere in the second layer with node coordinates (i, j, k) = (0, 0, 1) lays in the cavity formed by the spheres in the first layer with node coordinates (0, 0, 0), (1, 0, 0) and (0, 1, 0). To convert the lattice coordinates to Cartesian coordinates, the knowledge of the distance between the nodes in the three directions x, y, and z is needed. These distances

2.1. Fundamental assumptions The following four assumptions are at the basis of the new random walk model: [1] The particles move in a hexagonal close packed lattice (HCPL) from one node to another. This means that the trajectories of the particles are discretised. The HCPL is thus used as a framework to discretise particle motion in dense particulate systems.

Fig. 1. Hexagonal close packed lattice (HCPL).

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163

j,y

157

radial boundary confined random environment

k,z

i=0,j=1,k=0 j=constant i=0,j=0,k=0

i,x i=1,j=0,k=0

i=0,j=0,k=1

Fig. 2. Top view of the three-dimensional orthogonal lattice coordinate system defined in a hexagonal close packed lattice. pffiffi pffiffi are Δx ¼ Ds2 3, Δy ¼ Ds and Δz ¼ Ds3 6, where Ds is the diameter of the sphere (= the distance between two adjacent nodes).  pffiffi Each  node (i, j) in the second layer is translated with ðΔx; ΔyÞ ¼ Ds6 3 ; D2s with regard to the corresponding node in the first layer.

2.1.2. Real versus phantom particles To reduce calculation time and power, only a fraction, freal, of the total number of particles, Ntot, is assumed to actually move around in the lattice, implying that only a number Nreal ¼ intðf real Ntot Þ

ð2Þ

of particles is considered to act like so-called real particles. The number Nph ¼ Ntot −Nreal

ð3Þ

then represents the number of so-called phantom particles. All these phantom particles create the randomly changing environment. Phantom particles are not considered in calculations in terms of motion through the lattice, because they are given an initial position and then, at each step, they are randomly redistributed over the

non-occupied nodes of the HCPL. In this way, they indirectly affect the motion of the real particles without really moving. This will be discussed in Section 3. 2.1.3. Number of particles The total number of particles Ntot and the occupied node fraction fnode of the lattice nodes in the confined environment can be determined by   Ntot ¼ int ρp V CE f node ¼

ð4Þ

Ntot Nnode

ð5Þ

with ρp the particle number density (m − 3) in the confined environment with volume VCE and with Nnode the total number of lattice nodes in the confined environment, available for particles to occupy. The occupied node fraction, fnode, is defined as the fraction of nodes in the lattice that are occupied by particles. The value of Nnode of the lattice for a fixed value of VCE can be easily altered by decreasing or increasing the distance between the nodes of the lattice. By using Eqs. (2) and (3), the number of real and phantom particles can easily

Z

0

(a)

R

(b)

Fig. 3. Discretisation of the cylindrically confined environment into cylindrical shell control volumes: (a) side view, (b) top view.

158

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163

be determined. When the particles in the random walk model have a diameter that equals the distance between two adjacent nodes, the voidage εCE can be defined as the volume fraction of air in the confined environment. The number of particles in the confined environment is then determined by Ntot ¼

ð1−ε CE Þ⋅V CE Vp

ð6Þ

with Vp the volume of one particle. 2.2. Defining random motion in a hexagonal close packed lattice Random walks on hexagonal lattices have been studied in the past in a statistical way [26,27], but an actual model describing random motion in a hexagonal close packed lattice was not found in literature. The latter will be treated in this section. The lattice coordinates of the neighbouring nodes of a given node in the HCPL with coordinates (i, j, k) are given in Tables 1 and 2. It can be noticed that these coordinates depend on the parity of i and k, which is a consequence of the stacking of A- and B-layers. Pure random walk in a lattice comes with two requirements concerning symmetry, each on a different time-scale: 1. Scale of an elementary time-step: each possible jump from a certain node to a neighbouring node occurs with the same probability. 2. Scale of many elementary time-steps: in Tables 1 and 2, all possible changes in coordinates (for each parity of i or k) have to be symmetrically distributed. The first requirement concerns the symmetry at the level of one step. The second requirement concerns the symmetry at the level of many steps at which each possible jump listed in Tables 1 and 2 can be statistically considered to be a separate jump. There are 36 possible jumps when taking into account the parity of i and k: • 12 possible jumps to a neighbouring node in the same layer (Table 1: 6 possible jumps for each parity of i). • 12 possible jumps to a neighbouring node in the layer above (Table 2: 12 possible jumps with k + 1 for each parity of i or k) • 12 possible jumps to a neighbouring node in the layer below (Table 2: 12 possible jumps with k − 1 for each parity of i or k) To investigate the second requirement, it was analysed how all 36 possible jumps are distributed over all possible Δi, Δj and Δk. The number of jumps with a certain Δi, Δj and Δk were determined by summation in Tables 1 and 2. The result is shown in Table 3 (Δi and Δj) and in Table 4 ðΔkÞ. In order for all possible jumps to be symmetrically distributed, in Table 3, there should be symmetry of the numbers with respect to ðΔi; ΔjÞ ¼ ð0; 0Þ, Δi ¼ 0, Δj ¼ 0 and the diagonals. However, this is not the case. Jumps with ðΔi; ΔjÞ ¼ ð−1; −1Þ and ðΔi; ΔjÞ ¼ ðþ1; þ1Þ are privileged with respect to the jumps with ðΔi; ΔjÞ ¼ ð−1; þ1Þ and ðΔi; ΔjÞ ¼ ðþ1; −1Þ, respectively. Therefore, it can be concluded that, due to the definition of an orthogonal three-dimensional coordinate system in this study, it is not possible to define pure random motion in the lattice by the mechanism of elementary jumps to the twelve closest neighbouring nodes. Table 1 Lattice coordinates of the neighbouring nodes in the same layer as the node with coordinates (i,j,k). Node

i = even

i = odd

1 2 3 4 5 6

(i + 1,j,k) (i + 1,j-1,k) (i,j − 1,k) (i − 1,j − 1,k) (i − 1,j,k) (i,j + 1,k)

(i + 1,j + 1,k) (i + 1,j,k) (i,j − 1,k) (i − 1,j,k) (i − 1,j + 1,k) (i,j + 1,k)

Table 2 Lattice coordinates of the closest neighbours (white rows) and a second closest neighbour (grey row) in an adjacent layer of a node with coordinates (i,j,k).

k=odd

k=odd

k=even

k=even

i=even

i=odd

i=even

i=odd

(i+1,j,k±1) (i,j,k±1) (i,j+1,k±1) (i+1,j-1,k±1)

(i+1,j+1,k±1) (i,j,k±1) (i,j+1,k±1) (i+1,j,k±1)

(i,j,k±1) (i,j-1,k±1) (i-1,j-1,k±1) (i-1,j,k±1)

(i,j,k±1) (i,j-1,k±1) (i-1,j,k±1) (i-1,j+1,k±1)

Node 1 2 3 4

One way to get around this problem, is to allow a jump to a specific second closest neighbouring node in each adjacent layer. The lattice coordinates of this node are given in Table 2. This node is specific because it was chosen such that the asymmetry disappears (Table 3), hence the second requirement concerning symmetry is fulfilled. The first requirement for pure random motion is that each of the possible 14 elementary jumps (12 to the closest and 2 to the second closest neighbour) has the same probability to occur. Each jump can be modelled as a superposition of a vertical motion (kdirection) and a horizontal motion (i,j-direction). Concerning the vertical motion, for each elementary jump, there are four possibilities with Δk ¼ þ1, four with Δk ¼ −1, and six with Δk ¼ 0. To fulfil the first requirement, a uniform random number Nrand ∈ [0, 13] is needed. The following transformation is used: 0≤Nrand b6⇒Δk ¼ 0 6≤Nrand b10⇒Δk ¼ þ1 10≤N rand ≤13⇒Δk ¼ −1

Each jump now has a probability of 1/14 to occur, meaning that the first requirement for pure random motion is also fulfilled and random motion is defined in a hexagonal close packed lattice. For completeness, it should be noticed that introducing the idea of 14 elementary jumps (12 to the closest and 2 to the second closest neighbour) is only one way to solve the problem of symmetry. Another possible solution is keeping the 12 elementary jumps and appropriately biassing the jump probabilities. 3. Model implementation In the model description, only real particles move from one node to another while the ensemble of phantom particles constitutes the randomly changing environment by the randomly changing positions of the phantom particles. As will be detailed in this section, the environment can be described by a variable one-dimensional array of random numbers. The model implementation includes (1) the initialisation of the algorithm by allocating the particles (real and phantom) to the lattice nodes at time zero, (2) the implementation of random jumps and (3) the implementation of the randomly changing environment.

Table 3 Distribution of all jumps over Δi and Δj (white columns: 12 closest neighbours, grey columns: 12 closest + 2 second closest neighbours).

Δi = 0

Δi = −1 Δ j = −1 Δj=0 Δ j = +1

3 4 1

3 6 3

6 8 6

Δi = +1 6 8 6

1 4 3

3 6 3

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163 Table 4 Distribution of all the jumps over Δk (white columns: 12 closest neighbours, grey columns: 12 closest + 2 second closest neighbours).

Δk = −1 12

Δk = 0 16

12

Δk = +1 12

12

16

Table 5 Simulation identification number (sim _ id#) and number of real particles (Nreal) considered for a certain dimension Dcre of the confined environment and a certain real particle fraction freal. εCE

DCE (m)

0.40 0.60 0.80 0.40 0.60 0.80

0.025

freal = 0.004 Nreal

3.1. Implementing particle allocation at time zero To initialise the system, first, the particles needed to be allocated. First, a one-dimensional array A was generated, containing NA random numbers with a value between 0 and 1, generated with a pseudo-random number generator. Next, each lattice node in the confined environment was identified by a unique integer number q ∈ [1, NA]. The random number A[q] was then used to allocate the particles. For this, the occupied node fraction, fnode, was used. By convention, if A[q] ≤ fnode, the node was occupied by a particle (real or phantom). If A[q] > fnode, it was not. To determine whether the particle was real or phantom, the real particle fraction, freal, was used. If A[q] ≤ freal ⋅ fnode, the node was occupied by a real particle. If freal ⋅ fnode b A[q] ≤ fnode, a phantom particle was fixed at this location. 3.2. Implementing random jumps At each step and for each random jump of a real particle i ∈ [1, Nreal], it was first determined to which node the real particle was going to jump. Therefore, a one-dimensional array B was generated which contained Nreal random numbers with a value between 0 and 13, again generated with a pseudo-random number generator. Each random number B[i] was then used to determine to which node the real particle i would jump (see Section 2.2). Subsequently, it was checked whether this node was outside the boundary of the confined environment or already occupied by another real particle. If a real particle occupied this node, the jump did not take place and the particle i remained at its original position. If the node was not occupied by a real particle, A[q] was used to check whether this node was occupied by a phantom particle. If freal ⋅ fnode b A[q] ≤ fnode, the real particle remained at its position. If not, the real particle performed the jump to its new position. 3.3. Implementing a random environment

Ncv ¼

40,000

sim _ id# 1 2 3 4 5 6

tsim (days)

Nreal

2

1250

20

10,000

sim _ id#

tsim (days)

7 8 9 10 11 12

0.5

6

RCE W cv

ð7Þ

with Wcv the width of a shell and RCE the radius of the confined environment. The shell width Wcv in Eq. (7) was given a value of 0.0125 m. With a diameter, DCE, of the confined environment of

µ_x µ_y µ_z

2.5 2

Sim_10 Sim_3 Sim_1

1.5 1

Displacement (Ds)

To analyse the model, the particles in the confined environment were considered to be particles with a diameter equal to the distance between two adjacent nodes (Dp = Ds). As a case study, two cylindrical confined environments were considered with a diameter DCE of 0.025 m and 0.05 m, respectively. In both cases, the height HCE of the confined environment was assumed to be the same as the diameter. In this way, the particles could move through the lattice over a similar distance in both the radial and the axial direction so random walk in both directions could be compared. The voidage, εCE, was attributed values of 0.4, 0.6 and 0.8, respectively. For each dimension DCE and each value of εCE, two simulations were carried out with a real particle fraction, freal, of 0.004 and 0.001, respectively. This means that twelve simulations were carried out. The simulations were identified with a number, sim _ id#, as

0.05

5000

freal = 0.001

given in Table 5. The value of the distance between two adjacent nodes, Ds (= Dp for this case study), was given a value so that the number of real particles, Nreal (Eq. (2)), covered a broad range of values from 1250 to 40,000. For a certain value of the dimension DCE and the real particle fraction, freal, the number of real particles, Nreal, was kept constant (Table 5) for the different values of εCE by varying the distance, Ds, between two adjacent nodes. In this way, the confined particle motion could be statistically compared for different values of εCRE. The distance between two adjacent nodes, Ds, was 225 μm, 195 μm and 155 μm for a voidage, εCE, of 0.4, 0.6 and 0.8, respectively. A simulation lasted 2 × 10 4 elementary steps. During each step each real particle was subjected to a random jump. A log was kept of the displacement of each real particle in the x-, y- and z-direction with regard to its initial state. After completing each simulation, the parameters characterising the distribution of the displacements, namely the mean value Δr and the standard deviation, σ, were analysed. Then, it was considered how the dimension of the confined environment, the real particle fraction and the voidage affected the distribution. Finally, it was verified whether the algorithm was self-sustaining. Therefore, the cylindrical confined environment was divided in cylindrical shell volumes (or control volumes) as shown in Fig. 3. The number of control volumes, Ncv, was determined by

If the random numbers in array A would remain constant, the phantom particles would remain fixed at their initial position, meaning that the environment would not be random and certain preferential paths would be created which would prevent implementing genuine random walk. Therefore, at each step, a new array A was generated which virtually repositioned the phantom particles. 4. Model analysis

159

0.5 0

-0.5 -1 -1.5 -2 -2.5 0

2000

4000

6000

8000 10000 12000 14000 16000 18000 20000

Step Fig. 4. Displacement of the centre of gravity in the x-, y- and z-direction (multiples of Ds) versus the steps for sim _ id# 10, 3 and 1, respectively (see also Table 5).

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163

0.025 m and 0.05 m, two and four control volumes were considered, respectively. The number of particles in each control volume, Ntot, cv, was then calculated as Ntot;cv ¼ ρp V cv

ð8Þ

By checking whether the number of real particles in the control volumes (Nreal, cv) fluctuated around the initial value during a simulation, the self-sustaining nature of the algorithm was verified, i.e., no depletion or accumulation of real particles occurred during a simulation.

(a)

Sim_1 Sim_4 Sim_7 Sim_10

45 40

Standard Deviation (Ds)

160

35 30 25 20 15 10

5. Results and discussion 5

5.1. Simulation time, performance and stability

0

0

2000

4000

6000

8000 10000 12000 14000 16000 18000 20000

Step

(b)45

Sim_2 Sim_5 Sim_8 Sim_11

40 35

Standard Deviation (Ds)

The model was programmed in C++ in Microsoft Visual Studio 2005. The simulations were carried out on a personal computer with an Intel Core i7-2600 3.4 GHz processor. The times needed to complete each simulation, tsim, are listed in Table 5. As expected, tsim increased rapidly for an increasing number of real particles considered. This reveals that the algorithm was able to significantly reduce calculation time when the number of real particles was kept low. By performing three simulations with a significantly higher number of real particles (40,000) compared to all other simulations (≤10,000), the algorithm was also tested for its performance and stability. These three simulations each took 20 days during which no error notification or anomaly occurred, implying that the model was able to simulate random walk of a broad range of number of real particles.

30 25 20 15 10

5.2. Statistical analysis of the random walk

5 0 0

2000

4000

6000

8000 10000 12000 14000 16000 18000 20000

Step Sim_3 Sim_6 Sim_9 Sim_12

(c)

45 40

Standard Deviation (Ds)

In Table 6, for all simulations, the displacement of the centre of mass, Δμ r , of all real particles is shown after 5000, 10,000 and 20,000 steps, respectively. It can be noticed that the value of Δμ r was more or less of the order of magnitude of the particle diameter, which indicates that the particles did not drift. This is also shown in Fig. 4 where for three different simulations the displacement of the centre of mass in the x-, y-, and z-directions (Δμ x ; Δμ y ; Δμ z ) is plotted versus the number of steps. The displacements clearly fluctuated around zero with an amplitude more or less of the order of magnitude of the particle diameter. When looking at the standard deviation, σr, of the distribution of the displacements of the real particles in Table 7, it can immediately be observed that the standard deviation increased with increasing number of steps. This is also plotted in Fig. 5a, b, and c for all the simulations with voidage εCE equal to 0.4, 0.6, and 0.8, respectively. From each of these plots, it can be derived that the real particle fraction, freal, had a negligible effect on the standard deviation. However, it can be observed that the standard deviation increased when the dimensions of the confined environment increased. Indeed, in a larger confined environment the particles

35 30 25 20 15 10 5 0

0

2000

4000

6000

8000 10000 12000 14000 16000 18000 20000

Step Table 6 Displacement r (multiples of Ds) of the centre of mass of the real particles after 5000, 10,000 and 20,000 steps, respectively. εCE

0.4

freal

0.004

0.001

0.004

0.001

0.004

0.001

0.62 0.85 1.03 0.17 0.18 0.35

1.16 1.56 1.27 0.30 0.29 0.61

1.08 1.02 1.95 0.16 0.27 0.49

0.75 2.09 0.62 0.87 1.19 0.55

1.13 1.55 1.76 0.45 0.77 1.25

1.14 1.39 2.65 1.81 1.09 1.14

DCE (m)

0.025

0.05

0.6

0.8

Fig. 5. Standard deviation of the distribution of the particle displacements (multiples of Ds) versus the steps for simulations with (a) εCE = 0.4, (b) εCE = 0.6 and (c) εCE = 0.8 (see also Table 5).

Nstep

5000 10,000 20,000 5000 10,000 20,000

can travel over larger distances before encountering the boundaries of the confined environment. Looking at Fig. 5, it can be noticed that the relation between the standard deviation and the number of steps slightly changed for different values of the voidage εCRE. By fitting a power function to the plots in Fig. 5, it was analysed how the standard deviation of the simulated random walk in a confined

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163 Table 7 Standard deviation σr (multiples of Ds) of the distribution of the displacements of the real particles after 5000, 10,000 and 20,000 steps, respectively. εCE

0.4

freal

0.004

0.001

0.004

0.001

0.004

0.001

14.39 19.21 25.07 15.19 21.05 28.30

14.29 19.51 25.33 15.16 20.80 28.12

21.56 27.81 33.55 23.56 31.15 37.69

21.31 27.89 33.76 23.41 30.94 37.61

27.19 33.57 39.01 29.07 36.41 41.71

26.72 33.54 38.38 28.96 36.34 41.58

DCE (m)

0.025

0.05

0.6

0.8

5000 10,000 20,000 5000 10,000 20,000

environment related to the number of steps. For all simulations, this relation was of the form  b σ r ¼ a Nstep

ð9Þ

The values of the parameter a, the power coefficient b, and the coefficient of determination R 2 (i.e. a measure of how well the fitted curve (Eq. (9)) predicts the results of the model) are listed in Table 8. It can be noticed that the power coefficient b decreased and the parameter a increased with increasing εCE. Moreover, for all simulations, the power coefficient was smaller than 0.5, i.e. the value for a random walk in an infinite space (Eq. (1)). Since a value of 0.5 is valid for random walk in an infinite three-dimensional space, the value of the power coefficient can certainly not become larger than 0.5 for random walk in a finite dimensioned space. This can be explained by the finite dimensions of the confined environment, restricting the motion of the particles within the boundaries of the environment, hence suppressing the increase of σr when comparing to random walk in an infinite space. For the smallest value (εCE = 0.4), the standard deviation reached a value of 25Ds after approximately 17,000 steps, while for larger values of εCE, this value was already reached after 7000 steps. This can be easily explained because for lower voidages, hence larger particle densities, the particles will encounter more resistance to “walk” away from their initial position. The considered real particle fraction, freal, had no significant influence on the relation between the standard deviation and the number of steps. This is an indication that the approach in this model, namely replacing real particles by phantom particles that

Number of Real Particles in Control Volume

4000 CV1 CV2 CV1 CV2

3500 3000

Sim_3 Sim_3 Sim_12 Sim_12

2500 2000 1500 1000 500 0

0

2000

4000

6000

Table 8 Values of the parameter a and the power coefficient b describing the relation between the standard deviation of the displacement and the number of random steps. εCE

Nstep

8000 10000 12000 14000 16000 18000 20000

Step Fig. 6. Number of real particles in the first two control volumes versus the number of random steps, for simulations 3 and 12, respectively (see also Table 5).

161

0.40 0.60 0.80 0.40 0.60 0.80

DCE (m)

0.025

0.05

freal = 0.004 = 0.001

freal = 0.001

a

b

R2

a

b

R2

0.2920 0.5341 0.7261 0.2554 0.4727 0.7041

0.4539 0.4267 0.4137 0.4781 0.4515 0.4247

0.9985 0.9921 0.9859 0.9997 0.9944 0.9862

0.2821 0.5152 0.7438 0.2587 0.4714 0.7041

0.4579 0.4300 0.4104 0.4759 0.4514 0.4247

0.9990 0.9948 0.9854 0.9996 0.9945 0.9862

virtually change place every time step, is able to mimic random motion in the confined environment, regardless the value of freal. In Fig. 6, for simulations 3 and 12, the number of real particles in the first two control volumes (Nreal, cv) was plotted against the number of steps. It can be observed that the number of real particles in the control volumes fluctuated only slightly around the initial value, meaning that the algorithm was able to simulate a balanced random exchange of particles between the control volumes. Hence, the algorithm did not lead to a non-uniform particle distribution in the confined environment. The latter can also be seen in Fig. 7, where the particle distribution in a certain cross section of the confined environment was plotted after 50, 10,000 and 20,000 steps. The particle distribution remained uniform during the complete simulation. In Fig. 8, for simulation 4, a plot is shown of the initial and the final position of a few real particles. In Fig. 9, the simulated random walk of a real particle is plotted for simulation 9. From both figures it can be concluded that the real particles actually moved through the system and that they did not vibrate around their initial position. The above observations, together with the observation that the centre of mass fluctuated around zero, indicated that the consideration of only a fraction of the particles and the use of a pseudo-random number generator was already sufficient to model random walk since it did not lead to drift or non-uniform particle distribution. 6. Conclusions In this paper, a new discrete and statistical random walk model was presented that describes random walk of a particle in a confined and relatively dense environment defined on a hexagonal close packed lattice. The lattice is solely introduced in the model to discretise the random particle motion in the confined environment, not to describe or study particle motion on this particular lattice. The statistical approach enables the description of random walk without having to consider individual particle–particle collision events which in dense particle system cause the random walk of the particles. The model considered only a fraction of the total amount of particles (the so-called real particles) and hence reduced computational power significantly. The number of real particles had no significant influence on the random walk characteristics. The model was tested on a particulate system with a cylindrical geometry and proved to be reliable to simulate random walk. The algorithm, which makes use of a pseudo-random number generator, did not create drift or did not evolve to non-uniform particle distributions. This indicates that the pseudo-random numbers used in the algorithm were already sufficient to describe random motion in a multi-particle system. The relation between the standard deviation of the particle displacement distribution and the steps was characterised by a power coefficient that was smaller than 0.5, i.e. the value of the power coefficient for a random walk in an infinite space. This was due to the confinement of the random environment, which restricts the random motion of the particles within the boundaries of the system. Finally, it can be concluded that the modelling approach presented in this work to model random particle motion can be incorporated in the model

162

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163

y

y

y

x 50 random steps

x 10000 random steps

x 20000 random steps

Fig. 7. Particle distribution in a cross section of the confined environment after 50, 10,000 and 20,000 random steps, respectively.

that is currently being developed for fluid bed processes where particle–particle collisions occur frequently and for which a description of the individual particle motion is required. Nomenclature a parameter A one dimensional array (1xN) b power coefficient B one dimensional array (1xN) CE confined environment D diameter (m) DEM Discrete Element Method f fraction H height (m) HCPL hexagonal close packed lattice int real to integer conversion function N number r distance (m) R radius (m) sim_id# simulation identity number t time (days) V volume (m³) W width (m) μ centre of mass (Ds) σ standard deviation (Ds)

ε ρ Δ

voidage density (m − 3) displacement (distance) or interval (time)

Superscripts and subscripts

A cv lp node p ph r real rsp s sim steps tot x

array control volume lattice packing node particle phantom particles radius (polar coordinate) real particles random spherical packing sphere simulated/simulation steps total x-direction

0,025

0,05 0,020

0,04

Z

0,015

Z

0,03

0,010

0,02 0,005

0,01 0, 0,01

0,00

-0,005

0,00 -0,01

0,00

X

-0,01 0,01

0,02

-0,02

0,000

Y

Fig. 8. Initial and final positions of real particles for sim _ id# = 4 (see also Table 5).

0,000

X

0,005

0,010

0, 0,005 0,000 -0,005 -0,010 Y

Fig. 9. Simulated random walk of a real particle in a confined environment for sim _ id# = 9 (see also Table 5).

M. Vanderroost et al. / Powder Technology 221 (2012) 155–163

y z

y-direction z-direction

Acknowledgements The authors acknowledge the financial support of the Special Research Fund (BOF) of the Ghent University. References [1] F. Ronsse, J.G. Pieters, K. Dewettinck, Combined population balance and thermodynamic modelling of the batch top-spray fluidised bed coating process. Part I— model development and validation, Journal of Food Engineering 78 (2007) 296–307. [2] F. Ronsse, J.G. Pieters, K. Dewettinck, Combined population balance and thermodynamic modelling of the batch top-spray fluidised bed coating process. Part II — model and process analysis, Journal of Food Engineering 78 (2007) 308–322. [3] F. Ronsse, J.G. Pieters, K. Dewettinck, Modelling heat and mass transfer in batch, topspray fluidised bed coating processes, Powder Technology 190 (2009) 170–175. [4] M. Vanderroost, F. Ronsse, K. Dewettinck, J.G. Pieters, Modelling coating quality in fluidised bed coating: spray sub-model, Journal of Food Engineering 106 (2011) 220–227. [5] R. Luger, Exact non-parametric tests for a random walk with unknown drift under conditional heteroscedasticity, Journal of Econometrics 115 (2003) 259–276. [6] T. Chuluun, C.S. Eun, R. Kiliç, Investment intensity of currencies and the random walk hypothesis: cross-currency evidence, Journal of Banking & Finance 35 (2011) 372–387. [7] F.R. Hallett, P.A. Speight, R.H. Stinson, Introductory Biophysics, Methuen publications, Carswell Company Limited, USA, 1977. [8] P.G. Doyle, J.L. Snell, Random Walks and Electric Networks, Mathematical Association of America, Washington D.C, 1984. [9] Z. Burda, J. Duda, J.M. Luck, B. Waclaw, Localization of the maximal entropy random walk, Physical Review Letters 102 (2009) 160602. [10] O. Vitrac, M. Hayert, Effect of the distribution of sorption sites on transport diffusivities: a contribution to the transport of medium-weight-molecules in polymeric materials, Chemical Engineering Science 62 (2007) 2503–2521. [11] R. Caulkin, A. Ahmad, M. Fairweather, X. Jia, R.A. Williams, Digital predictions of complex cylinder packed columns, Computers and Chemical Engineering 33 (2009) 10–21.

163

[12] H. Nakamura, E. Abe, H. Yamada, Coating mass distributions of seed particles in a tumbling fluidized bed coater Part II. A Monte Carlo simulation of particle coating, Powder Technology 99 (1998) 140–146. [13] S. Abe, S. Thurner, Anomalous diffusion in view of Einstein's 1905 theory of Brownian motion, Physical A: Statistical Mechanics and its Applications 356 (2005) 403–407. [14] N. Ikeda, Network formed by traces of random walks, Physica A 379 (2007) 701–713. [15] L. da Fontoura Costa, G. Travieso, Exploring complex networks through random walks, Physical Review. E, Statistical, Nonlinear, and Soft Matter Physics 75 (2007) 16102-1-7. [16] G. Wei, New approaches to shapes of arbitrary random walks, Physica A 222 (1995) 155–160. [17] W. Duangkhamchan, F. Ronsse, F. Depypere, K. Dewettinck, J.G. Pieters, Comparison and evaluation of interphase momentum exchange models for simulation of the solids volume fraction in tapered fluidised beds, Chemical Engineering Science 65 (2010) 3100–3112. [18] W. Duangkhamchan, F. Ronsse, F. Depypere, K. Dewettinck, J.G. Pieters, CFD study of solids concentration in a fluidised-bed coater with variation of atomisation air pressure, Powder Technology 212 (2011) 103–114. [19] B.P.B. Hoomans, J.A.M. Kuipers, W.J. Briels, W.P.M. Van Swaaij, Discrete particle simulation of bubble and slug formation in a two-dimensional gas-fluidised bed: a hard-sphere approach, Chemical Engineering Science 51 (1) (1996) 99–118. [20] Y. Wang, M.T. Mason, Two-dimensional rigid body collisions with friction, Journal of Applied Mathematics and Mechanics 59 (1992) 635–642. [21] J. Schäfer, S. Dippel, D.E. Wolf, Force schemes in simulations of granular materials, Journal de Physique I 6 (1) (1996) 5–20. [22] K.D. Kafui, C. Thornton, M.J. Adams, Discrete particle-continuum fluid modelling of gas–solid fluidised beds, Chemical Engineering Science 57 (2002) 2395–2410. [23] R. Beetstra, M.A. van der Hoef, J.A.M. Kuipers, Numerical study of segregation using a new drag force correlation for polydisperse systems derived from Lattice-Boltzmann simulations, Chemical Engineering Science 62 (2007) 246–255. [24] R. Begum, M.A. Basit, Lattice Boltzmann Method and its applications to fluid flow problems, European Journal of Scientific Research (ISSN: 1450-216X) 22 (2) (2008) 216–231. [25] Johannes Kepler, Strena seu de nive sexangula (The six-cornered snowflake), MR0927925, 978-1589880535, 1611. [26] M. Koiwa, S. Ishioka, Random walk properties of lattices and correlation factors for diffusion via the vacancy mechanism in crystals, Journal of Statistical Physics 30 (1983) 477–485. [27] M. Sahimi, B.D. Hughes, L.E. Scriven, H.T. Davis, On Pólya random walks, lattice Green functions, and the bond percolation threshold, Journal of Physics A: Mathematical and General 16 (1983) L67–L72.