Analysis and optimization of material flow inside the system of rotary coolers and intake pipeline via discrete element method modelling

energies Article Analysis and Optimization of Material Flow inside the System of Rotary Coolers and Intake Pipeline via Discrete Element Method Modelling Jakub Hlosta 1,*ID , David Žurovec 1, JiˇríRozbroj 1,Álvaro Ramírez-Gómez 2, Jan Neˇcas 1and JiˇríZegzulka 1 1ENET Centre, Bulk Solids Centre-Czech Republic, VŠB-Technical University of Ostrava, 17. listopadu 15/2172, 70833 Ostrava-Poruba, Czech Republic; david.zur[email protected] (D.Z.); jiri.rozbr[email protected] (J.R.); [email protected] (J.N.); [email protected] (J.Z.) 2Department of Mechanical Engineering, Chemistry and Industrial Design, Technical University of Madrid, Ronda de Valencia 3, 28012 Madrid, Spain; alvaro.ramir[email protected] *Correspondence: [email protected]; Tel.: +420-597-329-371 Received: 7 June 2018; Accepted: 12 July 2018; Published: 14 July 2018   Abstract: There is hardly any industry that does not use transport, storage, and processing of particulate solids in its production process. In the past, all device designs were based on empirical relationships or the designer’s experience. In the field of particulate solids, however, the discrete element method (DEM) has been increasingly used in recent years. This study shows how this simulation tool can be used in practice. More specifically, in dealing with operating problems with a rotary cooler which ensures the transport and cooling of the hot fly ash generated by combustion in fluidized bed boilers. For the given operating conditions, an analysis of the current cooling design was carried out, consisting of a non-standard intake pipeline, which divides and supplies the material to two rotary coolers. The study revealed shortcomings in both the pipeline design and the cooler design. The material was unevenly dispensed between the two coolers, which combined with the limited transport capacity of the coolers, led to overflowing and congestion of the whole system. Therefore, after visualization of the material flow and export of the necessary data using DEM design measures to mitigate these unwanted phenomena were carried out. Keywords: discrete element method (DEM) modelling; design; equipment optimization; flow prediction; rotary cooler 1. Introduction Market pressure forces companies to develop increasingly optimized products with regard to transport efficiency, environmental impacts, and material and labor savings. Nowadays, to propose innovative processes, industries require the use of sophisticated simulation techniques. Some of these are based in the discrete element method (DEM), which use virtual models to simulate the dynamic flow of particulate solids [ 1 ] in transport, handling, storage, and process systems [ 2 ]. At present, DEM is drawing the attention of many designers and manufacturers of equipment as well as process engineers. It has become a useful tool to analyze industrial processes [ 3 ]. It makes it possible to test the basic functionality of devices [ 4 ] and their structural configurations before prototyping or to optimize existing equipment, resulting in considerable economic, time, labor, and material savings. One of the challenges of modeling the processes of particulate solids is the creation of a mathematical model for the bulk material itself, whose behavior is influenced by a great number of factors. The behavior of powders and bulk materials depends on external and internal factors that influence their bulk Energies 2018,11, 1849; doi:10.3390/en11071849 www.mdpi.com/journal/energies Energies 2018,11, 1849 2 of 18 density, internal and external angle of friction, angle of repose, cohesiveness, particle morphology and granulometry, flowability, and many other parameters [5]. Since the current requirements go beyond the simple transportation of materials from A to B but rather require a combination of the transportation with processes (i.e., drying, cooling, sorting, etc.), DEM opens up new possibilities in design and process optimization avoiding design errors before implementation. Recently, there has been a tendency to link the DEM method to other numerical simulation methods such ascomputational fluid dynamics (CFD) [ 6 ], finite element method (FEM) [ 7 ] or multi-body dynamics (MBD) to solve more complex problems [ 8 ]. Development and improvement of any structure is a time-consuming effort, but it is essential for business competitiveness and product sales. The using of numerical methods for study and optimization is very common and desirable today. These numerical techniques extend the possibilities of optimizing the construction elements of the conveyors as well as in the area of critical loads (stress) studies on the transport equipment [ 9 ]. In case of material flow through pneumatic transport systems, numerical methods are applied in the areas of new applications and designs of complex multiphase flow situations [ 10 ]. New theoretical and numerical approaches can be used to predict and evaluate possible safety risks during transport [ 11 ]. In the field of storage, it is possible to use numerical methods for solving complex tasks of flow of particulates from silos with respect to acting wall pressures, speed profiles, and particle shape [ 12 ]. It can be very challenging for some companies to design the correct pipeline for their equipment. This research motivates industrial practitioners and engineers to adopt DEM as an important tool for designing and troubleshooting of equipment. There are many processes and handling operations with particulates in the energy industry. The DEM can be used in wide scale of applications, i.e., coal and fly ash handling, grinding processes, storage bins design, etc. Bai and He (2015) published a study on the effect of coal loads on kinetic energy of balls for ball mills where the kinetic energy of grinding balls was optimized [13]. 2. Theoretical Background The basic principle of the DEM modeling procedure is based on a particle contact model, which is defined by parameters such as particle diameter, Poisson’s ratio, interaction coefficients and density. The particle simulation model defined in this way is calibrated based on experimental data. Finally, the data from the experimental and simulation model of the application (e.g., conveyor) is verified. The diagram of the modeling procedure is shown in Figure 1. The basic material optimization with Hertz-Mindlin contact model through the angle of repose has been performed in this study. In the next step, simulation with the original state of coolers was performed. This simulation corresponded with the behavior of the real equipment. Therefore, some geometric adjustments has been designed and tested by simulation again. Mass flow rates and particle trajectories in different parts of the cooling system were understood thanks to this study. The goal of this study is to present the principle and optimization process via DEM modelling in industrial scale equipment, where there are tonnes of material transported or handled per hour. Originally published by Cundall and Strack in 1979, the discrete element method has become over time a commonly used numerical tool for particle modeling [ 14 ]. In DEM models, the particulate material is modeled as an assembly of particles interacting with each other or with parts of the equipment or facilities where present (moving paddles of cooling conveyors, walls of hoppers, silos, and pipelines, etc.). The macroscopic behavior of particles is due to the microscopic interactions between the particles themselves and between the particles and boundary conditions. The determination of trajectory and velocity of particles is performed in discrete time steps. Much of the information about this behavior, such as collision frequency or interaction forces, is based on these calculations. Theoretical particle movement through other particles shows flow mechanisms at the scale and level of detail that would be otherwise very difficult to achieve by experimental means. Energies 2018,11, 1849 3 of 18 Energies 2018, 11, x FOR PEER REVIEW 3 of 19 Figure 1. General principle of discrete element method (DEM) modeling of every single process simulation. Originally published by Cundall and Strack in 1979, the discrete element method has become over time a commonly used numerical tool for particle modeling [14]. In DEM models, the particulate material is modeled as an assembly of particles interacting with each other or with parts of the equipment or facilities where present (moving paddles of cooling conveyors, walls of hoppers, silos, and pipelines, etc.). The macroscopic behavior of particles is due to the microscopic interactions between the particles themselves and between the particles and boundary conditions. The determination of trajectory and velocity of particles is performed in discrete time steps. Much of the information about this behavior, such as collision frequency or interaction forces, is based on these calculations. Theoretical particle movement through other particles shows flow mechanisms at the scale and level of detail that would be otherwise very difficult to achieve by experimental means. Reliability and the nature of discrete models depend entirely on simplifying many physical properties of particulate solids. These properties need to be described in DEM using several simple interaction parameters. Such simplifications are necessary and are widely used to solve complex problems within reasonable timeframes. In the literature, there has been little work done in verifying the flow of particulate solids in the global scale of the entire technology. If DEM is to become an equally important design tool, such as FEM or CFD numerical modeling, it is absolutely necessary to verify the capabilities of DEM by simulations, i.e., to provide insight into the interaction mechanisms in installations where it is difficult to monitor flow, obtain measurements or quantify its performance [15]. A number of problems with gravity tables are solved within gravity separation processes, e.g., a slope of the table or air intake [16]. Additive technologies address geometric optimization to improve the quality of the applied layer [17]. In the sorting areas, the optimal setting of the operating device is solved [18]. Mixing processes need to be investigated and compared with the properties of mixed materials [19]. It is necessary and possible to optimize new production technologies and equipment [20] in industry and engineering with regard to service life [21]. For example, there are problems with the flow of particulate solids when storing materials in silos [22] or during pneumatic transportation in pipelines [23]. In the pharmaceutical and chemical industry, material segregation or inappropriate distribution of the active substance may occur during processing. In all these industries, DEM can provide significant information for problem-solving [24]. The study of bulk materials can be carried out either as a continuum or as a particle conglomerate. As for the latter approach, a great deal of work has been published in recent years on the issue of physical interactions of individual particles [25]. A number of studies of micro-mechanical interactions between particles have been carried out, such as the stress state inside the material or the distribution of internal forces within the assembly. However, due to the lack of sufficiently sensitive Figure 1. General principle of discrete element method (DEM) modeling of every single process simulation. Reliability and the nature of discrete models depend entirely on simplifying many physical properties of particulate solids. These properties need to be described in DEM using several simple interaction parameters. Such simplifications are necessary and are widely used to solve complex problems within reasonable timeframes. In the literature, there has been little work done in verifying the flow of particulate solids in the global scale of the entire technology. If DEM is to become an equally important design tool, such as FEM or CFD numerical modeling, it is absolutely necessary to verify the capabilities of DEM by simulations, i.e., to provide insight into the interaction mechanisms in installations where it is difficult to monitor flow, obtain measurements or quantify its performance [ 15 ]. A number of problems with gravity tables are solved within gravity separation processes, e.g., a slope of the table or air intake [ 16 ]. Additive technologies address geometric optimization to improve the quality of the applied layer [ 17 ]. In the sorting areas, the optimal setting of the operating device is solved [ 18 ]. Mixing processes need to be investigated and compared with the properties of mixed materials [ 19 ]. It is necessary and possible to optimize new production technologies and equipment [ 20 ] in industry and engineering with regard to service life [ 21 ]. For example, there are problems with the flow of particulate solids when storing materials in silos [22] or during pneumatic transportation in pipelines [ 23 ]. In the pharmaceutical and chemical industry, material segregation or inappropriate distribution of the active substance may occur during processing. In all these industries, DEM can provide significant information for problem-solving [ 24 ]. The study of bulk materials can be carried out either as a continuum or as a particle conglomerate. As for the latter approach, a great deal of work has been published in recent years on the issue of physical interactions of individual particles [ 25 ]. A number of studies of micro-mechanical interactions between particles have been carried out, such as the stress state inside the material or the distribution of internal forces within the assembly. However, due to the lack of sufficiently sensitive measuring experimental equipment, mathematical models or computational simulations, approaches to the extensive problems of bulk mechanics have not been explored or extended [26]. As the world recedes gradually out off coal as a source of energy, cleaner and more competitive energy sources, such as renewable energy sources, will emerge. From biomass and waste, 13% of the world’s total energy consumption is projected to be from coal in 2040 [ 27 ]. Coal has experienced the greatest development as an energy source in 2000–2012. In addition, the use of oil and natural gas has experienced growth, followed by hydropower and renewable energy sources. The use of nuclear energy decreased in recent years because of nuclear disasters risk. Therefore, it is necessary to construct all energy devices with maximum reliability in terms of their number and importance. This study presents innovative ways of using the DEM method in optimizing a rotary cooler. Such cooler is Energies 2018,11, 1849 4 of 18 installed under a fluidized bed boiler of a power plant and a heating plant and serves for cooling and transporting of hot fly ash from the boiler to waste disposal. A comprehensive solution of material flow prediction and similar technological (cooling) processes includes the determination of flow, mechanical, physical, granulometric, and thermal properties of the conveyed particulate solids [ 28 ]. Various methods of determining the thermal conductivity of powders [ 29 ], which can be generally used in designing coolers, in calculation of the heat exchange surface, and the material retention time in the cooler, especially in the combination of CFD–DEM numerical methods [ 30 – 32 ], can be found in the literature. This study focuses on the shape optimization of the intake pipeline and the active parts of the cooling system. To solve the malfunction of one of two devices that are identical, an optimized solution has been found to ensure continuous operation of the device via DEM modeling in a comparative way. The DEM algorithm generally starts with placing the particles in the system. For the calculation of contact force, EDEM performs detection of contacts to determine whether the individual particles interact with each other in the current time step. Contact detection is time-consuming and requires substantial computing power. After contact detection is completed, the total force applied to each particle is determined and therefore the translational, as well as rotational motion, can be determined by integrating the Newton Equation (1) and its equivalent (2). Between each time step of the calculation process, the particles are moved based on the speed and acceleration calculated in the previous time step. Contacts are then used to calculate the reaction forces acting on each particle. This force determines its speed and acceleration in the next time interval. The resulting force and torque represent the sum of all the forces and torques applied to each particle, including gravity, fluid flow, magnetic fields, electrostatic fields, external forces and so on. In this study forces and torques caused by gravity, contacts, static, and rolling friction are considered. Therefore, the Newton equations of motion of the particle i in contact with the particle j the following form: midvi dt =∑Fn ij+Ft ij+mig (1) Iidωi dt =∑RVi×Ft ij −τr ij(2) where m i , I i , v i and ωi are the mass, the moment of inertia, the speed of movement and the angular velocity of the rotation of the particle i. F ijn and F ijt is the normal and tangential force induced by contact of particle i with particle j in the current time step. RV i represents the reaction vector between the center of the particle i and the contact point where the force Fijtis acting. The contact model is necessary to evaluate the mutual forces between the particles themselves but also between the particles and other solid bodies in their contact with each other. There are a number of different contact models, but there is no general rule as to which contact model is the best as these models are not equivalent and address different particle contact properties. For example, models based on Hertz theory consider particles in which the contacts lead to nonlinear elastic deformation, whereas the model of the linear spring damping assumes that this deformation is rather viscoelastic. All these models depend on the parameters whose values must be included in the input settings. The basic input values are the radius of the particles R, particle density ρ , Young’s modulus E (shear modulus G), and interaction coefficients. Tsuji et al. (1992) [ 33 ] proposed a nonlinear contact model based on the adaptation of the original model proposed by Cundall and Strack. This currently widely used Hertz-Mindlin model is described by Equations (3) and (4). This model included a non-linear member based on Hertz theory for normal contact. For tangential contact, Tsuji suggested to include viscous dissipation. The modified contact model, therefore, is as follows: Fn ij =−knδn3 2 ij −ηnvsh ij ·nijnij(3) Energies 2018,11, 1849 5 of 18 Ft ij =−ktδt ij −ηtvt ij(4) Slip speed v ijsh is replaced with relative tangent speed v ijt . Sufficiently high tangential forces cause the particles to slide relative to one another or to the surfaces with which they come in contact. EDEM software offers several different contact models and, if necessary, also the creation of user-defined models. Modified contact model Hertz-Mindlin represents the EDEM default contact model and is based on the model proposed by Tsuji et al. In the Hertz-Mindlin contact model, during the collision of the particle i with particle j the normal force Fijnacting on each particle is given by: Fn ij =−4 3E∗√R∗δn3 2 ij −2r5 6ψpknm∗vn ij (5) where E* is equivalent to Young’s module of two colliding particles, R* is an equivalent contact radius, δijn is normal displacement of particles by normal force, m* is the equivalent mass of particles, v ijn is the normal component of relative velocity, and the normal contact stiffness is k n = 2E*(R δn ) 0.5 . Damping coefficient ψ is a function of the coefficient of restitution e and ranges from 0 (absolutely viscous) to 1 (absolutely elastic). F ijt depends on tangential displacement δijt , relative tangential velocity v ijt as well as on tangential stiffness k t = 8G*(R δt ) 0.5 . The tangential force is limited by the Coulomb’s law of friction. Ft ij =−ktδt ij −2r5 6ψpktmvt ij (6) Damping coefficient ψis a function of the coefficient of restitution e and is defined as: ψ=−ln e q(ln e)2+π2(7) In the default EDEM Hertz-Mindlin contact model, the coefficient of restitution for the given case is constant with respect to the impact velocity, provided that other parameters are also constant. The two particles i and j have the equivalent contact radius R*, equivalent mass m*, equivalent Young’s Module E* and shear module G* defined as follows: R∗= 1 Ri+1 Rj!−1 =RiRj Ri+Rj(8) m∗= 1 mi+1 mj!−1 =mimj mi+mj(9) E∗= 1−υi Ei+1−υj Ej!(10) 1 G∗= 2−υi Gi+2−υj Gj!(11) E i and E j is Young’s modulus, G i and G j represent a shear modulus, υi and υj are Poisson’s ratios of particle i and particle j. If the simple Hertz-Mindlin contact model defined in Equations (5) and (6) is used, the following contact parameters must be defined for the simulation model in EDEM, for both particle–particle contacts and particle–geometry contacts: •Static friction coefficient µs •Coefficient of rolling friction µr •Coefficient of restitution e Energies 2018,11, 1849 6 of 18 At present, there are a number of sophisticated contact models that can simulate particular systems with a number of different properties. For example, various forms of powdery materials exhibiting plastic or elastic properties, cohesive soils, wet and sticky materials, or even rocks, using solid connections between the particles. Another modern trend is to connect DEM with other numerical methods and software (SW) such as CFD, FEM, MBD, etc. This interconnection improves numerical modeling to such an extent that it is possible to predict pipe wear during pneumatic transport of abrasive materials or to perform stress analysis of transport and mining machines. This, along with ever-increasing computational power opens up unlimited possibilities. 3. Material and Methods The rotary coolers in this study are used to cool fly ash leaving the boiler. The mechanical–physical properties of the fly ash sample were determined. The objectives of the operational optimization included the assessment of the suitability of the rotary cooler for the process, the assessment of the suitability of the flow path from the fluidized bed boiler through the inlet pipe (Figure 2) into coolers 1 and 2 (Figure 2), and assessment of the possibility of overloading the coolers under the transport capacity of about 13 t/h. The diagram of the whole assembly including the dimensional sketch is shown in Figure 2. Energies 2018, 11, x FOR PEER REVIEW 6 of 19 Ei and Ej is Young's modulus, Gi and Gj represent a shear modulus, υi and υj are Poisson's ratios of particle i and particle j. If the simple Hertz-Mindlin contact model defined in Equations (5) and (6) is used, the following contact parameters must be defined for the simulation model in EDEM, for both particle– particle contacts and particle–geometry contacts:  Static friction coefficient μs  Coefficient of rolling friction μr  Coefficient of restitution e At present, there are a number of sophisticated contact models that can simulate particular systems with a number of different properties. For example, various forms of powdery materials exhibiting plastic or elastic properties, cohesive soils, wet and sticky materials, or even rocks, using solid connections between the particles. Another modern trend is to connect DEM with other numerical methods and software (SW) such as CFD, FEM, MBD, etc. This interconnection improves numerical modeling to such an extent that it is possible to predict pipe wear during pneumatic transport of abrasive materials or to perform stress analysis of transport and mining machines. This, along with ever-increasing computational power opens up unlimited possibilities. 3. Material and Methods The rotary coolers in this study are used to cool fly ash leaving the boiler. The mechanical– physical properties of the fly ash sample were determined. The objectives of the operational optimization included the assessment of the suitability of the rotary cooler for the process, the assessment of the suitability of the flow path from the fluidized bed boiler through the inlet pipe (Figure 2) into coolers 1 and 2 (Figure 2), and assessment of the possibility of overloading the coolers under the transport capacity of about 13 t/h. The diagram of the whole assembly including the dimensional sketch is shown in Figure 2. Figure 2. Computer-aided design (CAD) data assembly with all active geometries which have been used for the discrete element method (DEM) simulation. Figure 2 shows all the active elements of the cooler geometry that come into contact with the conveyed material. The 3D computer-aided design (CAD) model represents a simplified form of the real device geometry in real operation. The assembly includes an unusual chute that feeds and distributes the material between the two coolers. Two coolers with an inner diameter of 0.866 m, a length of 6 m with an inclination of 5°, are directly connected to the chute, and the material is transported upwards. The outer shell is cooled by running water. The axial distance between the coolers is 1 m. 3.1. Fly Ash Characterization Before examining any conveyor, storage or process equipment working with bulk materials, it is first necessary to characterize the material to understand its behavior and specifics. Moreover, characterization of the material is also necessary for the creation of virtual material within the DEM Figure 2. Computer-aided design (CAD) data assembly with all active geometries which have been used for the discrete element method (DEM) simulation. Figure 2shows all the active elements of the cooler geometry that come into contact with the conveyed material. The 3D computer-aided design (CAD) model represents a simplified form of the real device geometry in real operation. The assembly includes an unusual chute that feeds and distributes the material between the two coolers. Two coolers with an inner diameter of 0.866 m, a length of 6 m with an inclination of 5 ◦ , are directly connected to the chute, and the material is transported upwards. The outer shell is cooled by running water. The axial distance between the coolers is 1 m. 3.1. Fly Ash Characterization Before examining any conveyor, storage or process equipment working with bulk materials, it is first necessary to characterize the material to understand its behavior and specifics. Moreover, characterization of the material is also necessary for the creation of virtual material within the DEM simulation. The fly ash was tested primarily for flow and shear properties that have a major effect on material flow through the device. The Schulze ring shear tester RST-01 was used to determine the angle of internal and wall friction, flowability, and cohesion. In addition, the angle of repose was determined by the method of “piling” where the material is released on a steel bowl of 9 cm in diameter. The material is dispensed from a hopper by a vibratory feeder which end is located above the center of the bowl at a height of 10 cm. The material is dosed as long as the cone is not formed on the entire base of the bowl and until the feed material starts to slip out of the bowl. The angle of repose is sometimes also used to characterize the flowability of powders. Materials with the angle of repose below 30 ◦ are classified as free-flowing, 30 ◦ –45 ◦ as easy-flowing, 45 ◦ –55 ◦ as cohesive and Energies 2018,11, 1849 7 of 18 above 55 ◦ as very cohesive to non-flowing materials [ 34 ]. According to the results, it is a non-cohesive, easy-flowing material that should not pose any problems during transportation and storage. 3.2. Virtual Material, Geometry Definition, and Contact Model The DEM model used in this study uses the “soft-sphere” (SSDEM) [ 35 ] method originally developed by Cundall and Strack [ 13 ]. With this method, the particles in contact are able to withstand small deformations, and these deformations are used to calculate the forces acting between the particles. Due to computational limits, there are some limitations in current DEM models. Because of the lack of standardized procedures, the determination of the input parameters necessary for simulation can introduce inaccuracies in the results obtained from the numerical models. One of the key issues is still the representation of actual particle shapes and the limited number of particles that can be modeled with a reasonable computational cost. To achieve an acceptable computational cost, most of the DEM simulations include only a limited number of particles with a diameter in the order of millimeters and a maximum particle count of hundreds of thousands. A simple spherical shape and a small number of particles are usually not comparable to real industrial systems, which consist of billions of small irregular particles with wide distributions of particle size. Therefore, a simulation of millions of irregular particles in a complex system cannot reasonably be performed with a single processor. In real industrial operations, some of the observed phenomena (such as mixing, drying, cooling, etc.) may take a long time in minutes and hours. Typical times that can currently be simulated using DEM at this scale are usually shorter. Parallel computing clusters can provide the required level of computational power required for larger DEM simulations. However, as computer hardware continually improves, it should be possible in the future to model complex geometries with a large number of non-spherical particles using readily available multiprocessor units. For instance, EDEM software since its 2017 version supports not only parallel computing on standard central processing unit (CPU) processors but also takes advantage of graphics processing unit (GPU) processors of graphic cards, which has become quite powerful due to computer games. The challenge for DEM modelers is to create simplified models that can be processed at reasonable time frames and are still relevant to real systems that work with millions of irregularly shaped particles with a wide size distribution. The main objective of DEM simulations is to simulate non-spherical particles efficiently. Various complex shape modeling methods have been published, for example, complex particles can be represented in DEM codes as a group of mutually permeable or non-protruding spheres of the same or different radius, ellipsoids or polygons. Models of non-spherical shapes that represent real particles can be complicated because more sophisticated models of particle interaction are needed rather than a simple contact spherical model. For different types of polygons, contact detection and subsequent computation of forces and torque can be complicated and computationally demanding. Using methods with particles created from a sphere, it is possible to achieve the computation efficiency to detect point contacts and calculation of reaction forces. The disadvantage of this method, however, is that accurate representation of a complex shape requires the clustering of many spherical particles, which in turn requires more computational memory. Another fundamental problem in the development of DEM for real industrial applications is the determination of the input parameters of the particulate material. The simulation input parameters are often not measured, and the values are estimated without proper justification because there is basically no stable method for their determination. To correctly interpret the simulation results, the correct input parameters should be selected, and the simulations should also be validated using experimental tests [36]. Parameters of the virtual material were chosen with respect to the results of the tests performed using the collected samples and by the empirical Equation (12), which indicates the influence of the individual parameters on the angle of repose, according to which the material was subsequently calibrated [ 37 ]. To optimize the computation time, a monodisperse virtual material of spherical particles with a diameter of D = 15 mm was chosen. Some different radiuses were tested, but the final particle radius has been chosen for a good compromise between computational time and particle/equipment Energies 2018,11, 1849 8 of 18 size ratio. The individual parameters for the virtual material model are shown in Table 1. The validation results via angle of repose test are shown in Table 2. The difference between experiment and DEM model is 1.4% in this test. Normally, the total uncertainty of validated simulations is about 5–10% if the quantifiable parameters are compared. Generally, many different tests can be used for validation of model in practice. Unfortunately, in the presented study a calibration at full-scale would not be feasible due to the very high temperature of material. The decisions taken during the whole validation and optimization process have been agreed with experienced practitioners with the intake and the cooling system. AoR =68.61·µs,pp0.27·µs,pw0.22·µr,pp0.06·µr,pw0.12·D−0.2 (12) Table 1. Discrete element method (DEM) input parameters; material properties and interaction coefficients for fly ash particles. Material Density Poisson’s Ratio Young’s Modulus Coeff. of Restitution Static Friction Rolling Friction (kg·m−3)(-) (MPa) (-) (-) (◦) Fly ash 7850 0.25 11 0.65 0.8 1.4 Steel 870 0.30 210 0.70 0.9 1.3 Table 2. Validation of virtual material via angle of repose test. Angle of Repose (◦)—Measurement No. Avg. 12345678 Experiment 38.5 39.2 38.7 38.4 37.0 38.1 37.1 37.3 38.1 DEM 35.0 43.5 38.8 36.5 37.0 37.0 40.1 40.9 38.6 The basic calibration of the material through angle of repose was performed, and the simulation of the un-modified cooler shows the same process behavior and malfunctions as well as real cooling equipment. By this mean of comparing, the material calibration and geometrical optimization were done for the industrial scale of cooler. Thus, the process is validated for the material, and then a comparative analysis is performed between the original and the optimized solution. The imported 3D CAD model is shown in Figure 3. This figure also shows the location of particle generation (particle factory) as well as seven virtual sensors that export the individual data obtained during the simulation. Sensor 1 (Intake) monitors the input of the material into the chute. The material is then separated and the pipeline continues to the cooler 1 (blue) and cooler 2 (red). Each of the coolers has its own sensor at the inlet, outlet and in the middle part. The analysis of material movement inside the chute and the cooler is described in Section 4. Energies 2018, 11, x FOR PEER REVIEW 9 of 19 Figure 3. Imported 3D CAD model and virtual sensors placing. 4. Results and Discussion In this study, an analysis of existing coolers and chute was first performed using DEM modeling. Once the causes and consequences of operational problems were identified, design measures were proposed, and a simulated operation under the same operating conditions was performed again. These conditions were based on actual operation. First, the material was fed to the chute from the particle factory at a mass flow rate of 6.5 t/h, and the filling of the interior of the cooler and the transport of material divided between the two coolers were monitored. In time T = 150 s, the mass of generated material was increased to the rate of 13 t/h, which caused problems in the actual operation. In all cases, the total simulation time was 330 s. Since this study is focused on the mechanics of the particles, i.e., flow of material throughout the system, the main monitored parameters included mass flow rates of the material in individual parts of the system, the velocity of the particle movement in the various axes of the coordinate system, the number of particle contacts with the heat exchange surface of the cooler as well as the time the material remained in the cooler. 4.1. Analysis of the Current State of Operation The basic parameter of each conveying device is its transport capacity, i.e., the volume/mass of the material per unit of time. The material flow values obtained from all seven sensors (Figure 3) are shown in Figure 4a. In the first stage of the study, a material flow simulation was performed in the existing facility. This helped create a reference simulation that shows the behavior of the material throughout the system. Figure 4 shows that at the generated mass flow rate of 6.5 t/h, the material fails to separate equally. Instead, 4.6 t/h was fed to the cooler 1 and 1.9 t/h to the cooler 2. The inclination of the 5° of the coolers clearly causes the accumulation of material in their rear part. Cooler 1 operates at the maximum level when the total flow rate is 6.5 t/h. In time T = 150 s, the volume of material generated in the feed pipeline was increased to 13 t/h, resulting in the immediate overloading of cooler 1. After the whole volume of the chute 1 is filled with the material, cooler 2 after a certain period of time, is overloaded as well. In time, T = 330 s, both coolers already operate at their maximum transport capacity of about 6 t/h. Based on the simulation of the process, it is obvious that the inclination of the coolers to some extent negates the paddle pitch on the cooler shafts. The purpose of these paddles is to move the material inside the cooler in the direction of transport and to mix and spread the material over the cooler body. When monitoring a few selected particles entering the cooler from the chute, their movement, time they remained in the cooler, and number of contacts with the cooling shell were analyzed. The average particle speed in the transport direction was 0.02 m∙s-1, which at the cooler length of 6 m means the average time spent in the cooler (cooling time) was 300 s. Analysis of transport and particle movement at the current state of the fly ash cooler system Figure 3. Imported 3D CAD model and virtual sensors placing. Energies 2018,11, 1849 9 of 18 4. Results and Discussion In this study, an analysis of existing coolers and chute was first performed using DEM modeling. Once the causes and consequences of operational problems were identified, design measures were proposed, and a simulated operation under the same operating conditions was performed again. These conditions were based on actual operation. First, the material was fed to the chute from the particle factory at a mass flow rate of 6.5 t/h, and the filling of the interior of the cooler and the transport of material divided between the two coolers were monitored. In time T = 150 s, the mass of generated material was increased to the rate of 13 t/h, which caused problems in the actual operation. In all cases, the total simulation time was 330 s. Since this study is focused on the mechanics of the particles, i.e., flow of material throughout the system, the main monitored parameters included mass flow rates of the material in individual parts of the system, the velocity of the particle movement in the various axes of the coordinate system, the number of particle contacts with the heat exchange surface of the cooler as well as the time the material remained in the cooler. 4.1. Analysis of the Current State of Operation The basic parameter of each conveying device is its transport capacity, i.e., the volume/mass of the material per unit of time. The material flow values obtained from all seven sensors (Figure 3) are shown in Figure 4a. In the first stage of the study, a material flow simulation was performed in the existing facility. This helped create a reference simulation that shows the behavior of the material throughout the system. Figure 4shows that at the generated mass flow rate of 6.5 t/h, the material fails to separate equally. Instead, 4.6 t/h was fed to the cooler 1 and 1.9 t/h to the cooler 2. The inclination of the 5 ◦ of the coolers clearly causes the accumulation of material in their rear part. Cooler 1 operates at the maximum level when the total flow rate is 6.5 t/h. In time T = 150 s, the volume of material generated in the feed pipeline was increased to 13 t/h, resulting in the immediate overloading of cooler 1. After the whole volume of the chute 1 is filled with the material, cooler 2 after a certain period of time, is overloaded as well. In time, T = 330 s, both coolers already operate at their maximum transport capacity of about 6 t/h. Based on the simulation of the process, it is obvious that the inclination of the coolers to some extent negates the paddle pitch on the cooler shafts. The purpose of these paddles is to move the material inside the cooler in the direction of transport and to mix and spread the material over the cooler body. When monitoring a few selected particles entering the cooler from the chute, their movement, time they remained in the cooler, and number of contacts with the cooling shell were analyzed. The average particle speed in the transport direction was 0.02 m · s −1 , which at the cooler length of 6 m means the average time spent in the cooler (cooling time) was 300 s. Analysis of transport and particle movement at the current state of the fly ash cooler system revealed problems in both the chute, which dispense the material between both coolers unevenly and in the conveying capacity of the coolers. The overloading of both coolers is shown in Figure 4b. The transport capacity limit of the existing equipment and the inlet drain was determined (Figure 4a). 4.2. Inlet Pipeline Design Optimization In the first phase of the optimization process of the entire system, it was necessary to design adjustments in the chute, to achieve the most uniform distribution of material between the two coolers in all circumstances. At the same time, the positions and inclinations of the inlet and outlet pipes should be respected. This section presents an analysis of the current state of the chute which does not fulfill the required function, which is to equally distribute the input materials into individual coolers. Given the current design of the chute, the volume and rate in which the material is fed into the chute are the most important factors. For the study of the division of the material in the chute, the chute was removed from the overall assembly. The “intake”, “chute 1” and “chute 2” virtual sensors were used to monitor the process of distribution of the input material into both outputs. Energies 2018,11, 1849 16 of 18 Nomenclature Symbols used D diameter (m) E Young modulus (Pa) e coefficient of restitution (-) F force (N) G shear modulus (Pa) I moment of inertia (kg·m2) k damping (N·m−1) m mass (kg) R radius (m) RV reaction vector (-) T time (s) v velocity (m·s−1) Greek letters δdisplacement (m) µfriction coefficient (-) νPoisson ratio (-) ρdensity (kg·m−3) ψdamping ratio (-) ωangular velocity (rad·s−1) Subscripts * equivalent i particle i j particle j n normal pp particle–particle pw particle–wall r rolling s static sh shear t tangential Abbreviations CAD Computer-Aided Design CFD Computational Fluid Dynamics CPU Central Processing Unit DEM Discrete Element Method FEM Finite Element Method GPU Graphics Processing Unit MBD Multibody Dynamics SW Software References 1. Wait, R. Discrete element models of particle flows. Math Model Anal. 2001,6, 156–164. [CrossRef] 2. Cleary, P.W. Discrete element modelling of industrial granular flow applications. Task. q.-Sci. Bull. 1998 ,2, 385–416. 3. Antonyuk, S.; Palis, S.; Heinrich, S. Breakage behaviour of agglomerates and crystals by static loading and impact. Powder Technol. 2011,206, 88–98. [CrossRef] 4. Owen, P.J.; Cleary, P.W. Prediction of screw conveyor performance using the Discrete Element Method (DEM). Powder Technol. 2009,193, 274–288. [CrossRef] 5. Feda, J. Mechanics of Particulate Materials—The Principles, 2nd ed.; Academia Prague: Prague, Czech Republic, 1982; ISBN 0-444-99712-X. Energies 2018,11, 1849 17 of 18 6. Hendrik, O.; Kerst, K.; Roloff, C.; Janiga, G.; Katterfeld, A. CFD-DEM simulation and experimental investigation of the flow behavior of lunar regolith JSC-1A. Particuology 2018, in press. [CrossRef] 7. Forsström, D.; Pär, J. Calibration and validation of a large scale abrasive wear model by coupling DEM-FEM: Local failure prediction from abrasive wear of tipper bodies during unloading of granular material. Eng. Fail. Anal. 2016,66, 274–283. [CrossRef] 8. Barrios, G.K.; Tavares, L.M. A preliminary model of high pressure roll grinding using the discrete element method and multi-body dynamics coupling. Int. J. Miner. Process. 2016,156, 32–42. [CrossRef] 9. Fedorko, G.; Ivanˇco, V. Analysis of Force Ratios in Conveyor Belt of Classic Belt Conveyor. Procedia. Eng. 2012,48, 123–128. [CrossRef] 10. Manjula, E.V.P.J.; Ariyaratne, W.K.H.; Ratnayake, C.; Melaaen, M.C. A review of CFD modelling studies on pneumatic conveying and challenges in modelling offshore drill cuttings transport. Powder Technol. 2017 , 305, 782–793. [CrossRef] 11. Chen, W.; Roberts, A.; Katterfeld, A.; Wheeler, C. Modelling the stability of iron ore bulk cargoes during marine transport. Powder Technol. 2018,326, 255–264. [CrossRef] 12. Markauskas, D.; Ramírez-Gómez, Á.; Kaˇcianauskas, R.; Zdanceviˇcius, E. Maize grain shape approaches for DEM modelling. Comput. Electron. Agric. 2015,118, 247–258. [CrossRef] 13. Bai, Y.; He, F. Modeling on the Effect of Coal Loads on kinetic Energy of Balls for Ball Mills. Energies 2015 ,8, 6859–6880. [CrossRef] 14. Cundall, P.A.; Strack, O.D.L. A discrete numerical model for granular assemblies. Geotech. 1979 ,29, 47–65. [CrossRef] 15. Marigo, M. Discrete Element Method Modelling of Complex Granular Motion in Mixing Vessels: Evaluation and Validation. Ph.D. Thesis, The University of Birmingham, Birmingham, UK, 2012. 16. Kannan, A.S.; Jareteg, K.; Lassen, N.C.K.; Carstensen, J.M.; Hansen, M.A.E.; Dam, F. Design and performance optimization of gravity tables using a combined CFD-DEM framework. Powder Technol. 2017 ,318, 423–440. [CrossRef] 17. Haeri, S. Optimisation of paddle type spreaders for powder bed preparation in Additive Manufacturing using DEM simulations. Powder Technol. 2017,321, 94–104. [CrossRef] 18. Zhao, L.; Zhao, Y.; Bao, C.; Hou, Q.; Yu, A. Optimisation of a circularly vibrating screen based on DEM simulation and Taguchi orthogonal experimental design. Powder Technol. 2017,310, 307–317. [CrossRef] 19. Pantaleev, S.; Yordanova, S.; Janda, A.; Marigo, M.; Ooi, J.Y. An experimentally validated DEM study of powder mixing in a paddle paddle mixer. Powder Technol. 2017,311, 287–302. [CrossRef] 20. Wang, D.; Servin, M.; Mickelsson, K.O. Outlet design optimization based on large-scale nonsmooth DEM simulation. Powder Technol. 2014,253, 438–443. [CrossRef] 21. Kulka, J.; Mantic, M.; Faltinova, E.; Molnar, V.; Fedorko, G. Failure analysis of the foundry crane to increase its working parameters. Eng. Fail. Anal. 2018,88, 25–34. [CrossRef] 22. Roberts, A.W.; Wensrich, C.M. Flow dynamics or ‘quaking’in gravity discharge from silos. Chem. Eng. Sci. 2002,57, 295–305. [CrossRef] 23. Klinzing, G.E. A review of pneumatic conveying status, advances and projections. Powder Technol. 2018 ,333, 78–90. [CrossRef] 24. Neˇcas, J.; Hlosta, J.; Žurovec, D.; Žídek, M.; Rozbroj, J.; Zegzulka, J. Feeder Type Optimisation for the Plain Flow Discharge Process of an Underground Hopper by Discrete Element Modelling. Adv. Sci. Technol. Res. J. 2017,11, 246–252. [CrossRef] 25. Roberts, A.W. Particle Technology—Reflections and Horizons. Chem. Eng. Res. Des. 1998 ,76, 775–796. [CrossRef] 26. Drescher, A.; Josselin, G. Photoelastic verification of a mechanical model for the flow of a granular material. J. Mech. Phys. Solids. 1972,20, 337–340. [CrossRef] 27. Raclavská, H.; Corsaro, A.; Juchelková, D.; Sassmanová, V.; Frantík, J. Effect of temperature on the enrichment and volatility of 18 elements during pyrolysis of biomass, coal, and tires. Fuel Process. Technol. 2015 ,131, 330–337. [CrossRef] 28. Belohlav, Z.; Brenková, L.; Hanika, J.; Durdil, P.; Rapek, P.; Tomasek, V. Effect of Drug Active Substance Particles on Wet Granulation Process. Chem. Eng. Res. Des. 2017,85, 974–980. [CrossRef] 29. Hlosta, J.; Zurovec, D.; Gelnar, D.; Zegzulka, J.; Necas, J. Thermal Conductivity Measurement of Bulk Solids. Chem. Eng. Technol. 2017,40, 1876–1884. [CrossRef] Energies 2018,11, 1849 18 of 18 30. Baniasadi, M.; Baniasadi, M.; Peters, B. Coupled CFD-DEM with heat and mass transfer to investigate the melting of a granular packed bed. Chem. Eng. Sci. 2018,178, 136–145. [CrossRef] 31. Bellan, S.; Matsubara, K.; Cho, H.S.; Gokon, N.; Kodama, T. A CFD-DEM study of hydrodynamics with heat transfer in a gas-solid fluidized bed reactor for solar thermal applications. Int. J. Heat Mass Transf. 2018 ,116, 377–392. [CrossRef] 32. Wu, H.; Gui, N.; Yang, X.; Tu, J.; Jiang, S. A smoothed void fraction method for CFD-DEM simulation of packed pebble beds with particle thermal radiation. Int. J. Heat Mass Transf. 2018,118, 275–288. [CrossRef] 33. 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. [CrossRef] 34. Raymus, G.J. Handling of Bulk Solids. In Chemical Engineers’ Handbook, 6th ed.; Perry, R.H., Green, D., Eds.; McGraw-Hill: New York, NY, USA, 1985; ISBN 978-0070494794. 35. Schwartz, S.R.; Richardson, D.C.; Michel, P. An implementation of the soft-sphere discrete element method in a high-performance parallel gravity tree-code. Granul. Matter 2012,14, 363–380. [CrossRef] 36. Hlosta, J.; Zurovec, D.; Rozbroj, J.; Ramírez-Gómez, Á.; Necas, J.; Zegzulka, J. Experimental determination of particle-particle restitution coefficient via double pendulum method. Chem. Eng. Res. Des. 2018 , in press. [CrossRef] 37. Beakawi Al-Hashemi, H.M.; Baghabra Al-Amoudi, O.S. A review on the angle of repose of granular materials. Powder Technol. 2018,330, 397–417. [CrossRef] 38. Bian, X.; Wang, G.; Wang, H.; Wang, S.; Lv, W. Effect of lifters and mill speed on particle behaviour, torque, and power consumption of a tumbling ball mill: Experimental study and DEM simulation. Miner. Eng. 2017 , 105, 22–35. [CrossRef] © 2018 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).