NUMERICAL AND EXPERIMENTAL INVESTIGATION OF MACROSCALE MIXING APPLIED TO PULP FIBRE SUSPENSIONS by Clara Gómez B.A.Sc. (Hons), Universidad Pontificia Bolivariana, Colombia, 2002 M.A.Sc., University of British Columbia, 2004 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in The Faculty of Graduate Studies (Chemical and Biological Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) June 2009 ©Clara Gómez, 2009 ii Abstract Effective mixing of pulp fibre suspensions is essential for many pulp and paper manufacturing processes. Pulp suspensions display non-Newtonian rheology and possess a yield stress, which complicates mixing. In order to improve our understanding of pulp mixing in agitated vessels, a series of studies was undertaken to assess the suitability of using computational fluid dynamics (CFD) to model these systems. CFD simulations for laboratory-scale and industrial-scale mixing chests were developed with the pulp suspensions treated as Bingham fluids. The computed flow fields were used to determine the dynamic response of the virtual mixers, which was then compared with experimental measurements providing very good agreement under conditions of moderate agitation. Comparison between calculations and measurements of torque and axial force was also good (relative average error of 12% to 24%). The simulation results provided insight into the mixing flow occurring within the systems, showing the formation of caverns around the impeller(s), the location of stagnation regions and the presence of channeling. However, the accuracy of these predictions was limited by the Bingham model used to describe the suspension's rheology and the uncertainty to which the suspension's yield stress could be measured. To assess the degree to which the approximated rheology contributed to the CFD results, the mixing of a model fluid having a well-defined rheology (Newtonian glycerin solution) was extensively investigated in a laboratory-scale vessel using a typical industrial geometry (rectangular chest, side-entering axial flow impeller). The flow fields were measured using particle image velocimetry (PIV) and compared with CFD computations for identical operating conditions. Good agreement was found (avg. 13.1% RMS deviation of local axial velocities) iii confirming that the approach used in the CFD model was adequate to calculate the complex mixing flow fields for a Newtonian fluid. These results encouraged further research on extending the combined CFD and PIV application on a system more closely representative of pulp suspension agitation. Calculations were then performed to select a suitable non-Newtonian model fluid and appropriate operating conditions to model pulp suspension mixing in the laboratory-scale vessel, with the dimensionless cavern diameter as the mixing criterion for conservation between the two systems. iv Table of Contents Abstract............................................................................................................................................ii Table of Contents........................................................................................................................... iv List of Tables ................................................................................................................................. ix List of Figures ................................................................................................................................. x Nomenclature.............................................................................................................................. xvii Acknowledgements ...................................................................................................................... xx Co-Authorship Statement............................................................................................................. xxi 1 Introduction............................................................................................................... 1 1.1 Background ........................................................................................................... 1 1.2 Literature Review .................................................................................................. 3 1.2.1 Mixing Process Considerations....................................................................... 3 1.2.2 Impeller Characteristics .................................................................................. 4 1.2.3 Power and Power Number .............................................................................. 5 1.2.4 Newtonian vs. Non-Newtonian Mixing.......................................................... 7 1.2.5 Pulp Suspension Mixing ................................................................................. 9 1.2.6 Characterization of Pulp Suspensions........................................................... 10 1.2.7 Design of Agitated Stock Chests .................................................................. 12 1.2.8 Mixing Dynamics in Agitated Chests ........................................................... 15 1.2.9 Computational Fluid Mixing......................................................................... 16 1.2.10 Stirred Tank Modeling................................................................................. 17 1.2.11 CFD Modeling of Pulp Suspension Mixing................................................. 19 1.2.12 Particle Image Velocimetry Principles and Application.............................. 20 v 1.3 Research Objectives ............................................................................................ 24 1.4 Thesis Outline...................................................................................................... 26 1.4.1 Modeling a Pilot Scale Pulp Mixing Chest using CFD ................................ 26 1.4.2 Computational Modeling of Industrial Pulp Stock Chests ........................... 27 1.4.3 Investigation of the Flow Field in a Rectangular Vessel Equipped with a Side Entering Agitator .......................................................................................... 27 1.4.4 Characterization of Axial Flow Impellers on Pulp Suspensions .................. 28 1.4.5 Carbopol as a Model Fluid for Studying Mixing of Pulp Fibre Suspensions... ....................................................................................................................................29 2 3 1.5 Contributions of this Work.................................................................................. 29 1.6 References ........................................................................................................... 32 Modeling a Pilot-Scale Pulp Mixing Chest using CFD ............................ 38 2.1 Introduction ......................................................................................................... 38 2.2 Experimental Techniques .................................................................................... 40 2.3 Results and Discussion........................................................................................ 44 2.4 Summary and Conclusions.................................................................................. 54 2.5 Nomenclature ...................................................................................................... 54 2.6 References ........................................................................................................... 56 Computational Modeling of Industrial Pulp Stock Chests .................... 58 3.1 Introduction ......................................................................................................... 58 3.2 Experimental ....................................................................................................... 64 3.2.1 Chest Geometry and Operating Conditions .................................................. 64 vi 3.2.2 Industrial Dynamic Tests .............................................................................. 66 3.3 Modeling Approach............................................................................................. 67 3.4 Results and Discussions ...................................................................................... 72 3.4.1 Machine Chest (MC-A) ................................................................................ 72 3.4.2 4 Latency Chest (LC-B)................................................................................... 80 3.5 Conclusions ......................................................................................................... 85 3.6 Nomenclature ...................................................................................................... 86 3.7 References ........................................................................................................... 88 Investigation of the Flow field in a Rectangular Vessel equipped with a Side-entering Agitator..................................................................................... 90 5 4.1 Introduction ......................................................................................................... 90 4.2 Experimental Methods ........................................................................................ 93 4.2.1 Experimental Setup ....................................................................................... 93 4.2.2 Experimental Procedures .............................................................................. 97 4.3 Computational Techniques................................................................................ 105 4.4 Results and Discussion...................................................................................... 109 4.5 Summary and Conclusions................................................................................ 120 4.6 Nomenclature .................................................................................................... 120 4.7 References ......................................................................................................... 123 Characterization of Axial Flow Impellers in Pulp Fibre Suspensions.... .....................................................................................................................................127 5.1 Introduction ....................................................................................................... 127 vii 6 5.2 Experimental .................................................................................................. ...130 5.3 Results and Discussion...................................................................................... 131 5.3.1 Impeller Characterization............................................................................ 131 5.3.2 Impeller Characterization using CFD ......................................................... 134 5.4 Summary and Conclusions................................................................................ 138 5.5 Nomenclature .................................................................................................... 138 5.6 References ......................................................................................................... 140 Carbopol as a Model Fluid for Studying Mixing of Pulp Fibre Suspensions ................................................................................................................ 142 6.1 Introduction ....................................................................................................... 142 6.1.1 Rheology of Pulp Fibre Suspensions and Carbopol Solutions ................... 144 6.2 Experimental ..................................................................................................... 147 6.3 Results and Discussion...................................................................................... 149 6.3.1 Selection of Scaling Criteria ....................................................................... 149 6.3.2 Pulp Suspension and Carbopol Rheology................................................... 151 6.3.3 7 Establishing Similarity Between Pulp and Carbopol Mixing Studies ........ 157 6.4 Summary and Conclusions................................................................................ 165 6.5 Nomenclature .................................................................................................... 165 6.6 References ......................................................................................................... 168 Conclusions ........................................................................................................... 172 7.1 Conclusions ...................................................................................................... 172 7.2 Recommendations for Future Work .................................................................. 177 7.3 References ......................................................................................................... 180 viii Appendices ............................................................................................................................. 181 Appendix A......................................................................................................................181 A.1 Details of the Pilot Mixing Chest CFD Model................................................. 181 A.2 References ........................................................................................................ 188 Appendix B......................................................................................................................189 B.1 Procedure for Estimation of Discretization Error in CFD Simulations ............ 189 B.2 References......................................................................................................... 191 Appendix C......................................................................................................................192 C.1 Extended Details on Impeller Characterization Using CFD............................. 192 C.2 Power and Force Number: CFD vs Experimental ............................................ 196 C.3 References......................................................................................................... 200 ix List of Tables Table 2.1 Comparison of model parameters estimated from experimental and CFD data. Cm = 3.3% FBK at Q = 7.9 L/min, 21.1 L/min and 37.1 L/min for selected impeller speeds (N)……………………………………………………………....47 Table 3.1 Industrial chest specifications……………………………………………………66 Table 3.2 Comparison of experimental and CFD simulated parameters for the stock chests……………………………………………………………………………..75 Table 3.3 Experimental data for LC-B……………………………………………………...81 Table 4.1. Calculations of discretization error using the GCI method (refer to Appendix B for a detailed description)……………………………………………………………………..109 Table 6.1 Comparison of industrial-scale pulp mixing chests with the pilot-scale chest used in scaling. Evaluated for a 3% FBK softwood pulp suspension………………..150 Table 6.2 Yield stress of pulp suspensions and carbopol solutions……………………….154 Table 6.3 Process parameters for pilot- and laboratory-scale mixing chests upon scale-down. Evaluated for a 3% FBK softwood pulp suspension…………………………...162 Table 6.4 Operating parameters for pilot-scale mixing chests with 3 wt. % pulp vs. laboratory-scale mixing chest with aqueous carbopol solutions………………..164 x List of Figures Figure 1.1 Dynamic model of the agitated chest (a) Ein-Mozaffari, 2002 and (b) Soltanzadeh et al., 2009………………………………………………………………..………16 Figure 1.2 Operation of the PIV system……………………………………………………..25 Figure 2.1 Boundary conditions specified in the computational domain. The projection shown corresponds to the z-y plane using the coordinate system chosen for the model……………………………………………………………………………..41 Figure 2.2 Block diagram of the dynamic model used to analyze the experimental and CFD generated data…………………………………………………………………....44 Figure 2.3 Velocity vectors (arrows) in three planes parallel to the impeller flow: in the plane x = 0 (the impeller mid-plane) and the planes x = ±9 cm (close to the impeller tips). Cm = 3.3% FBK suspension at Q = 7.9 L/min and N = 1132 rpm. Vectors are colored by velocity magnitude (m/s). Some detail in the impeller zone has been omitted for clarity………………………………………………………………..45 Figure 2.4 Input (upper)/Output (lower) signals for numerical (dashed lines) and experimental (solid lines) measurements for the scale-model chest. Cm = 3.3% FBK, Q = 7.9 L/min, N = 1132 rpm……………………………………………..46 Figure 2.5 The effect of impeller speed on channeling (f) as a function of the pulp flow rate (after Ein-Mozzafri et al., 2004)…………………………………………………48 Figure 2.6 The effect of the impeller speed on channeling at different flow rates. Comparison of experimental measurements and CFD predictions……………………………49 Figure 2.7 Input (upper)/Output (lower) signals for numerical (dashed lines) and experimental (solid lines) measurements for the scale-model chest. Cm = 3.3% FBK at Q = 37.1 L/min: (a) N = 1456 rpm; (b) N = 1203 rpm…………………..52 xi Figure 3.1 Grey box model for flow through pulp stock chests (after Soltanzadeh et al., 2009)……………………………………………………………………………..59 Figure 3.2 Machine chest (MC-A): (a) Multiblock structure used for grid generation; (b) outer surface of mesh used in the initial grid…………………………………….65 Figure 3.3 Latency chest (LC-B): (a) Multiblock structure used for grid generation; (b) outer surface of mesh used in the initial grid…………………………………………..69 Figure 3.4 Test for grid independence for MC-A: (a) magnitude of suspension velocity 8 cm in front of the impeller (along the y-axis); (b) magnitude of the suspension velocity 3m away from the impeller. (c) Effect of yield stress and yielding viscosity on computed velocity along the impeller axis…………………………72 Figure 3.5 Raw experimental bump test data for MC-A. Note problems with consistency meter calibration and compression of the data (particularly the exit data)………74 Figure. 3.6 Cavern formation in MC-A: (a) Iso-surface where the shear-stress equals the suspension yield stress; (b) Iso-surface where the velocity equals 0.01 VTip. (c) Trajectories of particles released at main stock inlet (Inlet 1)…………………...78 Figure. 3.7 MC-A: (a) input and (b) output signals for experimental and computational data. Consistency scale on the right axis corresponds to the computational output signals with adjusted gain and offset. Test conditions and model parameters as given in Table 3.2………………………………………………………………..79 Figure. 3.8 Raw experimental bump test data for LC-B……………………………………..81 Figure. 3.9 Cavern formation in LC-B: (a) Iso-surface where the shear-stress equals the yield stress; (b) Iso-surface where the velocity equals 0.01 VTip. (c) Trajectories of particles released at the stock inlet……………………………………………….83 xii Figure. 3.10 LC-B: (a) input and (b) output signals for experimental and computational data. Consistency scale on the right axis corresponds to the computational output signals with adjusted gain and offset. Test conditions and model parameters as given in Table 3.2………………………………………………………………..85 Figure 4.1 (a) Schematic of the experimental setup showing the front view of the mixing vessel, the PIV acquisition system and triggering device used. (b) Location and size of the measuring areas. Dashed lines illustrate the locations at which the velocity profiles were evaluated. The horizontal line is located at impeller height and the vertical line is located at a distance d cm (varied from 3−10cm) in front of the impeller. (c) PIV Illumination planes. (d) Photograph showing the flow patterns observed in the PIV illumination plane aligned with the impeller tip………………………………………………………………………………...95 Figure 4.2 Profile of the mean (a) axial and (b) vertical velocity at a line located at x*= 0.37; d = 6cm (see Fig. 4.1(b)) in the measurement plane aligned with the impeller tip (z = -5.25 cm), at N = 118 rpm, Re = 36……………………………………….101 Figure 4.3 Effect of ∆t on the mean value of the axial velocity measured at a line located at x*= 0.33; d = 4.5cm (see Fig. 4.1(b)) in the measuring plane aligned with the impeller tip. (a) N = 118 rpm, Re = 36 (b) N = 267 rpm, Re = 83……………...104 Figure 4.4 Effect of the phase angle θ on the mean value of the (a) axial and (b) vertical velocity measured at a line located at x*= 0.29; d = 3cm (see Fig. 4.1(b)), in the measuring plane aligned with the impeller tip. N = 267 rpm, Re = 83………...105 Figure 4.5 Computed vertical profiles of axial velocity in a plane aligned with the blade tip (a) for a vertical line at x* = 0.29; d = 3cm (see Fig. 4.1(b)) and (b) for a vertical line at x* = 0.48; d = 10 cm (see Fig. 4.1(b)). N = 327 rpm, Re = 101………...108 Figure 4.6 Mean planar velocity fields for Re =101 in the three vertical planes of measurement: (a) mid plane (z = 0 cm); (b) plane aligned with the impeller tip (z = -5.25 cm) and (c) plane located 2 cm away from the vessel wall xiii (z = -10.25cm). Velocity vectors colored by velocity magnitude (m/s)……….110 Figure 4.7 Experimental (left) and CFD (right) planar velocity fields in the center plane of measurement for (a) Re = 36, (b) Re = 82 and (c) Re = 101. Velocity vectors colored by velocity magnitude (m/s)……………………...……………………113 Figure 4.8 CFD and PIV results for horizontal profiles of normalized axial and vertical velocity components for a horizontal line located at the impeller center (see Fig. 4.1(b)), and covering the full length of the PIV image, in a plane aligned with the blade tip (z = -5.25cm). Results are shown for N = 386 rpm , Re = 120 and N = 118 rpm, Re = 36. The error bars indicate the standard deviation for a series of three measurements……………………………………………………………..115 Figure 4.9 CFD and PIV results for vertical profiles of normalized axial and vertical velocity components for a vertical line at x*= 0.34; d = 5cm (see Fig 4.1(b)), spanning the PIV image height, in a plane aligned with the blade tip (z = -5.25 cm). Results are shown for N = 118 rpm, Re = 36 and N = 267 rpm, Re = 82. The error bars indicate the standard deviation for a series of three measurements…………….117 Figure 4.10 CFD and PIV results for vertical profiles of normalized axial and vertical velocity components for a vertical line at x* = 0.48; d = 10 cm (see Fig. 4.1(b)), spanning the PIV image height, in a plane aligned with the blade tip (z =-5.25 cm). Results are shown for N = 118 rpm, Re = 36 and N = 267 rpm, Re = 82. The error bars indicate the standard deviation for a series of three measurements…………….118 Figure 4.11 CFD and PIV results for vertical profiles of normalized axial and vertical velocity components at Re = 82; (a) for line at x* = 0.34; d = 5cm (see Fig. 4.1(b)) in the mid plane of measurement (z = 0 cm), and (b) for line at x* = 0.48; d = 10 cm (see Fig. 4.1(b)) in the measurement plane located 2 cm away from the vessel wall (z = -10.25 cm). The error bars indicate the standard deviation for a series of three measurements.…………………………………………………………………..119 xiv Figure 4.12 Comparison of NP vs. Re for experimental and computational tests. The error bars are based on the accuracy with which the torque was measured (±0.002Nm)……………………………………………………………………120 Figure 5.1 Mixing vessel (T1) and impellers studied (Maxflo on the left and the marine impeller on the right)…………………………………………………………...131 Figure 5.2 NP and Nf for Maxflo impeller (D = 10.2 cm) in softwood pulp suspensions….132 Figure 5.3 NP and Nf for marine impeller (D = 9.0 cm) in softwood pulp suspensions……132 Figure 5.4 NP and Nf for Maxflo impeller (D = 10.2 cm) in hardwood and softwood pulp suspensions……………………………………………………………………..133 Figure 5.5 NP and Nf for Maxflo impellers as a function of D/T in hardwood pulp suspensions……………………………………………………………………..134 Figure 5.6 Velocity vectors colored by velocity magnitude (m/s) at plane z = 0. Hardwood pulp at Cm = 3%...................................................................................................136 Figure 5.7 Comparison of NP and Nf for experimental and computational tests for cylindrical vessel T1. D = 10.2 cm Maxflo impellers in a Cm = 3% hardwood pulp suspensions as a function of rotational speed…………………………………..137 Figure 6.1 Flow curves measured in six tests made of a 3 wt.% FBK pulp fibre suspension. The average flow curve is given by the solid line……………………………...151 Figure 6.2 Shear stress vs. shear rate. (a) Carbopol solutions of 0.1, 0.15 and 0.20 wt.%; (b) FBK pulp fibre suspensions of 0.75, 2.0, 2.5 and 3 wt.%. Vane in large cup geometry: Vane diameter = 2.5 cm, cup diameter = 10 cm, gap size = 3.75 cm; T = 23oC……………………………………………………………...................153 xv Figure 6.3 Instantaneous viscosity vs. shear stress. (a) Carbopol solutions of 0.1, 0.15 and 0.20 wt.%; (b) FBK pulp fibre suspensions of 0.75, 2.0, 2.5 and 3 wt.%. Vane in large cup geometry: Vane diameter = 2.5 cm, cup diameter = 10 cm, gap size = 3.75 cm; T = 23oC………………………………………………………………155 Figure 6.4 Predictions of shear stress when using the Bingham approximation against the averaged flow curve measured for a 3 wt.% pulp fibre suspension…………….156 Figure 6.5 Characterization of Maxflo II impeller (D =10.2 cm, D/T = 0.39) in softwood pulp suspensions and carbopol solutions. (a) NP vs. Re’; (b) NF vs. Re’……………..160 Figure A.1 Schematic of the experimental setup of the pilot scale pulp mixing chest. Inputoutput configuration simulated in chapter 2……………………………………181 Figure A.2 Three dimensional view of Fig. 2.1. (The impeller is centered on the back wall 12 cm above the floor and has an offset of 7cm on the z direction)……………….183 Figure A.3 Calculated and measured power input vs. impeller rotational speed. Cm = 3.3% FBK suspension at Q = 37.1 L/min. (The representative error bars are based on the accuracy of the torque meter (± 0.05 Nm) used to measure the experimental data)……………………………………………………………………………..187 Figure C.1 Isometric view of the grids used for simulations in the square (S1) and cylindrical baffled (T1) with the Maxflo impeller D=10.1cm. Final three-dimensional meshes of T1 and S1 had 463742 and 340575 computational cells respectively……………………………………………………………………..192 Figure C.2 Test for grid independence. The figures show the magnitude of velocity at lines located 3 cm in front of the impeller (along the x-axis) in the square vessel with hardwood pulp at Cm = 3%, N= 302 rpm .(a) line located at the plane z = 0 cm.(b) line located at the plane z = 5 cm……………………………………………….193 xvi Figure C.3 Velocity vectors colored by velocity magnitude (m/s) at planes z = 0 cm , z = -5 cm and z = 5 cm. Hardwood pulp at Cm = 3%, N = 404 rpm…………………..194 Figure C.4. Velocity vectors colored by velocity magnitude (m/s) at plane z = 0 cm. Hardwood pulp at Cm= 3%..................................................................................195 Figure C.5 Pathlines of massless particles released from a circular surface of D = 10.1cm located in front of the impeller at y = 1.2 cm (coordinate system is centered at the impeller) Hardwood pulp at Cm = 3%, N = 404 rpm…………………………...196 Figure C.6 Comparison of NP from experimental and computational tests for vessel S1. Maxflo impeller D = 10.2cm in a hardwood pulp at Cm= 3%.............................197 Figure C.7 Comparison of Nf from experimental and computational tests for vessel S1. Maxflo impeller D = 10.2cm in a hardwood pulp at Cm= 3%.............................198 xvii Nomenclature Bi Bingham number, dimensionless C impeller C factor, dimensionless Cm fibre mass concentration, % d particle displacement, mm or pixels D impeller diameter, m or cm f percentage of channeling, fraction k consistency factor (kg/ms) g gravitational acceleration, m/s2 Fa axial thrust generated by the impeller, N → F external body forces, N Gchest transfer function of the chest I unit tensor Ks Metzner- Otto impeller constant L mixing chest length, m LMo level momentum, m2/s2 m suspension property, m/s p static pressure, Pa. M torque, N.m Mo impeller momentum number, m4/s2 N impeller rotational speed, s-1 or rpm Fr Froude number, dimensionless NP power number, dimensionless xviii NQ impeller flow number, dimensionless Re Reynolds number, dimensionless P power, W Q impeller pumping rate (in chapter 1), m3/s R percentage of recirculation, fraction Re Reynolds number, dimensionless r position vector in the rotating frame t time, s T tank diameter, m or cm Td single time delay representing plug flow through chest T1 time delay for the channeling zone, s T2 time delay for the mixing zone, s u velocity, m/s → υ → absolute velocity, m/s υr relative velocity, m/s V suspension volume inside the chest, m3 VFM fully mixed volume, m3 W chest width, m Z suspension height, m Greek letters ∆t inter image time, µs Ω angular velocity, rpm γ shear rate, s-1 xix γavg average shear rate in the mixing vessel, s-1 µp plastic viscosity, kg/m.s µ fluid viscosity, kg/m.s µa apparent viscosity, kg/m.s µeff near impeller average or effective viscosity, kg/m.s µ0 yielding viscosity, kg/m.s ρ fluid density, kg/m3 = τ stress tensor, N/m2 τ shear stress, N/m2 τ1 time constant for the channeling zone, s τ2 time constant for the mixing zone, s τy yield stress, N/m2 xx Acknowledgements I offer my enduring gratitude and appreciation to my supervisors Prof. Chad Bennington and Prof. Fariborz Taghipour. I am forever indebted for their guidance and understanding through many difficult times and for giving me insightful scientific advice. Working in the Dynamic Mixing Group has been a privilege. I would like to thank all members, both past and present, for their hard work that has inspired so many ideas. In particular, I would like to thank Farhad, Sujit, Manish, Ali and Babak for our collaborative efforts during this work. Special thanks to Prof. Guy Dumont and Prof. Carl Ollivier-Gooch who have always been willing to help me through out these years. Thanks to Prof. James Olson for listening and offering his advice at the ending stages of this journey. I also thank Prof. Suzanne Kresta for her kindness and support. I am grateful for the financial support of the University of British Columbia UGF and NSERC. Special thanks to Tim Paterson, Ken Wong, Lisa Hudson, Helsa Leong, Richard Ryoo, David Roberts and Doug Yuen for their invaluable assistance. From the bottom of my heart, I thank my good friends for their unconditional support: Esther, Jens, Barbara, Joerg, Matt, Sean, Meggan, Pete, Tara, Austin, Louise, Chris, Desirae, Andres, Laura, Esteban and Gustavo. Their love and kindness never goes unnoticed. The most profound and loving thanks to my family; in particular, to my parents who have given me the greatest love and support anyone could ever expect. I am blessed to have them and I have missed them deeply during these years. They are my true inspiration. xxi Co-Authorship Statement Chapter 2 has been published: Ford† C., Bennington C.P.J and Taghipour F., 2007. Modeling a pilot scale pulp mixing chest using CFD. J. Pulp Pap. Sci., 33, 115−120. I conducted all the simulations presented. The written work was a collaborative effort between my supervisors Prof. Chad Bennington, Prof. Fariborz Taghipour and me. Chapter 3 has been submitted for publication: Bhattacharya S., Gomez† C., Sotanzadeh A., Taghipour F., Bennington C.P.J. and Dumont G.A., 2008. Computational modeling of industrial pulp stock chests. I conducted all the CFD simulations for the Machine chest A and geometry discretization for Latency chest B. Simulations for Latency chest B were conducted by Dr. Sujit Bhattacharya. Mill data was collected by Dr. Sujit Bhattacharya, Ali Sotanzadeh, Prof. Chad Bennington and I. Identification of dynamic parameters was conducted by Ali Sotanzadeh. The written work was a collaborative effort between Dr. Bhattacharya, Prof. Chad Bennington and me. Chapter 5 has been published. Bhole M., Ford† C. and Bennington C.P.J., 2009. Characterization of axial flow impellers in pulp fibre suspensions. Chem. Eng. Res. Des, 87, 4, 648−653. Dr. Manish Bhole conducted the experimental characterization of the impellers in pulp. I conducted all the CFD simulations presented. The written work was a collaborative effort between Dr. Bhole, Prof. Chad Bennington and me. xxii † Chapter 6 has been submitted for publication: Gomez C., Derakhshandeh B., Hatzikiriakos S.G. and Bennington. C.P.J., 2009. Carbopol as a model fluid for studying mixing of pulp fibre suspensions. Babak Derakhshandeh conducted the rheological experiments and parameter fitting for the pulp suspensions used. I conducted the impeller characterization tests, rheological tests and parameter fitting for the Carbopol solutions. I carried out all the similarity studies and scaling calculations presented. Prof. Savvas. Hatzikiriakos revised the manuscript and provided insightful advices on the rheological aspects addressed. The written work was a collaborative effort between my supervisor Prof. Chad Bennington and me. A version of Appendix A has been published. Ford† C., Ein-Mozaffari F., Bennington C.P.J. and Taghipour F., 2006. Simulation of Mixing Dynamics in Agitated Pulp Stock Chest Using CFD, AIChE Journal, 52, 10, 3562−3569. I conducted all the simulations presented. Prof. Farhad EinMozzafari carried out the dynamic experiments used for model verification. The written work was a collaborative effort between my supervisors Prof. Chad Bennington, Prof. Fariborz Taghipour and me. † C. Gomez, formerly C. Ford Chapter 1 1 Chapter 1 Introduction 1.1 Background Over the centuries, paper has been one of the drivers of world development, providing the means for people to record and communicate ideas, news and works of art. Despite the increased usage of plastics and the emergence of new information technologies, paper remains an essential communication medium in the modern world. Pulp and paper products span from hygiene items through special high quality print papers, to diverse packaging materials. The different pulps produced from mechanical and/or chemical processes are blended in order to manufacture paper products with the desired characteristics for their final use. Mixing and blending are essential operations throughout the pulp and paper manufacture, where there exists an extensive range of mixing applications and scales. In the macro scale, the agitation of pulp stock for homogenization and blending and the mixing of pulp suspensions with a variety of wet-end chemicals are the most common mixing operations in the production process. Pulp suspensions must be uniformly mixed with other pulp streams with different attributes, as well as with additives, fillers and chemicals prior to sheet forming. This mixing is readily accomplished in large mechanically stirred vessels known as agitated stock chests. However, the unique rheology of pulp fibre suspensions (which is non-Newtonian and displays a significant yield stress) creates significant challenges during mixing. Therefore, understanding the suspension behaviour is critical for the design of pulp mixing systems. References start on page 33 Chapter 1 2 The current design of these chests is based on limited published information or proprietary criteria held by mixer vendors based on experience with previous installations. The quality of the mixing produced by these designs is rarely measured. In many cases, agitated stock chests are designed to reduce fast or high frequency disturbances in pulp properties (mixture composition, fibre mass concentration, freenes, etc.) by acting as low pass filters. In this capacity, they complement the action of control loops designed to attenuate variability below their cut-off frequency between 0.05 rad/s and 0.005 rad/s (Bialkowski, 1992). An indication of poor performance of existing installations was obtained by Ein Mozaffari et al. (2004) from dynamic measurements made in various mills and confirmed in a pilot scale stock chest, showing that the dynamic response of these chests is far from ideal and their inability to attenuate variability in high frequency ranges (0.01 rad/s to 0.1 rad/s) allows disturbances to pass through to further stages in the papermaking process affecting paper quality and machine run-ability. A significant improvement in the design of agitated stock chests can be obtained by developing models which take into account the actual flow field inside the vessel. Computational fluid dynamics (CFD) makes the development of these models possible. CFD has been applied in a limited number of instances to pulp suspension mixing (Bakker, 1993; Bhattacharya et al., 2008; Bhole et al., 2008; Ford et al., 2007; Ford et al., 2006; Wikstrom and Rasmuson, 1998), but the agreements between the model and available experimental evidence has not been fully established. Part of this is due to the difficulties developing a suitable model for the pulp suspension rheology and obtaining meaningful data with which to validate model predictions. The inherent complexity of mixing flows, aggravated by the complex rheology of the suspension, brings high uncertainty to the CFD predictions. Thus, it is imperative to verify the References start on page 33 Chapter 1 3 simulations results with experimental measurements. Particle image velocimetry (PIV) is an advanced optical flow measurement technique that allows a non intrusive measurement of the instantaneous velocity field over a global domain. Since pulp suspensions are opaque, PIV can not be directly implemented. However using suitable model fluids at the laboratory scale enables the application of this technique thus providing the needed data for verification of CFD simulations. 1.2 1.2.1 Literature Review Mixing Process Considerations Fluid mixing is of vital importance in most production systems in the chemical process industries. There are many different ways to provide mixing action in a vessel, but we are primarily concerned with impeller-type agitators in mixing vessels. The mixing process essentially consists of a combination of the mechanical movement of an agitator and the resultant induced deformation and flow. The mechanisms of mixing can be summarized as diffusion, convection and bulk movement. A mixing device tries to incorporate all these mechanisms together in the most economical manner. Leaving aside capital expenditure, this means achieving the desired result in the shortest time using the least power. This in turn, means achieving the right balance between diffusion, convection and bulk movement, keeping in mind the properties of the materials to be mixed and the required final product (Sweeney, 1978). References start on page 33 Chapter 1 4 To design an effective stirred tank, an efficient impeller should be chosen for the process duty. The impeller power characteristics relate the effects of fluid properties, impeller diameter, speed, and geometry to power consumption. Power consumption is independent of process performance whether the mixing process objectives are being satisfied or not (Oldshue, 1983). Once the properties of the mixing fluid are determined and impeller type and system geometry have been selected, there are three variables remaining: power, speed and impeller diameter. However, power (P), is a function of rotational speed (N), and impeller diameter (D); so there are only two independent selections. Process correlations involve any two of these variables and most frequently involve power and impeller diameter. The mixer design is finalized with the mechanical design of the shaft, impeller blade thickness, baffle thickness and supports, inlet/outlet nozzles, bearing, seals, gearbox and support structures. A determination of design, installation and cost relationships for various combinations of power, speed and diameter must include the many mechanical alternatives possible. 1.2.2 Impeller Characteristics There are hundreds of impeller types in commercial use. Selection of the most effective impeller is based on the understanding of process requirements and physical properties. Impellers can be grouped as turbines for low to medium viscosity fluids and close-clearance impellers for high viscosity fluids. Thurbine impellers are further characterized based on the flow patterns, as radial flow and axial flow. A radial-flow impeller discharges flow along the impeller radius in distinct patterns. An axial flow impeller is one in which the principal locus of flow occurs along the axis of the impeller. Axial flow impellers (marine type impeller, pitched blade turbine, References start on page 33 Chapter 1 5 Maxflo, etc.) produce high flow with low turbulence. They also produce more flow per unit power than radial impellers and are more cost effective in flow-controlled operations such as solid suspension and blending (Dickey, 2001). The greatest percentage of all pulp agitation employs this type of impeller with unique proprietary designs (Yackel, 1990). 1.2.3 Power and Power Number Early investigators used dimensional analysis to derive an equation for agitator power. Buckingham's pi theorem gives the following general dimensionless equation for the relationship of variables (Uhl and Gray, 1966). ⎡⎛ ρ D 2 N f ⎢⎜ ⎣⎢⎝ µ ⎞ ⎛ DN 2 ⎟,⎜ ⎠ ⎝ g ⎞ ⎛ P ⎞ ⎛ D ⎞ ⎛ D ⎞ ⎛ D ⎞ ⎛ D ⎞ ⎛ D ⎞ ⎛ n2 ⎞ ⎤ , ⎟,⎜ ⎟ , ⎜ ⎟ , ⎜ ⎟ , ⎜ ⎟ , ⎜ ⎟ , ⎜ ⎟⎥ = 0 5 3 ⎟ ⎜ ⎠ ⎝ ρ D N ⎠ ⎝ T ⎠ ⎝ Z ⎠ ⎝ p ⎠ ⎝ w ⎠ ⎝ l ⎠ ⎝ n1 ⎠ ⎦⎥ (1.1) where D, T, Z, c, w, p, n, l, ρ , µ , P, N and g are impeller diameter, tank diameter, liquid depth, clearance of impeller off vessel bottom, blade width, pitch of blades, number of blades, blade length, fluid density, fluid viscosity, power, impeller rotational speed, and gravitational acceleration, respectively. For geometrically similar systems Eqn. (1.1) can be reduced to: ⎡⎛ ρ D 2 N ⎞ ⎛ DN 2 ⎞ ⎛ P ⎞ ⎤ f ⎢⎜ =0 ⎟,⎜ ⎟,⎜ 5 3 ⎟⎥ ⎣⎝ µ ⎠ ⎝ g ⎠ ⎝ ρ D N ⎠ ⎦ (1.2) Equation (1.2) may be written in the following form, which defines the power number (NP): N P = K ( Re) a ( Fr )b (1.3) where Re is the impeller Reynolds number (ratio of inertial forces to viscous forces): References start on page 33 Chapter 1 Re = 6 ρ D2 N µ (1.4) and Fr is Froude number (ratio of applied forces to gravitational forces): Fr = DN 2 g (1.5) All of the power a mixer supplies to a fluid produces flow (Q) and head (H) (or shear), given by Q α ND 3 (1.6) H α N 2 D2 (1.7) Therefore, the power consumed by a mixer can be obtained by multiplying the pumping flow and the head, and is given by P = N p ρ N 3 D5 (1.8) NP decreases as inverse of Reynolds number (Re -1) through the laminar flow regime, with the slope of the descent depending upon K (in eqn. 1.3), which is a function of the system geometry and fluid properties. At high Reynolds numbers, greater than about 105, the flow is turbulent and the power number is essentially constant, i.e. NP = B. In between the laminar and the turbulent region there exists a gradual transition zone in which no simple mathematical relationship exists between the power number and the Reynolds number (Harnby et al., 1985). High flow and low turbulence is produced in agitated pulp stock chests equipped with axial flow impellers, thus they generally operate in the laminar and transition regimes. 1.2.4 Newtonian vs. Non-Newtonian Mixing The viscosity of a Newtonian liquid depends only on temperature and pressure. At constant temperature and pressure, Newtonian behaviour is characterized as follows (Brown et al., 2004): References start on page 33 Chapter 1 7 • The shear stress is the only stress generated in simple shear flow. • The shear viscosity is independent of the shear rate. • The viscosity does not vary with the time of shearing and the stress in the liquid falls to zero immediately after the shearing is stopped. • The viscosities measured in different types of deformations (e.g., uniaxial extensional flow and simple shear flow) are always in simple proportion to one another. All liquids showing deviations from the behaviour above are non-Newtonian. The viscosity of a non-Newtonian fluid is not a coefficient of the shear rate but a function of it and is called the apparent viscosity (µa): µa = τ (1.9) . γ One of the most common types of non-Newtonian fluids found in mixing applications is timeindependent fluids, in which the shear stress at any point is dependent only on the instantaneous shear rate at that point. Other non-Newtonian fluids include time dependent and viscoelastic fluids. The mixing Reynolds number (Re) for a non-Newtonian fluid is commonly determined in the laminar flow regime by the method introduced by Metzner and Otto (1957). This method assumes that a representative near-impeller average shear rate ( γ avg ) can be calculated as a function solely of the impeller rotational speed i γ avg = K s N (1.10) References start on page 33 Chapter 1 8 where N is the impeller rotational speed (measured in rotations per second) and Ks was originally noted to be a constant that depends only on impeller geometry. This shear rate is then used to estimate a near-impeller average or effective viscosity µeff = µa = µ ( γ ) = τ ≠ constant γ (1.11) For the disk flat-blade turbine used in their study, Ks = 13 was determined and their approach was found to be a workable one for pseudoplastic fluids. Calderbank and Moo-Young (1959) studied various sizes and types of impellers in Bingham, pseudoplastic and dilatant fluids, and found Ks to be 10 instead of 13 for pseudoplastic fluids. They also indicated that this method was least effective for axial flow impellers operating in the transitional flow regime. The power number in this situation was over-predicted by 20-30%, especially for highly shear-thinning fluids. In the transitional flow regime, the impeller lift force becomes significant, and the viscous forces become less important. This flow regime is distinctly different from the laminar-flow situation (Re < 130), for which the Metzner and Otto method was primarily intended. The limitations of the Metzner and Otto Method have been documented by several authors (e.g. Posarac and Watkinson, 2000). Most shear thinning fluids exhibit a yield stress and/or have zero shear viscosity. Bertrand and Tanguy (1996) modified the Metzner and Otto approach, by including the fluid yield stress in the estimation of the impeller Re and introduced the so called Bingham number given in Eqn. (1.12) which determines whether the fluid yield stress ( τ y ) would contribute significantly to the power consumption of the impeller. References start on page 33 Chapter 1 Bi = 9 τy N µeff (1.12) They (Bertrand and Tanguy, 1996) found that the effect of the fluid yield stress on Ks and i consequently γ avg through Eqn.(1.10) is negligible for Bi < 7500. Multiple studies have extended the experiments performed by Metzner and Otto (1957) to the transitional regime for different agitators (e.g. Nagata et al., 1971; Wichterle et al., 1984; Wang and Yu, 1989) and they have all concluded that as Re increases, the dependence of the shear rate on speed is more complicated. Generally, they found that the shear rate in the transitional regime is proportional to the square root of the power input per unit mass, indicating that turbulence is contributing to the generation of shear. 1.2.5 Pulp Suspension Mixing Fibre suspensions are agitated in large vessels to blend various furnish components, create uniformity in the suspension composition and to modify fibre properties. Creation of suspension uniformity is critical for downstream processing, including the runnability of paper machines and the uniformity of pulp bleaching. Pulp suspensions do not conform to any of the usual criteria of solid suspensions. There is no measurable settling velocity and there is dewatering at the surfaces. A unique measurement of the viscosity of the slurry is impossible to obtain, but there is a “pseudo” viscosity which increases tremendously with concentration. Fibre suspensions are flow sensitive and, as such, References start on page 33 Chapter 1 10 the mechanism for agitation is the same as for fluid blending. This requires an impeller with high flow characteristics and low shear (Yackel, 1990). Pulp and paper mills use both cylindrical and rectangular mixing chests. Yackel (1990) suggested that the optimum shape for an agitated stock chest in order to achieve complete motion of the suspension with minimum power input, was a cube defined by L/W = 1.0, Z/W = 1.0; where L, W and Z are chest length, width and stock height, respectively. Since all chests can not be perfect cubes, general rules of thumb exist for their design, one of which is that the L/W ratio should not exceed 1.5. Chests with L/W ratio greater than 1.5 to 2.0 require the use of two agitators. The two units would be sized as though each were in a separate chest, with W/2 used as the effective W for each mixing zone, with each mixing zone designed for an L/W = 1.0 and Z/W = 1.0 (Reed, 1995). Oldshue (1971) compared the performance of side-entering agitators and top entering agitators in stock chests. He showed that side entering agitators in a chest with Z/T greater than 0.4, are more economical than top-entering agitators for the same total mixed volume (where T is the chest characteristic length, i.e. tank diameter in cylindrical chests and tank width in rectangular chests). The performance of top entering agitators in chests with Z/T ratios between 0.4 and 0.25 was very similar to that of side entering agitators. However, this Z/T ratio is not appropriate in terms of required power for complete motion. References start on page 33 Chapter 1 1.2.6 11 Characterization of Pulp Suspensions Pulp suspensions consist primarily of water and wood fibres. Independent, single fibres exist only as free moving particles at extremely low concentrations or intense shear fields. Different incoherent loose fibre bundles form at slightly higher fibre concentrations. These are usually unstable and can easily disentangle and reform when sheared. As fibre population is increased, coherent networks are formed. The forces at the fibre contact points create fibre aggregates or flocs that give mechanical strength to the network. Before motion can be initiated within a suspension, external forces sufficient to overcome the network forces must be applied (Bennington, 1988). Fluid and fibre-floc-network interactions are more complex than in particulate slurry systems. The mechanisms of shear depend on flow rate, fibre concentration, and the state of the suspension, which in turn depends on the fibre population, fibre type and the pulping process used to make the fibres. Fibres suspensions are highly structured and inhomogeneous. They are unlike conventional non-Newtonian fluids because the mechanistic sliding-layer, laminar-flow concept never occurs throughout the suspension at concentrations above the sedimentation concentration, although it may occur locally at fibre suspension/solid interfaces (Duffy, 2003). Many authors have reported different measurements of pulp network strength. Despite the differences among these results, they all show a marked dependency on consistency that can be represented by τ = aC mb (1.13) where τ is the strength of the pulp network in Pascal, Cm the mass concentration and a and b fitted constants. These parameters depend on many factors, including test procedure, the pulp References start on page 33 Chapter 1 12 type and degree of pulp treatment (beating) prior to testing. Values of the constant a between 1.18 and 24.5 and values of the exponent b between 1.25 and 3.02 have been reported in the literature (for Cm units of wt.% in eqn. 1.13)(Bennington et al., 1990). In suspension pipe flow, the yield stress results from the minor component of inhomogeneous strained fibres and flocs pressing against the pipe wall. Thus τ y depends strongly on surface roughness, fibre concentration and fibre type, and it is not a true total suspension property as conventionally described in non-Newtonian yield stress slurries. For these reasons, plus the fact that the suspension is not in laminar shear but it is rather in plug flow; it is unlikely that Bingham plastic models or similar homogeneous non-Newtonian models can adequately describe the highly structured mechanisms of flow of pulp fibre suspensions in pipes (Duffy, 2003). 1.2.7 Design of Agitated Stock Chests The design of stock chests is concerned largely with selecting the power necessary to ensure complete motion throughout the chest volume which can be very difficult to achieve particularly at higher mass concentrations. Suspension yield stress must be considered in mixer design. The impeller power for suspension agitation has been correlated with the apparent viscosity of the suspension, assuming Bingham plastic behaviour (Blasinski and Rzyski, 1972; Gibbon and Attwod, 1962): . µa = τ y + µp γ (1.14) . γ . where τ y is the yield stress of the suspension, µp the plastic viscosity, and γ the shear stress. References start on page 33 Chapter 1 13 It has been common practice in the pulp and paper industry to use the viscosity of water when dimensioning equipment where the suspension is exposed to high shear rates (Blasinski and Rzyski, 1972; Gibbon and Attwood, 1962; Head and Durst, 1957; Wikstrom and Rasmuson, 1998). Consequently, by applying the Metzner and Otto correlation, and assuming that the plastic viscosity is close to that of water, Eqn. (1.14) has been implemented as: µa ≈ τy (1.15) 10 N which can be used to calculate the impeller Reynolds number (Blasinski and Rzyski, 1972) ρ D 210 N 2 ⎛ ND ⎞ Re = = 10 ⎜ ⎟ τy ⎝ m ⎠ 2 (1.16) where m is a commonly used parameter that defines the resistance of fibre suspensions to shear (Head and Durst, 1957; Gibbon and Attwood, 1962). Considering suspension yield stress in the mixer design is crucial especially at higher mass concentrations. Yield stress fluids create caverns of effective motion exclusively around the impeller and stagnation elsewhere in the mixing vessel. (Amanullah et al., 1998; Elson et al., 1986; Solomon et al., 1981; Wichterle and Wein, 1975; Wichterle and Wein, 1981). The cavern boundary has been defined using different approaches. Following the force balance approach introduced by Wichterle and Wein (1981), the shear stress along the cavern boundary is equal to the yield stress of the fluid. A more complete force balance approach was introduced by Amanullah et al. (1998) which included both tangential and axial shear components in the determination of the total momentum imparted by an axial impeller. These forces are transported to the cavern boundary by the pumping action of the impeller. The model proposed for the size of the cavern is given by References start on page 33 Chapter 1 14 2 Re y ⎛ 2 4 N P ⎞ ⎛ Dc ⎞ ⎜ D ⎟ = π ⎜ N f + 3π ⎟ ⎝ ⎠ ⎝ ⎠ where N f = Fa ρN D 2 4 (1.17) 2 2 is the dimensionless axial force number and Rey = ρ N D τ y is called the yield stress Reynolds number. Cavern formation in pulp suspensions using side entering axial flow impellers has been addressed in a cylindrical vessel by Hui et al. (2009) and a modification to the model of Amanullah et al. (1998) was provided in order to account for the interaction between the cavern boundary and the vessel walls. Current design of industrial chests is based on proprietary criteria developed from experience. One common design method widely used in the industry and summarized by Yackel (1990) matches the momentum flux generated by an impeller with that needed to provide complete surface motion in the vessel, which has been correlated with a range of variables, including the chest geometry, impeller location, type of fibre and concentration and desired retention time. By determining the required momentum flux for a particular installation using the proposed graphical correlations, a suitable impeller may be selected provided the impeller momentum flux for the impeller is known. Scaling of agitated stock chests was addressed by various authors (Gibbon and Attwood, 1962; Yackel, 1990). One of the criteria presented is known as “level momentum”, defined as LMo = Mo V 2 3 = CN 2 D 4 V 2 3 (1.18) References start on page 33 Chapter 1 15 where V refers to the fluid volume, and C is a constant that depends on the impeller geometry and type. Based on this criterion, level momentum must be equal for the model and the prototype. 1.2.8 Mixing Dynamics in Agitated Chests The dynamics of pulp chests have been modeled assuming ideal mixing and first order behaviour (Walker and Cholette, 1958). This allows selecting the chest volume based on anticipated disturbances and degree of attenuation required. However, these studies did not consider nonideal suspension flow that can create short-circuiting and dead zones in the chest. Ein-Mozaffari (2002) studied the dynamics of pulp agitation using a 1:11 scale-model Plexiglas chest equipped with a side-entering impeller. It was found that even when complete surface motion was attained, considerable channeling, recirculation and stagnation existed below the suspension surface. The extent of non-ideal flow can be significant. At high input flow rates and low impeller speeds, the percentage of channeling was found to be as high as 90% of the pulp feed. Even when the scale chest was operated under conditions where the entire suspension in the tank was in motion, the percentage of channeling and dead volume was still 26% and 16%, respectively. Dynamic studies made in industrial chests confirmed these findings (Ein-Mozaffari et al., 2004). A grey-box dynamic model was developed by Ein-Mozaffari (2002) based on the observed nonideal behaviour of industrial pulp chests (Fig. 1.1(a)). The model included a well mixed zone (represented by a first order transfer function having a time constant τ2 and time delay T2) as well References start on page 33 Chapter 1 16 as the possibility of stock recirculation, R, within that mixing zone. A fraction of the entering flow, f, could by-pass the primary mixing zone and pass directly to the exit, although it could be involved in some limited mixing. Thus, a first order transfer function with time constant and delay (τ1 and T1) was included for this portion of the flow. Later on, this model was modified by using a single time delay (Td) to represent the contribution of the plug flow through the chest, and the recirculation effect was included in the mixing time constant (Fig. 1.1(b)). In this manner, the model became simpler yet more robust in providing an easier identification of the mixing parameters (Sotanzadeh et al., 2009). The system was excited in order to estimate the model parameters. Since the choice of excitation signal has a substantial influence on the accuracy with which the parameters could be estimated, a frequency modulated random binary signal was designed to provide excitation of the system at the critical chest frequencies (EinMozaffari, 2002). f u(s) G 1′ ( s ) = 1 1 + τ 1s y(s) e −Td s Mixing Zone 1-f (a) Figure 1.1. G ′2 ( s ) = 1 1 + τ 2s (b) Dynamic model of the agitated chest (a) Ein-Mozaffari, 2002 and (b) Soltanzadeh et al., 2009 1.2.9 Computational Fluid Mixing Modeling a stirred tank using computational fluid dynamics (CFD) requires consideration of many aspects of the process. First, the computational model requires that the volume occupied by the fluid in the vessel (domain), be described by a computational grid. In these cells the References start on page 33 Chapter 1 17 problem variables are computed and stored. The governing equations of the flow occurring in the mixing tank (continuity, momentum, energy, and species) need to be discretized and then solved in the computational grid which must fit the contours of complex geometries such as pitched impeller blades of different shapes depending on the particular impeller used. Detailed information regarding the discretization of the computational domain and the governing equations, as well as the numerical solution approach is available in multiple CFD references that have been summarized in my previous work (Ford, 2004). 1.2.10 Stirred Tank Modeling The problem of modeling stirred tanks has been dealt with by excluding the impeller region from the computational domain and prescribing, on the basis of simple models and/or experimental information, the flow field near the impeller as a stationary boundary condition (Brucato et al., 1998; Fokema et al., 1994; Ranade et al., 1989). These methods are called black box methods. As an alternative, the impeller has been modeled as a distributed source of momentum and turbulence quantities (Hutchings et al., 1989). The main limitation of the above approaches is their lack of generality. The near-field flow produced by a given impeller may depend on the overall tank geometry and the fluid properties, so that the impeller can not be uniquely characterized. Also, different choices of impeller boundary conditions lead to significantly different predicted flow fields. A number of approaches have been proposed which allow the explicit simulation of the whole flow field without any recourse to empirical data. The most popular ones are the rotating frame References start on page 33 Chapter 1 18 and multiple reference frames models (Luo et al., 1994) which are steady state and the sliding mesh model which is time dependent (Murthy et al., 1994). The rotating frame model solves the momentum equations for the entire domain in a rotating frame. Typically, the angular velocity of the primary rotating component is used as the angular velocity of the frame. The frame is assumed to rotate with the impeller (the impeller is at rest in the rotating frame). This simple steady-state model can not handle motion of elements such as baffles into or through the fluid. A modification of the rotating frame model is the multiple reference frames (MRF) model (Luo et al., 1994), which has been used in this work. More than one rotating or non rotating reference frame can be used in the simulation allowing for the modeling of baffled tanks and tanks with inflow and outflow ports. The impeller is at rest in the rotating frame and the tank walls and baffles are at rest in the stationary frame. The grid used for an MRF solution must have a perfect surface of revolution surrounding each rotating frame (Bakker and Marshall, 2004). In the implementation of the MRF approach, the calculation domain is divided into subdomains, each of which may be rotating with respect to the inertial frame. The governing equations in each subdomain are written with respect to that subdomain’s reference frame. Thus, the flow in stationary subdomains is governed by the continuity and momentum equations as follows: → ∂ρ + ∇ ⋅ (ρ υ ) = 0 ∂t (1.19) → →→ = → → ∂ ( ρ υ ) + ∇ ⋅ ( ρ υ υ ) = −∇p + ∇ ⋅ (τ ) + ρ g + F ∂t (1.20) References start on page 33 Chapter 1 19 = → → where p is the static pressure, τ is the stress tensor (described below), and ρ g and F are the gravitational force and external body forces respectively. The stress tensor is given by = ⎡⎛ →T → ⎞ 2 ⎠ 3 → ⎤ τ = µ ⎢⎜ ∇ υ + ∇ υ ⎟ − ∇ ⋅ υ I ⎥ ⎣⎝ (1.21) ⎦ where µ is the molecular viscosity, I is the unit tensor, and the second term on the right hand side is the effect of volume dilation. The flow in the rotating subdomain is governed by the “modified”continuity and conservation of momentum equations (additional terms include the acceleration of the rotating frame) as follows: → ∂ρ + ∇ ⋅ ( ρυr ) = 0 ∂t (1.22) → → → → → = → → ∂ ( ρ υ ) + ∇ ⋅ ( ρυr υ ) + ρ (Ω × υ ) = −∇p + ∇ ⋅ (τ ) + ρ g + F ∂t (1.23) → υr refers to the relative velocity and can be related to the absolute velocity by the following equation: → → → → υ r = υ − (Ω × r ) (1.24) → → Here, Ω is the angular velocity vector (that is, the angular velocity of the rotating frame) and r is the position vector in the rotating frame. → The continuity of the absolute velocity, υ , is enforced at the boundaries between subdomains. In this way, flow variables from one subdomain can be used in the calculation of fluxes at the boundary of the adjacent subdomain. References start on page 33 Chapter 1 20 The sliding mesh model is a time-dependent approach in which the grid surrounding the rotating components physically moves during the solution. The velocity of the impeller relative to the moving mesh region is zero, as is the velocity of the tank and other internals in the stationary mesh region. The grid moves in small discrete steps after which the conservation equations are solved until convergence is reached, and then the grid moves again. Information is passed through the interface from the rotating to the stationary regions. This is the most rigorous and computationally expensive solution method. Transient simulations using this model can capture low frequency oscillations in the flow field in addition to the ones that result from impeller baffle interaction (Bakker and Marshall, 2004). 1.2.11 CFD Modeling of Pulp Suspension Mixing Bakker and Fasano (1993) were the first to model the pulp flow in a stock chest with a side entering impeller using CFD. The impeller was simulated by defining a flat velocity profile at the impeller location and a combination of turbulent and laminar flow regimes were used in the tank. The rheological model used for the pulp suspension was not given, and the simulation results were only compared to visual observations of the flow field, with which the agreement was judged to be satisfactory. Roberg (1997) simulated the flow of paper pulp in a cylindrical stirred tank equipped with a sideentering mixer assuming laminar flow and power law rheology for the pulp suspension. The impeller was simulated as source of momentum. The results of the simulations made sense qualitatively but there were not compared against experimental measurements. References start on page 33 Chapter 1 21 The agitation of pulp suspension with a jet nozzle agitator was modeled in CFD by Wikstrom and Rasmuson (1998). The impeller was modeled by imposing boundary conditions on the external flow, the pulp suspension was described as Bingham fluid, and flow inside the chest modeled as being laminar. The flow field obtained from the CFD calculations increasingly deviated from the velocities measured as the distance from the impeller increased. Including a shear thinning parameter in the suspension rheology did not present significant improvement over the Bingham model and they suggested that the models were too rigid to describe the real behaviour of the suspension. In my previous work (Ford, 2004), a CFD model of a rectangular pulp-stock mixing chest was formulated with the rheology of the pulp suspension approximated using the Bingham plastic model. A specific user defined procedure was developed to determine the dynamic response of the chest from the numerical solution of the chest’s velocity fields and compare it with dynamic data obtained experimentally under identical conditions. Using this technique, it was shown that the numerical model captured the mixing dynamics of the scale-model chest fairly well, although it over-estimated the exponential extent of mixing in the chest. This was particularly noticeable in flow situations containing significant bypassing where the model over-estimated the extent of mixing in the bypassing flow as defined in the dynamic model shown in Fig 1.1(a). Departure of the simulated and measured dynamics was attributed to the rheology of the suspension which is not fully described by the Bingham model. The inherent complexity of mixing flows coupled with the complicated geometry of mixing systems and the intricate suspension rheology, bring uncertainty into the CFD predictions. Thus there exists high interest in experimental verification of the simulations. References start on page 33 Chapter 1 22 1.2.12 Particle Image Velocimetry Principles and Application Particle image velocimetry (PIV) is an optical method used to obtain instantaneous velocity measurements over a global domain. The use of PIV in stirred tanks is very recent though it has been used by various researchers to study flow patterns in different planes of the stirring vessel and validate CFD simulations (e.g. Bakker et al., 1996; Brown et al., 2004, Khopar et al., 2004). PIV is a velocimetry technique based on images of traces “seeding” particles suspended in the flow. Ideally, these particles should be perfect flow tracers, homogeneously distributed and their presence should not alter the flow properties. The local fluid velocity can be obtained by measuring the fluid displacement from multiple particle images and dividing that displacement by the time interval between the exposures. To get an accurate instantaneous flow velocity, the time between exposures should be small compared to the time scales in the flow; and the spatial resolution of the PIV sensor should be small compared to the length scales in the flow. The principal layout of a modern PIV system is shown in Fig. 1.2. The PIV measurement includes illuminating a cross-section of the flow field by a pulsing light sheet, recording multiple images of the seeding particles in the flow using a camera located perpendicular to the light sheet, and analyzing the images for displacement information (Raffel et al., 1998). Commonly a pulsed Nd:YAG laser (Neodymium Yttrium Aluminum Garnet) is used as the light source because of its high light intensity. Pulsed lasers need some time to build up energy before they can deliver a new pulse and the two images in a PIV image pair must be taken within a quite short period of time. Thus, it is common to use a laser with two cavities. The laser pulses have a duration time of 5-10 ns and the energy in one pulse can be up to 400mJ. Nd:YAG lasers emit References start on page 33 Chapter 1 23 light with a wavelength of 1064nm which is in the infrared range. For PIV purposes light with this wavelength is not very useful since most cameras have their maximum sensitivity in the blue-green part of the spectrum. Thus, the wavelength of the Nd:YAG laser is halved using a device called an harmonic generator, so that it becomes 532nm. PIV puts special demands on the camera to be used. For this purpose a camera with progressive scan architecture is utilized. The detailed description of the architecture and functioning of this charged couple device (CCD) can be found elsewhere (Raffel et al., 1998; Dantec dynamics Inc., 2003). The two exposures can be taken either as a double exposure of one image or as two different images. The double exposure technique is based on applying auto-correlation techniques to obtain the velocity information, resulting in a directional ambiguity (since it is impossible to tell whether an imaged particle was in the first or second exposure), thus the flow can only be in one direction when this method is used. Hence, the method with two images and cross-correlation analysis is more commonly used. When the cross-correlation of images has been performed, a measure of the displacement is found by detecting the location of the highest correlation peak. This measurement provides both the magnitude and direction of the flow with no ambiguity (Raffel et al., 1998; Adrian, 2005). The cross-correlation function is not calculated on the whole images but on smaller parts of these called interrogation areas. Most commonly used interrogation area dimensions are 64 x 64 pixels for autocorrelation, and 32 x 32 pixels for cross-correlation analysis. The maximum particle displacement between light pulses should be about one fourth of these dimensions, in References start on page 33 Chapter 1 24 order to achieve reasonable accuracy and dynamic range for PIV measurements (Raffel et al., 1998). Additionally, the sampled particle images should satisfy the Nyquist sampling criteria in order to measure the particle displacements with subpixel accuracy and thus enhance the dynamic range of the PIV measurements (Adrian, 2005). The aim of the cross correlation is to find the distance that the particle pattern has moved during the inter image time and translate this into a velocity measurement. The relation between velocities (u) and particle displacements (d) is simply u= d M ∆t (1.25) where M is the magnification and ∆t is the inter image time. The cross correlation calculated on one interrogation area results in one velocity vector. The cross-correlation of images can be seen as finding which relative displacement of the interrogation areas gives the best pattern match. This displacement should be proportional to the average velocity in the interrogation areas. Fast Fourier transform (FFT) techniques are used for the calculation of the correlation functions. Since the images are digitized, the correlation values are found for integral pixel values, with an uncertainty of ± 0.5 pixel. Different techniques such as Gaussian and parabolic fits have been used to estimate the location of the correlation peak with increased accuracy (e.g. ± 0.01 pixel in 8 bit digital cameras). Digital windowing and filtering techniques can be used prior to calculation of the correlation peaks in order to improve the results (Dantec dynamics Inc, 2003). References start on page 33 Chapter 1 25 Double pulsed laser Laser sheet Velocity field Field of view Seeded flow Imaging optics CCD Correlation ∆t IA Image frame 1 Image frame 2 IA Particle Images Figure 1.2. Operation of the PIV system The majority of PIV systems today are equipped with processing software to perform crosscorrelation of the two images, data validation techniques and calculation of statistical properties which include the mean of each velocity component, standard deviation of the mean and the covariance coefficient. References start on page 33 Chapter 1 1.3 26 Research Objectives To assess the capabilities and address the limitations of using CFD to model macro-scale mixing of pulp fiber suspensions. The research was divided in three areas with specific objectives as follows: I. CFD modeling of pulp suspension mixing using a Bingham plastic approximation and experimental dynamic tests for model evaluation: (i) To expand the range of operating conditions simulated with CFD for the pilot-scale mixing chest (Ford, 2004) and assess the model performance based on comparison with experimental dynamic tests (Ein-Mozzafari, 2002) (addressed in Chapter 2). (ii) To extend this CFD modeling approach to industrial mixing chest configurations and assess the accuracy of the simulation results using industrial data collected at the process operating conditions (addressed in Chapter 3). II. CFD modeling and PIV measurements of mixing of a viscous-Newtonian fluid in a geometrically-scaled mixing vessel: (iii) To measure fluid flow in a laboratory-scale mixing chest using PIV and a viscous Newtonian model fluid with the velocity measurements used to verify the accuracy of CFD simulations made under identical conditions (addressed in Chapter 4). III. Implementation of non-Newtonian model fluids for studying mixing of pulp fibre suspensions in the laboratory scale: References start on page 33 Chapter 1 27 (iv) To characterize the axial impellers used in pulp suspension agitation (experimentally and computationally) to determine the dependence of scaling parameters (NP and NF) on pulp type, suspension mass concentration, impeller rotational speed and vessel geometry (addressed in Chapter 5). (v) To identify appropriate scaling criteria for use in selection of a model non-Newtonian fluid and impeller operating conditions for reproduction of important mixing phenomena identified in industrial and pilot-scale pulp mixing chests (addressed in Chapter 6). 1.4 Thesis Outline Computational Fluid Dynamics (CFD) is frequently used together with experimental techniques to understand mixing performance. However, the accuracy of CFD solutions in pulp suspension mixing operations has not been established, primarily due to the difficulty obtaining sufficient measurements with which to verify the computational results. In this thesis, issues limiting the accuracy of CFD as applied to pulp suspension mixing operations are investigated. In the following, (sections 1.4.1−1.4.5) a brief outline of the manuscripts presented in chapters 2 to 6 of this thesis is presented. 1.4.1 Modeling a Pilot Scale Pulp Mixing Chest using CFD Extensive tests made on a 1:11 scale-model pulp chest demonstrated that the dynamic performance of agitated pulp stock chests was far from ideal, with a propensity for significant channeling, recirculation and stagnation flows. In particular, pulp suspension channeling was found to increase dramatically and suddenly for small changes in mixing conditions (including References start on page 33 Chapter 1 28 the suspension consistency, impeller speed or pulp production rate). As such minor changes in operating conditions potentially have dramatic consequences for industrial processes; we attempted to model the onset of this non-ideal behavior using computational fluid dynamics (CFD). Previous work had successfully modeled the mixing dynamics of the scale-model mixing chest using a commercial CFD software package (Fluent) with the pulp suspension approximated as a Bingham plastic (Appendix A). In this chapter, the ability of the simulation to predict the observed abrupt transition to severe channeling was investigated. 1.4.2 Computational Modeling of Industrial Pulp Stock Chests Agitated pulp stock chests are the most widely used mixers in pulp and paper manufacture. In the work presented here, the approach taken in chapter 2 was applied to model two industrial pulp stock chests. The first chest is rectangular, agitated using a single side-entering impeller, and feeds a mixture of chemical pulps at 3.5% mass concentration (Cm) to a paper machine. The second chest has rectangular geometry, with a mid-feather wall used to direct suspension flow through a U-shaped trajectory past four side-entering impellers. This chest is used to remove latency from a Cm = 3.5% thermo mechanical pulp suspension ahead of stock screening. Steady-state simulations were made corresponding to process conditions during mill tests. The calculated steady-state flows were then used to determine the dynamic response of the virtual chests, compared with experimental measurements and found to agree reasonably well. This chapter discusses difficulties encountered in characterizing the mixing (both experimentally and computationally) and the limitations of the industrial data. References start on page 33 Chapter 1 1.4.3 29 Investigation of the Flow Field in a Rectangular Vessel Equipped with a Side- Entering Agitator CFD simulations are commonly verified against force measurements made on the impeller or by using the dynamic response of the system (as used in the previous chapters). However, these techniques cannot fully validate the simulations. A number of factors contribute to errors in the computational results, and are related to a number of parameters including the geometrical complexity and/or numerical approaches implemented in the model, and the definition of the pulp suspension flow properties. In this chapter, the accuracy and efficiency of CFD simulations made on a scaled-down agitated stock chest when operated with a fluid of known rheology (a Newtonian viscous solution of glycerin), was assessed. Particle image velocimetry (PIV) was used to obtain velocity measurements over representative regions of the vessel and compared with the computational results. Very good agreement (an average of 13.1% RMS deviation between computed and measured local axial velocities) was found which indicated that the limiting factors in the proficient application of CFD for pulp suspension mixing are mostly related to rheological complexities (i.e. proper non-Newtonian viscosity definition and the implications of non-Newtonian rheology models on the computations) 1.4.4 Characterization of Axial Flow Impellers on Pulp Suspensions Two axial flow impellers commonly used in pulp agitation applications were characterized in hardwood and softwood low-consistency pulp fibre suspensions. The impellers operated in the laminar and transition-to-turbulence regimes with the suspension mass concentration significantly affecting both power and axial thrust numbers. Axial force numbers could be collapsed to a single operating curve using the yield stress Reynolds number, but the power References start on page 33 Chapter 1 30 number remained a function of suspension properties. CFD was used to model impeller flow using a Bingham approximation to describe the suspension rheology. Model agreement with the experimental measurements was generally good, with the values of NP and Nf determined to within 24% and 16%, respectively. The error was attributed to the uncertainty in the rheological characterization of the pulp suspension, particularly the yield stress. 1.4.5 Carbopol as a Model Fluid for Studying Mixing of Pulp Fibre Suspensions The selection of operating conditions for a geometrically-scaled laboratory pulp stock chest in which aqueous carbopol solutions were chosen as the model fluid is detailed here. These transparent solutions will permit the use of PIV to measure velocity profiles throughout the vessel and thus enable validation of computational models. To ensure that the carbopol solutions mix in a similar manner pulp fibre suspensions at the new scale it is necessary to select a mixing criterion for conservation. The dimensionless cavern size was chosen as the conserved quantity, although the implications of scaling based on other criteria were also examined. The pulp fibre suspensions and carbopol solutions were modeled as Herschel-Bulkley fluids, based on laboratory measurements. 1.5 • Contributions of this Work The limitations of the CFD model developed for the pilot-scale mixing chest were identified. Despite the good agreement found between experimental and the predicted system dynamics under moderate agitation and flow rates through the chest, the CFD solutions failed to predict the measured dynamic behaviour over the full spectrum of operating conditions. This analysis has again raised the awareness for extensive References start on page 33 Chapter 1 31 verification of CFD simulations prior to their implementation for process design and optimization, thus determining the direction of the present work in this area. • The simulation of the dynamics of industrial pulp mixing chests from CFD results was developed. This provides a valuable tool for evaluation of CFD simulations when velocity measurements are not available for validation purposes. Moreover, linking CFD to the dynamics of the mixing process has provided important information about the internal flows that lead to a particular system response following a change in input conditions. This work has directed us to re-visit generally adopted mixer design methods that rely on an ideal dynamic response of the system. • To the author’s knowledge, the simultaneous application of PIV and CFD to study mixing in side-entering, axially agitated rectangular vessels was presented for the first time. This has expanded the number of industrially relevant configurations that can be effectively addressed using these advance numerical and experimental tools. The tuning of the PIV measuring parameters described as part of this study provides a significant contribution to effective implementation of PIV in the study of three dimensional mixing flows. • The adequate characterization of impellers used in complex fluids has a significant effect in the assessment of mixer performance. Carrying out characterization experiments for each different impeller configuration in a complex fluid or suspension can be expensive and time consuming. Therefore, the presented implementation of CFD simulations for the characterization of impellers provides a powerful modeling tool accessible to most engineers involved with mixing studies. This work constitutes fundamental research that References start on page 33 Chapter 1 32 is fully transferable to many other industrial mixing applications that involve fluids with complicated rheology. In particular, accounting for the proper functionality of NP and NF on Re and fluid properties significantly improves the predictions of the existing cavern correlations for agitation of shear thinning and/or yield bearing fluids. • The use of simulant fluids is common practice in the mixing laboratory. Nevertheless, the adequate selection of the operating conditions of laboratory-mixers employing model fluids has not been addressed in the open literature. Thus, the protocol provided on the proper use of carbopol as a model fluid to study pulp suspension mixing represents a significant contribution to the mixing field. The method presented is not exclusive to pulp suspension agitation, and can be extended to a wide range of industrial applications. References start on page 33 Chapter 1 1.6 33 References Adrian, R.J., 2005. Twenty years of particle image velocimetry. Exp. Fluids, 39, 159−169. Amanullah, A., Hjorth, S.A. and Nienow, A.W., 1998. A new mathematical model to predict cavern diameters in highly shear thinning, power law liquids using axial flow impellers. Chem. Eng. Sci., 53, 3, 455−469. Bakker, A., Fasano J.B., 1993. A computational study of the flow pattern in an industrial paper pulp stock chest with a side-entering impeller. AIChE Symposium Series, 89 293, 118−124. Bakker, A., Myers, J., Ward, R.W. and Lee, C.K., 1996. The laminar and turbulent flow pattern of a pitched blade turbine. Trans. Inst. Chem. Eng., 74, 485−491. Bennington, C.P.J., 1988. Mixing pulp suspensions. PhD dissertation, University of British Columbia, Vancouver, Canada. Bennington, C.P.J., Kerekes, R.J. and Grace, J.R., 1990. The yield stress of fibre suspensions. Can. J. Chem. Eng., 68, 10, 748−757. Bertrand, F. and Tanguy, P., 1996. A new perspective for the mixing of yield stress fluids with anchor impellers. J. Chem. Eng. Jpn., 29, 1, 51−58. Bhattacharya, S., Gomez, C., Soltanzadeh, A., Taghipour, F., Bennington, C.P.J. and Dumont, G.A., 2008. Computational modelling of industrial pulp stock chests. Can. J. Chem. Eng., submitted for publication. Bhole, M., Ford, C. and Bennington, C.P.J., 2009. Characterization of axial flow impellers in pulp fibre suspensions. Che. Eng. Res. Des., 87, 4, 648−653. Bialkowski, W.L., 1992. Newsprint variability and its impact on competitive position. Pulp Pap. Can., 93, 11, T299−T306. References start on page 33 Chapter 1 34 Blasinski, H. and Rzyski, E., 1972. Mixing of non-newtonian liquids. Power consumption for fibrous suspensions. Inz. Chem, 2, 1, 169−182. Brucato, A., Ciofalo, M., Grisfi, F. and Micale, G., 1998. Numerical prediction of the flow fields in baffled stirred vessels: A comparison of alternative modeling approaches. Chem. Eng. Sci., 53, 21, 3653−3684. Calderbank, P.H. and Moo-Young, M.B., 1959. The prediction of power consumption in the agitation of non-newtonian fluids. Trans. Inst. Chem. Eng., 37, 26−33. Dantec Dynamics Inc., 2003. Dantec dynamics flow map PIV installation and user’s guide: General PIV reference. Duffy, G., 2003. The significance of mechanistic-based model in fibre suspension flow. Nord. Pulp Paper Res., 18, 1, 74−80. Ein-Mozaffari, F., 2002. Macroscale mixing and dynamic behaviour of agitated pulp stock chests. PhD dissertation, University of British Columbia, Vancouver, Canada. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2004. Dynamic mixing in agitated industrial pulp chests. Pulp Pap. Can., 105, 5, T 118−T 122. Elson, T.P., Cheesman, D.J. and Nienow, A.W., 1986. X-ray studies of cavern sizes and mixing performance with fluids possessing a yield stress. Chem. Eng. Sci., 41, 2555−2562. Fokema, M.D., Kresta, S.M. and Wood, P.E., 1994. Importance of using the correct impeller boundary conditions for cfd simulations of stirred tanks. Can. J. Chem. Eng., 72, 177−183. Ford, C., 2004. CFD simulation of mixing dynamics in agitated pulp stock chests, M.A.Sc. thesis, University of British Columbia, Vancouver, Canada. References start on page 33 Chapter 1 35 Ford, C., Bennington, C.P.J. and Taghipour, F., 2007. Modelling a pilot-scale mixing pulp mixing chest using CFD. J. Pulp Pap. Sci., 33, 115−120. Ford, C., Ein-Mozaffari, F., Bennington, C.P.J. and Taghipour, F., 2006. Simulation of mixing dynamics in agitated pulp stock chests using CFD. AIChE J., 52, 10, 3562−3569. Gibbon, J.D. and Attwood, D., 1962. Prediction of power requirements in the agitation of fibre suspensions. Trans. Inst. Chem. Eng., 40, 75−82. Head, V.P. and Durst, R.E., 1957. Stock slurry hydraulics. Tappi J., 40, 12, 931−936. Hui, L.K, Bennington, C.P.J., and Dumont, G.A., 2009. Cavern formation in pulp suspensions using side-entering axial flow impellers. Chem. Eng. Sci., 64, 509−519. Hutchings, B.J., Weetman, R.J. and Patel, B.R., 1989. Computation of flow fields in mixing tanks with experimental verification. Paper TN - 4810, ASME Annual Meeting, San Francisco. Khopar, A., Aubin, J., Ruibo-Atoche, C., Xuereb, C., Sauze, N.L., Bertrand, J. and V, V., 2004. Flow generated by radial flow impellers: PIV measurements and CFD simulations. Int. J. Chem. Reactor Eng., 2, A18, 1−16. Luo, J.Y., Gosman, A.D. and Issa, R.I., 1994. Prediction of impeller induced flow in mixing vessels using multiple frames of reference. Inst. Chem. Eng. Symposium Series, 136, 549−556. Metzner, A.B. and Otto, R.E., 1957. Agitation of non-newtonian fluids. AIChE J., 3, 1, 3−10. Murthy, J.Y., Mathur, S.R. and Choudhury, D., 1994. CFD simulations of flows in stirred reactors using a sliding mesh technique Inst. Chem. Eng. Symposium Series, 136, 341−348. Nagata, S., M. Nishikawa, H. Tada and S. Gotoh, 1971. Power consumption of mixing impellers in pseudo-plastic liquids, J. Chem. Eng. Jpn., 4, 72−76. References start on page 33 Chapter 1 36 Oldshue, J.Y., 1971. Agitation in groundwood blending. Pulp Pap. Can., 72, 12, T375−T380. Posarac, D. and Watkinson, P., 2000. Mixing of a lignin based slurry fuel. Can. J. Chem. Eng., 78, 265−270, Raffel, M., Willert, C. and Kompehans, J., 1998. Particle image velocimetry. A practical guide, Adrian, R.J., et al., Springer, Berlin, 1998. Ranade, V., Joshi, J.B. and Marathe, A.G., 1989. Flow generated by pitched blade turbines 2: Simulation using κ − ε model. Chem. Eng. Commun., 81, 197−224. Reed, C.S., 1995. Selecting the right equipment for agitation and blending, part 1. Tappi J., 78, 6, 252−254. Roberg, K.E., 1997. Numerical simulations of the flow of fibre supsension in a side-entering mixer. Récents Progrès en Génie de Procédés, 11, 51, 203−210. Solomon, J., Elson, T.P., Nienow, A.W. and Pace, G.W., 1981. Cavern size in agitated fluids with a yield stress. Chem. Eng. Commun., 11, 143−164. Soltanzadeh, A., G.A., Dumont, S., Bhattacharya and C.P.J., Bennington, 2009, Estimation of residence time and channeling in pulp mixers, in press Nordic J. Pulp and Paper Sci. Wang, K. and S. Yu, 1989. Heat transfer and power consumption of non-Newtonian fluids in agitated vessels, Chem. Eng. Sci., 44, 33−40. Walker, O.J. and Cholette, A., 1958. Determination of the optimum size and efficiency of stock chests. Part 1 : The ideal chest. Pulp Pap. Can., 59, 113−117. Wichterle, K. and Wein, O., 1975. Agitation of concentrated suspensions. 5th Int. Congr. Chem. Eng. Proc., Chem. Equip. Des. Autom., B4.6. References start on page 33 Chapter 1 37 Wichterle, K., M. Kadlec, L. Zak and P. Mitchka, 1984. Shear rates on turbine impeller blades, Chem. Eng. Commun., 26, 25−32. Wichterle, K. and Wein, O., 1981. Threshold of mixing of non-newtonian liquids. Int. Chem. Eng., 21, 116−120. Wikstrom, T. and Rasmuson, A., 1998. The agitation of pulp suspension with a jet nozzle agitator. Nord. Pulp Paper Res., 13, 2, 88−94. Yackel, D.C., 1990. Pulp and paper agitation: The history, mechanics and process, TAPPI press, Atlanta, 1990. References start on page 33 38 Chapter 2 Modeling a Pilot Scale Pulp Mixing Chest Using CFD 2.1 Introduction Effective mixing is critical in many pulp and paper process operations, including blending of pulp stocks, addition of wet-end chemicals and consistency control prior to papermaking, and steam addition and chemical contacting in pulp bleaching operations. Agitated pulp stock chests function as low-pass filters to reduce high-frequency variability in pulp properties, such as the suspension mass concentration and freeness, ahead of many paper making operations. This complements the action of control loops which can only control low-frequency variability below the loop cut-off frequency (Bialkowski, 1992). In spite of the importance of these applications, the problem of designing and scaling up agitated stock chests has been tackled mainly by means of semi-empirical methods. One common design method has been summarized by Yackel (1990) and is based on matching the momentum flux generated by an impeller with that needed to provide complete motion in the chest. Yackel uses the term ‘complete motion’ to indicate when the pulp suspension is in motion over the entire top surface of the mixing vessel (what we refer to as ‘complete surface motion’). Yackel’s definition agrees with that used for previous design of paper stock chests (Oldshue and Gretton, 1956). Ein-Mozzafari et al.(2003) studied the dynamics of pulp agitation using a 1:11 scale-model Plexiglas chest equipped with a side-entering impeller. They found that even when complete A version of this chapter has been published. Ford† C., Bennington C.P.J. and Taghipour F., 2007. Modeling a pilot scale pulp mixing chest using CFD, Journal of Pulp and Paper Science 33: 115-120 . †C. Gomez formerly C. Ford Chapter 2 39 surface motion was attained that considerable deviation from ideal mixing existed, including channeling and suspension stagnation below the pulp surface. Additional power was required to reduce or eliminate these regions of poor mixing. Dynamic studies made in industrial chests confirmed these findings (Ein-Mozaffari et al., 2004). Pulp suspensions display non-Newtonian rheology. They are fibrous networks that possess structure and strength that result from the interaction between neighbouring fibres. In order to create motion throughout the suspension, the shear stress imposed on the suspension must be greater than the network strength, which is measured as the yield stress. Pulp suspensions are often treated as plastic materials; they show little deformation up to the yield stress and flow once the yield stress is exceeded. Fibre suspensions also exhibit non-linear viscoplastic properties, such as shear thinning behaviour (Bennington, 1988; Wikstrom and Rasmuson, 1998), although, after the yield point, pulp rheology is not well characterized (Bennington, 1988). A significant improvement in stock chest design could be obtained by developing robust models that take into account the flow field inside the vessel. Computational fluid dynamics (CFD) makes the development of these models possible by solving the momentum and mass conservation equations that govern flow in the chest. In past work (Ford et al., 2006) commercial CFD software (Fluent v.6.1) was used to simulate pulp mixing in the laboratoryscale chest with the suspension rheology approximated using a Bingham plastic model. It was shown that the numerical model captured the mixing dynamics of the scale-model chest fairly well, although it over-estimated the extent of mixing in the chest. References start on page 56 Chapter 2 40 In laboratory work, Ein-Mozzafari (Ein-Mozaffari, 2002; Ein-Mozzafari et al., 2004) found that an abrupt transition from a well-mixed chest to one with significant channeling occurred for very small changes in chest operation, including small changes in suspension mass concentration, impeller rotational speed or chest production rate. In this paper, the CFD model (Ford, 2004) was extended to such conditions in order to examine the validity of the predictions calculated and to establish limits for the model. 2.2 Experimental Techniques The laboratory 1:11 scale-model chest used by Ein-Mozzafari (2002) was modeled by Ford et al. (2004, 2006). Pulp suspension was pumped from a feed tank and through the stock chest to a discharge tank. The stock chest (illustrated in Fig. 2.1) is rectangular, with dimensions: width: W = 40 cm, length: L = 53 cm and height: Z = 70 cm (the pulp suspension level was maintained at Z = 50 cm for all tests reported here). The impeller used was a scaled version of an industrial impeller (the Maxflo; Chemineer Inc., Dayton, OH) having a diameter of 16.5 cm. In the dynamic tests, the system was excited by the injection of saline solution into the pulp feed through a computer controlled solenoid valve. The conductivity of the input and output streams was measured using flow-through conductivity sensors and the data used to determine the extent of non-ideal flow in the system. Tests were made with a 3.3% softwood FBK (fully bleached kraft) pulp suspension over a range of impeller speeds (N = 1000 to 1730 rpm) and suspension flow rates through the chest (Q = 7.9 to 37.1 L/min). Detailed design information is available in earlier references (Ein-Mozaffari, 2002; Ein-Mozaffari et al., 2003; Ein-Mozzafari et al., 2004) and a schematic of the experimental set up is given in appendix A. References start on page 56 Chapter 2 41 The commercial CFD software package, Fluent v. 6.1 (Fluent Inc., Lebanon, NH) was used to simulate the suspension flow. The computational domain within the mixing chest and around the three-dimensional geometry of the impeller was discretized with mutiblock mesh structures using Gambit software (v. 2.1.6, Fluent, Lebanon, NH). Originally, the 3-D mesh of the model had 1.1 x 105 computational cells, consisting of tetrahedral and hexahedral elements. The mesh was refined to have 2.02 x 105 cells by means of an adaptive refinement of those cells with large velocity gradients based on the solution acquired in the coarser mesh, as outlined in Ford et al. (2006). A multi-reference frame solution approach (MRF) was utilized to obtain a steady-state solution of the governing equations in the computational domain, with coupling between reference frames made using a velocity transformation. Figure 1 shows the solution domain for the mixing chest (viewed along the y-z axis) and the locations of the boundary conditions specified in Fluent. Details of the model development are given in Ford et al (2006) and are summarized in appendix A of the this thesis Figure 2.1. Boundary conditions specified in the computational domain. The projection shown corresponds to the z-y plane using the coordinate system chosen for the model. References start on page 56 Chapter 2 42 The pulp suspension rheology was approximated using the modified Bingham model available in Fluent. Model parameters were based on measurements made on the pulp suspension used in the experimental tests. Below the yield stress, the suspension is modeled as an extremely viscous fluid with viscosity µo . Thus the apparent viscosity is given by µ = µo for τ ≤ τ y (2.1a) This is done to avoid computational instability, but it results in the CFD solution giving the suspension very-low velocities in locations where, in reality, it would be stationary. Above the yield stress, τ y , the suspension is described by µ= ⎡ ⎛τ τ y + k ⎢⎢γ n − ⎜⎜ y ⎝ µo ⎣⎢ γ ⎞ ⎟ ⎟ ⎠ n⎤ ⎥ ⎥ ⎦⎥ for τ > τ y (2.1b) with the consistency index, k, set to 0.001 Pa·s (the viscosity of water) and the power-law index, n, set to 1.0 (to approximate a Bingham plastic). The yield stress of the Cm = 3.3% FBK suspension was set to τy = 475 ± 62 Pa as measured by Ein-Mozzafari et al. (2005) and the yielding viscosity was estimated based on experimental data collected by Bennington (1988) as µo = 100 ± 43 Pa·s. The suspension flow observed in the experiments was primarily laminar (The non-Newtonian Reynolds number calculated for the mixing vessel was below 103). Consequently, the mass and momentum conservation equations were solved in the laminar regime using a second order upwind scheme. Iterations were continued until the scaled root mean square (RMS) residuals for each transport equation were below 1×10-5. Spatial grid independence was verified by demonstrating that additional grid lines near the impeller surface where velocity gradients and imbalance residuals were highest, did not change the calculated velocity magnitude measured at References start on page 56 Chapter 2 43 a line located 8 cm in front of the impeller by more than 3% (based on the integral of the velocity through the line defined by z = 15 cm and y = 12 cm across the chest width)(Ford et al., 2006). The dynamic tests conducted experimentally by Ein-Mozzafari (2002) were simulated numerically using virtual tracer addition (ρ = 1000 kg/m3) at the chest inlet and the transient solution of the species conservation equation. The steady-state velocity field computed for the test conditions was used in the simulation along with the same pseudo random binary (PRBS) input signal used in the experiment (details of this protocol are given in appendix A). The model output was directly compared with the experimental measurements. Additionally, the parameters of a black-box mixing model were estimated from both the experimental and numerical inputoutput data using the method described by Kammer et al. (2005) and compared. This dynamic model includes a well mixed zone (represented by a first order transfer function having a time constant τ2 and time delay T2) as well as the possibility of stock recirculation, R, back to it. A fraction of the entering flow, f, could by-pass the mixing zone and pass directly to the exit, although it could be involved in some limited mixing. Thus, a first order transfer function with time constant and delay (τ1 and T1) was included in the model for this portion of the flow. In the analysis that follows the extent of channeling is used as the mixing quality index. Mixing is improved (approaches ideal mixing) when f approaches 0. A block diagram of the model is given in Fig. 2.2. References start on page 56 Chapter 2 44 Figure 2.2. Block diagram of the dynamic model used to analyze the experimental and CFD generated data. 2.3 Results and Discussion The flow fields and impeller power consumption calculated by the CFD simulations agreed with visual observation of the flow (made of surface flows through the Plexiglas chest) and measured power use as outlined in previous work(Ford et al., 2006) . The velocity profile predicted by CFD for a typical set of conditions (a Cm = 3.3% FBK suspension with a flow rate through the chest of Q = 7.9 L/min and an impeller speed of N = 1132 rpm) is shown in Fig. 2.3. Under these conditions, suspension motion occurs throughout the chest although analysis shows it to be subject to significant channeling (f = 0.40). The flow field computed here (and observed visually) shows that regions of poor mixing (those with very low flow velocities) exist at the bottom of the chest towards the back wall (opposite the impeller). Here, pulp moves significantly slower than in the regions of bulk motion. In the upper region of the vessel, while flow is more vigorous, the location of the suspension outlet allows pulp to bypass the mixing zone contributing to the high degree of channeling measured. References start on page 56 Chapter 2 45 Figure 2.3. Velocity vectors (arrows) in three planes parallel to the impeller flow: in the plane x = 0 (the impeller mid-plane) and the planes x = ±9 cm (close to the impeller tips). Cm = 3.3% FBK suspension at Q = 7.9 L/min and N = 1132 rpm. Vectors are colored by velocity magnitude (m/s). Some detail in the impeller zone has been omitted for clarity. Figure 2.4 plots the suspension conductivity for the input and output signals as measured experimentally (the solid lines) and as calculated from a CFD simulation of the dynamic test (dashed lines) for the conditions noted above (Fig. 2.3). The excitation signal used in the CFD simulations did not include noise present in the experimental measurements and thus the simulated numerical output lacks a high-frequency noise component. Under these operating conditions the computed response using the CFD data agrees well with the measured data. The model parameters for the dynamic model (Fig.2.2) were calculated using data generated by the CFD simulations and compared with those obtained using experimental data. This data, given in Table 1 for a number of tests, shows that parameters calculated using the CFD flowfield are generally in the same range as the experimentally determined ones. An exception to References start on page 56 Chapter 2 46 this are the time constants determined for the mixing zones, particularly that in the bypass flow, τ1, which is found to be typically five times larger than the experimental value. When this occurs, the extent of channeling determined by CFD is increased over that found experimentally (as significantly greater mixing is estimated to occur in the bypass flow). 10 input conductivity (mS/cm) output Expt CFD 9 8 7 6 5 4 0 1000 2000 3000 4000 5000 t (s) Figure 2.4. Input (upper gray lines)/Output (lower black lines) signals for numerical (dashed lines) and experimental (solid lines) measurements for the scale-model chest. Cm = 3.3% FBK, Q = 7.9 L/min, N = 1132 rpm. Figure 2.5 plots the fraction of pulp bypassing the mixing zone (f) versus impeller speed for a range of suspension flow rates through the mixing chest as reported by Ein-Mozzafari et al. (2004). At a fixed flow rate, the extent of channeling increases as the impeller speed decreases; at a fixed impeller speed, the extent of channeling increases as the suspension flow rate increases. At low flow rates (Q ≤ 21.1 L/min), the increase in channeling with decreasing impeller speed occurs fairly smoothly as shown in Fig. 2.5. This figure also shows the sudden onset of significant channeling that can occur under certain operating conditions, typically at higher suspension flow rates and lower impeller speeds. In particular, for flow rates above References start on page 56 Chapter 2 47 28.6 L/min, a step change in channeling occurs for a very small (approximately 50 rpm) impeller speed change near a critical value. Table 2.1. Comparison of model parameters estimated from experimental and CFD data. Cm = 3.3% FBK at Q = 7.9 L/min, 21.1 L/min and 37.1 L/min for selected impeller speeds (N) Dynamic Parameters Q = 7.9 L/min f T1 (sec) τ 1 (sec) T2 (sec) τ 2 (sec) Q = 21.1 L/min f T1 (sec) τ 1 (sec) T2 (sec) τ 2 (sec) Q = 37.1 L/min f T1 (sec) τ 1 (sec) T2 (sec) τ 2 (sec) Experiment CFD N = 1258 rpm 0.35 0.50 45 60 13.3 75.5 140 295 966 1413 N = 1031 rpm 0.57 0.59 22 20 5.4 22.5 74 96 400 493 N = 1456 rpm 0.35 10 3.0 25 167 0.60 12 14.5 26 366 References start on page 56 Chapter 2 48 Figure 2.5. The effect of impeller speed on channeling (f) as a function of the pulp flow rate (after Ein-Mozzafri et al.,2004). Numerical simulations were performed using CFD at three different flow rates: Q = 7.9, 21.1 and 37.1 L/min. Figure 2.6 gives the predicted and measured values for the extent of channeling (f) as a function of N. Comparison of these data show that CFD models the flow correctly in two aspects: at a fixed impeller speed the extent of channeling increases as the flow through the chest rate is increased; and at a fixed flow rate the extent of channeling increases as the impeller speed is reduced. In the latter case, the trend predicted by CFD is fairly good at low flow rates (Q ≤ 21.1 L/min), particularly for lower impeller speeds (N ≤ 1300 rpm). However, at the highest flow rate (Q = 37.1 L/min), the CFD simulation was unable to predict the dramatic change in channeling that occurred around N = 1325 rpm. Here, values of f predicted by CFD were considerably below experimentally determined ones for N < 1300 rpm and above them for N > 1350 rpm. Indeed, the predicted extent of channeling did not change significantly with changes in impeller speed for any of the conditions examined. Further, the CFD model was References start on page 56 Chapter 2 49 unable to predict the elimination of channeling (f → 0) which occurred as N was increased above 1700 rpm. The inability of CFD to handle these two flow situations is examined below. Figure 2.6. The effect of the impeller speed on channeling at different flow rates. Comparison of experimental measurements and CFD predictions. The dynamic response calculated using the computed flow fields followed the experimental response in many circumstances. For example, for Q = 37.1 L/min and N = 1456 rpm, the response to the PRBS signal obtained shows good agreement with the experimental data (see Fig. 2.7(a)) despite the value of f being so different (see Table 2.1). However, the CFD simulation was not able to predict the sudden shift to channeling that occurred as the impeller speed was reduced. At an impeller speed of N = 1203 rpm (Fig. 2.7(b)), the modeled dynamic response (the dashed black line) follows the measured exit conductivity profile only for the first PRBS perturbation. Note that for the first step change, both the experimental and CFD data show some channeling in addition to a significant degree of mixing. However, for the second and all subsequent step changes, the experimental output signal changed – to one showing almost References start on page 56 Chapter 2 50 exclusive channeling (This is readily seen by comparing the experimental input and output signals in Fig. 2.7(b). When channeling is significant, the output signal closely follows the input signal with only a small time delay). This represents a dynamic change that cannot be captured by any steady-state simulation. A plausible explanation is that two flow situations exist for this set of experimental conditions, and that some disturbance triggered a transition between them. CFD was unable to predict the sharp transition, and only calculated one flow field. As impeller speed was reduced (in the CFD simulations) the extent of channeling (f) increased, but only gradually. No discontinuity was ever found. It was noted that in addition to the sudden change in mixing behaviour observed in the experimental tests, that similar and sudden onsets of channeling have been observed in industrial situations. This led us to try to improve the CFD simulation to see if its occurrence could be modeled. It was initially thought that the simulation was being affected by the boundary conditions used for the pulp inlet to the chest. Visual observation of experimental tests showed that the inlet flow could remain on the surface of the chest and be conveyed to the exit without entering the impeller mixing zone. By changing the momentum with which the simulated (CFD) flow entered the chest (by reducing its velocity) it was hoped to model the observed phenomena (by preventing penetration of the added flow to the mixed cavern). The inlet velocity was reduced in steps by increasing the diameter of the computed feed-pipe to the chest. These increased the time delays estimated for both the mixing and bypass zones but did not stop the flow from entering the mixing zone. The degree of channeling decreased (opposite to the effect predicted) and the simulated flow field never induced a sudden transition to extensive channeling (despite extending the range of impeller speeds investigated to even lower values than experimentally tested). Changes to the rheological parameters used to characterize the References start on page 56 Chapter 2 51 suspension (yield stress and yielding viscosity) or to the rheological model (modeling the pulp suspension as a shear-thinning suspension using the modified Hershel-Buckley model in Fluent) were likewise unable to model the sudden transition to channeling observed on reducing the impeller speed (Ford, 2004). It is likely that this sudden transition may be related to a phenomenon that cannot be described using the Bingham rheological model employed in our CFD simulations. Observation shows that part of the suspension added to the chest is conveyed directly to the chest exit along the suspension surface. A slip layer may exist between this surface flow and the flow generated by the impeller (which forms a cavern). The equations used to determine the dimensions of the cavern formed around the impeller (e.g. Amanullah et al, 1998) equate the force applied by the moving fluid at the interface with the yield stress of the suspension. If the yield stress is exceeded, the interface moves until a balance exists, i.e. the cavern grows. However, in a pulp fibre suspension, the force applied at the interface disrupts the fibre network and rolls the fibres into flocs. Flocs move more readily and are not as efficient as transferring force to the interface. The interface can still move outwards if the net force applied by the flowing suspension adjacent to it exceeds the yield stress of the suspension. Thus cavern size is still determined by the yield stress of the suspension – but not exactly as described by the Bingham model, in part due to slip at the interface. References start on page 56 Chapter 2 52 3.5 Input 3.0 Conductivity (mS/cm) Output Expt. CFD 2.5 2.0 1.5 N = 1456 rpm 0 (a) 200 400 600 800 1000 t (s) 3.5 Input 3.0 Conductivity (mS/cm) Output Expt. CFD 2.5 2.0 1.5 (b) N = 1203 rpm 0 200 400 600 800 1000 t (s) Figure 2.7. Input (upper gray lines)/Output (lower black lines) signals for numerical (dashed lines) and experimental (solid lines) measurements for the scale-model chest. Cm = 3.3% FBK at Q = 37.1 L/min: (a) N = 1456 rpm; (b) N = 1203 rpm. As the position of the interface changes (and it may not do so in a continuous manner), competition between flow exiting the chest and that entering the impeller suction changes, affecting the amount of surface flow entrained by the impeller and consequently the degree of References start on page 56 Chapter 2 53 channeling (Ein-Mozaffari et al., 2005). Changes in impeller speed and/or suspension flow alter the ratio between flow exiting the chest and entering the impeller zone. Over certain operating ranges this ratio may change significantly which would then dramatically change the degree of channeling in the system. This competition explains why CFD simulations are unable to predict the sudden increase in channeling observed experimentally under certain conditions. It may also explain why the predicted degree of channeling never reaches the experimentally found value of zero for high impeller speeds. The CFD simulations compute changes in velocity found for a continuous Bingham plastic fluid using parameters measured for the pulp suspension. The changes in velocity predicted are therefore more gradual than likely exist in the suspension and in the real pulp chest. The result is that some surface flow is always predicted to be directed towards the exit, whereas in reality it can be drawn into the impeller flow. The work presented in the following chapters of this thesis is attempting to address the issues identified here. Experimental investigations using transparent model fluids will permit the flow field to be measured using PIV (particle image velocimetry). This experimental data can then be compared directly with CFD simulations using a fluid having measurable rheology over the complete flow range investigated. This will help further validate the CFD model used in our work, as well as deal with boundary conditions including the free surface in the mixing chest. Other experiments made with pulp suspensions, are comparing the prediction of cavern volume with measurements made using a non-invasive imaging technology (ERT – electrical resistance tomography) coupled with Ultrasonic Doppler velocimetry (UDV) measurements of the suspension’s velocity across the cavern boundary (Hui et al., 2009) References start on page 56 Chapter 2 2.4 54 Summary and Conclusions Experimental results have shown that the scale-model mixing chest is susceptible to high percentage of non-ideal flow and that changes in chest performance can occur abruptly for modest changes in operating conditions. This same behaviour has been reported for industrial chests. The CFD simulation developed using an approximated Bingham rheology was not able to predict this sudden transition to channeling, although model agreement was good under many other flow conditions. This work has identified the limitations of the numerical modeling of mixing pulp suspensions when treated as a continuum fluid and has determined the direction of present work in this area. 2.5 Nomenclature D impeller diameter, m Cm mass concentration of pulp suspension (consistency), % f percentage of channeling, fraction k consistency factor (kg/ms) N impeller rotational speed, rpm Q pulp flow rate through the chest, m3/s R percentage of recirculation, fraction T1 time delay for the channeling zone, s T2 time delay for the agitated zone, s Greek letters γ shear rate, s-1 References start on page 56 Chapter 2 µ fluid viscosity, kg/m.s µ0 yielding viscosity, kg/m.s ρ fluid density, kg/m3 τ1 time constant for the channeling zone, s τ2 time constant for the mixing zone, s τy yield stress, N/m2 55 References start on page 56 Chapter 2 2.6 56 References Amanullah, A., Hjorth, S.A. and Nienow, A.W., 1998. A new mathematical model to predict cavern diameters in highly shear thinning, power law liquids using axial flow impellers. Chem. Eng. Sci., 53, 3, 455−469. Bennington, C.P.J., 1988. Mixing pulp suspensions, PhD dissertation. University of British Columbia, Vancouver, Canada. Bialkowski, W.L., 1992. Newsprint variability and its impact on competitive position. Pulp Pap. Can., 93, 11, T299−T306. Ein-Mozaffari, F., 2002. Macroscale mixing and dynamic behaviour of agitated pulp stock chests, PhD dissertation. University of British Columbia, Vancouver, Canada. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2003. Performance and design of agitated pulp stock chests. Appita J., 56, 2, 127−133. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2004. Dynamic mixing in agitated industrial pulp chests. Pulp Pap. Can., 105, 5, T 118− T 122. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2005. Suspension yield stress and the dynamic response of agitated pulp chests. Chem. Eng. Sci., 60, 2399−2408. Ein-Mozzafari, F., Kammer, L.C., Dumont, G.A. and Bennington, C.P.J., 2004. The effect of operating conditions and design parameters on the dynamic behaviour of pulp stock chests. Can. J. Chem. Eng., 82, 154−161. Ford, C., 2004. CFD simulation of mixing dynamics in agitated pulp stock chests, M.A.Sc. Thesis. University of British Columbia, Vancouver, Canada. References start on page 56 Chapter 2 57 Ford, C., Ein-Mozaffari, F., Bennington, C.P.J. and Taghipour, F., 2006. Simulation of mixing dynamics in agitated pulp stock chests using CFD. AIChE J., 52, 10, 3562−3569. Hui, L.K, Bennington, C.P.J., and Dumont, G.A., 2009. Cavern formation in pulp suspensions using side-entering axial flow impellers. Chem. Eng. Sci., 64, 509−519. Kammer, L.C., Ein-Mozzafari, F., Dumont, G.A. and Bennington, C.P.J., 2005. Identification of channeling and recirculation parameters of agitated pulp stock chests. J. Process Control, 15, 31−38. Oldshue, J.Y. and Gretton, A.T., 1956. Performance and design of paper stock mixers. Tappi J., 39, 6, 378−390. Wikstrom, T. and Rasmuson, A., 1998. The agitation of pulp suspension with a jet nozzle agitator. Nord. Pulp Paper Res., 13, 2, 88−94. References start on page 56 58 Chapter 3 Computational Modeling of Industrial Pulp Stock Chests 3.1 Introduction Effective operation of agitated pulp stock chests is crucial to the success of many processes in the pulp and paper industry, including those requiring blending of pulp suspensions, chemical contacting and consistency (fibre mass concentration) control. One key feature of these chests is their ability to attenuate high frequency disturbances in pulp properties, acting as low-pass filters to complement the action of process control loops. However, the ability to attenuate these disturbances depends on the effectiveness of mixing achieved, and it is important to ensure that mixing chests are properly designed and operated to achieve the required mixing quality. There are many challenges in meeting this objective, and this paper presents results from computational fluid dynamic (CFD) simulations carried out on two industrial-scale pulp stock chests as a part of these efforts. Pulp suspensions have unique rheology: they are non-Newtonian, multi-phase suspensions and possess a yield stress. This rheology contributes to the creation of undesirable flows in the mixers, including channeling (bypassing) and dead zones (stagnant regions), which degrade mixing quality. Stock chests are typically designed using semi-empirical rules based largely on past experience. This is amply illustrated in a review of chest design practices summarized by Yackel (1990). The design procedure recommended is empirical in nature and is based on A version of this chapter has been submitted for publication. Battacharya S., Gomez† C., Sotanzadeh A., Taghipour F., Bennington C.P.J. and Dumont G.A., 2008. Computational modeling of industrial pulp stock chests. † C. Gomez formerly C.Ford Chapter 3 59 matching the momentum generated by the impeller(s) (a function of impeller design, diameter and rotational speed) with that required to achieve complete motion across the suspension surface of the chest. The criterion of ‘complete surface motion’ reflects the fact that experimental data were often obtained from operating plants where only the observation of the suspension surface was possible. The required momentum is given as a function of the chest geometry, suspension mass concentration (consistency, Cm) and fiber type − factors accounted for using graphical correlations and correction factors. Although Yackel’s method provides a systematic approach to chest design, Ein-Mozaffari (2003) showed that it does not ensure that ideally mixed conditions are achieved. f u(s) e G 1′ ( s ) = 1 1 + τ 1s y(s) −Td s 1 Mixing Zone 1-f G ′2 ( s ) = 1 1 + τ 2s Figure 3.1. Grey box model for flow through pulp stock chests (Soltanzadeh et al., 2009). The dynamic response of a typical industrial chest to an inlet disturbance is often significantly different than that expected from an ideally mixed vessel (Bialkowski, 1992; Ein-Mozaffari et al., 2003). Extensive tests were made on a geometrically-scaled, rectangular, laboratory chest with a side entering impeller (Ein-Mozaffari, 2002) to study the interrelationships between chest geometry, pulp properties, process operating parameters and the resulting mixing efficiency. A dynamic mixing model based on a physical interpretation of flow in the stock chest was References start on page 88 Chapter 3 60 developed (Soltanzadeh et al., 2009). In this model (Fig. 3.1) the chest dynamic has been characterized with four parameters f, τ1, τ2, T1. f is the fraction of incoming pulp which bypasses the main mixing zone and proceeds directly to the chest exit. While most mixing occurs within this zone surrounding the impeller (represented by a first order transfer function with a time constant, τ2) a limited extent of mixing also occurs in the bypass stream (time constant τ1, where τ1<<τ2). A single time delay, T1, was used in the model to represent the contribution of plug flow through the chest and in the pipe inlet and outlet pipes. This gives the following transfer function (1 + τ Z s ) e−T s GChest = (1 + τ1s ) (1 + τ 2 s ) 1 (3.1) where τZ = f τ2 + (1 – f)τ1. Simplification of Eqn. (3.1) is possible if there is negligible mixing in the bypass stream (i.e. τ1 ~ 0). This gives 1 + τ Z′ s ) e −T s ( GChest = (1 + τ 2 s ) 1 (3.1a) with τ Z′ = fτ 2 . Finally, when there is no bypassing Eqn. (3.1a) reduces to a simple first order transfer function with delay GChest = e −T1s (1 + τ 2 s ) (3.1b) Experimental data was obtained by measuring the temporal variation of a suitable suspension property at both the inlet and outlet of the chest (in the laboratory − suspension conductivity following tracer addition; in the mill − suspension mass concentration, Cm, following a change introduced at the inlet). Model parameters were then identified using established methods (Soltanzadeh et al., 2009). This model was found to be more stable than a similar model used in earlier work (Ein-Mozaffari et al., 2005) with very similar values for the key parameters of f and References start on page 88 Chapter 3 61 τ2 identified. Once the model parameters were identified the fraction of suspension in the chest (VT) that is fully mixed (VFM) could be calculated VFM τ 2 (1 − f ) Q = VT VT (3.2) Ein-Mozaffari et al.(2003, 2004) found that f and VFM/VT were strongly influenced by impeller speed, pulp consistency, suspension flow rate (production rate), fibre type and the location of the stock exit and entrance with respect to the impeller. An increase in consistency, decrease in impeller speed or increase in stock flow increased the extent of channeling in the chest and reduced its ability to attenuate incoming disturbances. These trends were explained by the rheology of the pulp suspensions. When a yield stress fluid is agitated, a ‘cavern’ (an actively mixed region of fluid motion) forms around the impeller with its size determined by a balance between the stress applied by the impeller and the suspension yield stress (Wichterle and Wein, 1981). Factors that determine the applied stress include the impeller design, diameter and operating speed. For example, an increase in suspension mass concentration increases suspension yield stress, which for a fixed set of impeller operating conditions results in a decrease in cavern volume and a reduction in mixing quality (Ein-Mozaffari et al., 2005). The close correlation between the flow characteristics obtained from the analysis of laboratory tests (extent of channeling, extent of chest involved in mixing) and the physical interpretation of suspension flow based on pulp rheology provide a good basis for the use of the grey box model described previously for analysis of mixing performance. The experimental data also showed that poorly mixed regions existed below the suspension surface even it was in complete motion. Indeed, Ein-Mozaffari (2002) found that up to three References start on page 88 Chapter 3 62 times more power was required to produce complete motion in the vessel than specified by Yackel for ‘complete surface motion’. Creation of ideal mixing required even greater power. These findings indicate that stock chests designed using traditional methods are unlikely to achieve their design objectives. Performance analysis of five industrial stock chests confirmed this (Ein-Mozaffari et al., 2004). The fully mixed volume in these mixing chests was only 27% to 87% of the total volume resulting in effective time constants significantly lower than design values. Although grey box modeling provides an attractive means of identifying and quantifying mixing performance, it does not identify why sub-optimal mixing occurs. To obtain information on flow in these chests and to find relationships between suspension flow and mixing dynamics we are examining the ability of CFD to supply this information. A number of studies have previously used CFD to calculate steady-state flow fields in pulp mixers. Bakker and Fasano (1993) modeled agitation in a pulp stock chest with a side-entering impeller by defining velocity profiles at the impeller and using a combination of laminar and turbulent flow regimes in different regions of the tank. The flow field calculated matched that expected, but was not compared with experimental data to confirm the accuracy of the predictions. Roberg (1997) simulated mixing of a low consistency pulp suspension in a cylindrical stirred vessel agitated with a side-entering impeller. The author assumed laminar flow and compared calculated flow fields for Newtonian and non-Newtonian rheological models to highlight the similarities and differences between them. The impellers were modeled as momentum sources (obtained using the distribution of the measured axial and tangential forces produced by the impeller over the swept volume). Wikstom and Rasmuson (1998) modeled mixing of low-consistency (Cm = 2.9 to 3.2%) pulp suspensions by a jet nozzle agitator in a 24 m3 pilot-scale chest using the CFD code CFX. The outflow from References start on page 88 Chapter 3 63 the nozzle was treated as a velocity boundary and the suspension as a Bingham fluid. The computational results were compared with ultrasonic Doppler velocimetry (UDV) measurements made at four locations in the chest and found to agree in regions near the nozzle when a yield stress of 150 to 188 Pa was used in the model. Cavern sizes calculated using the CFD simulations (where the suspension velocity fell to 1% of the initial nozzle velocity) were much higher than those predicted using correlations available in the literature. Ford et al. (2006, 2007) used CFD software (Fluent 6.1), modified Herschel Bukley rheology, laminar flow and a multiple reference frames (MRF) approach to compute the steady-state flow fields in the laboratory-chest used by Ein-Mozaffari (2002). The flow fields were then used to compute the dynamic response to a virtual tracer input which permitted comparison with the model parameters determined experimentally. In general, the CFD simulations resulted in good prediction of chest dynamics indicating that the simulations reproduced the essence of the flow fields well. However, the CFD simulations were unable to predict sudden transitions in mixing quality that were measured under certain process conditions (Ford et al., 2007; Ford et al., 2006). The objectives of the present work are to extend our CFD modeling to more industrial relevant mixer configurations, assess the accuracy of the simulations using measured data and use details of the simulated flow-fields to draw insight into the effects of design parameters on the chest performance. References start on page 88 Chapter 3 3.2 3.2.1 64 Experimental Chest Geometry and Operating Conditions The two industrial chests selected for analysis are illustrated in Figs 3.2 and 3.3. The first chest has a capacity of 306 m3 and feeds a mixture of chemical pulps (softwood, hardwood and broke) at consistency, Cm = ~3.5% to a paper machine. This rectangular chest (length, L = 7.44 m; width, W = 7.16 m and height, H = 5.74 m with a nominal stock level of 0.5H) is agitated by a single side-entry impeller located on the wall opposite stock addition. The suspension is added to the top of the chest through the main stock inlet (Inlet 1) and falls to the suspension surface. Mixed suspension is withdrawn from the chest outlet located on the wall behind the impeller. The exit pipe extends 0.5 m into the chest and is bent downwards to reduce short-circuiting. However, as both the suspension inlet and outlet are located on the same side of the chest, some channeling is expected. A small amount of the discharge stream (4.4%) is directed to an in-line freeness tester and added back to the chest through Inlet 2. The suspension is agitated using a D = 0.91 m, three bladed, marine-type, side-entering impeller rotating at N = 445 rpm. The impeller is centered on the wall 0.66 m above the chest floor with an offset of 0.45 m from the back wall. Other than its volume, this chest is of similar geometry to the scale-model chest used in our laboratory (Ein-Mozaffari, 2002). Dynamic data were obtained by measuring the inlet and outlet mass concentration to a step change in inlet consistency. Pulp consistency was measured using a shear-based consistency sensor (Valmet SP, Metso Automation, Edmonton, AB) with a measurable range of 2.5 to 5% Cm and a repeatability of 0.005% Cm at 3% Cm at the inlet. Outlet consistency was measured by a microwave-based sensor (Kajaani MCA, Metso Automation) with repeatability of 0.01% consistency. References start on page 88 Chapter 3 65 Inlet 1 Symmetry boundary at liquid surface Inlet 2 Rotating zone No slip boundaries at walls and floor (a) Pump Suction (b) Figure 3.2. Machine Chest (MC-A): (a) Multiblock structure used for grid generation; (b) outer surface of mesh used in the initial grid. The second chest has a capacity of 340 m3 and is used to remove latency1 from a Cm = 3.5% thermo-mechanical pulp (TMP) suspension ahead of stock screening. This U-shaped chest (overall dimensions L × W × H = 11.3m × 6.1m × 4.9m with a stock level of 0.75H) uses four side-entry Lightnin A-312 impellers (NP = 0.36; C = 0.4, Npump = 0.56) running at 310 rpm 1 Latency describes a property of thermomechanical pulps where part of the mechanical strength of paper sheets produced from the pulp is lost due to curl in the fibers. Recovery of this latent strength requires straightening of the pulp fibres, which is accomplished by agitating a dilute (Cm ~ 3%) suspension at a temperature above the lignin softening point for a specified length of time (Beath et al., 1966) References start on page 88 Chapter 3 66 located at the chest corners (see Fig. 3.3). The mid-feather wall directs stock flow through the chest and minimizes channelling. Both inlet and outlet suspension consistencies are measured using BTG MEK-2000 series consistency sensors (BTG Americas Inc., Montreal, QC). A summary of geometry and operating conditions for both chests is provided in Table 3.1. Table 3.1. Industrial chest specifications Chest shape Dimensions: L×W×H (m) Normal operating level (% of H) Pulp type(s) Pulp mass concentration (%) Operating temperature (oC) Pulp flow rate (m3/min) Inlet pipe location and diameter (m) Outlet pipe diameter and location (m) Number of impellers and installed power (kW) Impeller type and diameter (m) Impeller rotational speed (rpm) Consistency sensors 3.2.2 Machine Chest – A Rectangular 7.44 × 7.16 × 5.74 Latency Chest – B U-shaped 11.3 × 6.1 × 4.9 50 75 Mixture of hardwood/softwood/broke 3.5 40 6.9 (+ 0.32 for the recirculated stock) Vertically from chest roof, 0.3 (both) Stock added to top of suspension at an angle, 0.2 0.3 (0.18 m from impeller) 0.35 (0.35 m from impeller 4) 1 @ 125 hp 4 @ 30 hp Black Clawson marine-type impeller side-entry configuration, 0.91 445 Inlet: Valmet SP Exit: Kajaani MCA Lightnin A-31 (highefficiency axial impeller) side-entry configuration, 0.76 310 BTG MEK 2000 series on both Thermomechanical 3.5 63 8–9 Industrial Dynamic Tests Measurements of suspension mass concentration were made at the chest inlet and outlets during a change in incoming stock consistency to measure the mixing dynamics. The requirements for obtaining good data for the analysis include: References start on page 88 Chapter 3 • 67 Process control loops must be disengaged so that the open-loop response to the disturbance is measured. (Often multiple control loops (stock level, consistency etc.) are nested with other operational units and must be placed on manual control). • Changes (often a step change) in input consistency must be sufficiently large to give measurable changes with an acceptable signal to noise ratio. (A step change in inlet consistency of 10 to 15% of its initial value was found to be suitable). • Sensors should be accurately calibrated (however, calibration errors must often be dealt with following data acquisition). • Data collection should be three to four times the nominal residence time of the chest to properly estimate the model parameters. • Data must be collected without compression. (Data management software often compresses data for long-term storage which results in a loss of frequency content). • The sampling rate for the data must exceed the Nyquist frequency of the highest desired measurement frequency (sampling times of 1 to 5 s are typical for industrial chests). We attempted to meet these requirements during our tests, although we were limited by operational constraints. 3.3 Modeling Approach The computational domain for each chest was developed using Gambit software (v.2.2.30, Fluent Inc., Lebanon, NH). Multi-block structures were created for each vessel to provide efficient allocation of computational cells, with a higher cell density in regions of high velocity gradients (e.g. the impeller zone) and near wall boundaries. Adjoining blocks allowed for a gradual transition in cell size (less than a four-fold change of the characteristic length of a cell between References start on page 88 Chapter 3 68 adjacent zones) with only tetrahedral and/or pyramid shaped cells used in regions containing structurally complicated components (i.e. the impellers), and tetrahedral or hexahedral cells used in the bulk of the vessel. The initial multi-block grid structures were further refined until the solutions were determined to be grid-independent. The chest floor and side walls were defined as no-slip boundaries while the liquid surface was treated as a symmetry boundary with no flux across it. The angular speed of the impeller and axis of rotation were specified for the zones surrounding the impellers with the impeller blades defined as walls with zero relative velocity. Uniform velocity profiles were specified at >10din upstream of the chest to permit full flow development for the inlets. Similarly, the outlet boundaries were located a distance of >10dout from the chests. Suspension rheology was approximated using a modified Herschel-Bulkley model. At very low strain rates ( γ ≤ τ y / µ 0 ) the suspension was treated as an extremely viscous fluid with viscosity µo . This avoids instability during computation and results in some regions in the chest having very low velocities where in reality they are stationary. As the strain rate was increased and the yield stress,τ y , was exceeded, the suspension rheology was described using a power law relationship. Thus µ = µo ⎡ ⎛τ τ y + k ⎢γ n − ⎜⎜ y ⎢⎣ ⎝ µo µ= γ ⎞ ⎟⎟ ⎠ n ⎤ ⎥ ⎥⎦ for τ ≤ τ y (3.3) for τ > τ y (3.4) References start on page 88 Chapter 3 69 The consistency index, k, was set to 0.001 Pa·s (the viscosity of water) and the power-law index, n, to 1.0 following Ford et al. (2006, 2007). The yield stress of the suspensions was estimated using correlations given by Bennington et al.(1990) and the yielding viscosity estimated using experimental data collected by Bennington (Bennington, 1988) as µo = 100±43 Pa.s. Impeller 3 Inlet Impeller 2 Impeller 4 Pump Suction (a) Impeller 1 (b) Figure 3.3. Latency Chest (LC-B): (a) Multiblock structure used for grid generation; (b) Outer surface of mesh used in the initial grid. The shaded regions in Fig. 3(a) are solved in a rotating reference frame. Commercial CFD software, Fluent 6.3, was used to simulate the flow field in the two chests using the steady-state, multiple reference frames (MRF) approach to model the impeller (rotating) regions. Continuity and momentum conservation equations were solved in the laminar regime with a second-order upwind scheme. Iterations were continued until the RMS scaled residuals for each transport equation fell below 10-5 and a steady velocity value was reached at References start on page 88 Chapter 3 70 various sampled locations in the domain. The converged steady-state flow field was then used to determine the dynamic response of the chest to a virtual change in tracer concentration (having the same properties as the pulp suspension) made at the chest inlet (and equivalent to that made in the mill test) using Fluent’s user defined function capability. The transient species transport equations for the virtual tracer were solved with a specified time step, typically around ∆t = 1s (following tests to determine independency from time discretization). Simulations were carried out for a period at least as long as that of the experimental tests. The temporal variations in tracer concentration were then analyzed to determine the dynamic characteristics of the virtual stock chest using Eqn. (3.1b). The initial computational domain for MC-A contained 12 zones and was refined using the evolving flow-field to obtain a final gradient adapted grid with 0.3 million cells (Ford et. al, 2006) (Fig. 3.2). Velocity profiles perpendicular to and 0.08m (a) and 3.0 m (b) in front of the impeller are given in Figs 3.4(a) & (b). Grid independence was verified by demonstrating that further refinement of the cells did not change the velocity magnitude at these locations by more than 3.5%. The multi-block grid developed for LC-B is shown in Figure 3.3 and was similarly refined to obtain grid independence. In the industrial chest, the inlet pipe terminates at the chest wall with incoming stock falling as a jet to the suspensions surface. To simplify the CFD computations, the location, angle and velocity of jet impact with the suspension surface was calculated and modeled as a virtual feed pipe joining the pulp bed at the calculated location. References start on page 88 Chapter 3 71 Figure 3.4. Test for grid independence for MC-A: (a) magnitude of suspension velocity 8 cm in front of the impeller (along the y-axis); (b) magnitude of the suspension velocity 3m away from the impeller. (c) Effect of yield stress and yielding viscosity on computed velocity along the impeller axis. (a) (b) References start on page 88 Chapter 3 72 (c) 3.4 3.4.1 Results and Discussions Machine Chest (MC-A) Computational simulation of MC-A began with selection of appropriate values for τy and µ0 in the rheological model (Eqns. (3.3) and (3.4)). This was difficult as the suspension feed was a mixture of pulps (and we were given only nominal composition values from the mill: 50 to 100% hardwood kraft, up to 50% softwood kraft and up to 40% broke). For the simulations reported here, the yield stress was taken as τy = 113 Pa based on the aggregate data reported in Bennington et al. (1990). Most computations were made using this yield stress value, although some simulations were initially made to determine the effect of suspension rheology on the flowfield computed in the mixing chest, as shown in Fig. 3.4c. Here suspension velocity in the xdirection (along the impeller axis) is shown for fixed set of operating conditions (N = 445 rpm; References start on page 88 Chapter 3 73 Q = 7.22 m3/min) for τy = 50, 113 and 580 Pa with µ0 = 100 Pa.s and for τy = 580 Pa with µ0 = 300 Pa.s. The velocity profiles are similar at the lower yield stresses. Velocity magnitudes are lower at higher yield stresses, most noticeable further from the impeller (due to the cavern formation in the chest). If an accurate method to measure suspension velocity existed, differences between calculated and measured velocity profiles could be used to judge the quality of the simulation. The steady-state flow field was computed for process conditions corresponding to those during collection of experimental data. Once a converged solution was obtained, a dynamic simulation to model the experimental bump test was made. This simulation was made using a time interval ∆t = 1s, which was found to give results independent of time discretization (time intervals from 0.25 to 5s were examined). The experimental data collected from the test made on MC-A is shown in Figure 3.5. Here the suspension inlet consistency was reduced in a step fashion from Cm = 3.67% to 3.60%, only a 2% decrease. Data was sampled every 5s from both the inlet and outlet consistency meters. The quality of the inlet data is good, with little evidence of data compression. However, the exit stream data is compressed, resulting in a loss of frequency content which could introduce artificial low frequency dynamics into the signal and making identification of the model parameters more difficult. A problem with the calibration of the consistency meters is also apparent. While the inlet consistency changed by only 0.075% the measured discharge consistency changed by more than 1%. Also, the consistencies are not equal at the beginning of the test despite the fact the chest is operating at steady-state. Note that poor sensor calibration is References start on page 88 Chapter 3 74 not as severe a problem as signal compression because the dynamics are important for parameter identification and the sensor gain and offset can be adjusted as needed. 0.0385 Exit Consistency Consistency at Inlet ( - ) 0.038 Inlet Consistency 0.0375 0.037 0.0365 0.036 0.0355 0.035 0 1000 2000 3000 4000 Time (s) 5000 6000 7000 8000 Figure 3.5. Raw experimental bump test data for MC-A. Note problems with consistency meter calibration and compression of the data (particularly the exit data). Although the conditions under which the mill data were obtained were less than ideal, they demonstrate typical problems faced in the field. In this case, as the pulp leaving this chest goes directly to the papermachine, a more dramatic change in inlet consistency is not possible without creating operating issues there. Further experimental data were not available from this chest due to the closure of the mill. Consequently, our analysis is based on the best data available. Direct and qualitative comparison of the simulated and experimental data provides an indication of the ability of CFD to properly simulate the mixing chest. The real-time data were analyzed using the system identification toolbox in MATLAB version 7 (R14) and the results given in Table 3.2. Analysis shows that only 2.4% of the incoming pulp References start on page 88 Chapter 3 75 bypasses the chest, with negligible mixing in the bypass stream (τ1 = 0). The nominal time constant of the chest (τnom) was 1271s (at Q = 7.22 m3/min), whereas τ2 was identified as 780s. Thus only 61% of the chest was involved in active mixing of the suspension. As little pulp bypassed the mixing zone, the balance of pulp in the chest is either stagnant or very slowly moving. Table 3.2. Comparison of experimental and CFD simulated parameters for the stock chests MC-A Operating Conditions Q (m3/min) H (%) Cm (%) τy (Pa) Model Parameters T1 (s) τ1 (s) τ2 (s) f Mixing Parameters VFM/VT LC-B Expt. CFD1 Expt.† CFD 7.2 NM 3.7-3.6 NM 7.2 50 3.7 113 8.2-9.6* 60-67.5 2.9-2.7 NM 10 64 3.0 90.5 98 0 780 0.024 191 0 1097 0.01 189 0 743 0.05 799 0 728 0 0.61 0.85 0.49 0.51 NM = not measured; *maximum reading on flow meter; † Test 1 in Fig. 3.8 Interpretation of the computed steady-state flow field supports these findings. The single impeller induces suspension flow away from the impeller towards the rear wall of the chest. Flow is turned by this wall, and loops back along the side walls, suspension surface and chest floor before being drawn back into the impeller suction. Some flow exits the chest, with competition between flow induced into the impeller suction and flow exiting the stock line governing the channeling observed. The region immediately surrounding the impeller forms a cavern having high shear and velocity gradients, and is well mixed. A number of methods can be used to identify the cavern using computed flow data. The first method locates the surface where References start on page 88 Chapter 3 76 the shear stress imposed by the flow exceeds the suspension yield stress, one criterion adopted to define caverns in yield stress fluids (Amanullah et al., 1998). The iso-surface where this condition is satisfied was calculated from the simulated flow field for MC-A and is shown in Fig. 3.6(a) with τy = 113 Pa. A second criterion used to identify caverns involves locating the enveloping surface where the magnitude of the local tangential velocity falls to some low value, typically 1% of the impeller tip speed (Jaworski et al., 1994). Figure 3.6(b) shows the surface where the local velocity is 0.21m/s, 1% of the impeller tip speed of 21m/s. Figure 3.6 confirms that a large region of MC-A is involved in weak suspension flow. While suspension flow almost completely fills the vessel, fibre floc disruption and good mixing occur in the high-shear region that immediately surrounds the impeller (where τ > τy). The region beyond this active zone is in motion due to circulation caused by the impeller and inlet/exit generated flows. The location of the exit pipe is close to the cavern boundary which determines whether flow will be directed to the impeller region or will exit the chest. Cavern size would decrease if the suspension yield stress increased (for example, through addition of a higher percentage of softwood pulp to the chest). Under these circumstances, the location of the exit pipe could increase undesired bypassing. Good chest design uses inlet/exit locations that force the suspension to traverse the well-mixed impeller zone. By visualizing path-lines of mass-less particles released at the inlet to the chest (Fig. 3.6(c)) we see that most are drawn into the impeller suction, in agreement with the experimental data. References start on page 88 Chapter 3 77 Figure 3.6. Cavern formation in MC-A: (a) Iso-surface where the shear-stress equals the suspension yield stress; (b) Iso-surface where the velocity equals 0.01 VTip. (c) Trajectories of particles released at main stock inlet (Inlet 1). (a) (b) References start on page 88 Chapter 3 78 (c) The dynamic simulations for MC-A are compared to real-time data collected during the bump test in Fig. 3.7. The computational input signal was chosen to match the experimental signal generated at the mill (three linear input regions were used) as shown in Fig. 3.7(a). Note that the high frequency content of the incoming stream was not simulated. The computational results are compared to the experimental data in Figure 3.7(b) (the scale has been adjusted to account for offset and gain issues in the mill data) and shows that the mixing dynamics are simulated well. The model parameters identified using the simulated (CFD) bump tests are included in Table 3.2. In this case, the identified time constant was about 1097s which gave a mixed volume equal to 85% of the available chest volume. The simulated MC-A displays better mixing than it is really attained. This can be related to a number of factors; particularly the approximations introduced in the rheology model used to describe the pulp suspension. Additionally, errors in flow measurements in the plant data and the quality of the bump test data contribute to differences References start on page 88 Chapter 3 79 between the dynamic parameters identified. Despite these deviations, the CFD results confirm the experimental analysis that the chest is not ideally mixed. Mill Data Simulation Consistency at Inlet ( - ) 0.0381 0.0376 0.0371 0.0366 0.0361 0.0356 0 2000 4000 Time (s) 6000 (a) 8000 Mill Data 0.0381 Simulation 0.03682 0.0376 0.03662 0.0371 0.03642 0.0366 0.03622 0.0361 0.03602 0.0356 0.03582 0 1000 2000 3000 4000 Time (s) 5000 6000 7000 Consistency at exit - CFD ( - ) Consistency at exit - MC-A ( - ) 0.03702 8000 (b) Figure 3.7. MC-A: (a) input and (b) output signals for experimental and computational data. Consistency scale on the right axis corresponds to the computational output signals with adjusted gain and offset. Test conditions and model parameters as given in Table 3.2. References start on page 88 Chapter 3 3.4.2 80 Latency Chest (LC-B) Latency chest LC-B receives its thermo-mechanical pulp feed from a refined stock chest located upstream. During normal operation, water is added into the feed line at the exit of the refined stock chest to maintain a steady consistency entering LC-B. A similar control loop maintains a steady discharge consistency from LC-B. Thus, to generate the step change shown in Fig. 3.8, the consistency controllers must be turned off and the incoming consistency changed manually (by 8 − 9% in this case). Unfiltered data was collected at a sampling time of 2s from both consistency sensors. The signal quality is good, but as observed in other mills, the consistency meters were not properly calibrated. Also, there appears to be a continuous oscillatory component present in the signal from the discharge stream that could be due to equipment problems such as valve stiction. Again, the problems with meter calibration and abnormalities in the transmitted signal make identifying the mixing model difficult. In an attempt to determine the variability of mixing quality over a range of operating conditions, tests were carried out for various operating modes (shown in Fig. 3.8). The result of the first step change is shown in Table 3.2 while the dynamic model identified for the complete set of operating conditions of Fig. 3.8 is given in Table 3.3. Mixing performance can vary significantly depending on chest operating conditions such as inlet consistency, fill volume etc. For example, the data from Table 3.2 and Fig. 3.8 correspond to a relatively high flow rate and suspension chest level. Under these conditions the time constant identified (743s) is about half of nominal value (τnom ~1450s) which could be related to bypassing ( i.e. when the level of the chest increases a greater portion of the input stream moves along the suspension surface which reduces the time constant), although this was not identified by the gray box model (f = 0.05). At lower suspension flow and stock chest levels (Test-2 in References start on page 88 Chapter 3 81 Table 3) the time constant increased to about 80% of the nominal residence time. Under similar conditions (Test-4) the fraction of bypassing increases. This latter finding agrees with dynamic tests made in the lab-scaled chest (Ein-Mozzafari, 2002). Table 3.3. Experimental data for LC-B (Tests shown in Fig. 3.8) Operating Conditions Q (m3/min) Cm (%) H (%) Model Parameters T1 (s) τ1 (s) τ2 (s) f Mixing Parameters VFM/VT Test – 2 Test – 3 Test – 4 8.0 2.9 58 9.5 2.7 57 8.0 2.94 57 432 0 1181 0.001 203 0 1210 0 319 0 1126 0.16 0.80 0.98 0.66 Figure 3.8. Raw experimental bump test data for LC-B (Tests 1−4 indicated by the gray dotted lines) References start on page 88 Chapter 3 82 CFD simulations were carried out to simulate the first step change shown in Fig. 3.8. The yield stress of the thermo-mechanical pulp was estimated to be 90 Pa (Bennington et al, 1990) with µ0 set to 100 Pa.s as before. Steady-state simulations were first carried out to obtain the flow-field in the vessel followed by dynamic simulations. ∆t = 1s was used in the dynamic simulations for LC-B, since this gave solutions independent of time discretization (verified by initial simulations at 0.5 and 5s). Analysis of the CFD simulation showed that there was negligible bypassing in the chest (as found experimentally). The simulated chest had a time constant of 728s compared to 743s determined from experimental data showing that only 50% of the chest volume was mixed. An examination of the computed flow field explains this low mixing volume. Pulp enters the chest from the inlet positioned directly above an impeller (Fig. 3.2) and moves in a u-shaped trajectory through the chest. Along this path the suspension is alternately exposed to high-shear mixed regions near the impellers and plug flow low-shear regions between them. Figures 3.9(a) and (b) show the caverns produced in this chest using the two methods discussed earlier. It is seen that the regions of higher shear and intense mixing occupy only a fraction of the chest volume, which is also indicated by the low mixed volume determined by the dynamic mixing model. Path-lines of mass-less particles released at the inlet (Fig. 3.9(c)) show movement of pulp from the inlet through each of the four impellers and subsequent discharge from the tank. These path lines also illustrate the recirculation paths followed by the suspension near the inlet and the first impeller. References start on page 88 Chapter 3 83 (a) (b) (c) Figure 3.9. Cavern formation in LC-B: (a) Iso-surface where the shear-stress equals the yield stress; (b) Iso-surface where the velocity equals 0.01 VTip. (c) Trajectories of particles released at the stock inlet. References start on page 88 Chapter 3 84 The dynamic simulations for LC-B were compared to the real-time data collected during the bump test in Fig. 3.10. Again the computational input signal was chosen to match the experimental signal generated at the mill (Fig. 3.10(a)). The chest response to the change (both measured and predicted) is given in Fig. 3.10(b). Again it was necessary to adjust both the consistency offset and gain to compare agreement between the model and experimental data. The most notable difference between the computed and experimental output signals is the much longer delay obtained from the simulated data. This is most noticeable once the computed input and output signals are used to identify the model parameters. The value of T1 identified from the computations is about four times larger than the one identified from the measured signals. One possible explanation for this is that the grey-box model (Fig. 3.1) does not include a term for recirculation within the chest, which is important in this configuration as shown in Fig. 3.9(c). Also, the location of the consistency meters in the mill was not identical to that of the inlet and outlet boundaries in the CFD model because of the restrictions imposed by the selected boundary conditions (uniform velocity profile at the inlet and outflow boundary condition at the outlet, both needed to be specified at least at 10d away from the chest). The dynamic model also does not account for level changes or flow rate changes, and requires that the bump chest be achieved in steady-state conditions. These experimental challenges must also be dealt with for a more accurate comparison between computed and experimental data to be made. References start on page 88 Chapter 3 85 0.0305 Mill Consistency at inlet ( - ) 0.03 Simulation 0.0295 0.029 0.0285 0.028 0.0275 0.027 0.0265 0.026 2000 4000 6000 Time (s) 8000 Consistency at exit - Mill ( - ) 0.0305 Simulation 0.03 Mill 0.0295 10000 (a) 0.031 0.03 0.029 0.029 0.0285 0.028 0.028 0.027 0.0275 0.026 0.027 0.025 0.0265 Consistency at exit - CFD (-) 0 0.024 0.026 0 2000 4000 6000 Time (s) 8000 10000 (b) Figure 3.10. LC-B: (a) input and (b) output signals for experimental and computational data. Consistency scale on the right axis corresponds to the computational output signals with adjusted gain and offset. Test conditions and model parameters as given in Table 2. 3.5 Conclusions CFD simulations provide insight into mixing occurring in industrial pulp stock chests, allowing the location of regions of active mixing (caverns) and the bypassing created by interaction between chest and impeller geometry and by suspension properties and process conditions, to be References start on page 88 Chapter 3 86 determined. The agreement between the computational simulations and the experimental measurements was good, with experimental trends displayed well by the simulations. However, several challenges remain both for simulating pulp mixing (including specifying the suspension rheology and improving the specification of the boundary conditions) and measuring mixing quality in the mill experimentally (creating a good bump test, obtaining uncompressed data having a good signal-to-noise ratio, dealing with instrument calibration errors and identifying the parameters of a suitable mixing model for the chest). Indeed, the largest challenge of all may be obtaining data from an operating mill at steady-state. 3.6 Nomenclature Cm mass concentration of pulp suspension (consistency) din diameter of inlet feed pipe, m dout diameter of exit pipe, m f fraction of the incoming feed stock bypassing the well-mixed region GChest transfer function representing mixing in the pulp chest h height of the inlet feed pipe above the pulp bed, m k consistency index, Pas L length of chest, m W Width of the chest, m H nominal height of the chest, m N impeller rotational speed, rpm D impeller diameter, m NP impeller power number, dimensionless References start on page 88 Chapter 3 C 87 impeller clearance from wall, m Npump impeller pumping number, dimensionless VTip tip velocity, m/s Q flow rate of pulp through the chest, m3/s R fraction of flow recalculating in mixing zone T1 time delay in bypass region, s VT total volume of pulp suspension in chest, m3 VFM volume of pulp suspension in the chest that is fully mixed, m3 Greek letters • γ shear rate, s-1 τy yield stress, Pa µo yielding viscosity, Pas τ1 time constant identified for the well-mixed zone, s τ2 time constant for mixing in bypass zone, s τnom nominal time constant of vessel, τnom=VT/Q, s References start on page 88 Chapter 3 3.7 88 References Amanullah, A., S.A., Hjorth, and A.W., Nienow, 1998. A new mathematical model to predict cavern diameters in highly shear thinning power law liquids using axial flow impellers. Chem Eng Sci, 53, 455–469. Bakker, A., and J.B., Fasano, 1997. A computational study of flow pattern in an industrial paper pulp stock chest with a side-entering impeller. AIChE Symposium Series, 89 , 293, 118−294. Also available at http://www.bakker.org/cfmbook/cfmbook.htm Beath, L.R., M.T., Neil, and F.A., Masse, 1966. Latency in mechanical wood pulp. Pulp and Paper Mag. Canada, 67, 10, T423−T430. Bennington, C.P.J., 1988, Mixing Pulp Suspensions. PhD dissertation, University of British Columbia, Vancouver, Canada. Bennington, C.P.J., R.J., Kerekes, and J.R., Grace, 1990. The yield stress of fibre suspensions. Can J Chem Eng, 68, 748−757. Bialkowski, W.L., 1992. Newsprint variability and its impact on competitive position. Pulp and Paper Canada, 93, 11, T299 – T306. Ein-Mozaffari, F., 2002. Macroscale mixing and dynamic behaviour of agitated pulp stock chests, PhD. Thesis, University of British Columbia. Ein-Mozaffari, F., G.A., Dumont, and C.P.J., Bennington, 2003. Performance and design of agitated pulp stock chests, Appita J., 56, 2, 127–133. Ein-Mozaffari, F., C.P.J. Bennington and G.A. Dumont, 2004. Dynamic mixing in agitated industrial pulp chests. Pulp Paper Can. 105, 5, 41−45. References start on page 88 Chapter 3 89 Ein-Mozaffari, F., C.P.J., Bennington and G.A., Dumont, 2005. Suspension yield stress and the dynamic response of agitated pulp chests. Chem Eng Sci, 60, 2399–2408. Ford, C., 2004. CFD Simulation of mixing dynamics in agitated pulp stock chests. MASc thesis, University of British Columbia, Vancouver, Canada. Ford, C., F. Ein-Mozaffari, C.P.J. Bennington, and F. Taghipour, 2006. Simulation of mixing dynamics in agitated pulp stock chests using CFD. AIChE J, 52, 10, 3562–3569. Ford, C., C.P.J. Bennington, and F. Taghipour, 2007. Modelling a pilot-scale pulp mixing chest using CFD. J. Pulp Paper Sci, 33, 115−120. Jaworski. Z., A.W. Pacek and A.W. Nienow, 1994. On flow close to boundaries in yield stress fluids. Chem. Eng. Sci, 49, 3321−3324 Roberg, K.E., 1997. Numerical simulations of the flow of fibre suspension in a side-entering mixer. Récents Progrès en Génie des Procédés, 11 (51), 203−210. Soltanzadeh, A., G.A., Dumont, S., Bhattacharya and C.P.J., Bennington, 2009. Estimation of residence time and chanelling in pulp mixers. Nordic J. Pulp and Paper Sci, in press. Wichterle, K. and O., Wein, Threshold of mixing of non-Newtonian liquids, Int. Chem Eng, 21, 116–121. Wikstrom, T. and A. Rasmuson, 1998. The agitation of pulp suspensions with a jet nozzle agitator. Nordic Pulp Pap. Res. 13, 2, 88–94. Yackel, D.C., 1990, Pulp and paper agitation: The history, mechanics and process. TAPPI Press, Atlanta, GA. References start on page 88 90 Chapter 4 Investigation of the Flow field in a Rectangular Vessel equipped with a Sideentering Agitator 4.1 Introduction Fluid mixing is a central unit operation in the chemical and biological process industries with the stirred tank being one of the most used devices. Empirical rules that rely on global operating parameters such as the power number and flow number have been widely used to design mixing systems. However, the investigation and optimization of these processes necessitates fundamental understanding of the flow structures and velocity distributions found within the agitated vessel (Mavros, 2001). The quality of mixing required varies according to the specific process at hand. For the most part, bulk motion is required in order to guarantee good circulation within the mixing vessel. In addition, in suspension and dispersion operations, acceptable local velocities in regions remote from the impeller are required. Therefore the selection of an impeller for a specific duty is governed by the flow patterns and the velocity profiles that can be generated in a particular fluid. The recent increased availability of sophisticated experimental and efficient computational tools enables extensive analyses of the flow field in stirred vessels that enhance our understanding of the phenomena occurring within these systems. A wealth of literature exists on the A version of this chapter will be submitted for publication. C. Gomez†, F. Taghipour C.P.J. Bennington., 2009. Investigation of the flow in a rectangular vessel equipped with a side entering agitator. †C. Gomez formerly C.Ford Chapter 4 91 implementation of CFD to simulate the flow fields generated in cylindrical agitated tanks with radial and axial impellers, for example: Adams and Barigou (2007), Couerbe et al. (2008), Khopar et al.(2004), Ranade et al., (2001), Torré et al. (2007), Yoon et al.(2001) and Zalc et al.(2001). A number of studies have also been made in non-cylindrical geometries (e.g. Bakker, 1993; Bhattacharya et al., 2008; Ford et al., 2007; Ford et al., 2006). The computational results must be validated against experimental data as their calculations depend strongly on several aspects including modeling of the system’s geometry, discretization method, mesh structure, viscosity definition and rheology model in the non-Newtonian case, and turbulence model used if any, for the momentum conservation equations. The flow fields obtained in a stirred tank have been studied experimentally using multiple techniques including hot wire anemometry, laser doppler anemometry, laser doppler velocimetry, particle image velocimetry, tomography, ultrasonic doppler velocimetry, laser induced fluorescence, etc. as reviewed by Mavros (2001). The vast majority of these experimental works have been in cylindrical vessels (baffled and unbaffled) equipped with centered top-entering agitators (e.g. Couerbe et al., 2008; Escudié and Liné, 2003; Khopar et al., 2004; Ranade et al., 2001; Torré et al., 2007; Yoon et al., 2001, Wu and Patterson, 1993). Extensive measurements of the flow induced by Rushton turbines have been reported employing both single point techniques such as laser doppler velocimetry (LDV) (e.g. Lee and Yianneskis, 1998; Stoots and Calebrese, 1995; Wu and Patterson, 1989) and global techniques like particle image velocimetry (PIV) (e.g. Bugay et al., 2002; Couerbe et al., 2008; Escudié and Liné, 2003; Khopar et al., 2004; Ranade et al., 2001; Sharp and Adrian, 2001; Torré et al., 2007; Yoon et al., 2001). The flow induced by centered top entering axial flow impellers has also been measured by various investigators, using LDV (e.g. Jaworski et al., 1996; Mavros et al., 1996; Schäfer et References start on page 123 Chapter 4 92 al., 1998); and PIV in cylindrical vessels (e.g. Bugay et al., 2002) and in square tanks (e.g. Kilander and Rasmuson, 2005). To our knowledge, no work has been reported in the open literature on the combined PIV and CFD study of flow fields in mixing vessels with side entering agitators. Given the wide application of this configuration for blending and sediment suspension in the chemical process industries, it is important to advance our knowledge and understanding on the flow and mixing in mechanically agitated vessels of this kind. Particle image velocimetry was selected for the investigation presented here. PIV is one of the most advanced experimental methods available to visualize the main characteristics of a flow field, providing a nonintrusive measurement of the instantaneous planar velocity field over a global domain. As the name suggests, it records the position of small tracer particles over a defined time interval to extract the local fluid velocity. PIV requires four basic components: an optically transparent test section that contains the fluid seeded with small, fluorescent, flowfollowing particles, an illuminating light source (laser), a recording device (a charge-coupled device (CCD) camera), and a computer algorithm to process (cross-correlate) the recorded images. Due to the complexities associated with the highly three-dimensional flow in stirred tanks, additional care must be taken in the application of PIV techniques in these configurations. The present work focuses on the application of CFD coupled with PIV to analyze the macroscale hydrodynamics in an axially agitated rectangular tank. The geometrical configuration of this agitated vessel was derived from the scale down of industrial systems for agitation of paper pulp suspensions (Gómez et al., 2008). Similar applications (e.g. large rectangular or cylindrical tanks equipped with side entering mixers) are found in other industries including the petroleum industry (for crude oil and gasoline storing) and the foods industry (for blending fats and oils) References start on page 123 Chapter 4 93 (Oldshue, 1983). The fluids encountered in these applications are generally highly viscous or non-Newtonian with complex rheology. In this study we investigate the global flow field of a viscous Newtonian fluid (98 wt.% glycerin solution) in an agitated vessel of non conventional, yet industrially relevant, geometry. This work establishes the basis of future studies on the effects of non-Newtonian behavior in these mixing systems. 4.2 4.2.1 Experimental Methods Experimental Setup Figure 4.1(a) shows a schematic of the experimental setup used. The system consists of a plexiglass model of a 1:18 scaled mixing chest with length, L = 36.5 cm, width, W = 24.5 cm, and height, Htank = 40cm. The rectangular tank is equipped with a side-entering scaled version of the Maxflo Mark II impeller (Chemineer, Inc., Dayton, OH) used industrially for paper stock agitation (3 bladed; D = 10.16 cm; pitch ratio = 0.44), driven by a 0.33 HP motor controlled by an inverter to enable digital control of the motor speed. The impeller is centered on the z direction (corresponding to the vessel width). The impeller clearance from the side wall was C = 7.5 (C/W = 3.2) and clearance from the bottom wall was 7.4 cm. A shaft-mounted rotary torque transducer (Staiger Mohilo 0260DM) was used to obtain accurate measurements of the torque (± 0.002 N·m). The rotation speed was varied from N= 59−386 rpm with a given speed maintained to ± 4 rpm during measurements. The flow conditions were characterized by the Reynolds number, Re, which was computed using measured values of viscosity and density References start on page 123 Chapter 4 Re = 94 ρ ND 2 µ (4.1) The working fluid was a 98 wt% glycerin solution of density ρ = 1.256 kg/m3 at 22°C. The rectangular vessel was filled to a height of H = 30.6 cm and left to stand overnight to allow entrained air to escape. Due to the hygroscopic nature of glycerin, measurement of the density and viscosity of the solution prior to testing was necessary to gurantee an accurate definition of the material properties for CFD simulations. The solution’s viscosity, µ , of 0.68 Pa·s, was measured immediately after the experiments at the exact same temperature registered for the fluid in the mixing vessel, 22°C, using a stress controlled parallel-plate rheometer (Physica MCR 501). The flow in the stirred vessel was investigated using PIV. The Flow-Map 2D PIV system from Dantec dynamics Inc. (Ramsey, NJ) was used to map the velocity profiles. It consists of a dualhead Nd:YAG laser (New Wave Research: Model SoloIII-15 Hz) that emits 10-ns laser pulses with energies of up to 50 mJ/pulse. The collinear frequency doubled beam (1064 nm to 532 nm) is spread into a 2 mm thick sheet through a set of cylindrical lenses attached directly to the laser head. The particle images were recorded on a high-resolution progressive-scan interline CCD camera (Hamamatsu, model HiSense MkII) having a resolution of 1344 × 1024 pixels. A 514nm line filter was used to minimize the detection of stray light on the cells of the CCD chip. The camera and laser were connected to a PIV hub with a programmable synchronizer that controlled the timing of the laser illumination and image acquisition. The hub is connected to a personal computer (PC) where the experimental sequence is defined and the resulting data stream is captured and analyzed asynchronously. References start on page 123 Chapter 4 95 (a) (c) (b) (d) Figure 4.1. (a) Schematic of the experimental setup showing the front view of the mixing vessel, the PIV acquisition system and triggering device used. (b) Location and size of the measuring areas. Dashed lines illustrate the locations at which the velocity profiles were evaluated. The horizontal line is located at impeller height and the vertical line is located at a distance d cm (varied from 3−10cm) in front of the impeller. (c) PIV Illumination planes. (d) Photograph showing the flow patterns observed in the PIV illumination plane aligned with the impeller tip. References start on page 123 Chapter 4 96 For the PIV experiments, the glycerin was seeded with ~50 µm polyamide 12 seeding particles (PSP, Dantec dynamics Inc.) of density ρp = 1.03 g/cm3. The seeding was added to smaller aliquots of glycerin (∼50ml) and then inserted below the free surface at a low rotational speed (N = 80−120 rpm). The impeller speed was then increased to 1200−1500 rpm to provide a uniform distribution of the seeding particles in the fluid. Nearly all of the particles remained suspended and followed the flow field closely. A small amount of particles floating on the surface were removed using a size 400 sieve. The particle size was carefully chosen so that once combined with the image size desired in our experiments (i.e. macro-scale global flow field measurements), the f-number of the lens and the laser intensity; the average size of the reflected image of the seeding particles on the CCD chip was about 2 pixels in diameter (Adrian, 1997). Hence, the particles were slightly larger than seeding materials used in other PIV studies on smaller mixing systems (e.g. Bugay et al., 2002; Escudié and Liné, 2003; Torré et al., 2007; Zalc et al., 2001). The velocity measurements taken were phase locked with respect to the position of the impeller in order to avoid the detection of the periodic fluctuations induced by the motion of the blades. The triggering device used for this purpose consists of an inductive proximity switch (Wenglor IB40) and a thin aluminum rod (of 1mm diameter) mounted on the impeller shaft. Each time the piece of metal passed the sensor, the image acquisition was triggered. The aluminum rod was located at θ=0° (aligned with the y axis), and at θ= 90° (aligned with the z axis), thus two different impeller positions were studied. The accuracy with which the impeller position could be determined with this device was estimated based on the instrument delay. The inductive sensor’s delay is 740 µs, which corresponds to 0.3° of a full impeller rotation at the minimum rotational speed tested, N = 59 rpm, and to 1.8° at the maximum one, N = 386 rpm. References start on page 123 Chapter 4 4.2.2 97 Experimental Procedures To measure the flow profiles, a 2D slice of the flow field is illuminated using a pulsing laser light sheet (pulse duration is 5−10ns) that freezes the position of the particles. A camera is used to record the particles position in the illuminated plane over a short time interval (∆t between laser pulses). Thus, two images of the fluid flow are obtained with a known time interval between them. Comparison of the two images by means of cross correlation techniques allows the determination of the corresponding velocity vector map in the plane of the laser sheet (Raffel et al., 1998). Two sequential images are sub-sampled using square regions called interrogation areas (IA) that define the sampling resolution. Within these sampling areas, an average spatial shift of seeding particles may be observed from the first camera image to the next, provided there is flow in the illuminated plane. Common particles need to exist in the IAs which are then correlated by crosscorrelation analysis using fast Fourier transform (FFT) techniques. A high cross-correlation peak is observed when many particles match up with their corresponding spatially shifted partners and small cross-correlation peaks are observed when individual particles match up with other particles entering or leaving the IA (either moving in the plane or crossing the plane), contributing to random correlations or noise. This analysis enables the calculation of the most probable particle displacement in each interrogation area. The technique is repeated at each IA thereby building the instantaneous velocity field. A digital filter was applied to eliminate the DC component of the measurement (noise-floor) prior to the determination of the correlation peaks and subpixel interpolation. This filter removes the noise introduced by multiple particle reflections and background light. Additionally, a References start on page 123 Chapter 4 98 Gaussian window function was applied to the interrogation areas (i.e. the window width was 85% of the IA width) in order to suppress information from particles near the edges of the IA where phantom correlations are likely to dominate (Brown et al., 2004). When window functions are used, it is recommended to increase the size of the IA (Dantec dynamics Inc., 2003). Consequently, the size of the IA used was 32 pixels and a 25% overlap between adjacent IAs was applied to optimize the use of information in the image map. An adaptive multi-pass cross-correlation technique was used for all velocity vector calculations in order to cover a larger area and still capture the flow details (Adrian, 1986). When this method is employed, the overlap of particles at small displacements is eliminated, and particle images that left the interrogation area during the time between the two light pulses (in plane drop out) can be captured, thus increasing the dynamic range (Raffel et al., 1998). For each instantaneous measurement, the displacement vectors in the correlation plane for which the ratio of the signal peak to the noise peak was greater than 1.3 were rejected. No further filtering or validation was employed on the remaining valid vector map so as to truly capture any peculiarities of the flow. Only the valid vectors were accounted for in the statistical-averaging procedure. Tests concerning the number of images or sample sets needed to get accurate time-averaged measurements were performed by taking a set of 1800 instantaneous measurement pairs at Re= 36 and Re=120. Vector statistics of subsets of 1600, 1400, 1200, 1000 and 800 samples were compared against the full set of 1800 measurements. It was determined that a sample size References start on page 123 Chapter 4 99 of 1000 double-images produced a stable time-averaged result. Thereafter, all subsequent measurements were performed with 1000 measurements. The accuracy of the velocity measurements depends on the seeding concentration, the time between laser pulses, the size of the PIV measurement area and the spatial resolution. In order to determine the flow field in a large domain, the size of the measurement area should be as large as possible, however the spatial resolution must be small enough so that essential velocity structures can be resolved and measured. Therefore, given our objective to investigate the global flow field in the mixing vessel, the largest possible image size that satisfied the PIV measurement criteria described by Adrian (1997), Raffel et al., (1998) and Dantec (Dantec dynamics Inc., 2003); was chosen. Accordingly, the conditions chosen to minimize error were: • The particle seeding concentration was adjusted so that a minimum of 6−10 particles per IA was attained, giving a high signal to noise ratio. • The particle displacement was always kept below 25% of the lenght of the IA by selecting the appropriate time delays between laser flashes as a function of the local hydrodynamics and the parallax effects due to the through plane motion of the particles. • A slightly defocused camera coupled with the selected image and particle size (~50 µm) was used in order to obtain adequate particle images (around two pixels in diameter) Size of the PIV measurement area The PIV measurement locations chosen are shown in Figure 4.1(b). The size of the largest image in which the above mentioned criteria were fulfilled was 21.3 × 16.2 cm2 (Fig. 4.1(b)). Three z planes (the coordinate system orientation is shown in Fig. 4.1(c)) were investigated: one References start on page 123 Chapter 4 100 plane located at the center of the impeller (z = 0 cm), one plane aligned with the impeller tip (1.7 mm from the blade tip, z = -5.25 cm) and one plane located 2 cm away from the side wall (z = -10.25 cm) (these three planes are located ∼5cm from each other (Fig. 4.1c)). Tests with a smaller image size were performed to validate the accuracy of the global velocity measurements taken. Data were obtained with a measurement area of 9.3 × 7.1cm2 located 4.3 cm in front of the impeller and covering 7.1 cm on the negative y direction measuring from the impeller center as illustrated by the gray lines in Fig. 4.1(b). The small measurement area examined is large when compared with experiments performed by Yoon et al. (2001) to characterize the flow in the immediate vicinity of a Rushston turbine, but it is comparable in size to the PIV areas used in many other studies (e.g. Bugay et al., 2002; Bujalski et al., 2006; Couerbe et al., 2008; Escudié et al., 2004; Escudié and Liné, 2003; Svensson and Rasmuson, 2006; Torré et al., 2007). The large PIV measurement area covers about 60% of the tank length and 53% of the liquid height, and it is comparable to the large image PIV experiments done by Bugay et al., (2002) and the standard image size used by Collignon et al., (2008). For both cases, the measuring area was divided into interrogation areas of 32×32 pixels. The size of these IAs was 5 × 5 mm2 and 2.2 × 2.2 mm2 for the large and the small sampling regions respectively. The results obtained from these two experiments are compared using the profiles of the two components of the mean velocity vector. The axial and vertical velocities and the axial and vertical positions are non-dimensionalized as follows: References start on page 123 Chapter 4 U* = 101 U U tip V* = V U tip x* = x L y* = y H (4.2) The vertical profiles of the axial velocity, U*, and the vertical velocity, V*, at Re = 36, are plotted at x*= 0.16 (d = 6cm in Fig. 4.1(b)) in Figure 4.2(a) and 4.2(b) (a) IA size, mm2 5×5 2.2 × 2.2 IA size, mm2 5×5 2.2 × 2.2 (b) Figure 4.2. Profile of the mean (a) axial and (b) vertical velocity at a line located at x*= 0.37; d = 6cm (see Fig. 4.1(b)) in the measurement plane aligned with the impeller tip (z = -5.25 cm), at N = 118 rpm, Re = 36. Time interval between laser pulses One of the main difficulties in providing adequate velocity measurements in a large domain, like the one chosen here, lies in selection of the time interval between laser pulses because of in plane and out of plane velocity values possible in this particular configuration. For instance, the fluid velocities are much higher near the impeller but the magnitude of the velocity diminishes rapidly in the flow with most of the flow captured in the large PIV image moving at a very low velocity. A compromise must be reached when selecting the time between pulses in order to enhance the number of valid velocity vectors within the image. Thus, some accuracy in the measurements is lost both in the immediate vicinity of the impeller and far away from it. References start on page 123 Chapter 4 102 The time delay between laser flashes (∆t) depends on the local velocities in the system. In practice, pulse separation must be optimized to ensure that the maximum possible number of particle pairs is captured. A reliable way to achieve this is to define the particle displacement as one quarter of the IA side, LIA (Dantec dynamics Inc., 2003; Raffel et al., 1998). The separation ∆t is then given as a function of the object to image scale factor (S), LIA and the maximum velocity within the system, the impeller tip speed (Utip), as follows: ∆t = 0.25 LIA S U tip (4.3) However, the through-plane movement of particles affects the accuracy of the in-plane measurements. Thus, ∆t and the laser sheet thickness ( zA ) must be chosen in relation to the velocity perpendicular to the measurement plane (Uz), so that the following relationship is also satisfied (Dantec, 2003): U z ∆t ≤ 0.25 zA (4.4) In the immediate vicinity of the impeller, Uz also reaches the maximum velocity value given by Utip at specific angles of rotation determined by the blade pitch angle. At this location, ∆t given by Eqn. (4.4) is 40% of the ∆t value from Eqn. (4.3). Considering the high viscosity of the working fluid and that a significant portion of the domain covered in the PIV measuring area is away from the impeller, the upper limit of the velocity range given by Utip can be relaxed to allow for a larger ∆t thus increasing the dynamic range of the measurements (Adrian, 1997). Previous PIV and LDV measurements in the vicinity of top-entering axial impellers (Escudié and Liné, 2003; Myers et al., 1997) have shown that the actual fluid velocities are only a fraction of the tip speed reaching a maximum of 30% of the tip speed at about 5mm below the impeller. We References start on page 123 Chapter 4 103 investigated the sensitivity of the measurements with respect to ∆t by substituting Uz = 0.5 Utip , Uz = 0.3 Utip and Uz = 0.2 Utip in Eqn. (4.4). Three different ∆t values were tested. For example, at Re = 36, ∆t values of ∆t1 = 1680 µs, ∆t2 = 2670 µs , ∆t3 = 4006 µs were applied. Similarly, at Re = 83, the ∆t values used were ∆t1 = 780 µs, ∆t2 = 1180 µs, ∆t3 = 1770 µs. Figure 4.3 shows the vertical profiles of axial velocity measured at a line located at x * = 0.33 (d = 4.5cm in Fig. 4.1(b)) using the three different ∆t values, for Re = 36 and 83. The major differences between measurements are found around the impeller which is located at y*= 0.39. The variation between measurements was quantified in terms of the root-mean-square (RMS) deviation: 1/2 ⎡ 1 i=N 2⎤ ⎢ N ∑ (u∆t1 − u∆t 2 ) ⎥ ⎦ RMS = ⎣ i =1 1/2 i=N ⎡1 2⎤ ⎢ N ∑ (u∆t 2 ) ⎥ ⎣ i =1 ⎦ (4.5) The deviation between measurements decreased as ∆t increased. For instance, the RMS deviation between the measurements made with ∆t1 and ∆t2 is was17% for Re = 36 and 11.6% for Re = 83. This RMS deviation value decreased to 6.7% for both Re numbers when the measurements taken with ∆t2 were compared with ∆t3. Therefore, all the subsequent measurements were taken using ∆t calculated from Eqn. (4.4) with Uz = 0.3 Utip at the measurement planes near the impeller, and Uz = 0.2 Utip at measurement plane near the wall. References start on page 123 Chapter 4 (a) 104 ∆t, µs 1680 2670 4006 ∆t, µs 780 1180 1770 (b) Figure 4.3. Effect of ∆t on the mean value of the axial velocity measured at a line located at x*= 0.33; d = 4.5cm (see Fig. 4.1(b)) in the measuring plane aligned with the impeller tip. (a) N = 118 rpm, Re = 36 (b) N = 267 rpm, Re = 83. Phase averaging The PIV data collected at two different impeller positions θ = 0° and θ = 90° were compared to determine whether the periodic fluctuations induced by the motion of the blades were significant. A location close to the impeller blade tip was chosen for the comparison of the velocity measurements. Figure 4.4 compares the mean values of the axial and vertical velocities, measured at the two different angles, at a line located at x* = 0.29 (d = 3cm in Fig. 4.1(b)) in a plane aligned with the impeller tip. The velocity profiles are very close to each other. The RMS deviation between these two measurements was 3.8% for V* and 4.1% for U*. Hence, it was confirmed that phase averaging is not essential in the case of shaped impellers, such as an axial impeller, as suggested by Bugay et al. (2002) and Weetman and Oldshue (1988). Measurement reproducibility was verified by repeating the full measurement set of 1000 instantaneous velocity fields three times and comparing the mean value of the axial and vertical References start on page 123 Chapter 4 105 velocities by computing the RMS deviation. Since a 25% overlap was used, the measuring plane contained 2310 vectors, and the average RMS deviation calculated between identical experiments was 6.0% for U* and 4.5% for V*. (a) (b) Figure 4.4. Effect of the phase angle θ on the mean value of the (a) axial and (b) vertical velocity measured at a line located at x*= 0.29; d = 3cm (see Fig. 4.1(b)), in the measuring plane aligned with the impeller tip. N = 267 rpm, Re = 83. 4.3 Computational Techniques The numerical simulations of the flow field in the rectangular mixing vessel were performed using FLUENT 6.3 (Lebanon, NH). Aside from the impeller shaft and the small mounting supports of the impeller hub into the shaft, the computational model of the system represents the entire mixer geometry exactly. Re from 18 to 120 were simulated by changing the agitation rate accordingly (using the measured viscosity, µ = 0.68 Pa.s and density, ρ = 1.256 kg/m3). The geometry was decomposed into four blocks to allow for a more efficient hybrid mesh (generated in GAMBIT 2.4): one cylindrical volume containing the impeller; two semiReferences start on page 123 Chapter 4 106 rectangular volumes above and below the cylinder and one rectangular volume connected to the common left side of these three volumes. An unstructured grid of tetrahedral and prism elements, growing in size from the impeller outwards, was created within the cylindrical volume. The remaining blocks were meshed using a structured grid of hexahedral elements, with controlled size near vessel walls and surfaces around the cylindrical block, which permitted a significant reduction in the number of computational cells. After discretization of the computational domain, the mass and momentum conservation equations for an incompressible fluid were solved using a multiple reference frame (MRF) model. The rotating frame was defined by the cylindrical volume containing the impeller. Within this sub-domain, the mass and momentum conservation equations include additional terms that account for the acceleration of the coordinate system given by the impeller rotational speed. The remaining volumes were stationary. In the solution process, the continuity of the absolute velocity is enforced at the boundaries between sub-domains. The impeller was defined as a wall with zero velocity within the moving reference frame. A no-slip condition was set at all solid/liquid interfaces. The free surface remained flat during the experiments in the range of Re investigated, so a free-slip condition was used on the upper surface of the vessel. The conservation equations were solved in the laminar regime using a second order upwind scheme for the convection terms. The SIMPLE pressure-velocity coupling method was used with the gradients evaluated based on nodal values. The convergence of the axial velocity at three different point locations in the domain (one near the impeller, one halfway to the wall References start on page 123 Chapter 4 107 opposite to the impeller, and one far away from the impeller) was continuously monitored during the computations. Iterations were continued until the RMS scaled residuals for the three velocity components and the pressure fell below 1×10-5. At this point, the axial velocities monitored were examined to ensure that they had reached a steady value. If this was not the case, the iterations were continued until these velocities settled. The time needed to reach a converged solution in a single PC (Intel Core 2 CPU at 2.66 GHz and 4GB of RAM) varied from 4 hrs at low Re to 8−12 hrs. at the highest Re. An initial mesh (grid 1) of 133,767 volumes, with individual cells ranging from 1.33×10-10 m3 to 1.41×10-6 m3 was constructed with higher mesh density near the impeller and vessel walls. The representative average cell size of grid 1, h1 was h1 = 5.89 ×10-3 m. This mesh was refined systematically in a stepwise fashion. The next finer mesh (grid 2) had 346,134 volumes, h2 = 4.29 ×10-3 m and a global refinement factor, r12 = h1/h2, of 1.37. The finest mesh (grid 3) was comprised 735,121 volumes, with individual cells ranging from 3.5 ×10-11 m3 to 3.03x10-7 m3 and h3 = 3.33 ×10-3 m. The global grid refinement factor r23 was 1.29. Figure 4.5 shows the computed vertical profile of axial velocity in a plane aligned with the impeller tip (z =5cm), at x* = 0.29 and x* = 0.48 (d = 3 cm and d = 10 cm in Fig. 4.1(b)) for all grids. In order to provide an estimation of the discretization error, the grid convergence index (GCI) method described by Celik (2004) and summarized in appendix B, was implemented. The method was applied for the axial velocity at two points within the profiles shown in Fig. 4.5, located near the impeller center (y*=0.4); where the main differences are observed between coarse and fine grid solutions. The solution in the finest grid (grid 3) was selected for evaluation against the PIV measurements, thus we present the error estimates on this grid. Since the method References start on page 123 Chapter 4 108 was applied in non-uniform grids and for field variables, the local cell size was used in the GCI calculations. The calculated local order of accuracy (p) along with the approximate relative error (ea23) and extrapolated error (eext23) are summarized in table 4.1 (definitions of these error estimates are given in appendix B). The numerical uncertainty in the fine grid solution for the axial velocities near the impeller, given by the fine grid converge index (GCIfine23) was found to be 2.9% 3 cm away from the impeller and 2.3% 10 cm away from the impeller. (a) (b) Figure 4.5. Computed vertical profiles of axial velocity in a plane aligned with the blade tip (a) for a vertical line at x* = 0.29; d = 3cm (see Fig. 4.1(b)) and (b) for a vertical line at x* = 0.48; d = 10 cm (see Fig. 4.1(b)). N = 327 rpm, Re = 101. References start on page 123 Chapter 4 109 Table 4.1. Calculations of discretization error using the GCI method (refer to Appendix B for a detailed description). 4.4 U at x*= 0.08, U at x*= 0.28, y* = 0.4 y* = 0.4 r12 1.65 1.65 r23 1.36 1.36 p 2.5 2.3 ea23 2.6% 1.9% eext23 2.2% 1.9% GCIfine23 2.9% 2.3% Results and Discussion Instantaneous velocity fields obtained from single measurements were averaged over a large number of events to determine statistical values significant to the global flow phenomena prevailing in the system. Hence, the mean “time-averaged” velocity field taken over 1000 measurements was used to reveal the prevailing macro-scale flow structures inside the mixing vessel. Figure 4.6 gives the averaged velocity fields at the three different measurement planes (shown in Fig. 4.1(c)) for Re = 101. It should be noted that measurements at the mid plane exclude a portion of the image (impeller region and the area behind it) because the impeller blocks the path of the laser sheet (Fig. 4.6(a)). The length of the vectors and color coding is based on the planar velocity magnitude normalized with respect to the tip speed. References start on page 123 Chapter 4 110 (a) (b) (c) Figure 4.6. Mean planar velocity fields for Re =101 in the three vertical planes of measurement: (a) mid plane (z = 0 cm); (b) plane aligned with the impeller tip (z = -5.25 cm) and (c) plane located 2 cm away from the vessel wall (z = -10.25cm). Vectors are colored by velocity magnitude (m/s). References start on page 123 Chapter 4 111 The velocity fields measured at the mid-plane (Fig 4.6(a)) show four recirculation loops, two of them below the impeller region, which are created by the impingement of one of the jets, generated by the impeller, on the bottom of the vessel. The other two loops are located above the impeller, and they form as a result of the upward jet flow generated at the tip of the impeller blades at this location. The size and location of this measuring plane (Fig 4.6(a)) limits the visibility of the upper loops, and only their bottom half can be appreciated. These main recirculation patterns evidence the axial flow structure generated by the Maxflo Mark II impeller. For all the PIV planes studied, the velocities measured are only a fraction of the impeller tip speed (planar velocity magnitudes were below 35% of Utip) due in part to the significant resistance to flow of the glycerin solution (i.e. high viscosity). The highest planar velocity magnitudes measured are located in the plane aligned with the blade tip, but this magnitude diminishes quickly away from the impeller. At this plane (Fig.4.6 (b)), the jet generated by the impeller moving upwards towards the free surface is seen. The position and shape of the recirculation loops previously identified in the mid plane measurements, changed significantly from one plane to another, confirming the three dimensional nature of the flow field. The planar velocity field at the measurement plane located 2 cm away from the side wall (Fig. 4.6(c)) shows the evolution of the flow near the wall boundaries. The centrifugal force generated by the tangential flow of the fluid following the rotating impeller creates radial flow towards the vessel walls. Measurements at this plane show how this secondary flow becomes axial, especially near the bottom wall of the vessel, and then turns vertically inward near the impeller wall. The bifurcation of the flow captured at the lower end of the image (x*= 0.22 and y*<0.25) at the impeller axial location, provides an indication of the boundary between the two major References start on page 123 Chapter 4 112 recirculation loops below the impeller. The planar velocity magnitudes are much lower in this plane compared to the other two planes examined, and no vortices can be identified at this location. Figure 4.7 compares the planar velocity vectors measured at the center plane (z = 0 cm) for Re = 36, 82, and 101 with those obtained using CFD results calculated with the fine grid. The two-dimensional CFD velocity fields on the right side of Fig. 4.7 were created using the axial and vertical components of the computed velocity vectors. These values were then used to calculate the planar velocity magnitude. No interpolation was applied, so the vectors and the spacing between them were determined by the mesh spacing in the CFD simulations. A number of observations can be made using these 2D-velocity fields. Note that details in the immediate vicinity of the impeller are excluded from the PIV measurements so that comparison is based on the bulk flow patterns at common areas between the two velocity fields. Also, the relative position of the impeller blades in the CFD simulations is not exactly the same as the one in the PIV images. For instance, in the computational model of the system, the relative impeller position is such that the mid-plane (z = 0 cm) intersects the impeller blade above the hub, whereas in the PIV images the same plane intersects the blade below the hub. This however does not prevent the analysis and comparison between them, as it has been already demonstrated that phase averaging is not essential. References start on page 123 Chapter 4 113 (a) (b) (c) Figure 4.7. Experimental (left) and CFD (right) planar velocity fields in the center plane of measurement for (a) Re = 36, (b) Re = 82 and (c) Re = 101. Vectors are colored by velocity magnitude (m/s). References start on page 123 Chapter 4 114 For all Re studied, both computational and experimental results show four recirculation loops, two above and two below the impeller (Fig. 4.7). The stagnation flow on the upstream faces of the blades leading to vortex flow on the downstream faces of the impeller blades can be easily appreciated from the CFD velocity fields. Very good agreement is obtained for the locations of the centers of these recirculation regions and the distribution of the planar velocity magnitude along the measurement plane. The distance between the recirculation centers and the blades increased with Re. This trend was captured by both modeling and experiments. The enhanced strength of the recirculation loops at higher Re was also confirmed in the measured and in the computed planar flow patterns. The experimental and computational velocity fields can be quantitatively compared by plotting the profiles of the velocity components along sample lines (see dashed lines in Fig. 4.1(b)). Figure 4.8 shows the axial and vertical velocity profiles on the measurement plane aligned with the blade tip (z = -5.25 cm), along a horizontal line located at the impeller center, extending from x* = 0.05 to x* = 0.63 (the domain length of the PIV image). The velocity measurements show good repeatability and the associated error bars indicate the standard deviation for a series of three separate measurements. The CFD results show very good agreement with the experimental measurements. The local maximum and minimum normalized velocities are very similar to the measured ones. As Re increases, the maximum axial velocity increases from about 5% of the tip speed at Re = 36 to 8.2% of the tip speed at Re = 120. The computed local maximum for the vertical velocity is slightly lower than the measured one. A second local maxima for the vertical velocity at Re = 120 was computed in the CFD simulations, but not captured by the PIV measurements. This is likely due to the limited resolution of the PIV measurements in the vicinity of the impeller, where the velocity gradients are higher (Keane and Adrian, 1992). In References start on page 123 Chapter 4 115 fact, the same interrogation area size (5 × 5 mm2) was used in the entire image for the PIV measurements, whereas an increased mesh resolution was used in the vicinity of the impeller for the CFD simulations, with control volumes of characteristic length as low as 0.3 mm. Re = 36 Re = 36 Re = 120 Re = 120 Figure 4.8. CFD and PIV results for horizontal profiles of normalized axial and vertical velocity components for a horizontal line located at the impeller center (see Fig. 4.1(b)), and covering the full length of the PIV image, in a plane aligned with the blade tip (z = -5.25cm). Results are shown for N = 386 rpm , Re = 120 and N = 118 rpm, Re = 36. The error bars indicate the standard deviation for a series of three measurements. References start on page 123 Chapter 4 116 The evolution of the planar velocity along the length of the vessel can be examined by plotting the axial and vertical velocity profiles along vertical lines positioned at different axial distances from the impeller (i.e. varying d in Fig. 4.1(b)). Flow structures in the impeller region can not be resolved accurately because they are smaller than the IA length. The CFD results computed using a laminar viscous model are unable to capture velocity fluctuations near the impeller due to eddies formed as a result of turbulence. Thus comparisons between PIV and CFD are made at locations where the PIV and the CFD model resolutions are both adequate. Figure 4.9 shows the measured and computed vertical profiles of the axial and vertical velocity in the measurement plane aligned with the blade tip, at a line located at x* = 0.34, (d =5 cm in Fig. 4.1(b)) for Re = 36 and 82. The horizontal axis represents the vertical distance along the line and the vertical axis represents the normalized velocity components. In all the cases, very good agreement is attained between the computed and measured velocities. The calculated RMS deviations (measured vs. computed) were 13.8 % for U* and 15% for V* at Re = 36, and 15% for U* and 10.9% for V* at Re = 101. In the same way, the vertical velocity profiles of the two velocity components at a line located at x*=0.48 (d =10 cm in Fig. 4.1(b)) were examined and are shown in Figure 4.10. Once again, very good agreement exists between experimental and numerical results. The calculated RMS deviations were 12.4 % for U* and 16% for V* at Re = 36, and 7.2% for U* and 20.1% for V* at Re = 82. The differences rendering this high RMS deviation for V* at Re =82, were found above the impeller (y* ≥ 0.39), suggesting that at this Re the deformation of the free surface (not considered in CFD) has some influence on the flow. No changes were observed on the surface during the experiments; hence the assumption of a flat upper surface was justified. Nevertheless, References start on page 123 Chapter 4 117 it is likely that a weak interaction with the surface exists which can not be easily observed (i.e. the free surface may rise a few millimeters). Re =36 Re =36 Re =82 Re =82 Figure 4.9. CFD and PIV results for vertical profiles of normalized axial and vertical velocity components for a vertical line at x*= 0.34; d = 5cm (see Fig 4.1(b)), spanning the PIV image height, in a plane aligned with the blade tip (z = -5.25 cm). Results are shown for N = 118 rpm, Re = 36 and N = 267 rpm, Re = 82. The error bars indicate the standard deviation for a series of three measurements. References start on page 123 Chapter 4 Re =36 Re =82 118 Re =36 Re =82 Figure 4.10. CFD and PIV results for vertical profiles of normalized axial and vertical velocity components for a vertical line at x* = 0.48; d = 10 cm (see Fig. 4.1(b)), spanning the PIV image height, in a plane aligned with the blade tip (z =-5.25 cm). Results are shown for N = 118 rpm, Re = 36 and N = 267 rpm, Re = 82. The error bars indicate the standard deviation for a series of three measurements. Similar comparisons between experimental and numerical results were made at the other two measuring planes and are shown in Figure 4.11 for Re = 82. The calculated RMS deviations were 13 % for U* at x* = 0.34 (d = 5cm in Fig. 4.1(b)) at the mid-plane of measurement and 17.2% for U* at x* = 0.48 (d = 10 in Fig. 4.1(b)) at the plane of measurement located 2 cm away from the vessel wall. References start on page 123 Chapter 4 (a) 119 (b) Figure 4.11. CFD and PIV results for vertical profiles of normalized axial and vertical velocity components at Re = 82; (a) for line at x* = 0.34; d = 5cm (see Fig. 4.1(b)) in the mid plane of measurement (z = 0 cm), and (b) for line at x* = 0.48; d = 10 cm (see Fig. 4.1(b)) in the measurement plane located 2 cm away from the vessel wall (z = -10.25 cm). The error bars indicate the standard deviation for a series of three measurements. Global parameters like the power number (NP) are conventionally used to characterize the mixing process. To complement the CFD evaluation procedure, the torque was measured and computed as the axial moment on the impeller center at each Re. The power draw and NP was then calculated and plotted against Re. Figure 4.12 shows a comparison of the NP vs. Re curves. The agreement between these curves is very good over the entire range of Re, with NP being under predicted by an average of 11%. The linear decay of NP with increasing Re confirms the laminar flow assumption implemented in the CFD model. References start on page 123 Chapter 4 120 Figure 4.12. Comparison of NP for experimental and computational tests. The error bars are based on the accuracy with which the torque was measured (± 0.002 Nm). 4.5 Summary and Conclusions The flow in a rectangular agitated vessel equipped with a side entering axial impeller was investigated numerically and experimentally over a range of Re between 20 and 120. Experiments were conducted using PIV to characterize the flow field from time averaged 2D velocity maps at three different locations relative to the impeller. The PIV measurement parameters were tuned for optimum capture of the velocity fields in a representative portion of the vessel. The planar velocity fields revealed characteristic flow structures underlying the global recirculation flow patterns occurring in the vessel. The system was modeled using CFD software (Fluent) and the PIV measurements were used to evaluate the model predictions. Special attention was directed towards building efficient and size controlled mesh structures, enabling fast and reliable computations that could be performed in one single conventional PC. The CFD predicted velocity vectors at the measurement planes agreed very well with the PIV data for global flow patterns and planar velocity magnitude. The predicted velocity profiles were examined at different positions within the measurement planes for various agitation rates in order References start on page 123 Chapter 4 121 to present an extensive evaluation of the CFD model. The experimental and computational results showed very good agreement in all cases. The NP vs. Re curve for the system was constructed from torque and rotational speed measurements and compared with the computational predictions resulting in a very adequate agreement. The accurate definition of the fluid’s viscosity in the CFD simulations was crucial to reach satisfactory results. The overall results from this study indicate that CFD can effectively be used to estimate multiple parameters relevant for different mixing applications. For example, the spatial heterogeneity within the mixer can be determined with confidence by means of computation of the rate of deformation tensor and strain rate. Following sufficient evaluation of the simulation results, CFD models can be applied with increased confidence on the re-design and optimization of similar mixing systems. 4.6 Nomenclature D impeller diameter, cm L mixing vessel lenght, cm Htank mixing vessel height, cm H liquid level in the mixing vessel, cm d axial distance in front of the impeller, cm z radial position from the impeller center, cm N impeller rotational speed, rpm Re Reynolds number NP power number LIA length of the interrogation area side References start on page 123 Chapter 4 122 S object to image scale factor zA laser sheet thickness, mm Utip tip speed, m/s Uz velocity perpendicular to the measurement plane, m/s U* normalized axial velocity V* normalized vertical velocity x* normalized axial position y* normalized vertical position p apparent local order of accuracy hi characteristic cell size of grid i r21 grid refinement factor from grid 1 to grid 2 r32 grid refinement factor from grid 2 to grid 3 ea 23 approximate relative error between the two finest grids (grids 2 and 3) eext 23 extrapolated relative error in the fine grid (grid 3) GCI fine 23 fine grid converge index Greek letters µ fluid viscosity, kg/m.s ρ fluid density, kg/m3 θ angle that defines the relative position of the impeller blades within one rotation, ° ∆t time delay between laser flashes, s References start on page 123 Chapter 4 4.7 123 References Adams, L.W. and Barigou, M., 2007. CFD analyisis of caverns and pseudo-caverns developed during mixing of non-newtonian fluids. Che. Eng. Res. Des., 85, A5, 598−604. Adrian, R.J., 1986. An image shifting technique to resolve directional ambiguity of double pulsed laser velocimetry. Applied Optics, 23, 3855−3858. Adrian, R.J., 1997. Dynamic ranges of velocity and spatial resolution of particle image velocimetry. Meas. Sci. Technol., 8, 1393−1398. Bakker, A., Fasano J.B. , 1993. A computational study of the flow pattern in an industrial paper pulp stock chest with a side-entering impeller. AIChE Symposium Series, 89 293, 118−124. Bhattacharya, S., Gomez, C., Soltanzadeh, A., Taghipour, F., Bennington, C.P.J. and Dumont, G.A., 2008. Computational modelling of industrial pulp stock chests. Can. J. Chem. Eng., submitted for publication. Brown, D., Pip, J., Middlenton., J., Papadopoulos, G. and Arik, E., 2004. Experimental methods, in: Paul, E.L., Atiemo-Oberg V.A. and Kresta S.M. (Eds.), Handbook of industrial mixing. Science and practice, John Wiley & Sons, New Jersey, pp. 145−174. Bugay, S., Escudié, R. and Liné, A., 2002. Experimental analysis of hydrodynamics in axially agitated tank. AIChE J., 48, 3, 463−475. Bujalski, J.M., Yang, W., Nikolov, J., Solnordal, C.B. and Schwarz, M.P., 2006. Measurement and CFD simulation of single-phase flow in solvent extraction pulsed column. Chem. Eng. Sci., 61, 2930−2938. Celik, I. B., 2004. Procedure for estimation and reporting of discretization error in CFD applications. Fluids Eng. Division of ASME. References start on page 123 Chapter 4 124 Collignon, M.L., Crine, M., Verdin, E., Chaubard, J.F., Peeters, L., Dessoy, S. and Toye, D., 2008. A study of the mixing by PIV and PLIF in bioreactor of cells animals culture. ISMIP VI, Niagara Falls, Canada. Couerbe, G., Fletcher, D.F., Xuereb, C. and Poux, M., 2008. Impact of thixotrophy on flow patterns induced in a stirred tank: Numerical and experimental studies. Che. Eng. Res. Des., 86, 545−553. Dantec dynamics Inc., 2003. Dantec dynamics flow map PIV installation and user's guide: General PIV reference. Escudié, R., Bouyer, D. and Liné, A., 2004. Characterization of trailing vortices generated by a rhuston turbine. AIChE J., 50, 1, 75−86. Escudié, R. and Liné, A., 2003. Experimental analysis of hydrodynamics in a radiallyagitated tank. AIChE J., 49, 3, 585−603. Ford, C., Bennington, C.P.J. and Taghipour, F., 2007. Modelling a pilot-scale mixing pulp mixing chest using CFD. J. Pulp Pap. Sci., 33, 115−120. Ford, C., Ein-Mozaffari, F., Bennington, C.P.J. and Taghipour, F., 2006. Simulation of mixing dynamics in agitated pulp stock chests using CFD. AIChE J., 52, 10, 3562−3569. Gómez, C., Derakhshandeh, B., Hatzikiriakos S.G.and Bennington, C.P.J., 2008. Carbopol as a model fluid for studying mixing of pulp fibre suspensions. Chem. Eng. Sci., submitted for publication. Jaworski, Z., Nienow, A.W. and Dyster, K.N., 1996. An LDA study of the turbulent-flow field in a baffled vessel agitated by an axial, down-pumping hydrofoil impeller. Can. J. Chem. Eng., 74, 13−15. References start on page 123 Chapter 4 125 Keane, R.D. and R.J., A., 1992. Theory of cross-correlation analysis of PIV images. Appl. Sci. Res., 49, 191−215. Khopar, A., Aubin, J., Ruibo-Atoche, C., Xuereb, C., Sauze, N.L., Bertrand, J. and V, V., 2004a. Flow generated by radial flow impellers: PIV measurements and CFD simulations. Int. J. Chem. Reactor Eng., 2, A18, 1−16. Kilander, J. and Rasmuson, A., 2005. Energy dissipation and macro instabilities in a stirred square tank investigated using an LE-PIV approach and LDA measurements. Chem. Eng. Sci., 60, 6844−6856. Lee, K.C. and Yianneskis, M., 1998. Turbulence properties of the impeller stream of a Ruhston turbine. AIChE J., 44, 1, 13. Mavros, P., 2001. Flow visualization in stirred vessels. A review of experimental techniques. Trans. IChemE., 79, part A, 113−127. Mavros, P., Xuereb, C. and Bertrand, J., 1996. Determination of 3-D flow fields in agitated vessels by laser doppler velocimetry. Effect of impeller type and liquid viscosity on the liquid flow patterns. Chem. Eng. Res. Des, 74A, 658−668. Myers, K.J., Ward, R.W. and Bakker, A., 1997. A digital particle image velocimetry investigation of flow field instabiliteis of axial flow impellers. J. Fluids Eng., 3, 623−632. Oldshue, J.Y., 1983. Fluid mixing technology, McGraw-Hill Publications Co., New York, 1983. Raffel, M., Willert, C. and Kompehans, J., 1998. Particle image velocimetry. A practical guide, Adrian, R.J., et al. (Eds), Springer, Berlin, 1998. Ranade, V.V., Perrard, M., Sauze, N.L., Xuereb, C. and Bertrand, J., 2001. Trailing vortices of rushton turbine: PIV measurements and CFD simulations with snap shot approach. Trans. IChemE., 79, 3−11. References start on page 123 Chapter 4 126 Schäfer, M., Yanneskis, M., Wäxhter, P. and Durst, F., 1998. Trailing vortices around a 45o pitched-blade impeller. AIChE J., 44, 6, 1233−1246. Sharp, K.V. and Adrian, R.J., 2001. PIV study of small scale flow structure around a Ruhston turbine. AIChE J., 47, 766−778. Stoots, C. and Calebrese, R.V., 1995. Mean velocity field relative to a rhuston turbine blade. AIChE J., 41, 1−11. Svensson, F.J.E. and Rasmuson, A., 2006. PIV measurements in a liquid-liquid system at volume percentages up to 10% dipersed phase. Exp. Fluids, 41, 917−931. Torré, J.P., Fletcher, D.F., Lasuye, T. and Xuereb, C., 2007. Single and multiphase CFD approaches for modeling partially baffled stirred vessels: Comparsion of experimental data with numerical predictions. Chem. Eng. Sci., 62, 6246−6262. Weetman R.J. and Oldshue J.Y., Power, flow and shear characteristics of mixing impellers. Proc. Eur. Conf. on Mixing, Pavia, Italy, May 24–26 (1988). Wu, H. and Patterson, G.K., 1989. Laser doppler measurements of the turbulent flow in stirred vessel to establish scaling rules. Chem. Eng. Sci., 44, 2207−2221. Yoon, H.S., Sharp, K.V., Hill, D.F., Adrian, R.J., Balachandar, S., Ha, M.Y. and Kar, K., 2001. Integrated experimental and computational approach to simulation of flow in a stirred tank. Chem. Eng. Sci., 56, 6635−6649. Zalc, J.M., M., A.M. and Muzzio, F.J., 2001. Extensive validation of computed laminar flow in a stirred tank with three rusthon turbines. AIChE J., 47, 10, 2144−2154. References start on page 123 127 Chapter 5 Characterization of Axial Flow Impellers in Pulp Fibre Suspensions 5.1 Introduction Mixing is an essential operation in pulp and paper manufacture. Pulp fibre suspensions are nonNewtonian and possess a yield stress. They must be uniformly blended with other pulp streams and with additives, fillers and chemicals prior to papermaking. The required mixing is commonly accomplished in mechanically stirred vessels (stock chests) using axial flow impellers designed for high flow, such as marine impellers, the Maxflo impeller (Chemineer Inc., Dayton, USA) and the A312 impeller (Lightnin Inc., Rochester, USA). Stock chest design involves selecting the correct impeller and power input required to ensure complete motion through a specified chest volume. Due to suspension rheology this can be a challenge even at low suspension mass concentrations (Cm = 2 to 4%). Chest design is based on proprietary criteria but typically follows the method summarized by Yackel (1990) which matches the fluid momentum generated by an impeller with that needed to provide complete motion over the suspension surface in the chest. The later is determined from correlations developed for a wide range of variables, including chest geometry (typically rectangular), fibre type and suspension concentration. The impeller momentum number, Mo, is given by Mo = CN 2 D 4 (5.1) where C is a constant that depends on impeller geometry and N and D are the impeller speed and diameter, respectively. C is usually available from the impeller manufacturer, often based on A version of this chapter has been published. Bhole M., Ford† C. and Bennington C.P.J., 2009. Characterization of axial flow impellers in pulp fibre suspensions. Chem Eng Res Des, 87, 4, 648−653. † C. Gomez formerly C.Ford Chapter 5 128 tests made in water. Yackel (1990) tabulates Mo values for several impellers over a range of operating conditions for use with his design procedure. Axial flow impellers are most often characterized by their power number, NP, and more recently by their axial force number, Nf NP = P FA where P = 2πNM and N f = 3 5 ρN D ρN 2 D 4 (5.2) with P being the power drawn, M the shaft torque and FA the axial thrust. NP for mixing of Newtonian fluids in standard vessels is available in the literature as a function of Reynolds number. For non-Newtonian fluids, the methodology developed by Metzner and Otto can be used to determine the appropriate Reynolds number for these correlations (Bakker and Gates, 1995; Metzner and Otto, 1957; Shekhar and Jayanti, 2003). Nf has not been widely adopted for mixer design although Yackel’s method essentially employs it. However, Yackel’s procedure has been shown to underestimate the power required for effective mixing of pulp suspensions (Ein-Mozaffari et al., 2004). Indeed, even when complete surface motion is attained, stagnant zones exist below the suspension surface which reduces mixing efficiency. The mixing of yield stress fluids results in formation of a cavern around the impeller where fluid motion is localized. Outside the cavern, fluid stresses fall below the yield stress and the fluid is stagnant (Solomon et al., 1981; Wichterle and Wein, 1981). Cavern size can be estimated by balancing the force transferred from the impeller to the cavern boundary once a cavern shape has been specified. A number of researchers have done this (Amanullah et al., 1998; Elson et al., 1986; Solomon et al., 1981), and Amanullah et al. included both NP and Nf in their equation for References start on page 140 Chapter 5 129 determining the surface area, As, of an unbounded spherical cavern created by an axial flow impeller, As = Re y D 2 Nf 2 ⎛ 4N P ⎞ +⎜ ⎟ ⎝ 3π ⎠ 2 (5.3) where the yield stress Reynolds number, Rey, is given by ρN 2 D 2 Re y = τy (5.4) with τy the suspension yield stress. When the parameters (Nf and NP) chosen for Eqn. 5.3 correspond to the actual mixer operating conditions (rather than values determined under turbulent conditions) the predictions are more accurate (Amanullah et al., 1998 ). Since the rheological behavior of pulps suspensions is not well characterized the measurement of these parameters becomes imperative. We have done this for two common axial flow impellers under typical pulp suspension mixing conditions. Computational fluid dynamics (CFD) has emerged as a powerful tool for characterizing impeller behaviour. While most simulations have focused on mixing single-phase Newtonian fluids under turbulent conditions, some more recent approaches have modeled agitation in non-Newtonian fluids (Adams and Barigou, 2007; Arratia et al., 2006; Ford et al., 2006; Kelly and Gigas, 2003). Although these simulations provide information on velocity profiles and permit identification of regions of poor mixing, their accuracy depends significantly on the characterization of the fluid rheology (Ford et al., 2007). We have used a commercial CFD package to determine impeller parameters under relevant pulp mixing conditions and compared them with our experimental measurements to determine the accuracy of this approach. References start on page 140 Chapter 5 5.2 130 Experimental The primary vessel used for impeller characterization (T1) is a fully baffled (four baffles of width B = 2.5 cm), slightly tapered, cylindrical vessel with diameter varying linearly from 29 cm at the top of the vessel to 26 cm at the bottom as shown in Fig. 5.1. Suspension height was maintained at Z = 29 cm. Impellers were mounted on an instrumented shaft driven by a 0.25 kW variable speed motor with the impeller off-bottom clearance, C1, maintained at 10 cm. Two impeller types were studied: a standard marine impeller (3-bladed, D = 9.0 cm, pitch ratio = 1.5) and a series of geometrically scaled Maxflo impellers (Chemineer Inc., Dayton, OH) (3 bladed, D = 10.2, 14.0 and 16.5 cm; pitch ratio = 0.44). Photographs of the impellers are included in Fig. 5.1. Impeller rotational speed was measured to ±1 rpm using an optical tachometer and shaft torque was measured to ±0.002 N⋅m using an in-line strain gauge. The axial force was measured gravimetrically using a balance (Mettler Toledo SB32001) accurate to ±0.1g although variations in axial force due to shaft vibrations at high rotational speeds could be as high as ±20g. Tests were made in the down-pumping mode with impeller speed varied from 50 to 800 rpm. T2 B Z C1 D T1 Figure 5.1. Mixing vessel (T1) and impellers studied (Maxflo on the left and the marine impeller on the right) References start on page 140 Chapter 5 131 Two industrial pulp fibre suspensions were studied. A softwood bleached kraft pulp (SW) with a length-weighted fibre length, lf, of 2.96 mm and a hardwood bleached kraft (HW) with lf = 1.28 mm were obtained from Domtar Inc. (Windsor, QC). Suspension mass concentration was varied from Cm = 0 to 4% with the suspension yield stress measured using a Haake RV12 viscometer giving the following correlations τ y = 12.3C m 2.05 for the softwood pulp (5.5) τ y = 1.25C m 3.00 for the hardwood pulp (5.6) with the coefficient of variation in yield stress averaging ±22% over the mass concentration range, in line with past measurements (Bennington et al., 1990). 5.3 5.3.1 Results and Discussion Impeller Characterization NP and Nf for agitation of softwood pulp suspensions with a Maxflo impeller (D = 10.2 cm) and the marine impeller are shown in Figs 5.2 and 5.3 (plotted against impeller rotational speed, N). For both impellers, the dependence of Nf and NP on N indicate that flow is in the laminar and transition-to-turbulence regimes. For a fixed rotational speed, both Nf and NP increase with suspension consistency indicting a significant influence of suspension rheology on impeller performance. For a given impeller speed and suspension consistency, the Nf and NP for the marine propeller are higher than for the Maxflo impeller due to its higher pitch ratio. References start on page 140 Chapter 5 132 100 100 10 10 1 Nf (-) NP (-) 1000 Cm (%) 0.0 1.0 2.0 3.0 4.0 1 0.1 10 0.1 100 Cm (%) 0.0 1.0 2.0 3.0 4.0 0.01 10 1000 100 1000 N (rpm) N (rpm) 1000 100 100 10 10 1 1 0.1 10 Nf (-) NP (-) Figure 5.2. NP and Nf for Maxflo impeller (D = 10.2 cm) in softwood pulp suspensions Cm (%) 0.0 1.0 2.0 3.0 4.0 0.1 100 N (rpm) 1000 0.01 10 Cm (%) 0.0 1.0 2.0 3.0 4.0 100 1000 N (rpm) Figure 5.3. NP and Nf for marine impeller (D = 9.0 cm) in softwood pulp suspensions NP and Nf were measured in both hardwood and softwood pulp suspensions from Cm = 1 to 4% and are plotted against Rey in Fig. 5.4. For Nf all curves collapse to a single operating line. However, the same is not true for NP where suspension properties other than the suspension yield stress affect the power number. Similar results are reported for pulp suspensions by EinMozaffari et al. (2003) and for Bingham plastic fluids by Nagata (1975). References start on page 140 Chapter 5 133 100 100 Cm(%) 1.0 2.0 3.0 4.0 SW Cm(%) 1.0 2.0 3.0 4.0 10 Nf (-) NP (-) 10 HW HW SW 1 1 0.1 0.1 0.01 0.1 1 10 100 1000 10000 0.01 0.01 Rey (-) 0.1 1 10 100 1000 10000 Rey (-) Figure 5.4. NP and Nf for Maxflo impeller (D = 10.2 cm) in hardwood and softwood pulp suspensions The effect of D/T ratio was studied by measuring NP and Nf as a function of Rey for various D/T ratios in a Cm = 3% hardwood pulp suspension using the Maxflo impeller. Vessel T1 was used for most tests in addition to a larger cylindrical baffled vessel (T2) having T = 57 cm which allowed D/T to be varied from 0.18 to 0.60. Figure 5.5 shows that NP is almost independent of D/T over the range considered. In Newtonian fluids NP is a function of D/T for axial flow impellers because the flow generated interacts with the vessel walls (and changes as D/T changes (Chapple et al., 2002). In pulp suspensions, the cavern created around the impeller often does not reach the vessel walls and thus no significant effect of D/T is observed. Nf does depend on D/T and further study is needed to quantify the effect. So far we have presented results for cylindrical baffled tanks. Industrial stock chests are usually rectangular due to space constraints and the ease of construction. Consequently some tests were made in a square tank (S1) with T = Z = 24 cm. NP and Nf for the D = 10.2 cm Maxflo impeller in a Cm = 3 % hardwood pulp suspension did not vary significantly (less than a 10% change) References start on page 140 Chapter 5 134 from the values measured in the cylindrical vessel (T1). Further details on these additional tests can be found in Appendix C. 10 100 D/T 0.18 0.39 0.54 0.63 D/T 0.39 0.54 0.63 T (m) 0.57 0.26 0.26 0.26 Nf (-) NP (-) 10 D (m) 0.102 0.102 0.140 0.165 D (m) 0.102 0.140 0.165 T (m) 0.26 0.26 0.26 1 1 0.1 0.1 1 10 100 1000 0.1 0.1 Rey (-) 1 10 100 1000 Rey (-) Figure 5.5. NP and Nf for Maxflo impellers as a function of D/T in hardwood pulp suspensions 5.3.2 Impeller Characterization using CFD Characterization of impellers for different fluid/suspensions can be expensive and time consuming. In this regard, using CFD to model impeller performance is promising. We used the commercial CFD software package, Fluent V6.3 (Fluent Inc., Lebanon, NH) to simulate the experimental tests performed using the D = 10.2 cm Maxflo impeller in a Cm = 3% hardwood suspension in vessel T1. The computational domain within the vessels was discretized with multi-block mesh structures using Gambit software (v. 2.4.6, Fluent, Lebanon, NH) with cell density optimized to capture flow details without being excessive. Increased mesh resolution was used near all solid boundaries with increases in mesh size limited to a growth rate of 20% (Fluent, 2006). References start on page 140 Chapter 5 135 A modified Herschel-Bulkley model was used to describe suspension rheology (refer to chapter 2, Eqns 2.1a & 2.1b), with the consistency index, k, set to 0.001 Pa·s (the viscosity of water) and the power-law index, n, set to 1.0. The yield stress of the 3 % hardwood suspension was set to its measured value of τy = 30 Pa and the yielding viscosity to µo = 100 Pa·s (Bennington et al., 1990). The final three-dimensional mesh of the cylindrical vessel (T1) had 463,742 computational cells (further details given in appendix C; Fig. C 1). The multiple reference frame approach was used for computation, with coupling between reference frames made using a velocity transformation. A no-slip condition was set at all solid/suspension boundaries and a free slip condition was set at the air/suspension interface. Mass and momentum conservation equations were solved in the laminar regime using a second order upwind scheme to limit numerical diffusion. For all the conditions simulated, the nonNewtonian Reynolds numbers defined by Attwood and Gibbon (1963) and Blasinski and Rzyski (1972) were below 103 indicating laminar or transitional flow, in agreement with the measured behavior (Figs. 5.2 and 5.3). Iterations were continued until the root mean square scaled residuals for each transport equation fell below 1×10-6, with convergence of axial velocity at three different point locations (one near the impeller and two remote from the impeller) monitored to verify that a steady value was reached. Grid independence was verified by demonstrating that further cell refinement near the impeller did not change the computed torque at the impeller center by more than 2% and that the magnitude of velocity calculated at two lines located in the impeller discharge did not change by more than 3% (illustrated in Fig. C.2) References start on page 140 Chapter 5 136 The calculated velocity field is shown in Fig. 5.6 at the vessel mid-plane for geometry T1 at three different impeller speeds. Only vectors having magnitudes below 0.05 m/s are plotted. Regions having velocity >0.05 m/s indicate rapid motion in the impeller (cavern) region and are left blank for clarity. The regions in blue (low velocity) are essentially stagnant. As impeller speed increases, the size of the cavern increases. Note that even at high impeller speeds the cavern does not significantly interact with the vessel walls. Additional figures from the CFD solution are provided in appendix C (Figs C.3 to C.6). Figure 5.6. Velocity vectors colored by velocity magnitude (m/s) at plane z = 0. Hardwood pulp at Cm = 3% Measured NP vs. N and Nf vs. N curves are compared with those obtained from the CFD simulations in Fig.5 7 for a Cm = 3% hardwood suspension in cylindrical vessel T1. Computed values of the shaft torque and the axial thrust acting on the impeller were used in Eqn. 5.2 to calculate NP and Nf. The agreement between experimental and computed curves is generally very good over the entire impeller speed range. NP is under-predicted by 24% except at the lower impeller speeds (13% deviation). Nf is under-predicted by 16%. Computations made for References start on page 140 Chapter 5 137 the square vessel (S1) show the same behaviour. The computed values of Nf were very similar for both geometries (within 4%), in agreement with the experimental results 10 Experimental CFD Experiment CFD Nf (-) NP (-) 10 1 0.1 100 1 0.1 1000 N (rpm) Figure 5.7. Comparison of NP and Nf 100 1000 N (rpm) for experimental and computational tests for cylindrical vessel T1. D = 10.2 cm Maxflo impellers in a Cm = 3% hardwood pulp suspensions as a function of rotational speed The fluid forces and torque are calculated using the pressure and velocity fields in the impeller which depend on the boundary conditions and material properties specified for the simulation, particularly the suspension viscosity. In the impeller vicinity, shear rates are high and the viscosity of the suspension is calculated with Eqn. 5.8 using parameters τy , µo, k and n. The calculated viscosity is sensitive to τy which is only measured to an accuracy of ±17% for the Cm = 3% hardwood suspension, equivalent to ±5 Pa. Consequently, the estimated viscosity near the impeller is only accurate to ±17% when all other parameters are known with total confidence (which they aren’t). When this uncertainty is considered, the deviation between the calculated and measured impeller parameters can be explained. References start on page 140 Chapter 5 5.4 138 Summary and Conclusions Two axial flow impellers used for the agitation of pulp fibre suspensions were characterized in hardwood and softwood pulp suspensions at low mass concentrations. The impellers operated in the laminar and transition-to-turbulence regimes in the suspensions, with the suspension mass concentration significantly affecting both power and axial thrust numbers. The axial force numbers could be collapsed to a single operating curve using the yield stress Reynolds number, although the power number remained a function of suspension properties. Computational fluid dynamics (Fluent) was used to model impeller flow using a Bingham approximation to describe the suspension rheology. The agreement of the CFD computations was good, with the values of NP and Nf determined to within 24% and 16% of the experimentally measured values. The discrepancy was attributed primarily to the uncertainty in measuring the suspension yield stress, which was to ±17% for the hardwood suspension at Cm = 3%. 5.5 Nomenclature Cm mass concentration of pulp suspension (consistency), % Mo impeller momentum number, m4/s2 C impeller constant, dimensionless N impeller rotational speed, rpm D impeller diameter, cm NP impeller power number, dimensionless Nf impeller axial force number, dimensionless FA axial force, N As cavern surface area, cm2 References start on page 140 Chapter 5 Rey yield stress Reynolds number, dimensionless T vessel diameter, cm B baffle width, cm Z suspension height, cm C1 impeller off-bottom clearance, cm k consistency index, Pas lf fibre length, mm 139 Greek letters τy yield stress, Pa τ shear stress, Pa µo yielding viscosity, Pas . γ shear rate, s-1 References start on page 140 Chapter 5 5.6 140 References Adams, L.W. and Barigou, M., 2007. CFD analysis of caverns and pseudo-caverns developed during mixing of non-newtonian fluids. Che. Eng. Res. Des., 85, A5, 598−604. Amanullah, A., Hjorth, S.A. and Nienow, A.W., 1998. A new mathematical model to predict cavern diameters in highly shear thinning, power law liquids using axial flow impellers. Chem. Eng. Sci., 53, 3, 455−469. Arratia, P.E., Kukura, J., Lacombe, J. and Muzzio, F.J., 2006. Mixing of shear thinning fluids with yield stress in stirred tanks. AIChE J., 52, 2310−2322. Attwood, D. and Gibbon, J.D., 1963. The agitation of paper stock. Pap. Technol., 4, 54−62. Bakker, A. and Gates, L.E., 1995. Properly choose mechanical agitators for viscous liquids. Chem. Eng. Prog., 91, 25−34. Bennington, C.P.J., Kerekes, R.J. and Grace, J.R., 1990. The yield stress of fibre suspensions. Can. J. Chem. Eng., 68, 10, 748−757. Blasinski, H. and Rzyski, E., 1972. Mixing of non-newtonian liquids. Power consumption for fibrous suspensions. Inz. Chem, 2, 1, 169−182. Chapple, D., Kresta, S.M., Wall, A. and Afacan, A., 2002. The effect of impeller and tank geometry on power number for a pitched blade turbine. Chem. Eng. Res. Des, 80, 364−372. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2003. Performance and design of agitated pulp stock chests. Appita J., 56, 2, 127−133. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2004. Dynamic mixing in agitated industrial pulp chests. Pulp Pap. Can., 105, 5, T 118− T 122. References start on page 140 Chapter 5 141 Elson, T.P., Cheesman, D.J. and Nienow, A.W., 1986. X-ray studies of cavern sizes and mixing performance with fluids possessing a yield stress. Chem. Eng. Sci., 41, 2555−2562. Fluent, 2006. V6.3 user's guide documentation, Lebanon, NH. Ford, C., Bennington, C.P.J. and Taghipour, F., 2007. Modelling a pilot-scale mixing pulp mixing chest using CFD. J. Pulp Pap. Sci., 33, 115−120. Ford, C., Ein-Mozaffari, F., Bennington, C.P.J. and Taghipour, F., 2006. Simulation of mixing dynamics in agitated pulp stock chests using CFD. AIChE J., 52, 10, 3562−3569. Kelly, W. and Gigas, B., 2003. Using CFD to predict the behavior of power law fluids near axial-flow impellers operating in the transitional regime. Chem. Eng. Sci, 58, 2141−2152. Metzner, A.B. and Otto, R.E., 1957. Agitation of non-Newtonian fluids. AIChE J., 3, 1, 3−10. Nagata, S., 1975. Mixing: Principles and applications, Wiley, New York. Shekhar, S.M. and Jayanti, S., 2003. Mixing of power law fluids using anchors: Metzner-Otto concept revisited. AIChE J., 49, 30−40. Solomon, J., Elson, T.P., Nienow, A.W. and Pace, G.W., 1981. Cavern size in agitated fluids with a yield stress. Chem. Eng. Commun., 11, 143−164. Wichterle, K. and Wein, O., 1981. Threshold of mixing of non-Newtonian liquids. Int. Chem. Eng., 21, 116−120. Yackel, D.C., 1990. Pulp and paper agitation: The history, mechanics and process, TAPPI press, Atlanta. References start on page 140 142 Chapter 6 Carbopol as a Model Fluid for Studying Mixing of Pulp Fibre Suspensions 6.1 Introduction Many fluids encountered in industrial mixing processes display complex rheological behavior. Special challenges are encountered when mixing yield stress fluids (as they produce caverns around impellers) or when rheological properties change during processing. Pulp fibre suspensions possess yield strength, and suspension mass concentration changes throughout pulp manufacture. To help understand suspension mixing behaviour, computational fluid dynamics (CFD) has been increasingly used, along with experimental techniques. However, difficulty determining suspension rheology increases the difficulty of CFD implementation and the uncertainty in their predictions. A significant challenge is to obtain sufficient and accurate experimental data from the opaque pulp suspensions to validate computational results. CFD has been used to study agitation of pulp suspensions in stock chests with good success (Bakker and Fasano, 1993; Ford et al., 2006; Ford et al., 2007; Bhattacharya et al., 2008; Bhole et al., 2008). These models only perform as well as the approximations included in their formulation, including those made in the choice of rheological model used to describe the pulp rheology, generally a modified Herschel-Bulkley model. Verifications of simulations have been conducted using torque and force measurements (Wikstrom et al., 1998; Ford et al., 2007; Bhole et al., 2008), visual observations of surface flows (Roberg, 1997; Bakker and Fasano, 1993) and A version of this chapter has been submitted for publication. C. Gomez†, B. Derakhshandeh, S.G. Hatzikiriakos, and C.P.J. Bennington., 2009. Carbopol as a model fluid for studying mixing of pulp fibre suspensions. † C. Gomez formerly C.Ford Chapter 6 143 assessment of the dynamic response of laboratory and industrial chests (Ford et al., 2006; 2007; Bhattacharya et al., 2008). While these comparisons provide confirmatory information on mixing occurring in the chest, they do not provide sufficient detail to fully validate the CFD simulations. Ideally, the simulated flow fields should be compared with velocity measurements made at specific locations in the vessel and cover the full spectrum of shear rates present. To obtain these data, non-intrusive flow measurement techniques, such as laser Doppler anemometry (LDA), particle image velocimetry (PIV), particle tracking velocimetry (PTV) and laser speckle velocimetry (LSV) have been used. However, many industrially important fluids, including paper pulp suspensions, are opaque so that optically-based techniques cannot be implemented. Consequently transparent simulant fluids that exhibit rheology similar to the industrial fluid must be used, usually in reduce-scaled vessels. One benefit of using a simulant fluid is that the rheology can often be characterized more precisely than that of the industrial fluid, increasing the accuracy of the computational predictions. An appropriate simulant fluid for studying pulp suspension mixing must have the correct rheological properties for the scale at which the measurements are to be performed, and may not fully replicate the suspension properties (Brown et al., 2004). After geometrical similarity has been established between the industrial and laboratory-scale mixers, scaling criteria must be compared to determine the operating conditions at the laboratory-scale that will mimic the mixing occurring at the process scale. Implementation of a particular scaling criterion (i.e. equal power input per unit volume, equal Reynolds number, etc.) yields a system of equations that can be solved to determine the appropriate operating conditions once a model fluid has been chosen. References start on page 168 Chapter 6 144 In this work we investigate the use of aqueous carbopol solutions as an appropriate simulant fluid for pulp suspension mixing studies. We first present a rheological study to establish similarities between the behaviour of carbopol solutions and pulp fibre suspensions. This behaviour is then incorporated in scale-down calculations to determine the appropriate operating conditions to reproduce important mixing conditions in a laboratory mixer. 6.1.1 Rheology of Pulp Fibre Suspensions and Carbopol Solutions Pulp fibre suspensions are composed of cylindrical fibres having lengths of 1−3 mm and aspect ratios from 60−400. The fibres typically flocculate into mass concentrations (called flocs) which move as discrete entities in the flow. Often multi-phase systems, such as pulp suspensions, can be considered to behave macroscopically as a single-phase fluid. However, this assumption neglects the effect suspension structure has on rheology. Although the assumption may work at certain scales, pulp suspensions can not be regarded as a continuum when the size of the dispersed phase (i.e. pulp fibres and flocs) is comparable to the geometrical dimensions of the flow field (Coussot and Piau, 1995). This is particularly important during characterization with conventional rheometric devices. For example, wall slip (a classical problem in rheology) is prone to occur in particulate suspensions (Swerin, 1998) and the limitation of fibre rotation in the narrow gaps employed for rheological measurements leads to fibre alignment which can result in measuring the viscosity of the suspending medium (water) instead of the suspension (Wikstrom et al., 1998). Due to these obstacles, rheological characterization of pulp suspensions is challenging, requiring careful interpretation of the experimental data and often resulting in large uncertainties in measured parameters. References start on page 168 Chapter 6 145 Fibre suspensions posses a significant yield stress (Bennington et al., 1990), and have been described using yield stress fluid models in some cases with a shear thinning component (Roberg 1997; Wikstrom et al., 1998; Pettersson, 2004)(Pettersson, 2004; Roberg, 1997; Wikstrom et al., 1998). The Hershel-Buckley model . τ HB = τ y + k γ n (6.1) is often used. The yield stress, τy, is a strong function of the suspension mass concentration, and its determination is accompanied by a large coefficient of variation, often reaching 50% (Bennington et al., 1990). The reported flow (n) and consistency (k) indices range widely for a given mass concentration. For example from n = 0.25 and k = 4 Pa.s0.25 (Roberg, 1997) to n = 0.07 and k = 40 Pa.s0.07 (Pettersson, 2004) for a 4 wt.% suspension. These indices also depend on the flow and geometry chosen for their measurement limiting their usefulness in other geometries. A suitable model fluid for pulp mixing studies should behave similarly to a pulp suspension and be readily characterized (so that computational and experimental data can be directly compared). Other attributes include transparency (to allow use of optically based experimental techniques), lack of toxicity and ease of handling. Possible candidates from among common model fluids include carboxymethyl cellulose, carbopol, and natrosol solutions. In this study we chose carbopol solutions because they are remarkably clear and approximate idealized yield stress fluids without significant viscoelasticity (Curran et al., 2002). The rheology of carbopol solutions is described well using the Herschel-Bulkley model (Roberts and Barnes, 2001; Curran et al., 2002; Dubash and Frigaard, 2007). References start on page 168 Chapter 6 146 A comparison of pulp suspension and carbopol solution rheology was made under similar conditions. A wide-gap cup and vane geometry was used with the vane geometry treated as a cylinder of fluid having the same diameter as the vane (Barnes and Nguyen, 2001; Cullen et al., 2003). This allowed the use of concentric cylinder geometrical equations for calculating the shear stress−shear rate relationship from torque and angular velocity measurements. For a vane having height, h, and radius, R, rotating at a constant angular velocity, the measured torque, Tq, is related to the shear stress at the vane periphery, σ, by (Krieger and Woods, 1996) σ= Tq (6.2) 2π R 2 h In the wide-gap geometry, the shear rate is calculated using the angular velocity, Ω, and σ (from Eqn. (6.2)) • γ = (2Ω) d (ln Ω) d (ln σ ) (6.3) This equation is applicable in the large gap geometry for any fluid (including non-Newtonian fluids with yield stresses) provided the no-slip boundary condition is valid (Barnes and Nguyen, 2001). An additional challenge when characterizing certain yield stress fluids is thixotropy associated with their structure (Moller et al., 2006). Both pulp suspensions and carbopol solutions are thixotropic due to their ability to form entanglements and cross-linkages, although on different time and length scales. In carbopol solutions, cross-linkages are formed between individual polymer particles (on the molecular scale) which are gelled together and act as a concentrated dispersion even at low concentrations. The exterior of these particles is covered with the free ends of the gel strands which interact with adjacent microgel particles to produce very high References start on page 168 Chapter 6 147 viscosity at low shear stresses (Carnalli and Naer, 1992). In pulp suspensions, physical crosslinkages form between fibres forming flocs that are disrupted by applied shear stress. The rheological behaviour of both systems is affected by their structure and hence by the history of shear experienced. Indeed, marked differences can be measured for identical pulp systems that have been prepared using different protocols. To minimize thixotropic effects and to establish a common structure for the systems they were pre-sheared and allowed to relax prior to rheological characterization. 6.2 Experimental Sifted carbopol 940 powder (Noveon Inc., Waterloo, ON) was dispersed in deionized water using an axial flow impeller (Lightnin A100) at ~200 rpm for 5 to 6 hours allowing complete hydration of the polymer and the formation of an aqueous solution. The solution was left to stand overnight to allow entrained air to escape, and then neutralized with 18 wt.% NaOH solution to a pH of 6.2 to 7.0 (2.3 parts NaOH with 18% per part of polymer by weight according to the manufacture’s recommendations (Noveon, 2002) while mixing at a lower rate (~150 rpm). The sodium hydroxide ionizes the carbopol, generating negative charges along the polymer back bone which uncoils the molecules into an extended conformation, thickening the solution and producing a yield stress (between pH 5 and 9)(Dubash and Frigaard, 2007). Carbopol solutions at 0.1, 0.15 and 0.2 wt% were prepared. Pulp fibre suspensions were prepared from bleached softwood kraft (FBK) (Domtar Inc., Windsor, QC) by breaking up the dried pulp sheets and rehydrating the fibres with tap water References start on page 168 Chapter 6 148 using a disintegrator (TMI, Montreal, QC) for 30 minutes at 120 rpm. FBK suspensions of 0.75, 2.0, 2.5 and 3.0 wt.% were prepared. The rheological properties of the pulp suspensions and carbopol solutions were measured using a Bohlin C-VOR (Malvern, Worcestershire, UK) rheometer in controlled stress mode. A fourvaned bob of 38.5 mm height and 25 mm diameter was placed in a cylindrical cup having an internal diameter of 100 mm and height of 64 mm giving a gap size of 37.5 mm. A series of tests were conducted to determine the optimum pre-shearing conditions for each system. For the pulp suspensions, tests were carried out at each suspension concentration with different levels of shear stress applied for different times. The optimum pre-shearing conditions were identified as those for which variation in pulp viscosity vs. shear rate following treatment was similar for tests conducted on different aliquots of the same pulp sample. These conditions were then used in subsequent tests: for pulp suspensions of 0.75, 2.0, 2.5 and 3.0 wt.%, 300 s of pre-shearing at 15, 80, 120 and 180 Pa, respectively, followed by relaxation for 300 s. The carbopol solutions at 0.1, 0.15 and 0.2 wt.% were pre-sheared at 50, 100 and 120 Pa for 100 s, respectively, followed by relaxation for 100 s. (Less extensive pre-shearing was needed for the carbopol solutions as they display less pronounced thixotrophy than the pulp suspensions.) Following pre-shearing and relaxation, the steady-state flow curves of the pulp suspensions and carbopol solutions were obtained by applying a logarithmic shear stress ramp with the strain response measured after 5 s at each applied stress. The torque versus angular velocity data were transformed to shear stress versus shear rate data using Eqns. (6.2) and (6.3). References start on page 168 Chapter 6 149 A geometrically-scaled Maxflo MarkII impeller (Chemineer Inc., Dayton, OH) having D = 10.2 cm was selected for the scale-down study presented here. This impeller was recently characterized (power number, NP, and axial force number, NF) in pulp suspensions by Bhole et al. (2008). The same procedure was followed for the carbopol solutions. 6.3 6.3.1 Results and Discussion Selection of Scaling Criteria The aim of this study is to define a model system for studying macroscale mixing of pulp suspensions at the laboratory scale. Thus, it is important to show that the model system provides a good representation of the original mixing process, allowing important flow characteristics to be recreated, including cavern formation, suspension channeling, recirculation and stagnation (Ein-Mozaffari et al., 2005; Ford et al., 2007). As the quality of mixing achieved is largely governed by the cavern size created in the mixer (Elson et al., 1986; Amanullah et al. 1998; EinMozaffari et al., 2005; Adams and Barigou, 2007), scale-down was based on achieving an equal dimensionless cavern size in the pulp suspension and model carbopol systems. The effect of other scaling criteria was also examined. A comparison of design and operating parameters for two industrial chests and a pilot-scale chest is given in Table 6.1. The cavern volume was not measured for the industrial chests (due to measurement limitations) although an estimate of the effective mixed volume fraction (EinMozaffari et al., 2005), VFM/VT, is given based on dynamic tests. The pilot-scale chest is geometrically similar to a typical single-impeller industrial chest, and the conditions shown are for operation when the effective mixed volume fraction is similar to industrial tests. References start on page 168 Chapter 6 150 Table 6.1. Comparison of industrial-scale pulp mixing chests with the pilot-scale chest used in scaling. Evaluated for a 3% FBK softwood pulp suspension. Fluid Properties τy (Pa) ρ (kg/m3) Chest Dimensions V, m3 D/T = D/W T = W (cm) Z (cm) L(cm) D (cm) Operating Parameters N (rpm) P (W) Tq (N.m) ηe (Pa.s) Dimensionless Parameters Fr NP NF Re' Dc/T VFM/VT Scaling Parameters Tq/V (N.m/m3) P/V (W/m3) St (m/s) Dc/D Industrial-scale 1-impeller Industrial-scale 4–impellers Pilot-Scale Chest (Scale 1) 113 1030 90.5 1030 154.7 1031 153 0.13 716 287 744 91 193 0.26 610 280 1130 76 0.11 0.41 40 50 55 16.5 445 74600 1601 1.5 310 71600 2206 1.75 1031 231 2.14 3.23 5.10 0.29 -4150 -0.61 2.07 0.49 --1750 -0.66 4.9 0.36 0.08 149 0.61 0.7 10.5 488 21.2 -- 11.4 371 12.3 -- 19.5 2100 8.91 1.73 References start on page 168 Chapter 6 6.3.2 151 Pulp Suspension and Carbopol Rheology Flow curves for six replicate tests made with the 3 wt.% pulp suspension (pre-sheared and relaxed as discussed above) are shown in Fig. 6.1. Each curve shows a yield stress and shear thinning behaviour. The averaged data is included as a solid line in Fig. 6.1. For pulp suspensions, the shear rate was determined with accuracy ranging from ±4.1 s-1 (0.75 wt.%) to ±46.8 s-1 (3.0 wt.%). The error for the carbopol solutions was considerably smaller, reaching a maximum value ±0.2s-1 at 0.2 wt.%. The qualitative similarity between the pulp suspensions and carbopol solutions is readily apparent as shown in Fig. 6.2. At low shear rates both systems exhibit the behavior of extremely viscous Newtonian fluids, followed by a shear-thinning region Shear stress, σ (Pa) with viscosity decreasing as shear rate increases. 10 3 10 2 10 1 10 3 wt.% FBK Flow curves -4 -3 10 -2 10 10 -1 10 0 1 10 2 10 10 3 Shear rate,γ (s ) -1 Figure 6.1. Flow curves measured in six tests made of a 3 wt.% FBK pulp fibre suspension. The average flow curve is given by the solid line. Multiple definitions of yield stress exist (Barnes and Nguyen, 2001) and different values are determined depending on the definition chosen. In this study, yield stress is defined as the value References start on page 168 Chapter 6 152 above which the strain response begins to increase significantly as a function of applied shear stress in the logarithmic shear ramp test described earlier. To more accurately identify the yield stress, additional tests were performed. Starting at small shear stress values, the shear stress was ramped linearly at 1 Pa/s and the compliance (strain/stress ratio) was measured. The “instantaneous viscosity” (reciprocal of the time-derivative of compliance) was plotted versus shear stress, as shown in Fig. 6.3. In the absence of a yield stress the instantaneous viscosity would rise indefinitely. The yield stress is identified as the stress where the instantaneous viscosity begins to decrease. Repeated determinations of pulp suspension yield stresses were accompanied by large standard deviations (±16 Pa for 3 wt.%) compared with those measured for carbopol solutions (±3 Pa for 0.2 wt%). Successive tests with different aliquots of the same carbopol solution yielded the same flow curves each time, but for pulp suspensions several tests (three to six) were performed and averaged to obtain the curves shown in Figure 6.2(b). The discrepancy between successive tests increased as pulp suspension concentration increased which is attributed to the large structural scale of the fibre suspensions. In rheological experiments, the sample is usually assumed to be homogeneous with the local shear rate derived assuming that the macroscopically imposed shear is distributed uniformly throughout the sample (with wall slip avoided). However, the irregularities (flocculation) of the pulp suspension does not permit a uniform shear distribution, and since the structure of the fibre network (floc arrangement) is different between samples, the distribution of shear is not repeatable despite pre-shearing the suspensions. References start on page 168 Chapter 6 153 Shear stress, σ (Pa) 100 10 1 wt.% carbopol 0.1 0.15 0.2 0.1 (a) 10 -4 -3 10 10 -2 -1 0 10 10 10 1 2 10 Shear rate, γ (s ) -1 Shear stress, σ (Pa) 1000 (b) 100 wt.% FBK 0.75 2 2.5 3 10 1 -4 10 10 -3 10 -2 10 -1 10 0 10 1 10 2 3 10 10 4 -1 Shear rate, γ (s ) Figure 6.2. Shear stress vs. shear rate. (a) Carbopol solutions of 0.1, 0.15 and 0.20 wt.%; (b) FBK pulp fibre suspensions of 0.75, 2.0, 2.5 and 3 wt.%. Vane in large cup geometry: Vane diameter = 2.5 cm, cup diameter = 10 cm, gap size = 3.75 cm; T = 23oC. The values of yield stress measured for the carbopol solutions and pulp suspensions are summarized in Table 6.2. For the pulp suspensions, these values agree with those reported by References start on page 168 Chapter 6 154 Bhole et al. (2008) and their dependence with concentration is similar to that reported by Bennington et al. (1990). For the carbopol solutions, the values agree with those reported by Roberts and Barnes (2001). Note that we have tabulated the yield values of the various carbopol solutions with pulp suspensions having similar values. However, there are substantial quantitative differences in their behavior despite having similar yield stresses. For example, the minimum level of shear required to induce a measurable deformation (proportional to the measured shear rate) in the carbopol solutions is considerably lower than that required for pulp suspensions. Also, the same imposed stress results in a much higher shear rate for the pulp suspensions. For instance, at a constant imposed shear stress of 5 Pa, the shear rate measured for the 0.75 wt.% pulp suspension is about one order of magnitude higher than that measured for the 0.1 wt.% carbopol solution. Thus, the “instantaneous viscosity” of a pulp suspension before yielding is much lower than that of the carbopol solution with matching yield stress. This attests to the severe interaction of carbopol microgel particles at low shear stresses, which is significantly stronger than that of fibre entanglement at the same shear stress. Table 6.2. Yield stress of pulp suspensions and carbopol solutions. FBK Pulp Suspensions Wt.% Yield stress (Pa) 0.75 12 ± 2 2.0 57 ± 7 2.5 104 ± 16 3.0 155 ± 16 Carbopol 940 Solutions wt.% Yield stress (Pa) 0.10 10.5 ± 0.6 0.15 61 ± 2.5 0.20 91 ± 3.0 - When fit to the Hershel-Buckley equation, the 3.0 wt.% pulp suspension gives . 0.21 τ HB = 155 + 135 γ R 2 = 0.6 3.0 wt% FBK (6.4) References start on page 168 Chapter 6 155 Instantaneous viscosity, η (Pa.s) 10000 1000 100 wt.% 10 carbopol 0.1 0.15 0.2 1 (a) σy (Pa) 10.5±0.6 61.4±2.5 91.5±3.0 1 10 100 Shear stress, σ (Pa) Instantaneous viscosity, η (Pa.s) 10000 1000 (b) 100 10 wt.% FBK 1 0.75 2 2.5 3 1 σy(Pa) 12±2 57±7 104±16 155±13 10 100 1000 Shear stress, σ (Pa) Figure 6.3. Instantaneous viscosity vs. shear stress. (a) Carbopol solutions of 0.1, 0.15 and 0.20 wt.%; (b) FBK pulp fibre suspensions of 0.75, 2.0, 2.5 and 3 wt.%. Vane in large cup geometry: Vane diameter = 2.5 cm, cup diameter = 10 cm, gap size = 3.75 cm; T = 23oC. The low R2 is expected, as shown by the flow curves in Fig. 6.1. To make this fit, the yield stress was determined (as described earlier and tabulated in Table 6.2). Non-linear regression is then used to determine k and n, with the accuracy of the fit affected by the yield stress used. The References start on page 168 Chapter 6 156 parameters in Eqn. (6.4) are measured to: τHB = 155 ± 16 Pa, k = 135 ± 40 Pa.sn and n =0.21 ± 0.03. Similar accuracy was obtained for suspensions at other concentrations: τ HB = 12.8 + 7.80 γ 0.25 R 2 = 0.81 0.75 wt.% FBK (6.5) τ HB = 57.2 + 13.7 γ 0.33 R 2 = 0.71 2.0 wt.% FBK (6.6) 2.5 wt.% FBK (6.7) τ HB = 104 + 14.8 γ 0.28 R 2 = 0.66 A comparison between the 3.0 wt.% suspension after yielding and the Bingham plastic approximation often used for pulp suspensions (Gibbon and Attwood, 1962; Blasinski and Rzyski, 1972; Wikstrom et al., 1998; Roberts and Barnes, 2001; Ford et al., 2006) is made in Fig. 6.4. The Bingham approximation is not an accurate depiction of suspension behaviour over the range of shear rates present in a pulp mixer. Shear stress, σ (Pa) 1000 100 10 -2 10 3.0 wt.% FBK measured Bingham k =0.001 Bingham k =0.01 Bingham k =0.1 Bingham k =1.0 10 -1 0 10 10 1 2 10 10 3 -1 Shear rate, γ (s ) Figure 6.4. Predictions of shear stress when using the Bingham approximation against the averaged flow curve measured for a 3 wt.% pulp fibre suspension. References start on page 168 Chapter 6 157 When the carbopol data is fit to the Hershel-Buckley model a much better fit is found τ HB = 10.5 + 2.67 γ 0.48 R 2 = 0.99 0.10 wt.% (6.8) τ HB = 61.4 + 46.8 γ 0.24 R 2 = 0.96 0.15 wt.% (6.9) τ HB = 91.4 + 46.0 γ 0.33 R 2 = 0.99 0.20 wt.% (6.10) The standard error of the estimated kHB and n was ± 0.06 and ±0.01 at 0.1 wt.%, ±1.4 and ±0.01 at 0.15 wt.% and ±2.4 and ±0.04 at 0.2 wt.% These parameters compare well with measurements reported by Robert and Barnes (2001). Handling the more concentrated carbopol solutions is difficult. At 0.2 wt.% the neutralized carbopol solution is a thick gel and removing gas bubbles becomes exceedingly difficult. Therefore the 0.2 wt.% solution is not a good choice for use in the lab-scale mixer as entrapped gas would limit the accuracy of optical based velocity measurements. 6.3.3 Establishing Similarity between Pulp and Carbopol Mixing Studies Industrial agitated pulp stock chests are often rectangular in shape, with chest length (L) to width (W) ratios (L/W) ranging from 1.0−2.0 for single impeller chests, to 1.5−2.0 for chests with 2 or 4 impellers. Chest height (H) to width ratios (H/W) range from 0.5−1.5. Optimum agitation and power efficiency is obtained for impeller diameter (D) to W ratios (D/W) of 0.4−0.5 (Gibbon and Attwood, 1962; Ein-Mozafarri, 2002), although in industrial configurations D/W ranges from 0.12 for single impeller chests to 0.3 for multiple impeller chests (Ein-Mozafaarri, 2002; EinMozaffari et al., 2004; Bhattacharya et al., 2008). This is due to the desire to minimize the initial capital cost of the mixers, but results in increased power requirements and reduced motion in the References start on page 168 Chapter 6 158 vessel. A rectangular chest having parameters within these ranges was used by Ein-Mozzafari (2002) in his pilot-scale mixing studies. The 18:1 laboratory-scale mixer chosen was designed to achieve geometrical similarity with the pilot-scale mixer and typical industrial mixers (Tables 6.1 and 6.3). A Maxflo Mark II impeller (D = 10.1 cm) was used with D/W = 0.4; L/W = 1.4; liquid level (Z) to width ratio, Z/W = 1.25 and impeller clearance (C) from wall to chest width ratio, C/W = 0.175. Although the D/W chosen was higher than found for industrial chests, industrial mixing conditions could be approximated by operating the impeller at lower speeds (Bhole et al., 2008). In the laminar to transition regimes at which industrial chests operate (Gibbon and Attwood, 1962; Ein-Mozaffari, 2002; Ford et al., 2006) the effect of D/T on NP vs. Re is not significant (which confirms the ability to model industrial-scale mixing at a smaller scale with higher D/W ratios). The volume of the laboratory-scale chest is 30% of the pilot-chest (i.e. 0.03 vs. 0.11 m3). To achieve dynamic similarity on the laboratory-scale, relevant dimensionless parameters must be identical to those of at the larger scale. However, it is impossible to achieve this for all relevant parameters and compromises must be made. In this study we matched the dimensionless cavern size and examined the effect of this choice on other relevant parameters. To characterize impeller behaviour in the pulp suspensions and carbopol solutions required specification of the non-Newtonian Reynolds number, Re' Re' = ρ ND 2 η (6.11) where the viscosity, η, for a Hershel-Buckley fluid is given by References start on page 168 Chapter 6 159 . η= τy + kγ n (6.12) . γ A representative shear rate was used to evaluate Re ' based on the Metzner-Otto relationship . ( γ e = K s N )(Metzner and Otto, 1957). Although the applicability of this relationship has been debated for non-Newtonian systems having a yield stress (Hemrajani and Tatterson, 2004), it has been effectively used when 0.15 < Re' < 620 (Adams and Barigou, 2007) as in the present study. This allows Eqn. (6.11) to be rewritten as . n ηe ≈ τy + kγ e . γe = τ y + k ( K s N )n (6.13) Ks N where ηe is the effective viscosity and Ks = 10 (Metzner and Otto, 1957)(reported for marine-type impellers). Re’ is then computed using Eqn. (6.10) with η = ηe. NF and NP are correlated with Re' as shown in Fig. 6.5. The NF vs. Re' and NP vs. Re' curves do not collapse to a single operational curve for either system. This behavior suggests that Ks is not independent of the rheological properties of each system. While the applicability of the Metzner-Otto relationship to viscoplastic fluids and the dependence of Ks on rheological parameters has been studied (Hirata and Aoshima, 1996; Anne-Archard et al., 2006) it is beyond the scope of the present work. For all the scaling calculations we have used the curves measured for each system. References start on page 168 Chapter 6 160 100 NP (-) 10 1 3.0 wt% FBK 2.0 wt% FBK 0.1 wt% carbopol 0.15 wt% carbopol 0.1 0.1 1 10 100 1000 (a) Re' 100 NF(-) 10 1 0.1 0.01 0.1 3.0 wt% 2.0 wt% 0.1 wt% 0.15 wt% 1 FBK FBK carbopol carbopol (b) 10 100 1000 Re' Figure 6.5. Characterization of Maxflo II impeller (D =10.2 cm, D/T = 0.39) in softwood pulp suspensions and carbopol solutions. (a) NP vs. Re’; (b) NF vs. Re’. There are a number of equations used for predicting cavern size in yield stress fluids (Elson et al., 1986; Amanullah et al., 1998). The force model proposed by Ammanullah et al. (1998) uses NP, NF and a Reynolds number to predict the size of unbounded spherical caverns created by axial flow impellers and was chosen for this work. References start on page 168 Chapter 6 161 2 ρ N 2 D2 ⎛ Dc ⎞ = ⎜ ⎟ πτ y ⎝ D⎠ NF 2 ⎛ 4NP ⎞ +⎜ ⎟ ⎝ 3π ⎠ 2 (6.14) where Dc is the cavern diameter and (ρN2D2)/τy is termed the yield stress Reynolds number. Scale-down is made from the pilot-scale mixer used by Ein-Mozaffari (typical operating conditions are provided in Table 6.1). If pulp is used as the working fluid in the laboratory-scale chest, scaling based on equal power per unit volume, Reynolds number and dimensionless cavern diameter are given in Table 6.3. For P/V scaling, impeller speed increases but Dc/D decreases from 2.1 to 1.7. Tq/V and St also decrease. For Re' scaling, impeller speed increases and P/V increases significantly. Dc/D increases due to variation in viscosity with shear rate at the different scales. Scaling based on maintaining Dc/D yields a slightly lower Re' and a higher P/V (although lower than with Re' scaling). Tq/V and St remain relatively constant, supporting the constant torque per unit volume scaling recommended by some pulp mixer vendors. From this analysis, it appears that scaling to constant Dc/D gives the closest results to the larger scale, but at a higher P/V. The latter is not an issue when chemical reactions are not occurring. The operating conditions for the laboratory-scale chest must be changed to maintain similar mixing conditions when a model fluid is used. Table 6.4 shows the consequences of using 0.1 and 0.15 wt.% carbopol solutions in place of the 3 wt.% pulp suspension. Two scaling criteria are investigated: equal Re’ and equal Dc/D. Calculations begin by assuming an impeller rotational speed, computing the effective shear rate using the Metzner-Otto correlation and calculating ηe using Eqn. (6.13). The other parameters are then computed (calculation of Re’ allows NP and NF to be determined using the data in Fig. 6.5) with N changed iteratively until the scaling criterion is met. References start on page 168 Chapter 6 162 Table 6.3. Process parameters for pilot- and laboratory-scale mixing chests upon scaledown. Evaluated for a 3% FBK softwood pulp suspension. Pilot-scale Chest Fluid Properties τy (Pa) k(Pa.s) n ρ (kg/m3) Chest Dimensions V, m3 T = W (cm) Z (cm) L(cm) D (cm) Zc (cm) C (cm) D/T = D/W Z/W L/W C/W Zc/W Scaling Criterion Operating Parameters N (rpm) P (W) Tq (N.m) ηe (Pa.s) Dimensionless Parameters Fr NP NF Re' Dc/T Scaling Parameters St (m/s) Tq/V (N.m/m3) P/V (W/m3) Dc/D Laboratory-scale Chest 155 135.9 0.21 1031 0.11 40 50 55 16.5 12 7 0.4 1.25 1.4 0.17 0.30 Equal P/V 0.03 25 31 36 10.2 7 4 0.4 1.25 1.4 0.17 0.30 Equal Re’ Equal Dc/D 1031 231 2.14 3.23 1368 56.2 0.39 2.54 1736 99.5 0.55 2.08 1659 88.3 0.51 2.16 4.9 0.36 0.08 149 0.61 5.41 0.42 0.08 96.2 0.50 8.71 0.36 0.09 149 0.67 7.9 0.37 0.12 137 0.61 8.91 19.5 2100 1.73 7.31 14.7 2100 1.56 9.27 20.5 3721 1.80 8.86 19.2 3340 1.73 References start on page 168 Chapter 6 163 For the 0.1 wt.% carbopol, equal Re' scaling results in lower N and P/V and higher Dc/D. This is due to differences in the rheology of the two systems. At the same shear rate, the carbopol solution displays a lower viscosity than pulp suspension (slope of the shear stress vs. shear rate curve of Fig. 6.2). Thus, the impeller speed must be reduced to satisfy the scaling criteria. Note that the diameter of the cavern formed is larger than the one calculated for the pulp suspension (at an even larger N) and is larger than the width of the vessel. Equation (6.14) is no longer applicable once the cavern reaches the vessel wall, and the large cavern size means that mixing quality will not be similar at the two scales. Constant Dc/D scaling results in lower N, Re' and P/V. Despite these differences, the fully mixed volume fraction of this system, as defined by Ein-Mozaffari et al., (2002; 2005) is equivalent to that in the pilot-scale mixer and to that estimated for typical industrial chests (Ein-Mozaffari et al., 2004; Bhattacharya et al., 2008). This means that the percentage of carbopol yielded and in active motion is the same as in the pulp system. Consequently, mixing quality is expected to be similar. Re’ is lower due to the higher apparent viscosity of the carbopol solution at this impeller speed, and thus the model system will operate more towards the laminar regime than the pulp suspension. The same calculations were performed using the 0.15 wt.% carbopol solution. The higher yield stress, consistency index and flow index at this concentration gives a significantly higher effective viscosity at low impeller speeds, and N must be increased to satisfy the scaling criteria. Both Re' and Dc/D scaling give similar values of cavern size and Reynolds number, but the higher operating impeller speed and the difficulty characterizing the impeller at this concentration (rod climbing was observed in these tests) are harbingers of difficulties. We noted References start on page 168 Chapter 6 164 earlier that the 0.2 wt.% carbopol solution was unacceptably difficult to work with. For these reasons the 0.1 wt.% carbopol solution was chosen as the most appropriate model fluid for our studies. Table 6.4. Operating parameters for pilot-scale mixing chest with 3 wt.%pulp suspension vs. laboratory-scale mixing chest with aqueous carbopol solutions. Pilot-mixer Pulp suspension 3 wt.% Dimensionless Parameters Fr NP NF Re ' Dc/T Scaling parameters St (m/s) Tq/V (N.m/m3) P/V (W/m3) Dc/D 0.1 wt.% Equal Re’ Scaling Criterion Fluid Properties τy (Pa) k(Pa.s) n ρ (kg/m3) Operating Parameters N (rpm) P (W) Tq (N.m) ηe (Pa.s) Aqueous carbopol solutions in laboratory-scale mixer 155 136 0.21 1031 0.15 wt.% Equal Dc/D Equal Re’ 10.47 2.65 0.47 998 Equal Dc/D 61 46.67 0.24 998 1031 231 2.14 3.23 389 2.35 0.058 94 0.37 0.04 1074 21.6 0.192 1002 19.1 0.182 0.45 1.28 1.25 1.32 4.9 0.36 0.08 149 0.61 0.44 0.78 0.19 149 1.05 0.03 8.65 1.08 12.7 0.61 3.34 0.34 0.11 149 0.65 2.90 0.37 0.12 131 0.61 8.91 19.5 2100 1.73 2.08 2.15 87.8 2.25 0.50 1.41 13.9 1.73 5.74 7.19 808 1.77 5.35 6.80 713 1.73 References start on page 168 Chapter 6 6.4 165 Summary and Conclusions The mixing and handling of non-Newtonian fluids and suspensions occurs across a large range of industries, although their rheology is not fully understood. In many cases, their opacity makes studying them difficult. We have focused on one particular example, which is the macro-scale mixing of pulp fibre suspensions. To study these systems both experimentally and computationally requires the use of model fluid systems. In the present work, we overcome experimental difficulties using aqueous carbopol solutions that display qualitatively similar flow behaviour to that of the pulp suspensions. However, careful examination of the rheology of both systems is necessary to choose the appropriate conditions to use. A step by step approach leading to a successful scale-down of pulp mixing in the laboratory is presented. Scaling based on relative cavern size was found to be most appropriate to mimic the mixing quality attained with a 3 wt.% pulp suspension. Different concentrations of carbopol solutions were tested, with the 0.1 wt.% solution chosen as it permitted achieving the cavern size scaling criterion without significantly departing from the flow regime present in pulp mixing at the larger scale. 6.5 Nomenclature h vane height, cm R vane radius, cm Tq torque, N.m n flow index k consistency index, Pa.s References start on page 168 Chapter 6 166 Ks Metzner and Otto impeller constant D impeller diameter, cm T characteristic diameter of cylindrical mixing vessel or characteristic width of rectangular mixing vessel, cm W mixer width, cm H mixer height, cm L mixer length, cm Z liquid level in the chest, cm C impeller clearance from the wall at which the impeller is located, cm Zc impeller clearance from the chest bottom, cm V mixer volume, m3 VFM fully mixed volume, m3 Dc cavern diameter, cm P power, watt St tip speed, m/s Re Reynolds number Re' non-Newtonian Reynolds number NP power number NF axial force number Fr Froude number Greek letters σ shear stress, Pa Ω angular velocity, rad/s References start on page 168 Chapter 6 τy yield stress, Pa τ HB shear stress predicted by the Herschel-Bulkley model, Pa . γe . 167 shear rate, s-1 γe effective shear rate, s-1 η non-Newtonian viscosity, Pa.s ηe effective viscosity; non-Newtonian viscosity at the effective shear rate, Pa.s References start on page 168 Chapter 6 6.6 168 References Adams, L.W. and Barigou, M., 2007. CFD analysis of caverns and pseudo-caverns developed during mixing of non-Newtonian fluids. Chem. Eng. Res. Des. 85, A5, 598−604. Amanullah, A., Hjorth, S.A. and Nienow, A.W., 1998. A new mathematical model to predict cavern diameters in highly shear thinning, power law liquids using axial flow impellers. Chem. Eng. Sci. 53, 3, 455−469. Anne-Archard, D., Marouche, M. and Boisson, H.C., 2006. Hydrodynamics and Metzner-Otto correlation in stirred vessels for yield stress fluids. Chem. Eng. J. 125, 15−24. Bakker, A. and Fasano, J.B., 1993. A computational study of the flow pattern in an industrial paper pulp stock chest with a side-entering impeller. AIChE Symposium Series 89 v. 293, 118−124. Barnes, H.A. and Nguyen, Q.D., 2001. Rotating vane rheometry - a review. J. Non-Newtonian Fluid Mechanics 98, 1, 1−14. Bennington, C.P.J., Kerekes, R.J. and Grace, J.R., 1990. The yield stress of fibre suspensions, Can. J. Chem. Eng. 68, 10, 748−757. Bhattacharya, S., Gomez, C., Soltanzadeh, A., Taghipour, F., Bennington, C.P.J. and Dumont, G.A., 2008. Computational modelling of industrial pulp stock chests. Can. J. Chem. Eng. (submitted for publication). Bhole, M., Ford, C. and Bennington, C.P.J., 2009. Characterization of axial flow impellers in pulp fibre suspensions, Chem. Eng. Res. Des., 87, 4, 648−653. Blasinski, H. and Rzyski, E., 1972. Mixing of non-Newtonian liquids. Power consumption for fibrous suspensions. Inz. Chem. 2, 1, 169−182. References start on page 168 Chapter 6 169 Brown, D., Pip, J., Middlenton, J., Papadopoulos, G. and Arik, E., 2004. Experimental methods, in Handbook of industrial mixing: Science and practice (Paul, E.L., Atiemo-Obeng, V.A. and Kresta, S.M. Eds.). John Wiley & Sons, New Jersey, pp. 145−174. Carnalli, J.O. and Naer, M.S., 1992. The use of dilute solution viscometry to characterize the network properties of carbopol microgels. Colloid Polymer Sci. 270, 2, 183−193. Coussot, P. and Piau, J.M., 1995. A large-scale field coaxial cylinder rheometer for the study of the rheology of natural coarse suspensions. J. Rheol. 39, 1, 105−124. Cullen, P.J., O'Donnell, C.P. and Houska, M., 2003. Rotational rheometry using complex geometries - a review. J. Texture Studies 34, 1, 1−20. Curran, S.J., Hayes, R.E., Afacan, A., Williams, M.C. and Tanguy, P.A., 2002. Properties of carbopol solutions as models for yield-stress fluids. J. Food Scie. 67, 1, 176−180. Dubash, N. and Frigaard, I.A., 2007. Propagation and stopping of air bubbles in carbopol solutions, J. Non-Newtonian Fluid Mechanics 142, 1−3, 123−134. Ein-Mozaffari, F., 2002. Macroscale mixing and dynamic behaviour of agitated pulp stock chests. PhD dissertation, University of British Columbia, Vancouver, Canada. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2004. Dynamic mixing in agitated industrial pulp chests. Pulp Paper Can. 105, 5 T118−122. Ein-Mozaffari, F., Bennington, C.P.J. and Dumont, G.A., 2005. Suspension yield stress and the dynamic response of agitated pulp chests. Chem. Eng. Sci. 60, 2399−2408. Elson T.P., D.J., C. and Nienow, A.W., 1986. X-ray studies of cavern sizes and mixing performance with fluids possessing a yield stress. Chem. Eng. Sci. 41, 2555−2562. References start on page 168 Chapter 6 170 Ford, C., Bennington, C.P.J. and Taghipour, F., 2007. Modelling a pilot-scale mixing pulp mixing chest using CFD, J. Pulp Pap. Sci. 33, 115−120. Ford, C., Ein-Mozaffari, F., Bennington, C.P.J. and Taghipour, F., 2006. Simulation of mixing dynamics in agitated pulp stock chests using CFD. AIChE J. 52, 10, 3562−3569. Gibbon, J.D. and Attwood, D., 1962. Prediction of power requirements in the agitation of fibre suspensions, Trans. I. Chem. Eng. 40, 75−82. Hemrajani, R.R. and Tatterson, G.B., 2004. Mechanically stirred vessels, in Handbook of industrial mixing: Science and practice (Paul, E.L., Atiemo-Obeng, V.A. and Kresta, S.M. Eds.). John Wiley & Sons, New Jersey, pp. 345-390. Hirata, Y. and Aoshima, Y., 1996. Formation and growth of cavern in yield stress fluids agitated under baffled and non-baffled conditions. Chem. Eng. Res. Des. 74, 4, 438−444. Krieger, I.M. and Woods, M.E., 1966. Direct determination of flow curves of non-Newtonian fluids. 4. Parallel-plane rotational viscometer. J. App. Phys. 37, 13, 4703−4704. Metzner, A.B. and Otto, R.E., 1957. Agitation of non-Newtonian fluids. AIChE J. 3, 1, 3−10. Moller, P.C.F., Mewis, J. and Bonn, D., 2006. Yield stress and thixotrophy: On the difficulty of measuring the yield stress in practice. Soft Matter 2, 274−283. Noveon, 2002. www.pharma.noveonic.com. Bulletin 10: Neutralization procedures. Pettersson, J.A., 2004. Flow and mixing of pulp suspensions. PhD dissertation. Chalmers Technical University, Gotenburg, Sweden. Roberg, K.E., 1997. Numerical simulations of the flow of fibre suspension in a side-entering mixer. Récents Progrès en Génie de Procédés 11, 51, 203−210. References start on page 168 Chapter 6 171 Roberts, G.P. and Barnes, H.A., 2001. New measurements of the flow-curves for carbopol dispersions without slip artifacts. Rheologica Acta 40, 499−503. Swerin, A., 1998. Rheological properties of cellulosic fibre suspensions flocculated by cationic polyacrylamides. Colloids and Surfaces A − Physicochemical and Engineering Aspects 133, 3, 279−294. Wikstrom, T., Defibrator, S. and Rasmuson, A., 1998. Yield stress of pulp suspensions: The influence of fibre properties and processing conditions. Nordic Pulp Pap. Res. J. 13, 3, 243−250. References start on page 168 181 APPENDICES Appendix A: A.1 Details of the Pilot Mixing Chest CFD Model The laboratory 1:11 scale-model chest used by Ein-Mozaffari (Ein-Mozaffari, 2002) in his experimental studies was modeled using the commercial CFD software package, Fluent V6.1 (Fluent Inc., Lebanon, NH) .The computational domain within the chest was generated and discretized using Gambit software (v. 2.1.6, Fluent, Lebanon, NH). An schematic of the pilot experimental set up is shown in Fig A.1. Feed Discharge Tracer (NaCl solution) Figure A.1. Schematic of the experimental setup of the pilot scale pulp mixing chest. Inputoutput configuration simulated in chapter 2. The flow in the mixing chest is unsteady in an inertial frame of reference because the rotor/impeller blades sweep the domain periodically. However, in the absence of stators or A version of this appendix has been published. Ford† C., Ein-Mozaffari F., Bennington C.P.J. and Taghipour F., 2006. Simulation of Mixing Dynamics in Agitated Pulp Stock Chest Using CFD, AIChE Journal, 52, 10, 3562-3569 †C. Gomez formerly C. Ford Appendix A 182 baffles, it is possible to perform calculations in a domain that moves with the rotating impeller. In this case, the flow is steady relative to the rotating (non-inertial) frame and the velocity of the coordinate system is included in the equations of motion describing the flow (Luo et al., 1994). To implement this multiple reference frame (MRF) numerical approach, the calculation domain is divided into subdomains, each of which may be rotating with respect to the inertial frame. The governing equations in each subdomain are written with respect to the subdomain’s reference frame and the continuity of the absolute velocity is enforced at the boundaries between subdomains. Figure A.2 shows a three dimensional view of the solution domain for the mixing chest and the locations of the boundary conditions for the mixing vessel. For the rotating and stationary subdomains, the angular velocity and axis of rotation in each reference frame was specified. The rotation axis origin was located at the center of the impeller hub with its orientation in the positive z direction. The impeller speed was specified in the moving sub-domain and the impeller blades were treated as walls that moved with zero relative velocity. Likewise, the chest walls were given a velocity of zero in the absolute reference frame. The velocity inlet was defined as having a flat profile ten diameters upstream of the chest with a magnitude ν = 4Q π D2 and a negative y orientation. The exit boundary condition was specified as fully developed flow at the end of a pipe with an L/D ratio of ten. References start on page 188 Appendix A 183 inlet Stationary wall with zero shear Moving wall at zero velocity relative to impeller angular velocity 16.5 Outlet 20 cm 14.5 cm Stationary wall (zero absolute velocity) 20 cm MRF with N rpm angular velocity Stationary subdomain 50 cm 53 cm 40 cm Figure A.2. Three dimensional view of Fig. 2.1. (The impeller is centered on the back wall 12 cm above the floor and has an offset of 7cm on the z direction). In general, the viscosity for non-Newtonian fluids is a function of all three invariants of the rate of deformation tensor, D . However, for non-Newtonian liquids, µ can be considered to be a function of shear rate, γ , which is related to the second invariant of D and is defined as γ = (A.1) D:D with D given by ⎛ ∂u j ∂u i D=⎜ + ⎜ ∂x ∂x j ⎝ i ⎞ ⎟ ⎟ ⎠ i, j = 1,2,3 (A.2) The apparent viscosity of the suspension was approximated accordingly, as a function of the shear rate as defined in Eqn A.1, using the modified Herschel-Bulkley model available in Fluent. At very low strain rates ( γ ≤ τ y / µ 0 ) the suspension is modeled as an extremely viscous fluid with viscosity µo . As a result, regions in the suspension that are stationary in reality have a very References start on page 188 Appendix A 184 low velocity in the CFD solution. As the strain rate is increased and the yield stress, τ y , passed, suspension behaviour is described using a power law, as described by Eqns. 2.1a and 2.1b. Since suspension flow observed in the experiments is primarily laminar, the mass and momentum conservation equations were solved in the laminar regime using a second order upwind scheme. (Although flow near the impeller might be turbulent, the fluctuation velocities will diminish quickly in the areas outside of the impeller zone due to the fibrous structure of the suspension. Considering the flow as laminar in the entire domain will not impact the flow field significantly). Further, the non-Newtonian Reynolds number calculated for the impeller using the method of Gibbon and Attwood (1962) and Blasinski and Rzyski (1972) was below 103 and therefore in the laminar regime according to criteria used to characterize agitated vessels. To obtain converged solutions of the CFD model, calculations were started at low impeller speeds and the impeller speed gradually increased until the desired operating condition was reached. Iterations were continued until the scaled residuals for each transport equation were below 10-5. Spatial grid independence was verified by demonstrating that additional grid lines near the impeller surface did not significantly change the calculated velocity magnitude measured in a line 8 cm in front of the impeller centered on the impeller hub. This location was chosen because it is critical for the numerical solution as the highest velocity gradients and mass imbalance residuals were located in this area. The final dynamic adapted mesh used (which contained 202,010 cells) was further refined performing one and two additional velocity gradient adaptations (426,520 and 786,422 cells, respectively). Increasing the number of cells by a factor of four did not change the velocity magnitude in front of the impeller by more than 3%. References start on page 188 Appendix A 185 Computations were carried out using a 2.4 GHz Pentium 4 CPU with 1.0 GB of RAM. Second order convergence was typically achieved after 3−5 hours. After a mesh independent result was achieved, the possible physical modeling errors were evaluated by comparing the simulation results with measured dynamic data. This was accomplished by creating time-dependent solutions to the steady-state model results and comparing the model output directly with the experimental measurements, or by fitting both the simulation and experimental data sets to an independent grey-box mixing model and comparing model parameters. The transfer function of the chest ( Gchest) obtained from the grey box model shown in Fig 2.2 is given by Gchest = f e −T1s + 1 + τ 1s (1 − f )(1 − R ) e −T2 s 1+τ 2s (A.3) R eT2 s 1− 1 +τ 2s To estimate the model parameters, the system was excited. The choice of excitation signal has a substantial influence on the accuracy with which the parameters could be estimated, thus, a frequency modulated random binary signal was designed to provide excitation of the system at the critical chest frequencies (Ein-Mozaffari, 2002) The same pseudo random binary input used by Ein-Mozaffari(2002). for the experimental tests was simulated in Fluent using a user- defined function written in the C programming language. To obtain raw data that could be directly compared with the experimental measurements as shown in Figs. 2.4 and 2.7, the mass fraction input concentration of tracer was set to a value References start on page 188 Appendix A 186 analogous to the measured conductivity values (as conductivity is proportional to mass concentration) and the input/output data recorded. The discretization step size, ∆t , was initially set one order of magnitude smaller than the time constant estimated for the mixing system, and was then further reduced to ensure convergence after 20 − 30 iterations at each time step (Fluent, 2004) Dynamic data obtained using different step sizes at each flow rate were compared to verify discretization independency. For example, at Q = 37.1 L/min, the model-generated dynamic response was compared with sampling intervals of ∆t =1 and 0.5 seconds. The predicted responses overlapped, indicating temporal independence. Similar results were found for tests made at Q = 7.9 L/min for time step sizes of 5 and 2.5 seconds. It turns out that the sampling intervals chosen for obtaining the dynamic response of the virtual system were identical to those used by Ein-Mozaffari (2002) for acquiring the original experimental data. The dynamic parameters were estimated from the input-output data using the numerical method described by Kammer et al (2005). The power input to the impeller (P) was calculated from the CFD simulation data using P = 2πNM where M is the moment vector about the center of the impeller. The calculated power input matched the experimentally determined values well as shown in Fig.A.3, particularly at high impeller speeds where the deviation was as low as 1%. For all the conditions studied, the values were over-predicted by CFD although the dependence of power input as a function of rotation speed was very similar. This discrepancy may be due to a number of the assumptions made in the simulation, including: the rheological model chosen to describe suspension rheology, the use of the laminar-flow model in the region near the impeller, and/or the inability of the model to capture the small recirculation zones and swirls near the blade tips. References start on page 188 Appendix A 187 700 Measured Power Computed Power 600 P (W) 500 400 300 200 1000 1100 1200 1300 1400 1500 1600 N (rpm) Figure A.3. Calculated and measured power input vs. impeller rotational speed. Cm = 3.3% FBK suspension at Q = 37.1 L/min. (The representative error bars are based on the accuracy of the torque meter (± 0.05 Nm) used to measure the experimental data). References start on page 188 Appendix A 188 A.2 References Blasinski, H. and Rzyski, E., 1972. Mixing of non-Newtonian liquids. Power consumption for fibrous suspensions. Inz. Chem, 2, 1, 169−182. Ein-Mozaffari, F., 2002. Macroscale mixing and dynamic behaviour of agitated pulp stock chests, PhD dissertation, University of British Columbia. Fluent, 2004. V 6.2 user's guide documentation, Lebanon, NH. Gibbon, J.D. and Attwood, D., 1962. Prediction of power requirements in the agitation of fibre suspensions. Trans. Inst. Chem. Eng., 40, 75−82. Kammer, L.C., Ein-Mozzafari, F., Dumont, G.A. and Bennington, C.P.J., 2005. Identification of channeling and recirculation parameters of agitated pulp stock chests. J. Process Control, 15, 31−38. Luo, J.Y., Gosman, A.D. and Issa, R.I., 1994. Prediction of impeller induced flow in mixing vessels using multiple frames of reference. Inst. Chem. Eng. Symposium Series, 136, 549−556. References start on page 188 189 Appendix B B.1 Procedure for Estimation of Discretization Error in CFD Simulations† Step 1: Define a representative grid size h. For field variables, the local cell size can be used. ⎡1 h=⎢ ⎣N 1/3 ⎤ (∆Vi ) ⎥ ∑ i =1 ⎦ N (B.1) where Vi is the volume of the ith cell, and N is the total number of cells used Step 2: Select three significantly different set of grids, and run simulations to determine the value of a key variable ( φ ) important to the objective of the study. It is desirable that the grid refinement factor, r=hcoarse/hfine, be about or greater than 1.3 (Celik, 2004). Step 3: Let h3 < h2 < h1 and r12 = h1 / h2 , r23 = h2 / h3 , and calculate the apparent order, p, of the method using Eqns. (B.2) through (B.5) p= 1 ln ε12 / ε 23 + q( p ) ln(r23 ) ⎛r p −s⎞ q ( p ) = ln ⎜ 23 p ⎟ ⎝ r12 − s ⎠ s = 1⋅ sign(ε12 / ε 23 ) ε12 = φ1 − φ2 ε 23 = φ2 − φ3 φk denoting the solution on the kth grid. Step 4: Calculate the extrapolated values from †This appendix contains extended details on the calculations of discretization error presented in Chapter 4. (B.2) (B.3) (B.4) (B.5) Appendix B 190 p φext 23 = (r23 φ1 − φ2 ) (B.6) (r23 p −1) Step 5: Calculate and report the following error estimates, along with the apparent order p: Approximate relative error: ea 23 = φ3 − φ2 φ3 (B.7) Extrapolated relative error: eext 23 = φext 23 − φ1 φext 23 (B.8) The fine grid converge index: GCI 23 fine = 1.25ea 23 r23 p − 1 (B.9) Table 4.1 illustrates this calculation procedure for the three grids used in this study. Note that this does not account for modeling assumptions. The numerical uncertainty in the fine grid solution is reported as the GCIfine23 percentage References start on page 191 Appendix B 191 B.2 References Celik, I. B., 2004. Procedure for estimation and reporting of discretization error in CFD applications. Fluids Eng. Division of ASME. References start on page 191 192 Appendix C C.1 Extended Details on Impeller Characterization Using CFD Figure C.1. Isometric view of the grids used for simulations in the square (S1) and cylindrical baffled (T1) with the Maxflo impeller D=10.1cm. Final three-dimensional meshes of T1 and S1 had 463742 and 340575 computational cells respectively. This appendix contains extended details on the CFD simulations in chapter 5. Appendix C 193 0.040 340575 cells 434543 cells 527565 cells 0.035 Velocity Magnitude m/s 0.030 0.025 0.020 0.015 0.010 0.005 -0.15 -0.10 -0.05 0.00 0.05 0.10 0.15 X position (cm) (a) 0.040 340575 cells 434543 cells 527565 cells 0.035 Velocity Magnitude m/s 0.030 0.025 0.020 0.015 0.010 0.005 -0.15 -0.10 -0.05 0.00 0.05 0.10 0.15 X position (cm) (b) Figure C.2. Test for grid independence. The figures show the magnitude of velocity at lines (along the x-axis) located 3 cm in front of the impeller in vessel T1 with hardwood pulp at Cm= 3%, N= 302 rpm .(a) line located at the plane z = 0 cm.(b) line located at the plane z = 5 cm. References start on page 200 Appendix C 194 The calculated velocity field captures the circulation pattern of the suspension in the vessels adequately from a qualitative point of view. Velocity vectors plots are used to investigate in detail the computed velocity field in the vessels. Figure C.3 shows the planar velocity fields at planes z = 0 cm, z = 5 cm and z = -5cm for both the cylindrical baffled and the square vessel at N = 404 rpm. The velocity magnitude range covers values from 0 m/s to 0.02 m/s, the equivalent of 1% of the tip speed. Regions with velocities that are higher than 0.02 m/s, are displayed as empty spaces and represent the caverns created by the impeller motion within the vessel volumes when using the cavern boundary definition presented by Jaworski and Nienow (1994). Figure C.3. Velocity vectors colored by velocity magnitude (m/s) at planes z = 0 cm, z = -5 cm and z = 5 cm. Hardwood pulp at Cm = 3%, N = 404 rpm. Figure C.4 illustrates the velocity vectors at the mid plane (z = 0 cm) for three different impeller speeds (low, medium and high) for the square vessel. The same velocity magnitude range was References start on page 200 Appendix C 195 used for all plots (from 0 m/s to 0.005 m/s). Areas in the mid plane, at which the velocity exceeds 0.05m /s, are displayed as empty spaces in the vector plots. Figure C.4. Velocity vectors colored by velocity magnitude (m/s) at plane z = 0 cm. Hardwood pulp at Cm= 3%. Another way to visualize the velocity field is by illustrating the mixing patterns calculated in the numerical solution. Figure C.5 displays the flow patterns generated by tracking massless particles released from a circular surface with D = 10.1 cm, located right in front of the impeller. In the vicinity of the impeller, the particles follow the impeller motion closely. After a certain period of time, the particles evolve and leave this area following the mixing patterns shown, in agreement with the expected flow patterns produced by axial impellers (Oldshue, 1983). References start on page 200 Appendix C 196 Figure C.5. Pathlines of massless particles released from a circular surface of D = 10.1cm located in front of the impeller at y = 1.2cm (coordinate system is centered at the impeller) Hardwood pulp at Cm= 3%, N = 404 rpm C.2 Power and Force Number: CFD vs. Experimental The power number was computed according to Eqn. (5.2), using the calculated torque from the CFD solution at the selected range of impeller speeds. The NP vs N (rpm) curves calculated from CFD and experimental measurements are shown in figure C.6 for the square vessel. The shape of the power number curves is adequately followed by the simulation results. However, NP is under predicted by avg. 22-24% for most situations, regardless of the particular geometry, except at lower impeller speeds (N = 50, 100 & 160 rpm) where the deviation between experimental and CFD calculated values is about 8% for the square vessel and 13% for the cylindrical vessel. At a fixed impeller speed, both the measured and predicted torque was about avg. 12% higher in the square tank, which is explained by the presence of baffles in the cylindrical configuration. The presence of baffles aids in keeping the flow symmetric, which gives lower forces perpendicular to the shaft that result in lower torque variations and bending moments. References start on page 200 Appendix C 197 Experimental CFD NP (-) 10 1 0.1 100 1000 N rpm Figure C.6. Comparison of NP from experimental and computational tests for vessel S1. Maxflo impeller D = 10.2cm in a hardwood pulp at Cm= 3%. Several axial loads are imposed on the impeller shaft including the weight, the pressure forces and in the case of an axial flow impeller, the axial thrust upwards (Dickey and Fasano, 2004) The predicted axial thrust acting on the impeller was obtained from the computed vertical fluid force acting on the impeller surface. This was used to calculate the axial thrust number as defined by Eqn. (5.1) over the selected range of impeller speeds. NF vs. N (rpm) curves calculated from CFD and experimental measurements are shown in figure C.7 for the square vessel. The shape of the axial number vs. impeller speed curves is properly described by the CFD simulation results. NF is under predicted by avg. 16 % and 22% for the cylindrical fully baffled vessel and the square vessel respectively. At a fixed impeller speed, the measured thrust in the square vessel was avg. 7% higher than measured in the cylindrical vessel. This is small difference was not observed in the NF values computed from CFD which are very close to each other (within 4% ) for all circumstances and only at highest impeller speed NF computed from the square tank is 2% larger than NF from the cylindrical vessel. References start on page 200 Appendix C 198 10 Nf (-) Experimental CFD 1 0.1 100 1000 N rpm Figure C.7. Comparison of Nf from experimental and computational tests for vessel S1. Maxflo impeller D = 10.2cm in a hardwood pulp at Cm= 3%. The fluid forces and torque calculated in CFD are computed from the pressure and velocity fields calculated over the impeller region which are entirely dependant on the boundary conditions and material properties specified for the simulation, more particularly so, on the fluid viscosity. In the vicinity of the impeller, where shear rates are particularly high, the viscosity of the suspension is calculated by Eqn. (2.1b). This equation defines the non Newtonian viscosity η of . the yielded suspension as a function of γ , defined by the parameters τ , µo, k and η. The calculated value of η is especially sensitive to the value of the yield stress τy which carries a standard deviation of about 36% equivalent to ± 12 Pa, from the experimental measurements. Consequently, the estimated viscosity in the near impeller region is accurate within ± 36% of its actual value when all other parameters in Eqns. (2.1a) and (2.1b)are defined with full confidence (i.e. with out any other error added). In reality, the parameters n = 1.0, k = 0.001 Pa.s come from physical reasoning instead of an actual fitting of the suspension’s flow curves to the model and References start on page 200 Appendix C 199 thus they carry an amount of error that can not be quantitatively accounted for, because, the after yield behavior of pulp suspensions has not yet been well characterized (Bennington, 1988). Moreover, there exists a significant difficulty associated with measuring and modeling the actual rheology of fiber suspensions that is beyond the scope of the work presented in chapter 5. When the uncertainty associated with the viscosity definition of the suspension in the CFD simulations is considered, the deviation of the predicted power numbers and axial thrust numbers from the experimentally determined ones is justified. References start on page 200 Appendix C 200 C.3 References Bennington, C.P.J., 1988. Mixing pulp suspensions, PhD dissertation, University of British Columbia, Vancouver, Canada. Dickey, D.S. and Fasano, J.B., 2004. Mechanical design of mixing equipment, in: Paul, E.L., et al. (Eds.), Handbook of industrial mixing. Science and practice, John Wiley & Sons, New Jersey, pp. 145−174. Jaworski, Z., Pacek, A.W. and Nienow, A.W., 1994. On flow close to boundaries in yield stress fluids. Chem. Eng. Sci., 49, 3321−3324. Oldshue, J.Y., 1983. Fluid mixing technology, McGraw-Hill Publications Co., New York. References start on page 200
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Numerical and experimental investigation of macro-scale...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Numerical and experimental investigation of macro-scale mixing applied to pulp fibre suspensions Gómez, Clara 2009
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Numerical and experimental investigation of macro-scale mixing applied to pulp fibre suspensions |
Creator |
Gómez, Clara |
Publisher | University of British Columbia |
Date Issued | 2009 |
Description | Effective mixing of pulp fibre suspensions is essential for many pulp and paper manufacturing processes. Pulp suspensions display non-Newtonian rheology and possess a yield stress, which complicates mixing. In order to improve our understanding of pulp mixing in agitated vessels, a series of studies was undertaken to assess the suitability of using computational fluid dynamics (CFD) to model these systems. CFD simulations for laboratory-scale and industrial-scale mixing chests were developed with the pulp suspensions treated as Bingham fluids. The computed flow fields were used to determine the dynamic response of the virtual mixers, which was then compared with experimental measurements providing very good agreement under conditions of moderate agitation. Comparison between calculations and measurements of torque and axial force was also good (relative average error of 12% to 24%). The simulation results provided insight into the mixing flow occurring within the systems, showing the formation of caverns around the impeller(s), the location of stagnation regions and the presence of channeling. However, the accuracy of these predictions was limited by the Bingham model used to describe the suspension's rheology and the uncertainty to which the suspension's yield stress could be measured. To assess the degree to which the approximated rheology contributed to the CFD results, the mixing of a model fluid having a well-defined rheology (Newtonian glycerin solution) was extensively investigated in a laboratory-scale vessel using a typical industrial geometry (rectangular chest, side-entering axial flow impeller). The flow fields were measured using particle image velocimetry (PIV) and compared with CFD computations for identical operating conditions. Good agreement was found (avg. 13.1% RMS deviation of local axial velocities) confirming that the approach used in the CFD model was adequate to calculate the complex mixing flow fields for a Newtonian fluid. These results encouraged further research on extending the combined CFD and PIV application on a system more closely representative of pulp suspension agitation. Calculations were then performed to select a suitable non-Newtonian model fluid and appropriate operating conditions to model pulp suspension mixing in the laboratory-scale vessel, with the dimensionless cavern diameter as the mixing criterion for conservation between the two systems. |
Extent | 13564023 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-06-26 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NonCommercial-NoDerivatives 4.0 International |
DOI | 10.14288/1.0228853 |
URI | http://hdl.handle.net/2429/9707 |
Degree |
Doctor of Philosophy - PhD |
Program |
Chemical and Biological Engineering |
Affiliation |
Applied Science, Faculty of Chemical and Biological Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2009-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nc-nd/4.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2009_fall_gomez_clara.pdf [ 12.94MB ]
- Metadata
- JSON: 24-1.0228853.json
- JSON-LD: 24-1.0228853-ld.json
- RDF/XML (Pretty): 24-1.0228853-rdf.xml
- RDF/JSON: 24-1.0228853-rdf.json
- Turtle: 24-1.0228853-turtle.txt
- N-Triples: 24-1.0228853-rdf-ntriples.txt
- Original Record: 24-1.0228853-source.json
- Full Text
- 24-1.0228853-fulltext.txt
- Citation
- 24-1.0228853.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0228853/manifest