CFD-DEM Simulation of a Coating Process in a Fluidized Bed Rotor Granulator

Abstract

Coating of particles is a widely used technique in order to obtain the desired surface modification of the final product, e.g., specific color or taste. Especially in the pharmaceutical industry, rotor granulators are used to produce round, coated pellets. In this work, the coating process in a rotor granulator is investigated numerically using computational fluid dynamics (CFD) coupled with the discrete element method (DEM). The droplets are generated as a second particulate phase in DEM. A liquid bridge model is implemented in the DEM model to take the capillary and viscous forces during the wet contact of the particles into account. A coating model is developed, where the drying of the liquid layer on the particles, as well as the particle growth, is considered. The simulation results of the dry process compared to the simulations with liquid injection show an important influence of the liquid on the particle dynamics. The formation of liquid bridges and the viscous forces in the liquid layer lead to an increase of the average particle velocity and contact time. Changing the injection rate of water has an influence on the contact duration but no significant effect on the particle dynamics. In contrast, the aqueous binder solution has an important influence on the particle movement.

1. Introduction

For various products in the pharmaceutical, chemical, and food industries, the coating of solid particles with different materials is an important processing step [1]. A wide range of equipment is available for coating, such as drum coaters, bubbling and spouting fluid beds, and fluidized bed rotor granulators. In all these processes, it is crucial to consider the change of the particle dynamics due to an increasing amount of liquid in the apparatus. While agglomeration is encouraged in some processes, it should be avoided in others. Especially in the pharmaceutical industry, round, coated pellets are needed for oral drug delivery [2]. The fluidized bed rotor granulator was developed in order to be able to carry out the individual process steps of spheronization, coating, and drying in one apparatus. It consists of a rotating round base plate with an annular gap between the rotor and the wall. Air flows through the gap, which fluidizes and dries the pellets. The coating liquid is sprayed onto the pellets from a horizontal injection nozzle. Due to its specific fluid dynamics and mechanical conditions in the process chamber, pellets with high strength, smooth surface, and high sphericity could be obtained [3,4,5].

For the optimization and control of coating processes, it is essential to know and describe the particle kinetics and dynamics, as well as the growth of the coating layer on the particles. These data can be determined with numerical simulation. Continuous Euler–Euler flow models were often used to describe the gas phase and particles in the fluidized beds [6,7,8,9,10]. These approaches allow a relatively fast calculation but simplify the particle-particle interactions, and details of the individual particle movements are lost. Due to the increasing computing power, the more computationally intensive Eulerian–Lagrangian approaches are used increasingly to simulate multiphase flows. Thereby, the computational fluid dynamics (CFD), is coupled with the discrete element method (DEM) [11]. With CFD, the flow field of the gas in the processing apparatus can be determined, which treats the fluid phase as a continuum. In DEM, the interactions are calculated for every single particle based on models describing the physical properties of the particles, such as the adhesion, and their mechanical behavior during slow, fast, and repeated stressing [12,13].

Neuwirth et al. [3] investigated the particle dynamics in a fluidized bed rotor granulator under dry conditions experimentally and with CFD-DEM simulations. A good agreement of the simulation results with the experiments was shown [14,15]. However, only a few studies have been carried out to examine the wetting of particles caused by an additional spraying process. Fries et al. [16] simulated the moisture distribution in a fluidized bed using the CFD-DEM approach by introducing a liquid containing a spray domain. Suzzi et al. [17] investigated a tablet-coating process, applying the same method of a liquid containing spray domain to DEM simulations. In these works, no droplets are simulated, but the particles or tablets gain a certain amount of liquid by passing this spray domain. To consider sprayed droplets in the simulation, a second discrete phase needs to be added to the DEM model. Since in a conventional spray process a huge number of small droplets are generated, such simulation approaches become very computationally expensive. Sutkar et al. [18] presented a three-phase approach of gas, droplet, and particle flow for a 2D spout fluidized bed with draft plates. Börner et al. [19] simulated a spraying process in a 3D fluidized bed granulation. However, for reasonable time scales, a similarity model-based approach was used to scale the individual droplets. In our previous work [20], we simulated the two-phase spray nozzle with CFD-DEM to obtain the droplet distribution in the fluidized bed apparatus. In the next studies, this time-consuming simulation approach for the sprayed droplets was simplified with the spray domain [1,21] or homogenous wetting of all the particles in the bed [22]. These models are faster; however, they do not take into account the real droplet deposition processes in the spray zone.

In this work, the particle dynamics during a coating process in a fluidized bed rotor granulator was studied. The coating process was simulated with a 3D CFD-DEM coupling in order to obtain additional information about the influence of the injected liquid on the particle dynamics. The capillary forces due to the existence of liquid layers after droplet deposition were taken into account based on the model of Tsuazawa et al. [3], and the viscous forces were calculated according to models of Lian et al. [23] and Popov [24]. In addition, a new coating model is implemented into the DEM code. It describes the droplet deposition on the particles and the wall, as well as the drying of the liquid layer over time. The simulation results of the fluidization process without liquid injection are compared with two cases with different injection rates of water, as well as one case where an aqueous binder solution is sprayed in.

2. Model Description 2.1. CFD Modeling

In this study, a one-way approach is used, which considers the influence of the gas phase on the solid phase, but neglects the influence of the solid phase on the gas phase. With this choice of coupling approach, it was possible to simulate a real process time of about one minute in reasonable calculation time. For the CFD, the gas phase is considered as continuum and the geometry of the apparatus is discretized by mesh cells. In each volume element of the CFD mesh, the gas flow is calculated by solving the volume-averaged Navier–Stokes equations by applying the semi-implicit method for the pressure-linked equations (SIMPLE) algorithm [25]. The equations describing the heat, mass, and momentum transfer in the gas phase can be written as a general equation:

∂_ρgϕ∂t+∇(_{ρg ug}ϕ−_Γϕ∇ϕ)=_Sϕ

where_ρgis the gas density,ϕdescribes the scalar quantity of the gas (e.g., the velocity components, the temperature, turbulent kinetic energy dissipation rate, and turbulent kinetic energy),_ugis the gas velocity,_Γϕspecify the diffusion coefficient, and_Sϕrepresent the source term, which accounts for any sources or sinks that either create or destroy the scalar quantity of the gas.

2.2. DEM Modeling

In DEM, firstly introduced by Cundall and Strack [26], the motion of each single particle is calculated using motion equations for translation and rotation, according to Newton and Euler, respectively. The translational velocity of a particle_vp,iwith the mass_mp,iis determined by calculating the gravitational force and the sum of the contact forces_Fc,ijacting on the particle i due to interactions with other particles j or the apparatus wall:

_mp,id_vp,idt=_mp,ig+∑j=0k_Fc,ij

wheregis the gravitational acceleration.

For the used one-way coupling of the CFD with the DEM, the model of Cundall and Strack is extended by the drag and pressure gradient forces to take into account the influence of the gas phase on the particles:

_mp,id_vp,idt=−_Vp,i∇p+12_{CD ρg urel}2_Ai+_mp,ig+∑j=0k_Fc,ij

where_Vp,i,∇p,_urel, and_Aiare the particle volume, the pressure gradient, the gas velocity relative to the particle, and the frontal projection area of the particle, respectively. The contact forces are described in detail in the following section. The drag coefficient_CDdepends on the particle Reynolds numberR_ep:

R_ep=_urel⋅_ρg⋅_dp,iηvis

where_dp,iis the particle diameter and_ηvis the dynamic viscosity. The correlation used for different flow regimes is summarized in Table 1 [27,28].

For the description of the rotation, the angular velocity_ωp,iof each particle with the moment of inertia_Jp,iand radius_rc,ijis calculated by the sum of the torques_Mt,ijcaused by the tangential forces acting on the particle and_Mr,ijdue to the rolling friction if the particle rotates:

_Jp,id_ωp,idt=∑j=0k(_Mt,ij+_Mr,ij)

2.3. Contact Force

The contact force is calculated according to the viscoelastic soft-sphere model. The normal component of the contact force in Equation (6) is based on the theory by Hertz [29] and Tsuji [30]. The tangential component of the contact force is calculated according to a no-slip approximation of the model by Mindlin and Deresiewicz [31], as proposed by Tsuji (Equation (7)) [30]. In addition, to describe the wet contacts in the simulation of the coating process, viscous forces_Fvis,n, and capillary forces_Fcapare considered. Both forces are described in detail in the following two sections. The normal and tangential component of the contact force, index n, and t, respectively, are calculated as:

_Fc,ij,n=−_{kn δn3/2 nij}−_{ηn up,ij,n}+_Fvis,n+_Fcap

_Fc,ij,t={−_{kt δt tij}−_{ηt up,ij,t}+_Fvis,tif|_Ft,ij|≤μ|_Fn,ij|−μ|_Fn,ij|_tij+_Fvis,tif|_Ft,ij|>μ|_Fn,ij|

wherekdescribes the spring stiffness coefficient,δis the displacement, _nij, as well as_tijare normal and tangential unit vectors,_up,ij,trepresents the tangential component of the collision velocity vector of the particle with another particle or a wall andμis the sliding friction coefficient. Considering viscoelastic deformation behavior that causes energy dissipation, the damping factorηfor the normal and tangential direction is included:

_ηn=2α^m* _knδn1/4

_ηt=2α^m* _kt

whereαrepresents a function of the restitution coefficient, and^m* is the reduced mass of the contact partners. A detailed description of this contact model can be found in Antonyuk et al. [32] and Deen et al. [11].

2.3.1. Capillary Force

During a particle-particle or particle-wall contact, where one of the contact partners has a liquid layer in the contact area, a liquid bridge forms. Thus, capillary forces need to be taken into account. In the literature, several models are available to calculate the capillary forces for symmetric pendular bridges [23,33,34,35,36]. These models were derived on the basis of the total liquid-bridge energy, and the pressure difference across the liquid bridge. In this study, the capillary forces are modeled according to Israelachvili [36]. Fujihashi et al. [37] reported that the model proposed by Israelachvili is in closer accordance to the experimental data presented by Willet [38] than the model of Rabinovich [35]. This model by Israelachvili was already successfully applied to simulate a mixing process in DEM by Y. Tsunazawa et al. [39].

Based on the different geometric shapes, a distinction is made between particle-particle contact and particle-wall contact. Thus, the capillary force acting between the two particles, _Fcap,ppand the capillary force acting between a particle and a wall, _Fcap,pware calculated by the following expressions:

_Fcap,pp=−4πRγcos(θ)1+^{(1+2_VbπRh2−1)−1}

_Fcap,pw=−8πRγcos(θ)1+^{(1+2_VbπRh2−1)−1}

where_Vbis the volume of the liquid bridge, which is calculated as the sum of the liquid volume of each contact partner,hdescribes the distance between the particle centers or the particle center and the wall,θrepresents the wetting angle, andγis the surface tension.

Since the particle radius increases individually per particle during the coating process, the effective contact radiusRis used:

R=_{ri rjri}+_rj

where_riand_rjare the radii of the two contact partners.

A liquid bridge forms when the distance between the contact partners is less than the thickness of the liquid layers_h0on the particles or wall. It is assumed that the liquid always spreads ideally over the surface and a uniform thin film is formed. The thickness of the liquid layer on the particle_h0,pis calculated by subtracting the radius of the dry particle volume_rifrom the radius of the wetted particle determined by the volume of the liquid on the particle_Vp,i,liq, and its volume_Vp,i,dry:

_h0,p=3(_Vp,i,liq+_Vp,i,dry)4π3−_ri

The increase of the dry particle volume by drying is explained in detail in Section 2.4. For the calculation of the layer thickness on the wetted wall, the liquid volume of the wetted wall grid cell_Vcell,iis divided by the surface of the corresponding mesh cell of the wall_Scell:

_h0,w=_Vcell,iScell

A more detailed description of the allocation of the liquid during the simulated coating process can be found in the section on the coating model.

At a critical distance of the contact partners, the liquid bridge breaks. Lian et al. [40] proposed an approximation of the rupture distance_hrupobtained from the numerical solutions of the Young–Laplace equation, which depends on the wetting angle and the volume of the liquid bridge:

_hrup=(1+0.5θ)_Vb1/3

2.3.2. Viscous Force

In addition to the capillary force, a viscous force also acts due to the liquid. The dimensionless capillary number, which describes the relation between theviscous force and surface force, is used to assess the relevance of the viscous forces [41]:

Ca=_{ηvis up,ij}γ

where_up,ijdescribes the relative velocity of the colliding particles or particle with a wall.

From the Reynolds lubrication theory, considering the relative translational collision velocity of the particles in the normal direction, the viscous force in the normal direction_Fvis,n can be determined as [23,42]:

_Fvis,n=6π_ηvis ^R2 _up,ij,nh

where_up,ij,n represents the relative velocity of the colliding particles or the particle with a wall in the normal direction. The distance between the contact partners, at which the viscous force must be taken into account, is physically limited by the roughness of the respective surfaces. Therefore, a minimum distance is set for the calculation of the viscous. Similar treatments can be found in the literature [43,44,45].

In the tangential direction, the viscous force is calculated as proposed by Popov [24]:

_Fvis,t= 2π_ηvisR_up,ij,tln(1+R2h)

2.4. Coating Model

In order to consider a coating process of the particles, which is always a wetting and drying process, a new coating model is implemented in DEM. In the model, the particle volume consists of a dry particle volume_Vp,i,dryand a liquid volume on the particle_Vp,i,liq:

_Vp,i=_Vp,i,dry+_Vp,i,liq

Depending on this volume, the mass of the particle is recalculated in each time step. The droplets are generated as a second particulate phase in DEM, and their movement is simulated. In contrast to the solid particles, no contact forces are considered for the droplets. Depending on the type of the contact, different cases are distinguished for the droplet. If two spray droplets meet, this contact is ignored and the droplet trajectories are unaffected. If a contact with a particle occurs, the volume of the droplet_Vdropis transferred to the particle as a liquid volume of the particle_Vliq,i, and the droplet is deleted from the simulation. At each further contact of the particle with a droplet, the volume of the droplet is added to the liquid volume of the particle:

_Vp,i,liq=_Vliq,i+∑i=0k_Vdrop,i

In the case of a contact with an apparatus wall, which is discretized by a mesh, the volume of the droplet is assigned to the grid cell in the contact. Similar to the particles, each further droplet is added to the liquid volume of the grid cells:

_Vcell,i=_Vliq,i+∑i=0k_Vdrop,i

This creates the so-called “liquid layer” on the particle or grid cell of the wall. With the assumption of an ideal wetting over the entire particle surface and cell surface, the film thickness could be calculated as expressed in Equations (13) and (14).

Furthermore, liquid transfer during a contact, whereas at least one of the contact partners has a liquid layer, needs to be considered. As described in Section 2.3, the volume of the liquid bridge is equal to the sum of the liquid volume of each contact partner.

_Vb,ij,z+1=_Vliq,i,z+_Vliq,j,z

wherez is the number of time steps. When the particles rebound, half of the volume of the liquid bridge remains as liquid volume on each contact partner. This assumption is also used by other authors [46,47].

_Vb,ij,z+1/2=_Vliq,i,z+2=_Vliq,j,z+2

In Figure 1, the liquid transfer during the contact of the two particles with a liquid volume is shown.

The implemented drying process depends on the drying rate_φR, which specifies the amount of liquid mass that dries per time unit. In the case that a process with saturated air is simulated where no drying of the liquid volume takes place, the drying rate is zero. Otherwise, the drying coefficient is greater than zero and the liquid volume of the particle and the wall decreases in each time stepΔtuntil no liquid volume is left. At each time step of the simulation, the liquid volume_Vliq,i,z+1is recalculated:

_Vliq,i,z+1=_Vliq,i,z−_φRΔt_ρl

where_ρlis the density of the liquid.

In addition, a loading coefficient of the liquid α is introduced to consider the solid concentration in a solution used for coating. Therefore, it defines the ratio of the liquid volume converted to dry volume on the particle during the drying process. The assumption is made that the properties of the liquid do not change during the drying process, and thus, there is no change in concentration. With a loading coefficient of zero, no dry volume is formed, the liquid layer only decreases over time. If the loading coefficient is one, the decrease of the liquid volume corresponds to the increase of the dry volume per time step. The dry particle volume_Vp,i,dryconsists of the initial volume of particle_Vp,i,iniand a dried volume_Vdry,i:

_Vp,i,dry=_Vp,i,ini+_Vdry,i

Similar to the liquid volume, the dry volume is recalculated in each time step of the simulation:

_Vdry,i,z+1=_Vdry,i,z+α_φRΔt_ρl

This results in a total particle volume afterz+1time steps:

_Vp,i,z+1=_Vp,i,ini+_Vdry,i,z+α_φRΔt_ρl+_Vliq,i,z−_φRΔt_ρl+∑i=0k_Vdrop,i

The radius_riused for calculation of the thickness of the liquid layer in Equation (13), thus increases during the drying process:

_ri,z+1=3(_Vp,i,ini+_Vdry,i,z+α_φRΔt_ρl)4π3

The liquid volume of the grid cell afterz+1time steps is expressed as:

_Vcell,i,z+1=_Vliq,i,z−_φRΔt_ρl+∑i=0k_Vdrop,i

For the wall cells, the formation of a dry layer is not considered, only the drying process. In this model, no heat transfer between the gas phase and solid phase is taken into account since it is assumed that both phases have the same temperature. Drying of the particles becomes unimportant when the gas and liquid temperatures are low and the moisture content of the air is high. The coating model is shown schematically in Figure 2.

3. Simulation Setup 3.1. Geometry of the Rotor Granulator

In this study, the coating process in a rotor granulator with an inner diameter of 295 mm is simulated (Figure 3). The central element of the rotor granulator is the horizontally rotating unstructured plate with a diameter of 268 mm located in the middle of the apparatus, on which the starting cores or pelletized materials are located. The gas flows vertically into the particle bed via a two-millimeter-wide annular gap between the rotor plate and the apparatus wall, as can be seen in Figure 3a. The inflow surface for the gas into the process chamber is much smaller than that of the conventional fluidized beds, resulting in high velocity in the gap. A single-fluid nozzle with an outlet diameter of 1 mm is installed horizontally in the process chamber to inject the liquid for the coating of the particles. The nozzle is located 43.5 mm above the rotor plate.

3.2. Simulation Parameters

The geometry of the apparatus was discretized into hexagonal mesh cells. The mesh in the region where particles are located consists of approximately equally sized cells. The gas field was solved with an ANSYS Fluent 2019R2 (ANSYS, Inc., Canonsburg, PA, USA, 2019). For the CFD simulation, an inlet velocity of 1 m/s, and outlet pressure of 1.01325 bars were taken as boundary conditions. This resulted in an average gas velocity in the gap of 27 m/s. Air at 26 °C was used for the gas phase and it was considered as compressible gas in the simulation. The rotor disc rotated with 300 rpm. The CFD simulation parameters are listed in Table 2.

The particulate phases were calculated with EDEM 2020 (DEM Solution Ltd., Edinburgh, United Kingdom, 2019) The particles were generated above the rotor plate. Similar to the CFD simulation, the rotor rotated with 300 rpm. For comparison with future experiments in a second study, round particles consisting of a ceramic core and a shell of polyvinyl butyral are simulated. The particle properties needed for the DEM simulations were obtained by several experimental setups. The coefficient of the restitution was measured with a free-fall device, as described in the work of Grohn et al. [2]. The Young’s modulus was determined with a Nanoindenter (Hysitron TI Premier, Bruker Corporation, Billerica, MA, USA), and the friction coefficients with a Texture Analyser^® (TA.XTplus, Stable Micro Systems, London, UK). The particle properties needed for the DEM simulations are summarized in Table 3.

In this study, four different cases were analyzed to investigate the influence of the coating liquid on the particle behavior. In the first case, no liquid was injected, while in the second and third cases, different volumes of water per second were injected. To obtain a double liquid rate, the number of drops was kept constant, but the drop diameter was increased from 0.2 mm to 0.252 mm. In the fourth case, an aqueous binder solution, water with 12 mass percent of PHARMACOAT^® 606 (Hydroxypropyl Methylcellulose, Shin-Etsu Chemical Co., Ltd., Chiyoda-ku, Japan), was used, which has a dynamic viscosity more than 600 times higher than water. For the generation of the droplets, a virtual factory was created at the end of the nozzle. The droplets had an initial velocity of 35 m/s in the y-direction (Figure 3), and the spray angle was 60° (Figure 3b). The properties of the injected coating liquid were also measured. The surface tension was determined with a drop contour analysis system (OCA 15EC, DataPhysics Instruments, Filderstadt, Germany), the density of the liquid with a density meter (DMA 55, PhysicaMesstechnik GmbH, Ostfildern, Germany), and the dynamic viscosity with a rotational viscometer (RHEOTEST^® RN 4.1, RHEOTEST Medingen GmbH, Medingen, Germany). In Table 4, the parameters of the droplets for the different simulation cases are given. For the drying of the particles, the assumption is made that a layer of 20 μm dries in 300 s, resulting in a drying rate of 0.0015 mg/s.

4. Results 4.1. Influence of Liquid Injection on Particle Dynamics

The particle dynamics in the fluidized bed can change during a coating process. The sprayed liquid droplets are deposited on particles leading to an increased particle volume and generate the additional forces in particle contacts, such as capillary and viscous forces. The impact of these interactions depends on the amount and viscosity of the sprayed binder solution, drying kinetics, and the mixing behavior. To understand the characteristic particle dynamics in the rotor granulator, the particles colored by velocity magnitude, are shown in Figure 4 for the dry case 60 s after the process started. One can see that the fastest particles are located on the plate, and with increasing height, the velocity of the particles decreases. The rotating plate in the rotor granulator accelerates the particles in the toroidal direction upon contact. Moreover, the centrifugal force causes the particles to move outwards towards the wall of the apparatus. Since the gap with the incoming gas is also located there, the particles are additionally accelerated vertically upwards here. With increasing flight height, the particles slow down more and more due to the particle or wall impact, friction, and gravitational force. Subsequently, they are pushed towards the center by new particles accelerated by the upward flowing gas in the gap. Then, they land back on the plate, where the particles again move towards the wall due to the prevailing centrifugal forces.

In order to analyze the effect of the liquid injection on the particle dynamic, one-minute process time was simulated for all cases and the results were compared. The injection started after simulation of two seconds of the dry process, because than the steady state of the particle movement is reached. First, the dry case is compared with the wet case II, where 400,000 water droplets per second with a diameter of 200 μm are injected. During this time, a volume of 96.7 mL water was added to the system, corresponding to a liquid mass of 10.4% of the total mass of the particles. In Figure 5, the time-averaged solid volume fraction for the dry and wet case can be seen. For a better overview, and since all particles are close to the apparatus wall due to centrifugal forces, only the radial positions between 80 mm and 150 mm from the center are shown. The values are averaged over the last 5 s of the simulations. The highest solid concentration could be observed in a zone above the plate near to the gap for both cases (dashed lines in Figure 5). However, in the wet case, the concentration remains high, even directly above the gap. The attracting forces due to the liquid bridges cause the particles to stick together, which results in a higher packing density. In the dry case, particles are frequently localized even at a height of more than 40 mm near the wall, in contrast to the wet case, where the particle bed is more homogeneously distributed and particles do not reach a height above 40 mm. One can see that the particle bed is lower and denser in the wet case.

While the average contact force of the particles stays constant over time for the dry case, it increases when the injection starts (Figure 6a). After around 5 s of liquid injection, the mean value of the particle contact forces oscillates around a constant average value, which is 21.8% higher than in the dry case. Comparing the average velocities of the particles, one can see that the average particle velocity is also higher in the wet case (Figure 6b). Due to the additional adhesive forces caused by the liquid, the particles are in longer contact with the fast-rotating rotor plate, and the average velocity increases at the beginning of the injection, and stabilizes around an 8.7% higher average value.

During the coating process, the amount of water in the system increases linearly because of the constant injection volume flow and the constant drying rate. In Figure 7a, the constant increase of the total liquid volume of all particles in the process over time is shown. After 58 s of injection, 91 mL of liquid remains distributed on the particles, corresponding to 94.1% of the injected liquid. A total of 5.6 mL is on the walls of the apparatus. The low drying rate results in very low liquid reduction over time. Due to the wetting of the particles, liquid bridges are formed in the wet contacts. The number of formed liquid bridges increases strongly in the first three seconds. The rupture of the liquid bridges starts later when more volume of the liquid is deposited on the particles (Equation (15)). Over time, a balance between newly formed and destroyed bridges is reached, with an average number of 1581 of liquid bridges (Figure 7b). Agglomerates are formed from maximum six particles, but on average with three to four.

Due to the capillary and viscous forces caused by the liquid bridge, additional adhesive forces act on the particles, and as a consequence, the velocity of the particles also changes, accordingly. To analyze the effect of the liquid in more detail, Figure 8 compares the velocity distribution of the particles in the dry and wet cases.

In both cases, two peaks can be seen, one at a velocity of 0.2 m/s and the second at 0.6 m/s. In the dry simulation, over 30.1% of the particles have a velocity under 0.2 m/s, whereas in the wet case, it is only around 16.0% of the particles. The particles with a low velocity are the particles which have direct contact with the wall. These particles are slowed down by friction with the stationary wall. Most of the particles in the bed move at a velocity of 0.6 m/s. The highest velocity reaches the particles after contact with the rotor plate due to the high rotational speed of the plate. In the dry case, the highest peak is the first one, whereas in the wet case the second is the highest. Moreover, the distribution around the velocity of 0.6 m/s is more homogeneous in the wet case than in the dry case. While in the dry simulation there are no more particles with a velocity faster than 1.2 m/s, over 4.4% of the particles in the wet case have such a velocity. The additional adhesive forces caused by the capillary bridge result in a longer contact duration. The average contact duration increases from 0.058 ms in the dry case I to 0.127 ms in the wet case II (Table 5). As the contact duration of the particles becomes longer, the time interval in which kinetic energy is transferred becomes greater as well and more energy can be transferred from the rotating plate as well as faster particles to the slower ones.

In Figure 9, the distributions of the liquid thickness on the particles for the second case after 15 s, 30 s, 45 s, and 60 s of process time are shown. After the first 15 s, all particles are wetted, but some get liquid only by a transfer during wet contact with another wet particle or another wet wall and not directly by a droplet. After 30 s, all particles had at least one direct collision with a droplet and the average layer thickness was 86.3 μm. Over time, the average liquid layer on the particles increases continuously. After 60 s of simulation, corresponding to 58 s of liquid injection, an average liquid layer of 178.7 μm is reached.

4.2. Influence of the Water Spray Rate and Coating Solution In order to investigate the influence of the injected liquid rate on the particle behavior, in the simulation case III, the droplets were generated with a volume twice as large as in the first wet case (case II). In addition, a widely used coating substance, Pharmacoat, was simulated in the case IV, and the obtained particle behavior was compared with the other cases.

In Figure 10, the average contact force of the particle (Figure 10a) and the average particle velocity (Figure 10b) for the double water injection, as well as for the injection of Pharmacoat solution, are shown. Analogous to the results of the first wet case II (Figure 6), the forces and the velocity increase in the first five seconds after the injection started and then stabilized. Table 5 also shows the mean values of the force and velocity in the stabilized regimes. While the liquid rate does not seem to influence the average contact force and the average particle velocity, it changes for the simulation with Pharmacoat. Analyzing the average contact forces in the case of the Pharmacoat solution, one can see that they increase strongly in the first 10 s and are about twice as high as for the cases with the water injection. Due to the high dynamic viscosity of the Pharmacoat solution, the resulting viscous force in the contact increases (Equations (17) and (18)). However, after a certain amount of Pharmacoat solution in the process, the average contact forces then decrease again and are only 12% higher at the end of the simulation. Similar to the contact forces, the particles have a higher average velocity at the beginning of the liquid injection than in the cases with water injection. Then, the velocity decreases over time and is 6.7% lower at the end compared to case II. The increasing amount of Pharmacoat solution in the process results in higher energy dissipation during contact. The capillary number could be used here to describe the influence of the viscous forces on energy dissipation. The relative particle collision velocity of 0.155 m/s corresponds to the capillary number of around 0.002 for the case II with water injection. Whereas the dynamic viscosity of the solution with Pharmacoat is over 600 times higher than that of water, resulting in a capillary number of 1.3, which indicates a dominant influence of viscous forces.

Comparing the number of liquid bridges over time of the wet case III with the wet case II, no influence of the injected liquid rate could be found (Table 5). The average number of liquid bridges stabilizes even in the case of the double injection at around 1600 bridges. With Pharmacoat, about 22.7% less liquid bridges (in average 1222 liquid bridges) are formed than in case II.

In the wet case III, the average contact duration increases by about 23% (Table 5). The reason is the higher liquid volume on the particles. Since the rupture distance of the liquid bridge increases with increasing liquid volume, the liquid bridge is more stable over a longer distance between the particles after the rebound. As a result, the time until the liquid bridge is destroyed is also longer and, as a consequence, the contact duration. In the case of the simulated Pharmacoat solution, the average contact time increases slightly compared to case II. The higher acting viscous forces during contact result in greater energy dissipation, which prolongs the contact time. The fact that the difference is not greater is due to the dependence of the capillary force on the surface tension (Equations (10) and (11)). The surface tension of Pharmacoat is only slightly more than half the value of water.

No significant differences in the velocity distribution can be seen between single injection and double injection (Figure 11a compare to Figure 8b). The distribution in the case of the Pharmacoat solution differs more significantly from the wet case II (Figure 11b compared to Figure 8b). Particles that repeatedly stick to the wall for a short time and, therefore, have a velocity of 0 m/s are present. In contrast, there are particles, which are significantly faster (more than 1.6 m/s). These are particles that stick to the rotor, whereby high kinetic energy is transferred from the rotor plate. The particles approach the circumferential speed of the rotor, which has a maximum speed of 4.2 m/s.

In Figure 12, the average growth rate of the dry layer on the particles during the simulated process with the injection of the Pharmacoat solution is shown. As soon as a particle has a liquid volume of the Pharmacoat solution, the liquid begins to dry and a dry layer is formed. For this reason, there is no dry layer in the first few seconds. Since the drying rate is constant, the average dry layer also increases linearly over time. After 60 s, an average dry layer thickness of 0.05 μm on the particles is reached. It must be taken into account that the assumed drying rate at 26 °C is quite low. Furthermore, the coating process usually takes much longer than one minute to achieve the desired thickness of the dry layer.

5. Conclusions A coating process in a rotor granulator was simulated with CFD-DEM considering capillary and viscous forces between particles wetted due to the deposition of the liquid droplets or liquid transfer during collision with wet particles or wall. A new coating model was implemented to take into account a volume exchange of liquid during contact, as well as the drying of the liquid layer on the particle. In order to take the liquid injection into account, a single-fluid nozzle is also simulated. The droplets were generated as a second particle phase in DEM. The effect of the liquid injection on the particle behavior was analyzed with four different cases: one dry case, two cases with different rates of water injection, and one case of a solution of six mass percent of Pharmacoat 606 in water. The liquid on the particles led to an increase of the average contact forces of the particle and average particle velocity because of a better momentum transfer from the rotor plate to the wet particles, which stick significantly longer to the rotating plate than in the dry case. Moreover, the average contact duration increases considerably from the dry case to the wet cases. A double water injection led to more liquid volume on the particles, resulting in a greater rupture distance and, therefore, in an increase in the contact duration. In the first 40 s of the Pharmacoat injection, the average contact forces and the average velocity of the particles are higher than in the cases with water injection. After that, the average contact force and average velocity decreases over time and are lower at the end. The velocity distribution in the case of the Pharmacoat solution differs significantly from the wet cases since the particles adhere more often to the wall or the rotor, which results in a wider distribution. The average contact time increases slightly compared to the first wet case due to the higher energy dissipation by greater viscous forces caused by the higher dynamic viscosity of the Pharmacoat solution.

Figure 1. Liquid transfer between two colliding particles, where both particles have a liquid volume.

Figure 2. Schematic of the coating model: Droplet hits a particle (a) and droplet hits a cell of the wall grid (b).

Figure 3. Geometry of the simulated rotor granulator: Side view (a) and top view (b).

Figure 4. Particles colored by velocity magnitude for the dry case after 60 s.

Figure 5. Time-averaged solid volume fraction for the dry case (a) and for the wet case (b).

Figure 6. Average contact force acting on the particle (a) and average particle velocity (b) for the dry and wet case.

Figure 7. Total liquid volume of the particles over time (a) and number of liquid bridges over time (b).

Figure 8. Velocity distribution of the particles for the dry case I (a) and in the wet case II (b).

Figure 9. Cumulative distribution of the liquid layer thickness on the particles during liquid injection (case II).

Figure 10. Average particle velocity (a) and average contact force acting on the particle (b) for the double water injection case III, and the Pharmacoat 606 case IV.

Figure 11. Velocity distribution of the particles for the double water injection (case III) (a) and for the injection with the Pharmacoat 606 solution (case IV) (b).

Figure 12. Average dry layer of the particles coated with the Pharmacoat 606 solution.

Particle Reynold Number	Regime	Drag Coefficient
Rep≤0.25	Laminar	CD=24/Re
0.25<Rep≤1000	Transition	CD=24·(1+0.15·Re0.687)Re
Rep>1000	Turbulent	CD=0.44

Parameters of the Gas	Value
Number of grid cells (-)	2,743,200
Time (-)	steady
Turbulence model (-)	k-ε
Inlet velocity (m/s)	1
Outlet pressure (bar)	1.01325
Temperature (°C)	26
Rotation plate (rpm)	300
Specific heat (kJ/(kg·K))	1.006
Thermal conductivity (W/(m·K))	0.0242
Dynamic viscosity (kg/(m·s))	1.7894∙10−5
Gas density (-)	Ideal gas

Parameters of the Particle	Value
Number of particles (-)	8139
Bed mass (kg)	0.923
Diameter (mm)	4.2
Density (kg/m3)	2922
Young’s modulus (GPa)	1.68
Poisson’s ratio (-)	0.25
Static friction coefficient (-)	0.456
Rolling friction coefficient (-)	0.065
Restitution coefficient (-)	0.73

Parameters of the Droplets	Water Injection Case II	Water Double Injection Case III	Water with Pharmacoat 606 Case IV
Droplet generation rate (1/s)	400,000	400,000	400,000
Volume flow (mL/min)	100	200	100
Diameter (mm)	0.2	0.252	0.2
Density (kg/m3)	1000	1000	1026
Surface tension (mN/m)	75.75	75.75	42.91
Dynamic viscosity (mPas)	1.01	1.01	620.43
Wetting angle (°)	28.65	28.65	29.12

	Dry Case I	Water Injection Case II	Water Double Injection Case III	Water with Pharmacoat 606 Case IV
Contact force (N)	0.083 ± 0.007	0.101 ± 0.008	0.099 ± 0.006	0.148 ± 0.027
Particle velocity (m/s)	0.596 ± 0.003	0.647 ± 0.003	0.645 ± 0.003	0.647 ± 0.021
Contact time (ms)	0.058 ± 0.045	0.127 ± 0.007	0.165 ± 0.004	0.129 ± 0.005
Number of liquid bridges (-)	-	1581 ± 85	1605 ± 91	1222 ± 82

Author Contributions

Conceptualization, P.G.; Methodology, P.G.; Software, M.L. and P.G.; Validation, P.G.; Formal Analysis, P.G.; Investigation, P.G.; Particle production, T.O.; Measurement of Particle Properties: P.G. and S.A.; Data Curation, M.L. and P.G.; Writing-Original Draft Preparation, P.G.; Writing-Review & Editing, P.G.; Visualization, P.G.; Supervision, S.H. and S.A.; Project Administration, P.G.; Funding Acquisition, S.H. and S.A. All authors have read and agreed to the published version of the manuscript.

Funding

This study was conducted within the research project, which is funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation), Project-ID-AN 782/12-1, what the authors gratefully acknowledge.

Conflicts of Interest

The authors declare no conflict of interest.

References

1. Heinrich, S.; Dosta, M.; Antonyuk, S. Multiscale Analysis of a Coating Process in a Wurster Fluidized Bed Apparatus. Mesoscale Modeling in Chemical Engineering Part I; Elsevier: Amsterdam, The Netherlands, 2015; pp. 83-135. ISBN 9780128012475.

2. Grohn, P.; Weis, D.; Thommes, M.; Heinrich, S.; Antonyuk, S. Contact Behavior of Microcrystalline Cellulose Pellets Depending on their Water Content. Chem. Eng. Technol. 2020, 43, 887-895.

3. Neuwirth, J.; Antonyuk, S.; Heinrich, S.; Jacob, M. CFD-DEM study and direct measurement of the granular flow in a rotor granulator. Chem. Eng. Sci. 2013, 86, 151-163.

4. Korakianiti, E.S.; Rekkas, D.M.; Dallas, P.P.; Choulis, N.H. Sequential optimization of a pelletization process in a fluid bed rotor granulator. J. Drug Deliv. Sci. Technol. 2004, 14, 207-214.

5. Langner, M.; Kitzmann, I.; Ruppert, A.-L.; Wittich, I.; Wolf, B. In-Line particle size measurement and process influences on rotary fluidized bed agglomeration. Powder Technol. 2020, 364, 673-679.

6. Gidaspow, D. Multiphase Flow and Fluidization. Continuum and Kinetic Theory Descriptions; Elsevier Science: San Diego, CA, USA, 1994; ISBN 9780122824708.

7. Yang, N.; Wang, W.; Ge, W.; Li, J. CFD simulation of concurrent-up gas-solid flow in circulating fluidized beds with structure-dependent drag coefficient. Chem. Eng. J. 2003, 96, 71-80.

8. Armstrong, L.M.; Gu, S.; Luo, K.H. Study of wall-to-bed heat transfer in a bubbling fluidised bed using the kinetic theory of granular flow. Int. J. Heat Mass Transf. 2010, 53, 4949-4959.

9. Loha, C.; Chattopadhyay, H.; Chatterjee, P.K. Effect of coefficient of restitution in Euler-Euler CFD simulation of fluidized-bed hydrodynamics. Particuology 2014, 15, 170-177.

10. Adamczyk, W.P.; Klimanek, A.; Białecki, R.A.; Węcel, G.; Kozołub, P.; Czakiert, T. Comparison of the standard Euler-Euler and hybrid Euler-Lagrange approaches for modeling particle transport in a pilot-scale circulating fluidized bed. Particuology 2014, 15, 129-137.

11. Deen, N.G.; van Sint Annaland, M.; Van Der Hoef, M.A.; Kuipers, J.A.M. Review of discrete particle modeling of fluidized beds. Chem. Eng. Sci. 2007, 62, 28-44.

12. Salikov, V.; Heinrich, S.; Antonyuk, S.; Sutkar, V.S.; Deen, N.G.; Kuipers, J.A.M. Investigations on the spouting stability in a prismatic spouted bed and apparatus optimization. Adv. Powder Technol. 2015, 26, 718-733.

13. Breuninger, P.; Weis, D.; Behrendt, I.; Grohn, P.; Krull, F.; Antonyuk, S. CFD-DEM simulation of fine particles in a spouted bed apparatus with a Wurster tube. Particuology 2019, 42, 114-125.

14. Neuwirth, J.; Antonyuk, S.; Heinrich, S. Particle dynamics in the fluidized bed: Magnetic particle tracking and discrete particle modelling. In Proceedings of the Powders and Grains 2013 7th International Conference on Micromechanics of Granular Media (AIP 2013), Sydney, Australia, 8-12 July 2013; pp. 1098-1101.

15. Neuwirth, J. Charakterisierung und Diskrete-Partikel-Modellierung des Strömungs- und Dispersionsverhaltens im Rotorgranulator, 1st ed.; Cuvillier Verlag: Göttingen, Germany, 2017; ISBN 9783736994768.

16. Fries, L.; Dosta, M.; Antonyuk, S.; Heinrich, S.; Palzer, S. Moisture Distribution in Fluidized Beds with Liquid Injection. Chem. Eng. Technol. 2011, 34, 1076-1084.

17. Suzzi, D.; Toschkoff, G.; Radl, S.; Machold, D.; Fraser, S.D.; Glasser, B.J.; Khinast, J.G. DEM simulation of continuous tablet coating: Effects of tablet shape and fill level on inter-tablet coating variability. Chem. Eng. Sci. 2012, 69, 107-121.

18. Sutkar, V.S.; Deen, N.G.; Patil, A.V.; Salikov, V.; Antonyuk, S.; Heinrich, S.; Kuipers, J.A.M. CFD-DEM model for coupled heat and mass transfer in a spout fluidized bed with liquid injection. Chem. Eng. J. 2016, 288, 185-197.

19. Börner, M.; Bück, A.; Tsotsas, E. DEM-CFD investigation of particle residence time distribution in top-spray fluidised bed granulation. Chem. Eng. Sci. 2017, 161, 187-197.

20. Heine, M.; Antonyuk, S.; Fries, L.; Niederreiter, G.; Heinrich, S.; Palzer, S. Modeling of the Spray Zone for Particle Wetting in a Fluidized Bed. Chem. Ing. Tech. 2013, 85, 280-289.

21. Dosta, M.; Antonyuk, S.; Hartge, E.-U.; Heinrich, S.; Dosta, M. Parameter Estimation for the Flowsheet Simulation of Solids Processes. Chem. Ing. Tech. 2014, 86, 1073-1079.

22. Antonyuk, S.; Heinrich, S.; Smirnova, I. Discrete Element Study of Aerogel Particle Dynamics in a Spouted Bed Apparatus. Chem. Eng. Technol. 2012, 35, 1427-1434.

23. Lian, G.; Thornton, C.; Adams, M.J. Discrete particle simulation of agglomerate impact coalescence. Chem. Eng. Sci. 1998, 53, 3381-3391.

24. Contact Mechanics and Friction. Physical Principles and Applications, 1st ed.; Popov, V.L. (Ed.) Springer: Berlin, Germany, 2010; ISBN 978-3-642-10802-0.

25. Patankar, S.V. Numerical Heat Transfer and Fluid Flow; CRC Press: Boca Raton, FL, USA, 2018; ISBN 9781315275130.

26. Cundall, P.A.; Strack, O.D.L. A discrete numerical model for granular assemblies. Géotechnique 1979, 29, 47-65.

27. Oertel, H.; Laurien, E. Numerische Strömungsmechanik; Vieweg + Teubner Verlag: Wiesbaden, Germany, 2003; ISBN 978-3-528-03936-3.

28. Schiller, L.; Naumann, A. A drag coefficient correlation. Z. Ver. Dtsch. Ing. 1935, 77, 318-320.

29. Hertz, H. Über die Berührung fester elastischer Körper. J. Rein. Angew. Math. 1881, 92, 156-171.

30. Tsuji, Y.; Tanaka, T.; Ishida, T. Lagrangian numerical simulation of plug flow of cohesionless particles in a horizontal pipe. Powder Technol. 1992, 71, 239-250.

31. Mindlin, R.D.; Deresiewicz, H. Elastic Spheres in Contact under Varying Oblique Force. J. Appl. Mech. 1953, 20, 327-344.

32. Antonyuk, S.; Heinrich, S.; Tomas, J.; Deen, N.G.; Van Buijtenen, M.S.; Kuipers, J.A.M. Energy absorption during compression and impact of dry elastic-plastic spherical granules. Granul. Matter 2010, 12, 15-47.

33. Mikami, T.; Kamiya, H.; Horio, M. Numerical simulation of cohesive powder behavior in a fluidized bed. Chem. Eng. Sci. 1998, 53, 1927-1940.

34. Muguruma, Y.; Tanaka, T.; Tsuji, Y. Numerical simulation of particulate flow with liquid bridge between particles (simulation of centrifugal tumbling granulator). Powder Technol. 2000, 109, 49-57.

35. Rabinovich, Y.I.; Esayanur, M.S.; Moudgil, B.M. Capillary Forces between Two Spheres with a Fixed Volume Liquid Bridge: Theory and Experiment. Langmuir 2005, 21, 10992-10997.

36. Israelachvili, J.N. Intermolecular and Surface Forces; Elsevier: Amsterdam, The Netherlands, 2011; ISBN 9780123919274.

37. Fujihashi, D.; Tsunazawa, Y.; Tokoro, C.; Sakai, M. DEM Simulation of Particle Behavior in Pan-Type Pelletizer Considering the Effect of the Capillary Force. J. Soc. Powder Technol. Jpn. 2014, 51, 828-836.

38. Willett, C.D.; Adams, M.J.; Johnson, S.A.; Seville, J.P.K. Capillary Bridges between Two Spherical Bodies. Langmuir 2000, 16, 9396-9405.

39. Tsunazawa, Y.; Fujihashi, D.; Fukui, S.; Sakai, M.; Tokoro, C. Contact force model including the liquid-bridge force for wet-particle simulation using the discrete element method. Adv. Powder Technol. 2016, 27, 652-660.

40. Lian, G.; Thornton, C.; Adams, M.J. A Theoretical Study of the Liquid Bridge Forces between Two Rigid Spherical Bodies. J. Colloid Interface Sci. 1993, 161, 138-147.

41. Reynolds, O., IV. On the theory of lubrication and its application to Mr. Beauchamp tower's experiments, including an experimental determination of the viscosity of olive oil. Philos. Trans. R. Soc. 1886, 177, 157-234.

42. Adams, M.J.; Perchard, V. The cohesive forces between particles with interstitial liquid. Inst. Chem. Eng. Symp. 1985, 91, 147-160.

43. Nase, S.T.; Vargas, W.L.; Abatan, A.A.; McCarthy, J. Discrete characterization tools for cohesive granular material. Powder Technol. 2001, 116, 214-223.

44. Anand, A.; Curtis, J.S.; Wassgren, C.; Hancock, B.C.; Ketterhagen, W.R. Predicting discharge dynamics of wet cohesive particles from a rectangular hopper using the discrete element method (DEM). Chem. Eng. Sci. 2009, 64, 5268-5275.

45. Washino, K.; Miyazaki, K.; Tsuji, T.; Tanaka, T. A new contact liquid dispersion model for discrete particle simulation. Chem. Eng. Res. Des. 2016, 110, 123-130.

46. Schmelzle, S.; Nirschl, H. DEM simulations: Mixing of dry and wet granular material with different contact angles. Granul. Matter 2018, 20, 19.

47. Shi, D.; McCarthy, J. Numerical simulation of liquid transfer between particles. Powder Technol. 2008, 184, 64-75.