Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

On the fluid mechanics of viscoplastic particle suspension fractionation : understanding multilayer spiral… Al-Shibl, Mohammad 2016

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
24-ubc_2016_september_alshibl_mohammad.pdf [ 2.66MB ]
Metadata
JSON: 24-1.0305821.json
JSON-LD: 24-1.0305821-ld.json
RDF/XML (Pretty): 24-1.0305821-rdf.xml
RDF/JSON: 24-1.0305821-rdf.json
Turtle: 24-1.0305821-turtle.txt
N-Triples: 24-1.0305821-rdf-ntriples.txt
Original Record: 24-1.0305821-source.json
Full Text
24-1.0305821-fulltext.txt
Citation
24-1.0305821.ris

Full Text

On the Fluid Mechanics of ViscoplasticParticle Suspension Fractionation:Understanding Multilayer SpiralPoiseuille Flow in an AnnulusbyMohammad Al-ShiblB.Sc., Mechanical Engineering, King Saud University, 1993M.Sc., Mechanical Engineering, King Fahd University of Petroleum and Minerals,2006A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFDOCTOR OF PHILOSOPHYinThe Faculty of Graduate and Postdoctoral Studies(Chemical and Biological Engineering)THE UNIVERSITY OF BRITISH COLUMBIA(Vancouver)July 2016c© Mohammad Al-Shibl 2016AbstractThe focus of the present work is the study of laminar spiral multi-layer vis-coplastic flow in annular geometries. The motivation for the present workstems from an interest in utilizing such flow in fractionation of particle sus-pensions. The work is presented in four studies.In the first study we solve the fully developed condition of the flowanalytically. This solution is considered a reliable reference to validate theremaining numerical studies. Also it provides a means to test the stabilityof the flow.The second study is related to the fractionation of particle suspensionsutilizing the solution demonstrated in the first study. We develop fractiona-tion curves of particles of different sizes in fluids with different rheology. Wedevelop a code to simulate thousands of flow cases (a flow case has a uniquecombination of streams flowrates) of known fluids properties. We predict thefractionation operating window (some range of streams flowrates) needed forsuccessful fractionation.iiAbstractIn the third study, we examine the flow in the full annulus geometryincluding the entrance region of the flow. Here, we estimate the flow entrancelength in order to design the length of the mixing zone of the two streams inthe continuous fractionation device. We study the effect of Kelvin-Helmholtzand density current instabilities on the flow.In the fourth study, we attempt to design the continuous fractionationdevice in which we use the results of the first three studies together withthe analysis of the constraints imposed by the physical construction of thedevice. We explore the flow behavior in the exit region of the device andsuggest some guidelines to achieve successful fractionation accordingly.Results of this work show that spiral multi-layer Poisuille viscoplasticflow can be stable and the range of stability expands with increasing fluidsyield stress. Gravity current instability is evident in density unstable casesand the flow can be stabilized by increasing the fluids yield stress. KelvinHelmholtz instability was not found on conditions tested. Viscoplastic flowentrance length was found to be shorter than the equivalent Newtonian onefor the same range of Reynolds number.iiiPrefaceThree chapters of this thesis will be submitted for publication in refereedjournals. Following are the detail of the submission:1. Chapter three, Alshibl M. and Martinez D.M. Continuous fractiona-tion of particle suspensions in a viscoplastic fluid using spiral multi-layer Poiseuille flow. I solved the flow analytically, conducted the studyand wrote the paper under the supervision of Dr. Martinez.2. Chapter four, Alshibl M. and Martinez D.M. Study of entrance regionspiral multi-layer Poiseuille flow of viscoplastic fluid in annular geome-tries with density current effects. Dr. Martinez supervised the researchand I solved the flow numerically, validated the solution, planned andperformed the simulations.3. Chapter five, Alshibl M. and Martinez D.M. Study of the operationand efficiency of the continuous fractionation system of particle sus-pensions in a viscoplastic fluid using spiral multi-layer Poiseuille flow.ivPrefaceThis research has been supervised by Dr. Martinez. I designed the de-vice, solved the flow numerically, validated the solution, planned andperformed the simulations.vTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . xxiiiDedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xxiv1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52.1 The current pulp fractionation methods . . . . . . . . . . . . 52.2 Particles motion in Newtonian suspensions . . . . . . . . . . 8viTable of Contents2.3 Non-Newtonian fluids rheology . . . . . . . . . . . . . . . . . 142.4 Particle motion in viscoplastic fluids . . . . . . . . . . . . . . 172.5 The novel fractionation method . . . . . . . . . . . . . . . . 202.5.1 Fractionation principle . . . . . . . . . . . . . . . . . 202.5.2 Development of an efficient continuous fractionationprocess . . . . . . . . . . . . . . . . . . . . . . . . . . 232.5.3 Introduction to the academic problem . . . . . . . . . 253 The fully developed problem and particle fractionation . . 373.1 Method of solution . . . . . . . . . . . . . . . . . . . . . . . 443.2 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 453.3 Stability definition . . . . . . . . . . . . . . . . . . . . . . . . 463.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 474 The full annulus problem and density current effects . . . 564.1 Mathematical modeling . . . . . . . . . . . . . . . . . . . . . 574.2 Method of solution . . . . . . . . . . . . . . . . . . . . . . . 644.3 Grid independence test . . . . . . . . . . . . . . . . . . . . . 654.4 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 654.5 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 684.5.1 Entrance length . . . . . . . . . . . . . . . . . . . . . 68viiTable of Contents4.5.2 Tracking the interface location . . . . . . . . . . . . . 714.5.3 Axial velocity profile development . . . . . . . . . . 714.5.4 Kelvin Helmholtz instability . . . . . . . . . . . . . . 724.5.5 Density difference effects . . . . . . . . . . . . . . . . 765 The continuous fractionation device and fractionation effi-ciency . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 795.1 Mathematical modeling and method of solution . . . . . . . 805.2 Grid independence test . . . . . . . . . . . . . . . . . . . . . 815.3 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 825.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 835.4.1 Axial velocity profile development and entrance length 845.4.2 Outlet separator effects on interface . . . . . . . . . . 865.4.3 Flow rates design . . . . . . . . . . . . . . . . . . . . 885.4.4 Exit region effects on fractionation diagram . . . . . 896 Summary and conclusions . . . . . . . . . . . . . . . . . . . . 95Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97viiiTable of ContentsAppendicesA Fluid rheology and stability . . . . . . . . . . . . . . . . . . . 111B Fractionation computer program operation . . . . . . . . . 117C Parabolic equations versus elliptic partial differential equa-tions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119ixList of Tables3.1 A summary of the runs conditions simulated for each param-eter: inner and outer fluids yield stresses, plastic viscositiesand radius ratio. . . . . . . . . . . . . . . . . . . . . . . . . . 494.1 A summary of the runs conditions simulated for the entrancelength study. Simulations are for single fluid of Newtonian orviscoplastic behaviour. . . . . . . . . . . . . . . . . . . . . . 624.2 A summary of the runs conditions simulated for each study:Interface behaviour, Density current and Kelvin Helmoholtzinstabilities. . . . . . . . . . . . . . . . . . . . . . . . . . . . 635.1 A summary of the simulation runs conditions at: Rez =149.3, Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093. . . . . . . 825.2 fractionation case properties. . . . . . . . . . . . . . . . . . . 91xList of Figures1.1 A schematic of the geometry showing the fully developed re-gion, the annulus and the fractionation device. . . . . . . . . 12.1 Schematic of a pressure screen. . . . . . . . . . . . . . . . . . 62.2 Schematic of a hydrocyclone. . . . . . . . . . . . . . . . . . . 72.3 Schematic illustrating the yielded (region (1)-white) and un-yielded flow regions (region (2)-blue) around a spherical par-ticle in a yield stress fluid, according to the numerical resultsby Beris et al. (1985). This figure was reproduced from Putzet al. (2008). . . . . . . . . . . . . . . . . . . . . . . . . . . . 18xiList of Figures2.4 A demonstration of the fractionation of a bi-disperse suspen-sion of spherical particles. In (a), an image of the suspensionis given before the commencement of the centrifuge. (b) Thestate of the suspension after the application of the centrifugalforce. It should be noted that most of the darker particles areon the periphery of the centrifuge. This figure was reproducedfrom Madani et al. (2010a). . . . . . . . . . . . . . . . . . . . 212.5 A demonstration of the fractionation of a bidisperse suspen-sion of cylindrical particles. In (a) an image of the suspensionis given before the commencement of the centrifuge. (b) is thestate of the suspension after the application of the centrifugalforce. This figure was reproduced from Madani et al. (2010a). 222.6 A schematic cross section of the continuous fractionation device. 243.1 A schematic of the model geometry. . . . . . . . . . . . . . . 373.2 Dimensionless axial velocity profile comparison between presentwork and the work of Madani et al. (2013). . . . . . . . . . . 453.3 Multi-layer flow stability categories: (a) an unstable flow, (b)an inner stable flow and (c) an outer stable flow. The dashedlines are for the inner fluid velocity profiles. . . . . . . . . . . 48xiiList of Figures3.4 The solution of equations 3.22-3.23 for the conditions givenas series 1 in table 3.1. The colour map defines the stabilityof the flow state, as given in figure 3.3. The dashed linesrepresent the interface position and the solid lines representthe pressure drop rate. . . . . . . . . . . . . . . . . . . . . . . 523.5 Fractionation curves for spherical stainless steel particles indifferent yield stress fluids. The critical rotational speeds ωare denoted by the colored lines and the black lines for thecharacteristic speed Uc. Radius ratio, κ = 0.6 , outer radius,R = 0.127 m, plastic viscosity, µ = 1.0833 Pa.s, spheresdensity, ρs = 7800 kg/m3 and suspension density, ρf = 1000kg/m3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 533.6 Settling and residence rates in the mixing zone of the device. 54xiiiList of Figures3.7 Flow stability and operating diagram. Residence to settlingtime ratio in blue lines, interface locations in black dashedlines and dimensionless pressure drop rate in solid black lines.The annulus axial length L = 50 cm, inner fluid yield stressτ(1)y = 1 Pa, outer fluid yield stress τ(2)y = 0.5 Pa, radiusratio, κ = 0.6 , outer radius, R = 0.127 m, plastic viscosity,µ(1) = µ(2) = 0.75 Pa.s, spheres density, ρs = 7800 kg/m3and suspension density, ρf = 1000 kg/m3, ω = 37.2 rad/sand pressure drop rate of (50− 500) Pa.s. . . . . . . . . . . . 554.1 A schematic of the model geometry. . . . . . . . . . . . . . . 574.2 Full annulus model grid independence forRez=437.17, Reθ=8.7,At=0, B(1)=0.058 and B(2)=0.035 : (a) volume fraction dis-tribution: the red and blue colors denote the inner and outerfluids, respectively. (b) grid independence. . . . . . . . . . . . 66xivList of Figures4.3 Full annulus model validation with analytical solution: (a)volume fraction distribution for case 1, (b) volume fractiondistribution for case 2, (c) validation (case 1: Rez=437.17,Reθ=8.7, At=0, B(1)=0.058 and B(2)=0.035) and (d) vali-dation (case 2: for Rez=820.57, Reθ=8.7, At=0, B(1)=0.031and B(2)=0.019). The red and blue colors denote the innerand outer fluids, respectively. . . . . . . . . . . . . . . . . . . 674.4 Entrance length. a) Volume fraction distribution: the redand blue colors denote the inner and outer fluids, respectively.b) axial velocity contours with enlarged with enlarged viewsof velocity vectors near the wall and in mid-gap, c) Center-line velocity normalized with analytical fully developed valueversus axial location, d) Entrance length for Newtonian andviscoplastic fluids in an annulus of radius ratio of 0.8. . . . . 704.5 Axial velocity profile development: a) Axial velocity profilesat inlets, entrance region and fully developed region. b) Axialvelocity contours, Rez=437.17, Reθ=8.7, At=0, B(1)=0.058and B(2)=0.035. . . . . . . . . . . . . . . . . . . . . . . . . . 734.6 Axial velocity profiles, Rez=5248, Reθ=7.25, At=0.05, B(1)=B(2)=0and Ri=0.08. . . . . . . . . . . . . . . . . . . . . . . . . . . . 75xvList of Figures4.7 Density difference effect: Volume fraction distribution for At= 0.2 at different Reynolds and Bingham numbers. The redand blue colors denote the heavy and light fluids, respectively. 785.1 A schematic of the continuous fractionation device model ge-ometry. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 805.2 Fractionation device model grid independence forRez = 149.3,Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093 : (a) volume frac-tion distribution and (b) grid independence. . . . . . . . . . . 835.3 Fractionation device model validation with analytical solutionat Rez = 149.3, Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093.(a) volume fraction distribution. (b) validation. . . . . . . . . 845.4 Axial velocity profile development: (a) Volume fraction dis-tribution, (b) axial velocity profile development (before, in-side and after mixing zone). Simulated at Rez = 149.3,Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093. . . . . . . . . . . 865.5 Fractionation device volume fraction distribution showing theeffect of outlet separator on interface location, Rez = 149.3,Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093. . . . . . . . . . . 87xviList of Figures5.6 Interface behaviour at the exit region, Rez = 149.3, Reθ =7.88, B(1) = 0.155 and B(2) = 0.093 at: (a) ri = 0.84,Q1/Q2 = 2.98, (b) ri = 0.75, Q1/Q2 = 1.22, (c) ri = 0.72,Q1/Q2 = 0.96 and (d) ri = 0.69, Q1/Q2 = 0.75. . . . . . . . . 885.7 Settling and residence rates in the mixing zone of the device. 925.8 Fractionation diagram for the case detailed in table 5.2. Timeratios in blue lines, interface locations in black dashed linesand dimensionless pressure drop rate in solid black lines. . . . 935.9 Fractionation diagram after including the effect of outlet sep-arator. Case properties are per table 5.2. Time ratios in bluelines, interface locations in black dashed lines and dimension-less pressure drop rate in solid black lines. . . . . . . . . . . . 935.10 Fractionation diagram for an operating window criterion of40%. Case properties are per table 5.2. Time ratios in bluelines, interface locations in black dashed lines and dimension-less pressure drop rate in solid black lines. . . . . . . . . . . . 945.11 Fractionation diagram for an operating window criterion of30%. Case properties are per table 5.2. Time ratios in bluelines, interface locations in black dashed lines and dimension-less pressure drop rate in solid black lines. . . . . . . . . . . . 94xviiList of FiguresA.1 The solution of equations 3.22-3.23 for the conditions givenas series 2 in table 3.1. The colour map defines the stabilityof the flow state, as given in figure 3.3. The dashed linerepresents the interface position and the solid line representsthe pressure drop rate. . . . . . . . . . . . . . . . . . . . . . . 112A.2 The solution of equations 3.22-3.23 for the conditions givenas series 3 in table 3.1. The colour map defines the stabilityof the flow state, as given in figure 3.3. The dashed linerepresents the interface position and the solid line representsthe pressure drop rate. . . . . . . . . . . . . . . . . . . . . . . 113A.3 The solution of equations 3.22-3.23 for the conditions givenas series 4 in table 3.1. The colour map defines the stabilityof the flow state, as given in figure 3.3. The dashed linerepresents the interface position and the solid line representsthe pressure drop rate. . . . . . . . . . . . . . . . . . . . . . . 114A.4 The solution of equations 3.22-3.23 for the conditions givenas series 5 in table 3.1. The colour map defines the stabilityof the flow state, as given in figure 3.3. The dashed linerepresents the interface position and the solid line representsthe pressure drop rate. . . . . . . . . . . . . . . . . . . . . . . 115xviiiList of FiguresA.5 Maximum allowable flow rates Q in (l/s) for different stain-less steel spherical particles in a fluid having a yield stressτy = 1 Pa. The annulus axial length L = 50 cm, radius ra-tio, κ = 0.6 , outer radius, R = 0.127 m, plastic viscosity,µ = 1.0833 Pa.s, spheres density, ρs = 7800 kg/m3 and sus-pension density, ρf = 1000 kg/m3: a) fractionation curve, b)maximum allowable flow rates. . . . . . . . . . . . . . . . . . 116xixNomenclatureRoman lettersGˆ Pressure drop rate Pa/mpˆ Pressure PaQˆ Flow rate kg/m3Rˆ Annulus outer radius mrˆ Radial position mRˆi Annulus inner radius muˆ Velocity m/sUˆo Characteristic velocity m/sAt Atwood numberB Bingham numberxxNomenclatureDh Hydraulic diameter mF Applied force NLe Entrance length mRe Reynolds numberRi Richardson numberSubscriptsθ Tangentialc Criticalr Radialz AxialGreek lettersτ Shear stress Paτy Yield stress Paα Concentrationγ Shear rate s−1ωˆ Rotational speed rad/sxxiNomenclatureκ Annulus radius ratioµ Dynamic viscosity pa.sρ Density kg/m3xxiiAcknowledgementsAcknowledgment is due to the College of Engineering of King Saud Univer-sity for supporting this research.I would like to express my sincere gratitude to my advisor Prof. MarkMartinez for the continuous support of my Ph.D study and related research,for his patience, motivation, and immense knowledge. His guidance helpedme in all the time of research and writing of this thesis. I could not haveimagined having a better advisor for my Ph.D study.Besides my advisor, I would like to thank the rest of my thesis committee:Dr. Kendal Bushe, and Dr. Heather Trajano, for their insightful commentsand encouragement, but also for the hard question which incented me towiden my research from various perspectives.Last but not the least, I would like to thank my parents and my familyfor supporting me spiritually throughout writing this thesis and my life ingeneral.xxiiiDedicationxxivChapter 1IntroductionThe focus of the present work is the study of laminar multi-layer viscoplasticflow in annular geometries as shown in figure 1.1. This type of flow is usedin the novel fractionation process of particle suspensions that was advancedby Madani et al. (2010b). They proved this new fractionation methodologyon a batch wise basis. The purpose of this work is to extend the new methodapplicability to continuous fractionation for use in industrial applications.Figure 1.1: A schematic of the geometry showing the fully developed region,the annulus and the fractionation device.The motivation for the present work stems from an interest in the frac-1Chapter 1. Introductiontionation of particle suspensions especially with papermaking fibre suspen-sions. Fractionation is a vital part in papermaking process. It allows theprocessing of low quality fractions of pulp suspensions to improve paperquality and reduce the consumption of the required chemicals, energy orassets.Additionally, it enables papermakers to produce products over a widerange of properties that can be used in different applications. To elaboratemore, if we fractionate a mixture of long and short pulp fiber suspension, weare going to get a short fiber fraction that is known for having better print-ing capabilities and a long fiber fraction that has a higher tensile strength.In addition to that, we can blend these fractions at different ratios to getdifferent paper products with a variety of properties that can be used indifferent paper applications.We conducted four different but complementary studies to facilitate thismulti-layer flow in the continuous fractionation process, these studies arecarried out in order to design an industially applicable process. The firststudy deals with the fully developed condition of the flow. This conditionallows us to solve the problem analytically and hence provide a reliablereference to validate the remaining numerical studies. This solution alsoprovides a means to test the stability of the flow.2Chapter 1. IntroductionThe second study is related to the fractionation of particles in particlessuspensions utilizing the solution demonstrated in the first study. We de-velop particle fractionation curves of particles of different sizes in fluids withdifferent rheology. We develop a code to simulate thousands of flow cases(a flow case has a unique combination of streams flowrates) that we mayhave for fluids of known rheology and for particles of different sizes. We pre-dict the fractionation operating window (some range of streams flowrates)needed for successful fractionation.In the third study, we examine the flow in a full annulus geometry in-cluding the entrance region of the flow. Here, we estimate the flow entrancelength in order to design the length of the mixing zone of the two streams inthe continuous fractionation device. We study the effect of Kelvin-Helmholtzand density current instabilities on the flow.In the fourth study, we attempt to design the continuous fractionationdevice in which we use the results of the first three studies together withthe analysis of the constraints imposed by the physical construction of thedevice. We explore the flow behavior in the exit region of the device andsuggest some guidelines to achieve successful fractionation accordingly.This thesis is presented in 6 chapters. In chapter 1 we present the mo-tivation of this work. Chapter 2 gives a background on the fractionation3Chapter 1. Introductionmethods, fractionation efficiency, presents the new fractionation principleand the potential flow characteristics to be used in this new methodology.The fully developed flow problem and particle fractionation studies are pre-sented in chapter 3. We present the full annulus problem solution anddensity current analysis in chapter 4 and the continuous fractionation de-vice study in chapter 5. We end the work by conclusions in chapter 6 andbibliography.4Chapter 2BackgroundIn this chapter we will shed light on the available fractionation methods,their efficiency and the need for a new fractionation technique. We will dis-cuss in more detail particle motion with both Newtonian and non-Newtoniansuspensions as they are the core of the fractionation process. The proposedcontinuous fractionation device that employs the new technique is exploredin terms of identifying the nature of the flow to be used.2.1 The current pulp fractionation methodsCurrently, there are two methods to fractionate pulp fibers: pressure screens(figure 2.1) and hydro-cyclones (figure 2.2). In pressure screens, the fibersuspension enters an annular gap between a rotor and an outer cylinder withsmall openings (either slots or holes), separation is then achieved based onthe length of the fiber as small fibers pass through the openings. Pressurescreens usually are equipped with different types of rotors that continuously52.1. The current pulp fractionation methodsdisperse fibers to prevent them from accumulating on the screen surface. Anumber of studies discuss pressure screens design and operation, for example,Sloane (2000) and Julien Saint Amand and Perrin (1999) gave a generalreview of pressure screening systems, their design and efficiency.Figure 2.1: Schematic of a pressure screen.On the other hand, hydrocyclones separate fibers based on their specificsurface or surface area per gram, resulting in separation of large diameterthin walled fibers from small diameter thick walled ones. Hydrocyclonesare conical or partly cylindrical devices with no moving parts. Here thepulp suspension flows into the hydrocyclone with a tangent angle creatinga swirling flow. Due to centrifugal force the dense materials are pushed to-ward the outside and collected at the bottom of the hydrocyclone as rejects.62.1. The current pulp fractionation methodsFigure 2.2: Schematic of a hydrocyclone.Paavilainen (1992) studied experimentally the hydrocyclone operation andshowed that they are indeed able to fractionate fibers based on its specificsurface.Despite the fact that pressure screens and hydro-cyclones are consideredthe industrially approved methods, these methods have relatively low ef-ficiencies. The flow within these devices is complex due to the stochasticnature of turbulence, the presence of boundaries, and the long range hydro-dynamic interactions between suspended fibers, Madani et al. (2010b). All72.2. Particles motion in Newtonian suspensionsthese facts have a negative impact on the efficiency of both pressure screensand hydro-cyclones.To address the problem of low efficiency of these devices, we will explorefractionation efficiency in more detail. Fractionation efficiency is built uponan understanding of the motion of different classes of particles in a flowingNewtonian suspension. In the next section we briefly review a number ofstudies on particle motion for the purpose of understanding why the effi-ciency in these devices is low. Afterwards, we discuss particle behaviour inviscoplastic fluids and how this behaviour makes utilizing viscoplastic fluidsin fractionation a wise decision.The novel fractionation principle of Madani et al. (2010b) is then intro-duced with some preliminary results from batchwise experiments. Followingthat, the continuous fractionation method proposed by Madani et al. (2010b)using viscoplastic fluids is then introduced. And we finalize this chapter withan introduction to the academic problem.2.2 Particles motion in Newtonian suspensionsThe motion of particles in Newtonian suspensions is extremely complex.Studies by Batchelor (1972), and more recently by Mackaplow and Shaqfeh(1998) were concerned with the simplest case of settling of an isolated rod82.2. Particles motion in Newtonian suspensionsin an unbounded fluid experiencing stokes flow described by the followingequation.Vsed =∆ρd216µ[(ln 2r+0.193+O(ln 2r)−1)g+(ln 2r−1, 807+O(ln 2r)−1)(p.g)p](2.1)where ∆ρ is the difference in density between the fibre density ρ and thesurrounding fluid density ρf ; r is the aspect ratio of the fibre defined byl/d where d is the diameter of the fibre ; µ is the viscosity of the fluid ; gis the acceleration due to gravity ; and p is the unit vector that indicatesfibre orientation.They found that rods -based on their orientation- are subjected to lateralmovement while settling vertically under the effect of gravity with the driftvelocity strongly dependent on fiber orientation, this is not the case forspheres.On the other hand Jayaweera and Mason (1965) reported that orienta-tion does not affect settling at low Re, (Re < 0.01) and cylinders fall withthe same altitude that they were released with. We conclude that there isno unique terminal velocity for settling cylinders, and the terminal velocitydepends on the orientation and aspect ratio of the cylinder. This differencein terminal velocities may give a basis for separating particles mechanically.92.2. Particles motion in Newtonian suspensionsIn more realistic situations, i.e., higher concentration and Reynolds num-bers, long range disturbances are introduced to the system resulting in adistribution for the settling velocity. Disturbances are created by presenceof other particles and/or increasing the Reynolds number.Let us start looking at the effect of the ‘wake of a moving particle in thevicinity of another. In this regard, Happel and Brenner (1965) explainedthat the physics becomes more complex as each individual fibre settles androtates under the influence of the ‘wakes’ or ‘long-range hydrodynamic dis-turbances’ of the other settling particles. This in turn leads to inhomoge-neous settling rates and local floc formation. Floc formation was observedalso by Kumar and Ramarao (1991) in monodisperse glass fibre suspensions.In an interesting finding, Herzhaft and Guazzelli (1999) reported that inthe dilute regime, the ensemble-averaged settling velocity actually increaseswith concentration and may exceed the velocity of an isolated particle. Theyfound also that most of the fibers are aligned in the direction of gravity underthese flow conditions.Some numerical simulations conducted by (Mackaplow and Shaqfeh (1998);Butler and Shaqfeh (2002); Koch and Shaqfeh (1989)) for the flow un-der these conditions and in the limit of Re = 0, using slender body the-ory, suggested that the settling suspension should segregate into particle102.2. Particles motion in Newtonian suspensionsclumps. The simulations were in good agreement with experimental resultsof Herzhaft and Guazzelli (1999).For higher Reynolds numbers which are more relevant to the paper-making industry, Marton and Robie (1969) studied fibres settling in water.It has been shown that, unlike in the Stokes’ regime, at Reynolds num-bers Re ∼ O(1) isolated non-spherical particles tend to exhibit preferentialorientation during settling. The same conclusions were drawn by other in-vestigations (Jayaweera and Mason (1965); Feng et al. (1994); Jianzhonget al. (2003)) when studying flows at these Reynolds numbers, in the dilutelimit; they have shown that isolated non-spherical particles tend to exhibitpreferential orientation during settling.The torque induced on thin cylinders causes the body to rotate into astable position with its symmetry axis aligned horizontally. In a relativelyrecent work Holm et al. (2004) showed that at elevated concentrations, thelong range hydrodynamic interactions perturb the flow field to the pointwhere recirculation or swirling like structures are apparent.In a later work Salmela et al. (2007) conducted a more comprehensivestudy on the motion of settling particles as a function of concentration,aspect ratio, fluid viscosity, and fibre length for both monodisperse and bi-disperse suspensions and used a similar index of refraction matching tech-112.2. Particles motion in Newtonian suspensionsnique. They tracked the motion of tracer fibres in three-dimensions whichallowed them to study fibre fractionation in this equipment.They reported a non-monotonic behaviour in the average settling velocityas a function of concentration with the maximum velocity occurring at avolume fraction of 0.05 % similar to that reported by Herzhaft and Guazzelli(1999). Holm et al. (2004) also observed non-monotonic behaviour with themaximum in the initial settling speed occurring at a crowding number of N∼ 16 for papermaking fibres.In terms of orientation behaviour, Salmela et al. (2007) showed com-plex orientation behaviour for suspensions under these conditions. For lowconcentration they showed that the suspension fibers tend to orient in thehorizontal state and then tend to align in the vertical state as concentrationincreases.One vital conclusion can be drawn from some unpublished results ofSalmela and his coworkers that gives a rise to an understanding of the pos-sibility to fractionate based upon settling. Although not reported in Salmelaet al. (2007), they studied the motion of a bi-dispersed suspension comprisedof glass rods of two different aspect ratios r=23 and 50 and measured theensemble-averaged settling velocity of each class of particle.They measured the settling velocity as a function of the initial suspension122.2. Particles motion in Newtonian suspensionsconcentration. They found that under dilute conditions each fibre fractionsettles at statistically different velocities. This means that there is a possibil-ity for fractionation in this regime. However, with increasing concentration,the differences between the terminal velocities diminish until their differenceis indistinguishable and consequently deter the possibility of fractionationin this regime.In summary, from the above literature, we found that the settling velocityof isolated particles in an unbounded fluid is dependent on density, diameter,orientation and aspect ratio of the cylindrical particles. In the creeping flowregime (Re 1) separation of fibers based on difference in settling velocitiescan be achieved only in extreme dilute suspensions and on condition thatthe orientation distribution remains constant during descent.However, at higher concentration for all Reynolds number ranges separa-tion of particles cannot be achieved due to the long range hydrodynamic in-teractions between particle that leads to floc formation and chaotic swirlingstructures. Salmela et al. (2007) reported that there is no statistical differ-ence between the settling velocities of particles in bi-disperse suspensions inany reasonable concentrations.We conclude from above review that particles motion in suspensionswith concentration Reynolds number relevant to the industrial suspensions132.3. Non-Newtonian fluids rheologyis extremely complex. Chaotic behaviour is evident due to long range hy-drodynamic interactions between particles. This behaviour is observed inthe ideal lab situations, when it comes to the real life industrial equipmentlike cyclones and pressure screens the flow is event more complex beingthree-dimensional, time-dependent and turbulent in nature.Now before we move on to discuss particles motion in viscoplastic fluids(where the novel fractionation method is all about) it is a good idea tointroduce the difference between Newtonian and non-Newtonian fluids andwhere viscoplastic fluids stand.2.3 Non-Newtonian fluids rheologyFluids are divided to Newtonian and non-Newtonian based on their behaviorin the basic shear diagram (shear stress versus shear rate). The Newtonianfluids exhibit a linear relation between shear stress and shear rate with theline passing through the origin.τ = µγ (2.2)Where τ is the shear stress, µ is a constant value of dynamic viscosityand γ is the shear rate. A typical example of a Newtonian fluid is water. Allother types of fluids are non-Newtonian, meaning that either the relation142.3. Non-Newtonian fluids rheologybetween shear stress and shear rate is not linear or the line does not passthrough the origin or the material exhibits time-dependent behavior.Non-Newtonian fluids are divided into several types: shear-thinningfluids, shear-thickening fluids, yield stress (viscoplastic) fluids and time-dependent fluids.With the shear-thinning fluids, the shear stress - shear rate curve passesthrough the origin and is concave downwards, with an increase in shear rateresulting in a less than proportional increase in shear stress or in other wordsthe dynamic viscosity is not constant.For the shear-thickening fluids, the shear stress - shear rate curve alsopasses through the origin but is concave upwards, meaning that, an increasein shear rate results in a more than proportional increase in shear stresswith a variable dynamic viscosity.For time-dependent fluids, the shear stress versus shear rate data ob-tained in ascending order of shear rate are different from those obtained indescending order. If the latter is lower, the fluid is called time-dependentshear-thinning (thixotropic), while if the data obtained in descending orderis higher, the fluid is called time-dependent shear-thickening (antithixotropic).In yield stress fluids, the flow may not commence until a threshold valueof stress (τy, yield stress) is exceeded. A yield stress fluid with shear thinning152.3. Non-Newtonian fluids rheologybehavior is called Herschel-Bulkley fluid, while a yield stress fluid with alinear shear-stress shear rate relation is called Bingham plastic fluid.To be consistant with related previous studies in the literature, we uti-lize the Bingham plastic model to describe viscoplastic fluids rheology. Therelevant mathematical formulations of Bingham plastic fluids are the consti-tutive equation and the dimensionless Bingham number defined as follow:τ = τy + µpγ (2.3)B =τyLµpVo(2.4)where τ is the shear stress, γ the shear rate, µp the plastic viscosity, Bthe Bingham number, L the length scale and Vo the velocity scale. If τ < τythe Bingham plastic fluid behaves as a solid, otherwise it behaves as a fluid.If τy = 0, this model reduces to the Newtonian fluid.Now since the novel fractionation technique utilizes yield stress fluids asthe carrying medium, we need to know how particles behave or move in sucha medium. The next section will deal with the motion of particles in yieldstress fluids.162.4. Particle motion in viscoplastic fluids2.4 Particle motion in viscoplastic fluidsParticle motion in viscoplastic fluid is more complex than that in Newtonianfluid. If we consider the simplest case of one particle in an unboundedviscoplastic fluid, here the particle has to overcome the force created bythe yield stress to start moving. If we denote the net force applied on theparticle by Fa , then there should be a critical force value that is requiredto initiate motion (overcome the fluid yield stress resistive force).To elaborate more, let us consider an example of a particle immersed in amaterial having a yield stress τy. If Fa is not sufficient to overcome the yieldstress then the particle will be suspended indefinitely. Motion will commencewhen Fa > τyAe, where Ae is the surface area over which the force due tothe yield stress is applied. The exact value for Ae however remains an openquestion. The above relationship can be made dimensionless by scaling eachside by a characteristic area of the particle. In this case if we divide eachside of the equation by D2, where D is the diameter of thr particle, we candefine a dimensionless force ratio F which represents the criteria for motion:F =FaτyD2>AeD2(2.5)The bound where motion begins is defined as the critical force ratio Fc.172.4. Particle motion in viscoplastic fluidsThe particle movement depends solely on determination of Ae. Several stud-ies were conducted by researchers to determine Ae, Andres (1961) reportedthat for spherical particles settling will start when Fc = 4.8. Beris et al.(1985) on the other hand suggested that there are two yield surfaces aroundthe spherical particle; a surface with kidney shape away from the sphereand two triangular shapes attached to the leading and trailing edges of thesphere with critical force value of Fc ∼ 11 to start motion, figure 2.3 .Figure 2.3: Schematic illustrating the yielded (region (1)-white) and un-yielded flow regions (region (2)-blue) around a spherical particle in a yieldstress fluid, according to the numerical results by Beris et al. (1985). Thisfigure was reproduced from Putz et al. (2008).Experimentally, Tabuteau et al. (2007); Jossic and Magnin (2001); Lax-ton and Berg (2005) reported the critical force value for motion in the range16 < Fc < 25 and Madani et al. (2010a) in the range 13 < Fc < 29. In182.4. Particle motion in viscoplastic fluidsthe same work Madani et al. (2010a) reported Fc for cylinders and bentcylindrical rods oriented in parallel and perpendicul to the direction of theapplied centrifugal force in the range of 28 < Fc < 517 and 87 < Fc < 423respectively.Chhabra (2006) gave a comprehensive summary of settling and sedi-mentation in different fluids including viscoplastic fluids, he concentratedon terminal velocity and drag coefficients values that are used usually in de-sign purposes for engineering applications. Putz et al. (2008) experimentallyestimated the yield surface shape for settling spheres in carbopol solutionsby ovoid spheroid that has a major axis dimension of more than 5 times theradius of the sphere.To summarize, what is clear from this body of work is that during settlingthe flow is confined in the vicinity of the particle within an envelope the sizeof which is related to the yield stress of the material. For particles to settle,a critical force must be applied to overcome the resistance created by theyield stress. Relatively little is known about the shape of this surface andthe magnitude of the applied force to create motion. What can be said isthat the unyielded envelope is larger than the body itself but its shape hasyet to be determined rigorously.With this understanding of particle behaviour in suspensions we are192.5. The novel fractionation methodready now to introduce the novel fractionation technique developed by Madaniet al. (2010b) in the following section.2.5 The novel fractionation methodIn this section, we explain the new particle suspension fractionation methodthrough the work of Madani et al. (2010b). We explore the possible contin-uous fractionation equipment and the flow associated with it. We introducethe academic problem and state our aim and objectives in this project.2.5.1 Fractionation principleBy understanding motion of particles in yield stress fluids, we can see thatthe differences in particles static stability in yield stress fluids when sub-jected to one applied force make a novel criteria for separation of particles.It is informative at this point to present Madani et al. (2010a,b, 2011)’s re-sults in this regard as they considered two different spherical particles groupsof equal dimensions but of different densities.They rotated each particle in a centrifuge at exactly the same angularvelocity and at equal radial positions from the axis of rotation (both par-ticle types experienced the same acceleration field). Particles with greaterdensity, experienced a larger centrifugal force as they have larger mass.202.5. The novel fractionation methodSeparation of these two particles types occurs when the centrifugal forceapplied on type 1 particles, (F1) is greater than the critical force (τyAe)and at the same time, this critical force is greater than the centrifugal forceapplied on type 2 particles, (F2). Or in other words F1 > τyAe > F2.Figure 2.4 shows the images of suspensions before and after rotation ofthe centrifuge, the spheres with the higher density are shown in black andthat with the smaller density in red. Both sphere types are of same numberand distributed evenly throughout the centrifuge.Figure 2.4: A demonstration of the fractionation of a bi-disperse suspensionof spherical particles. In (a), an image of the suspension is given before thecommencement of the centrifuge. (b) The state of the suspension after theapplication of the centrifugal force. It should be noted that most of thedarker particles are on the periphery of the centrifuge. This figure wasreproduced from Madani et al. (2010a).212.5. The novel fractionation methodFigure 2.5: A demonstration of the fractionation of a bidisperse suspensionof cylindrical particles. In (a) an image of the suspension is given beforethe commencement of the centrifuge. (b) is the state of the suspension afterthe application of the centrifugal force. This figure was reproduced fromMadani et al. (2010a).They chose an appropriate rotation rate that led to separation (thelighter particles were statically stable in the suspension while the heavierparticles moved to the periphery of the centrifuge).Similar treatment was done to cylindrical rods as shown in figure 2.5.Here the heavier (type1) rods migrated towards the periphery of the cen-trifuge under the effect of the centrifugal force while lighter (type2) rodsstayed stable. Madani et al. (2011) performed more experiments on pa-permaking fibre suspensions demonstrating the utility of this approach andhas shown that indeed separation may proceed based upon fibre length or222.5. The novel fractionation methodcoarseness; their work was however conducted on a batchwise process. Thepurpose of this work is to advance this technology by demonstrating itspotential in an industrially relevant suspension, namely a papermaking sus-pension.2.5.2 Development of an efficient continuous fractionationprocessMadani et al. (2010b) has introduced a novel continuous fractionation methodbased upon the principles described earlier. They consider a flow field inwhich the particles translate in the axial direction which is perpendicular tothe motion induced by the centrifugal force. The axial motion is decoupledfrom the radial (centrifugal) motion. In addition, they emphasise that asignificant portion of the flow to be unyielded. To acheive this type of flow,they consider a pressure driven flow imposed onto solid body rotation.A schematic of the continuous device that uses this method is shown inFigure 2.6. Here particles, suspended in a viscoplastic fluid, are introducednear the inner cylinder (inlet 1) and are transported axially towards theexit; a similar viscoplastic fluid is introduced near the outer cylinder (inlet2) but no particles are present in this fluid. The walls of each cylinder rotateat the same angular speed and in the same direction. The size of the un-232.5. The novel fractionation methodyielded portion of the flow is controlled by the pressure drop only allowingthe operator to control separation and production rates independently. Thedevice has two exits to collect the sorted particles streams.Figure 2.6: A schematic cross section of the continuous fractionation device.The entire device rotates and creates a centrifugal field during transport.As described earlier, if an appropriate centrifugal force is selected, a certainfraction of particles will migrate towards the periphery and exit near theouter cylinder. Hence fractionation may occur. To date, the device has yetto be built as more information is required about the details of the flow field.242.5. The novel fractionation methodSo basically we have two viscoplastic layers driven by axial pressure dropand rotational angular velocity (solid body rotation), such a flow is referredto as a multi-layer spiral Poiseuille viscoplastic flow.It is the purpose of this work to develop the design criteria for a stablerotating multi-layer flow. A review of this flow stability will follow in thenext section.2.5.3 Introduction to the academic problemThe flow under consideration is a multi-layer spiral Poiseuille viscoplasticflow in annular geometry that is utilized in the novel fractionation equip-ment. The most important scientific question in designing such equipmentis its flow stability.In this work we study the stability of the flow using the Bingham plasticmodel in order to understand the instability bounds due to the presence ofdifferent fluids yield stresses, centrifugal forces resulting from rotation andpresence of rigid walls containing the flow.In the Bingham plastic constitutive model the flow may not commenceuntil a threshold value of stress (yield stress) is exceeded. Also, the fluidin this model has a linear shear-stress shear rate relation. So basically thismulti-layer flow may have both yielded and un-yielded regions depending on252.5. The novel fractionation methodthe local shear stress values compared to each fluid unique yield stress value.Unyielded regions behave like rigid bodies undergoing linear or rotationalmotion.We can think of three sources for instability in Poiseuille multi-layerviscoplastic flow. First, is the instability caused by a yielded interface in theevent of having yielded velocity profiles of both fluids sharing the interface.Second is the instability caused by Kelvin Helmholtz instability since wehave parallel multi-layer flow with different layer velocities.The third instability may be caused by solid body rotation in the case ofdensity unstable multi-layer flow (having a heavier inner fluid that is closerto the axis of rotation and a lighter outer fluid that is away from the axis ofrotation).To adress the questions associated with these instabilities and the be-haviour of the multi-layer annular flow in both the fully developed and entryregions, we review the related literature in four categories: i) non-Newtonianfluid flows in ducts and annular gaps, ii) entrance length, iii) interfacial sta-bility and iv) density current instability.Now let us start with the first catagory. Solutions of the steady non-Newtonian fluid flow in ducts and annular gaps do exist in the literature foriso-dense cases. Details can be gained from the work of Fredrickson and Bird262.5. The novel fractionation method(1958); Shul’Man (1970); Grinchik and Kim (1972); Anshus (1974); Hanks(1979); Papanastasiou (1987); Fordham et al. (1991) . Some researchers likeBittleston et al. (2002) and Pelipenko and Frigaard (2004) were interested indisplacement viscoplastic flows in annular gaps with emphasis on eccentricityand annulus inclination.When it comes to introducing wall rotation, Bittleston and Hassager(1992) were one of the first researchers to study viscoplastic fluid flow inducts with rotation. They investigated the flow in a concentric annulus whenthere is both axial and tangential flows. The tangent flow arises from therotation of the inner cylinder of the annulus. They used the Bingham plasticmodel in their analysis. They solved the flow analytically by considering theflow in a slot and solved the full problem numerically. They found the criticalrate at which the inner wall has to move to reduce the plug size to zero.One of the few research groups that performed an analytical analysis inaddition to Bittleston and Hassager (1992) were Liu and Zhu (2010). Theystudied the axial Couette Poiseuille flow of Bingham fluids through concen-tric annuli and presented eight different forms of velocity profiles dependingon values of three dimensionless parameters: the Bingham, axial Couettenumbers and the radius ratio.Some researchers used the Herschel-bulkley model in their analysis. For272.5. The novel fractionation methodexample Nouar (1998) in his work reported the critical wall movement rateto reduce the plug size to zero and Escudier et al. (2002) who used the finite-volume method to solve the flow for the fully developed laminar flow throughan eccentric annulus with inner cylinder rotation. Also, Meuric et al. (1998)numerically solved the laminar flow in vertical concentric and eccentric an-nuli for an imposed axial flow with inner wall rotation. They showed thatthe inclusion of rotational effects, for a fixed pressure gradient, is likely toincrease the axial volumetric flow rate over non-rotating situations in con-centric geometries. However, in eccentric annuli, the situation is reversedand the flowrate gradually decreases as the rotation rate is increased.Some of the researchers like Wang focused on tracking the yield surfacesand mobile plug zones. In his work,Wang (1997) used the finite elementmethod to study Bingham fluid flow in eccentric annuli and in an L-shapedducts. He found the plug regions and tracked the yield surfaces in thesegeometries. In a later work,Wang (1998) studied viscoplastic fluid flow in asquare duct using a finite-element method and tracked the yield surface andmobile plug zones.When it comes to studying the linear stability of the flow, insight can begained from the work of Peng and Zhu (2004) and Madani et al. (2013). Pengand Zhu (2004) studied Bingham-plastic fluid flow between two concentric282.5. The novel fractionation methodcylinders rotating independently and with axial sliding of the inner cylinder(spiral Couette flow). They found that islands of instability, which are foundin the spiral Couette flow of Newtonian fluids, may not exist owing to theeffect of yield stress. Also they reported that both the rotation of the outercylinder and a decrease of the gap between the cylinders have stabilizingeffects.Madani et al. (2013) investigated the linear stability of both Newtonianand Bingham fluids in spiral Poiseuille flow in the annular gap betweentwo co-rotating cylinders. They found that solid body rotation increasesthe margin of stability for both Newtonian and Bingham fluid flow cases.Additionally, the flow is linearly stable for all B > 0 where B is Binghamnumber.In terms of the shape of axial velocity profiles they found that the steady,fully-developed spiral Poiseuille flow consists of an unyielded region in thecenter of the channel, for finite B (Bingham number) , bounded by twoyielded regions. The position of the yield surfaces was found as part of thesolution methodology reported by Liu and Zhu (2010) and is dependent onlyon B; the swirl component does not affect the position of the plug.What is clear from this body of literature is that besides previous worksthere is few, if any, studies we could find in the literature addressing multi-292.5. The novel fractionation methodlayer flows in swirl, it is a relatively unexplored area. In related work, theliterature is both substantial and protracted for single-fluid viscoplastic flowin annular gaps or in displacement flow, i.e. one fluid contacting anotherfluid on a plane perpendicular to the axis of motion.We continue the discussion by considering the entry region solutions inthe literature (second category). Perhaps among the first related studies forNewtonian fluid flows are the work of Heaton et al. (1964) who reported thehydrodynamic entry length for different annular geometries. Shortly afterthat McComas (1967) reported the laminar hydrodynamic entrance lengthfor annular geometries. Shah (1978) explored more geometries includingcircular and noncircular ducts, parallel plates, rectangular, equilateral tri-angular, and concentric annular ducts. They proposed a correlation forlaminar hydrodynamic entry length for these geometries. After that, Feld-man et al. (1982) extended the research to eccentric annular geometries andestimated numerically the hydrodynamic entrance length.In the most recent studies on Newtonian flows, Durst et al. (2005) es-timated the entrance length for pipe flow of both laminar and creepingregimes. Poole and Ridley (2007) modified the Newtonian entrance lengthcorrelation proposed by Durst et al. (2005) to be applicable to power lawfluids. However, for annular geometries, entrance length was reported only302.5. The novel fractionation methodfor higher ranges of laminar regime by researchers like Maia and Gasparetto(2003).For viscoplastic fluid flows some researchers like (Chen et al. (1970);Shah and Soto (1975); Nowak and Gajdeczko (1983); Vradis et al. (1993);Min et al. (1997)) proposed correlations for the pipe flow enterance lengthas a constant number times Reynolds number or (Le/Dh = CRe) in thelaminar regime ignoring the lower range of Reynolds number flows wherethe entrance length collapse to zero based on their correlations. Here, Leand Dh are the entrance length and hydraulic diameter respectively. Thisprobelm was addressed by Ookawara et al. (2000) and more recently by Pooleand Chhabra (2010) who took into consideration the effect of the diffusiondominated (low Reynolds number) pipe flow in their correlations.The general outcome of the literature in this subject is that the entrancelength for annular geometry Newtonian flows is usually shorter than thecorresponding ones of pipe flows for the same values of Reynolds numbers.Also, there is little (if any) studies on entry region flows in lower rangesof laminar and creeping regimes of annular flows for both Newtonian andnon-Newtonian fluids.Let us now consider the third cateogry of this literature review, theinterfacial stability caused by a yielded velocity profile at the interface. A312.5. The novel fractionation methodconsiderable amount of literature studies the interfacial instability of multi-layer viscoelastic flows, in contrast, a few researchers like Frigaard (2001),Moyers-Gonzalez et al. (2004) and Huen et al. (2007) studied this instabilityin multi-layer viscoplastic flows. They reached to a conclusion that a stablemulti-layer viscoplastic flow is possible if we are able to preserve an unyieldedregion at the interface. To do that, the applied shear stress at the interfacehas to be less than that of the yield stress of one of the fluids that are sharingthis interface. So, effectively the velocity profiles are flat for at least one ofthe fluids sharing the interface.Now we turn our attention to the last category of the literature review,the problem of density difference between layers. A density current or some-times called gravity current is usually a horizontal flow in some gravitationalfield that is driven by a density difference. Due to its wide existence in na-ture, gravity currents were intensively investigated by researchers.Let us start with studies involving Newtonian fluids. Benjamin (1968)was one of the first people to study gravity currents mathematically. Hisstudies involved flows of moving heavy fluids with stagnant lighter fluids.Before that, Abbott (1961) and Barr (1967) reported experimental resultsabout intrusions and collisions between fluids of differing density. In aneffort to do a comprehensive review of the literature, Simpson published322.5. The novel fractionation methodsome work (Simpson (1982), Simpson (1999)) summarizing the research inthe field of gravity currents with examples in the atmosphere and in theocean.For non-Newtonian flows, the work is substantial and prolonged forpower law fluid flows. To mention some, Di Federico et al. (2012), Chowd-hury and Testik (2012), Di Federico et al. (2006b), Di Federico et al. (2006a),Longo et al. (2013), Vola et al. (2004), Gratton et al. (1999), Pascal (2003),Longo et al. (2015), Longo et al. (2013) and Ciriello et al. (2015). Theirwork was related to different studies like gravity currents in a porous layer,creeping gravity currents and propagation of free-surface channelized viscousgravity currents. They consider different geometries like horizontal rectan-gular tanks, horizontal rigid plane and horizontal rectilinear channels forboth steady and unsteady flows.It is clear that most of the gravity current literature is related to power-law fluids in spreading flows of one fluid over another mostly in horizontalplanes or rectilinear channels, etc. To the best of author’s knowledge, theproblem of density current of two viscoplastic fluids in an annulus with solidbody rotation has not been explored previously.To summarize, studies in the literature are few, if any, for: multi-layerviscoplastic flows with swirl, or annular flows entrance length in lower ranges332.5. The novel fractionation methodof laminar/creeping regimes for both Newtonian and non-Newtonian fluids,and, density current instability in bounded multi-layer flows especially an-nular flows. In related work the literature is substantial for: single fluidviscoplastic flows/displacement flows, higher ranges of laminar pipe and an-nular flows, and, density current studies of power law fluids spreading flowsin planes or rectilinear channels. For interfacial stability, a stable multi-layerviscoplastic flow is possible if we maintain at least one velocity profile of thefluids sharing the interface to be flat.Based on this summary of the literature, our aim is to investigate themulti-layer spiral Poiseuille viscoplastic fluid flow in annuli. We vary therheology of the fluid spatially in cross stream direction in order to allow usstudy the layered annular flows. We explore the entrance length for bothNewtonian and viscoplastic fluids in lower ranges of laminar and creepingregime. In addition we explore the effect of density difference between layers.At this point let us define the objectives of this research. Our main goalis to design a continuous fractionation equipment that employs the novelfractionation method. To acheive this goal, we need to fulfill the followingobjectives:• To understand the stability of a multi-layer annular flows in whichthe viscoplastic fluids are of different rheology. Here we will analyze342.5. The novel fractionation methodthe stability of two fluids in contact, in swirling Poiseuille flow. Wesolve analytically the fully developed problem of the Spiral Poiseuillemulti-layer flow for stability of the flow in the interface area.• We will examine computationally (in a labrotary device using 2D ax-isymmetric CFD simulation) the full annulus problem including esti-mating the entrance length of the flow to predict the required axiallength of the fractionation device. In this simulation we will examinethe effect of the design and elongational stresses in the entry region,on final interface position. In addition, we will explore the effect ofdensity current and Kelvin-Helmholtz instabilities on the flow stabil-ity. We will benchmark this computational solution with the fullydeveloped analytical solution.• Study computationally the continuous fractionation device includingthe inlet and exit conditions. We will simulate the flow field in alaboratory device using 2D axisymmetric CFD simulation. We willcompare the computational solution of the flow in the device with thefully developed analytical solution.There are some limitations to our analysis as this work is intended toprovide the theoretical basis to apply the novel fractionation principle on352.5. The novel fractionation methodcontinuous industrial applications and to provide a design for the fraction-ation device. Additional experimental work is needed to validate the op-eration and the efficiency of this new principle in a physical fractionationprocess.The first limitation is that we do not consider chemical reactions be-tween fluid layers. As an example the fractionation efficiency in fluids withdifferent pH concentration need to be explored experimentally. In this workwe assume no chemical reactions between fluid layers.The second limitation is that we do not consider diffusion between layersdue to particle movement between layers. As discussed earlier, the effectof solid particles on the flow is limited to a small envelope surroundingthe particle especially in the laminar regime under consideration which hastypically a very large Peclet number. For large Peclet number, the diffusiveeffects are usually expected to be limited over the simulations time scales.36Chapter 3The fully developed problemand particle fractionationWe consider the fully developed fluid flow in an annular space formed be-tween an inner and outer cylinders of radii rˆ ∈ [κRˆ, Rˆ], where κ is a constantdefined in the interval (0, 1), figure 3.1. In this work we define the hattednotation to represent a dimensional quantity.Figure 3.1: A schematic of the model geometry.The inner and outer cylinders rotate with the same angular speed ωˆ. The37Chapter 3. The fully developed problem and particle fractionationflow is subjected to a constant pressure gradient Gˆ, in the axial direction.We consider a multi-layer flow where two or more fluids flow with theirinterface aligned parallel to the direction of flow.Each individual fluid has its unique physical and rheological proper-ties. To be consistent with the work of Bittleston and Hassager (1992) andMadani et al. (2013) that considered a single viscoplastic fluid, we utilizethe Bingham plastic rheological model to describe multi-layer fluids rheologywith the apparent viscosity µˆ(k) and yield stress τˆy(k) defined for each fluid,where k = 1, 2, 3...etc. The fluids are incompressible with the same densityρˆ.The location of each layer is defined apriori by the position of its upperboundary rˆ(k)l . Hence each layer is defined in the region (rˆ(k−1)l , rˆ(k)l ). Forsimplicity in the analysis we will drop the superscript and assume that theseproperties vary as a function of radial position.To develop an appropriate model for the problem, we assume the flowto be laminar, steady, and axisymetric, the fluids to be Non-Newtonian(viscoplastic) with constant properties. In addition, we assume no chemicalreactions and the wall surfaces to be smooth and free from active impurities.The flow in the previously described geometry is governed by the basicpartial differential equations which result from the consideration of both the38Chapter 3. The fully developed problem and particle fractionationconservation of mass principle (continuity equation) and the conservation ofmomentum principle (Navier-Stokes equations). These equations are givenbelow for laminar, steady and incompressable flow.∇.uˆ = 0 (3.1)ρˆ(uˆ.∇) uˆ = −∇pˆ+∇.τˆ (3.2)where uˆ is the velocity, pˆ the pressure and τˆ the deviatoric stress tensor.we non-dimensionalize the Navier-Stokes equations using a length scale ofRˆ and a viscosity scale of µˆ(1). We define a velocity scale Uˆo, time scale tˆoand pressure scale τˆc of :Uˆo =√√√√(κωˆRˆ)2 +(GˆRˆ2µˆ(1))2(3.3)tˆo =ρˆRˆ2µˆ(1)(3.4)τˆc =µˆ(1)UˆoRˆ(3.5)Where Gˆ is the pressure drop per unit length. In the definition of Uˆo,the two terms reflect the effect of rotational speed ωˆ and pressure drop rate39Chapter 3. The fully developed problem and particle fractionationGˆ on the flow cheracteristics.Using these scalings, and omitting the hat notation for dimensionlessvariables, the scaled constitutive equations for the fluid are:τij =(µ(r) + B(r)γ˙)γ˙ij ⇔ τ > B(r)γ˙ij = 0⇔ τ ≤ B(r)(3.6)Where γ˙ and τ are the rate of strain and stress tensors, respectively and Bis Bingham number. These are defined by:γ˙ =(∑ 12γ˙ij γ˙ij) 12(3.7)τ =(∑ 12τijτij) 12(3.8)whereγ˙ij = ∇(uij + uji) (3.9)With these, we find that this flow is characterized by five dimensionlessgroups, the axial and tangential Reynolds numbers, Rez and Reθ, the Bing-ham number , B, the ratio of the swirl and axial velocities, ω, and the ratioof the radii of the two cylinders, κ defined as :40Chapter 3. The fully developed problem and particle fractionationRez =ρˆUˆoRˆµˆ(1)Reθ =ρˆωˆRˆ2µˆ(1)B(r) =τˆy(r)Rˆµˆ(1)Uˆoω =ReθRezκ =RˆiRˆwhere Rˆi is the inner radius of the annulus.If fully developed, the equations of motion reduce to:1rduθdθ+duzdz= 0 (3.10)ddr(r2τrθ) = 0 (3.11)1rddr(rτrz) = G (3.12)withγ˙ =1√2(γ˙2rθ + γ˙2rz) 12 (3.13)τ =1√2(τ2rθ + τ2rz) 12 (3.14)41Chapter 3. The fully developed problem and particle fractionationγ˙rθ = rddr(uθr)(3.15)γ˙rz =ddr(uz) (3.16)Integrating the governing equations 3.11 and 3.12 leads to:τrθ =C1r2(3.17)τrz =Gr2+C2r(3.18)To advance, we need to eliminate γ˙. To do so, we first square the stress-strain relationship. i.e.τ2ij =(µ(r) + B(r)γ˙)2γ˙2ij ⇔ τ > B(r)γ˙2ij = 0⇔ τ ≤ B(r)(3.19)And then summing the terms,∑τ2ij =(µ(r) + B(r)γ˙)2∑γ˙2ij ⇔ τ > B(r)∑γ˙2ij = 0⇔ τ ≤ B(r)(3.20)42Chapter 3. The fully developed problem and particle fractionationOrγ˙ = τ−B(r)µ(r) ⇔ τ > B(r)γ˙2 = 0⇔ τ ≤ B(r)(3.21)Substituting equations 3.15, 3.16, 3.17 and 3.18 into equation 3.21 yields thefinal system of equations:ddr(uθr)=C1r31µ(r)(1− B(r)τ)⇔ τ > B(r) (3.22)ddr(uz) =(Gr2+C2r)1µ(r)(1− B(r)τ)⇔ τ > B(r) (3.23)The governing Equations are subject to the following boundary conditionsincluding the walls non-slip conditions:• At interface:r = ri• A per unit length pressure drop across the annulus of:Gˆ433.1. Method of solution• At the inner wall:uz(κ) = 0uθ(κ) =κωˆRˆUˆo• At the outer wall:uz(1) = 0uθ(1) =ωˆRˆUˆo3.1 Method of solutionHere, the final system of equations together with the boundary conditionsforms a boundary value problem. The equations are integrated and solvedusing Matlab built-in program bvp4c. The program bvp4c uses a colloca-tion method that results in a system of nonlinear algebraic equations that issolved by a variant of Newtons method. This involves many partial deriva-tives of several kinds. To make solving BVPs as easy as possible, the programapproximates these partial derivatives with finite differences, Shampine et al.(2005).443.2. Validation3.2 ValidationWe compare the present model with the work of Madani et al. (2013) forsingle viscoplastic fluid in a concentric annulus of 0.8 radius ratio. Validationis carried out for two values of Bingham number, B = 0.5 and 0.3.Figure 3.2 shows the calculated dimensionless axial velocity profile com-pared to Madani et al. (2013) values. The agreement is very good with themean percentage difference falls below 1.5%.Figure 3.2: Dimensionless axial velocity profile comparison between presentwork and the work of Madani et al. (2013).453.3. Stability definition3.3 Stability definitionAs highlighted earlier, a stable multi-layer viscoplastic flow is possible if theapplied shear stress at the interface is less than that of the yield stress ofone of the fluids that are sharing this interface, (Frigaard (2001); Moyers-Gonzalez et al. (2004); Huen et al. (2007)). Having the applied shear stressbelow the yield stress means that the fluid will not deform and we will haveplug flow. This plug area is needed for the fractionated particles to penetratethrough in order to be separated from other particles in the suspension.Here the force required for the targeted particles (of certain geometryand density) to pass the plug region is known because we know the fluidyield stress value that we need to overcome. But in the yielded region ofthe flow on the other hand the applied stress on fluid is more than the yieldstress and all particles of different densities and geometries move with theflow.Based on the above definition of stability, the flow will be divided to thefollowing sub-categories with regard to stability:• An unstable flow where the applied shear stress on both fluids at theinterface is higher than the yield stress of each fluid. Here, both fluidsare yielded at the interface, figure 3.3-a.463.4. Results• A stable flow. In this type, the applied shear stress at one or moreof the two fluids in the interface is lower than its yield stress value.This situation makes the flow un-yielded on at least one side of theinterface. If the flow is un-yielded in the inner fluid region, we denotethe flow as inner stable flow, figure 3.3-b and if the flow is un-yieldedon the outer fluid region, we would call it an outer stable flow, figures3.3-c.3.4 ResultsWe study the effect of the yield stress, plastic viscosity and radius rationumbers defined respectively as follows:τ (r)y =τˆ(r)yτˆ(1)yµ(r) =µˆ(r)µˆ(1)κ =RˆiRˆwhere Rˆi is the inner radius of the annulus.The inner and outer fluids flow rates are defined respectively as:473.4. ResultsFigure 3.3: Multi-layer flow stability categories: (a) an unstable flow, (b)an inner stable flow and (c) an outer stable flow. The dashed lines are forthe inner fluid velocity profiles.Q1 =Qˆ1QˆoQ2 =Qˆ2QˆoWhere Qˆ1 and Qˆ2 are the dimensional inner and outer fluids flow rates.Qˆo is the reference flow rate defined using the annulus cross-sectional area483.4. ResultsSeries τ(1)y τ(2)y µ(1) µ(2) κ1 (0,1,2) 0.5 0.75 0.75 0.62 1 (0,0.5,2) 0.75 0.75 0.63 1 0.5 (0.375,0.75,1.5) 0.75 0.64 1 0.5 0.75 (0.375,0.75,1.5) 0.65 1 0.5 0.75 0.75 (0.6,0.8,0.9)Table 3.1: A summary of the runs conditions simulated for each parameter:inner and outer fluids yield stresses, plastic viscosities and radius ratio.Aˆc and the reference speed Uˆo:Qˆo = AˆcUˆoThe standard case has the following parameter values:τ (1)y = 1, τ(2)y = 0.5, µ(1) = 0.75, µ(2) = 0.75, κ = 0.6The detailed runs conditions are shown in Table 3.1.One representative case is shown ,i.e. series 1 as defined in the table 3.1,in figure 3.4. This figure outlines the flow state (contour) and defines thepressure drop and interface position at different inlet flow rate conditions.What is clear in this image is that stable flows are achievable and the region493.4. Resultsof stability grows with increasing yield stress.Similar results can be found for the other cases considered (which areshown in Appendix A). Our findings can be summarized into a number ofqualitative “rules of thumb”(a) Increasing yield stress increases the region of stability while increasingthe pressure drop required to maintain the desired flow rate.(b) The size of the stability window was fairly insensitive to the viscosityratio or channel size.At this point we turn our attention to attempting to use these resultsto estimate if separation is possible, given a stable flow field. We representthis methodology as a toy problem to outline the key features required tomake such an estimate. This discussion is by no-means a complete analysis.For separation to occur in the annulus, we need to know the differencein start criteria, as given in figure 3.5, as an example, and the settling rateof the particle introduced on the inner radius. The settling rate and the gapsize represent the settling time in the mixing zone which must be less thanthe residence time in the device, figure 3.6.As such we built a small algorithm to study the potential combination ofoperating conditions to separate two particles for a given rheology of fluids.Formally we seek a solution for the rotational rate and flowrates to separate503.4. Resultsparticles in a known viscoplastic fluid, subject to the interface being stableand the settling time being smaller than the residence time. The settlingrate of particles in a viscoplastic fluid where obtained from Derksen et al.(2011).The test case which we will discuss is similar to what we will conducton the laboratory device. Here, we will attempt to determine the flowratesto separate a 2.5 mm steel spheres from a dilute suspension of spheres withsmaller diameter.The spheres are suspended in a fluid with a yield stress of 1 Pa and islayered onto another fluid with a yield stress of 0.5 Pa. The first step in thisanalysis is to define the critical rotation rate to cause motion. We examinethe data from Madani and define the critical rotation rate, i.e.ωc >√FcτyD2∆ρV Rwhere Fc is the critical force as determined by Madani et al. (2010a) fordifferent particle sizes. ∆ρ is the difference between particle and suspensiondensities, V is the volume of the spherical particle, R is the distance betweenthe particle and axis of rotation, τy is the suspension yield stress and D isthe particle diameter.With this we solve equations 3.22-3.23, for the given fluid, and attempt513.4. ResultsFigure 3.4: The solution of equations 3.22-3.23 for the conditions given asseries 1 in table 3.1. The colour map defines the stability of the flow state,as given in figure 3.3. The dashed lines represent the interface position andthe solid lines represent the pressure drop rate: (a) τ(1)y = 0, (b) τ(1)y = 1,(c) τ(1)y = 2.523.4. ResultsFigure 3.5: Fractionation curves for spherical stainless steel particles in dif-ferent yield stress fluids. The critical rotational speeds ω are denoted bythe colored lines and the black lines for the characteristic speed Uc. Radiusratio, κ = 0.6 , outer radius, R = 0.127 m, plastic viscosity, µ = 1.0833Pa.s, spheres density, ρs = 7800 kg/m3 and suspension density, ρf = 1000kg/m3.to find the subset of cases in which the settling time is smaller than theresidence time in a stable operating window. The results are shown infigure 3.7. The blue lines in this case represent the ratio of residence tosettling times. When this ratio is greater than unity, we consider this to bean acceptable solution.A number of potential flow configurations are evident (stable, unstableand operating window). The operating window configuration is defined asa stable flow with at least 50% of the inner fluid region being unyielded.533.4. ResultsFigure 3.6: Settling and residence rates in the mixing zone of the device.If the inner fluid is mostly yielded, then there will be no enough space forthe targeted particles to penetrate through to the other layer and achieveseparation. In addition, the flow will be at risk of losing its stability.In this chapter we solve analytically the fully developed region of the flowof interest. We validate the solution against previous single fluid viscoplasticanalytical work. We categorize the stability of the flow based on the yieldingcondition of the fluids sharing the interface. We present the effect of fluidsreheology on the flow stability.In addition, we show that the novel fractionation method can be usedeffectively to fractionate particles suspensions in continuous mode. We pro-vide some useful curves and diagrams that characterise the flow with regardsto its stability and its ability to successfully fractionate particles suspensions.To do so, we define the fractionation operating window, we present the crit-ical rotational speed to fractionate particles of different sizes in fluids ofvariable yield stress values.543.4. ResultsFigure 3.7: Flow stability and operating diagram. Residence to settling timeratio in blue lines, interface locations in black dashed lines and dimensionlesspressure drop rate in solid black lines. The annulus axial length L = 50 cm,inner fluid yield stress τ(1)y = 1 Pa, outer fluid yield stress τ(2)y = 0.5 Pa,radius ratio, κ = 0.6 , outer radius, R = 0.127 m, plastic viscosity, µ(1) =µ(2) = 0.75 Pa.s, spheres density, ρs = 7800 kg/m3 and suspension density,ρf = 1000 kg/m3, ω = 37.2 rad/s and pressure drop rate of (50−500) Pa.s.Moreover, we develop a flow and stability diagram for a specific frac-tionation case (known inner and outer fluids rheology, annulus geometryand targeted particles size and density). In this diagram we show the areasof stability, instability, operating window and other useful information tothe operator.55Chapter 4The full annulus problemand density current effectsIn the previous chapter we examined the fully developed solution of theproblem in order to have the first base for designing the continuous device.The second base will be solving the entrance region of the multi-layer vis-coplastic flow inside an annulus. This intermediate stage is needed primarilyto design the mixing point of the continuous device. In addition to our treat-ment of the yielded interface instability in the previous chapter, we extendthe discussion to the Kelvin Helmholtz and density current instabilities inour effort to address all aspects of the flow stability.In this chapter, we aim at obtaining a numerical solution for the fullsystem of Navier stokes equations of the pressure driven laminar forcedmulti-layer flow problem of an annulus with rotating walls of a constantangular velocity. The layers are of viscoplastic fluids with variable density564.1. Mathematical modelingand are aligned in parallel to the annulus axis.4.1 Mathematical modelingWe consider the fluid flow in the annular space formed between an innerand outer cylinders of radii rˆ ∈ [κRˆ, Rˆ], where κ is a constant defined inthe interval (0, 1), figure 4.1. In this work we define the hatted notation torepresent a dimensional quantity.Figure 4.1: A schematic of the model geometry.The inner and outer cylinders rotate with the same angular speed ωˆ. Theflow is subjected to a constant pressure gradient Gˆ, in the axial direction.We consider a multi-layer flow where two fluids flow with their interfacealigned parallel to the direction of flow.We follow the same mathematical modeling as highlighted in chapterthree taking into consideration that k is limited to k = 1, 2 since we have574.1. Mathematical modelingonly two fluids. Here, 1 denotes the inner fluid and 2 denotes the outerfluid. Both the inner and outer fluids enter the annulus with a uniformvelocity profile. In developing our model we assume the same assumptionsas we did in the previous study with the exception of the fully developedflow assumption.We utilize the concentration α1 to model the change in concentrationbetween pure fluids 1 and 2 with α1 ∈ (0, 1). After scaling , the non-dimensional continuity and Navier-Stokes equations are∇.u = 0 (4.1)(1 + φAt)Rez(u.∇) u = −∇p+∇.τ (4.2)where u is the velocity, p the pressure and τ the deviatoric stress tensor.Here the function φ(α1) = 1 − 2α1 has the range of φ ∈ (−1, 1) forα1 ∈ (0, 1). The additional dimensionless parameter is the Atwood number,At which results from the variable density nature of the model. The Atwoodnumber is defined asAt =ρ1 − ρ2ρ1 + ρ2584.1. Mathematical modelingwhere, ρ1 and ρ2 are the densities of the heavier and lighter fluids, respec-tively.It is imperative at this point to define the Richardson number since wewill explore the Kelvin Helmholtz instability of the flow. The Richardsonnumber (Ri) is the dimensionless number that expresses the ratio of thebuoyancy term to the flow gradient term.Ri =buoyancy termflow gradient term=gρ∇ρ(∇u)2where g , ρ and u are acceleration of gravity, density and velocity, respec-tively.In the current numerical solution, all control volumes are filled with oneor more of the phases (fluids). The κth fluids volume fraction in a controlvolume is called ακ and we can have ακ = 0 if the control volume is emptyof the κth fluid, ακ = 1 if the control volume is full of the κth fluid and0 < ακ < 1 if the control volume contains an interface between the κth fluidand other fluids. Here, in each control volume, the volume fractions of allphases sum to unity.594.1. Mathematical modelingn∑κ=1ακ = 1 (4.3)Where n is the number of phases.All properties appearing in the continuity and momentum equations arevolume fraction averaged properties. For example, the volume fraction av-eraged density is computed as followsρ =n∑κ=1ρ(κ)ακ (4.4)The governing equations are subject to the following inlet, outlet andnon-slip boundary conditions:• At the inlet (z = 0):A uniform velocity profile is used for the z-component of the fluidsvelocityuz1 =uo1Uˆouz2 =uo2Uˆowhile the r and θ components of the velocity are assumed to be zero.• At the outlet (z = L):We impose zero axial gradients at the outlet (outflow condition influent).604.1. Mathematical modeling• At the inner wall, we impose the non-slip condition:uz(κ) = 0uθ(κ) =κωˆRˆUˆo• At the outer wall, we impose the non-slip condition:uz(1) = 0uθ(1) =ωˆRˆUˆoWe conduct different experiments to study the hydrodynamic entrancelength. Here we include reference to the Reynolds number based on thehydraulic diameter and average velocity, Reh, to compare the current resultswith previous literature work. Reh is defined asReh =ρˆ(1)VˆavDˆhµˆ(1)where Vˆav is the flow average velocity and Dˆh is the annulus hydraulicdiameter. Details of entrance length study runs conditions are shown in table4.1. Another set of runs are conducted to study the interface behaviour,density current instability and Kelvin Helmholtz instability. These runsconditions are detailed in table 4.2.614.1. Mathematical modelingStudy Reh Rez B(1)0.10 8.1 x 101 0.0001.00 8.1 x 102 0.000Newtonian 10.00 8.1 x 103 0.00072.00 5.8 x 104 0.000100.00 8.1 x 104 0.000200.00 1.6 x 105 0.0000.10 5.4 x 102 0.047Viscoplastic 1.00 1.3 x 103 0.01910.00 8.7 x 103 0.003Table 4.1: A summary of the runs conditions simulated for the entrancelength study. Simulations are for single fluid of Newtonian or viscoplasticbehaviour.624.1. Mathematical modelingStudy Rez Reθ At B(1) B(2)Interface behaviour 4.4 x 102 8.70 0.00 0.058 0.0358.2 x 102 8.70 0.00 0.031 0.0194 x 102 0.00 0.20 0.00 0.00Density current instability 4.2 x 102 166.00 0.20 0.00 0.002 x 104 166.00 0.20 0.062 0.0623 x 104 166.00 0.20 0.064 0.064Kelvin Helmholtz Instability 5.2 x 103 7.25 0.05 0.000 0.000Table 4.2: A summary of the runs conditions simulated for each study:Interface behaviour, Density current and Kelvin Helmoholtz instabilities.634.2. Method of solution4.2 Method of solutionFor this case we have to solve the flow in the entire annulus including theentry region. This means that we have to deal with the full form of thegoverning equations, there is no analytical means to solve the governingequations and hence they will be treated numerically.We solve the governing equations of the conservation of mass and mo-mentum using a finite difference technique. A numerical algorithm will beused to solve for the three velocity components and pressure of the solu-tion domain. The problem under consideration is an axi-symmetric steadylaminar problem; the governing equations are solved numerically using theFLUENT software. The numerical model uses a control-volume-based tech-nique to convert the governing equations to algebraic equations that can besolved simultaneously.This control volume technique is based on discretizing the domain intodiscrete control volumes using a computational grid. Then integrating thegoverning equations about each control volume, yielding discrete equationsthat conserve each quantity on a control-volume basis. After that the dis-cretized equations are linearized. The resultant linear equation system issolved to yield updated values of the dependent variables.The differencing schemes used are both formally second-order in ac-644.3. Grid independence testcuracy: central differencing is used for the diffusive terms and a second-order up-winding scheme for the convective terms. Coupling of the pressureand velocity was achieved using the well-known semi-implicit method forpressure-linked equations (SIMPLE) implementation of Patankar (1980).4.3 Grid independence testTo achieve grid independence, the solution was computed for different gridsizes to assess the effects of mesh refinement. The grid sizes used were100x266, 150x400 and 225x600 referred to as Grid-1, Grid-2, and Grid-3 respectively. Figure 4.2 shows the volume fraction distribution inside theannulus and the parameter used to check the grid independency. The param-eter is the fully developed axial velocity. Percentage difference between allgrids falls below 1%. Consequently, it was concluded that Grid-1 (100x266)is satisfactory and gives an accurate solution.4.4 ValidationTo check the adequacy of the present numerical solution, we compared thefully developed axial velocity profile of the numerical solution with the cor-654.4. ValidationFigure 4.2: Full annulus model grid independence for Rez=437.17,Reθ=8.7, At=0, B(1)=0.058 and B(2)=0.035 : (a) volume fraction distribu-tion: the red and blue colors denote the inner and outer fluids, respectively.(b) grid independence.responding analytical solution as seen in figure 4.3. Here, the annulus has aradius ratio of 0.8. The inner fluid enters the annulus in the lower quarterof its gap.Two cases (case1 and case2) were conducted to validate the numericalcode. For case 1, figure 4.3-c, the average percentage difference between thenumerical solution and the analytical solution lies within 2.8% differenceindicating good agreement.For the other case or case 2, we increase the two streams flow rates totest the code at different interface location. The comparison between the664.4. ValidationFigure 4.3: Full annulus model validation with analytical solution: (a)volume fraction distribution for case 1, (b) volume fraction distribution forcase 2, (c) validation (case 1: Rez=437.17, Reθ=8.7, At=0, B(1)=0.058 andB(2)=0.035) and (d) validation (case 2: for Rez=820.57, Reθ=8.7, At=0,B(1)=0.031 and B(2)=0.019). The red and blue colors denote the inner andouter fluids, respectively.674.5. Resultsfully developed axial velocity profile of the numerical solution with the cor-responding analytical solution for this case is shown in figure 4.3-d. Theaverage percentage difference between the numerical solution and the ana-lytical solution lies within 3% difference that reflects again a good agreementbetween the solutions.4.5 ResultsIn this subsection we estimate the entrance length for the flow under con-sideration, track the interface location, look at the axial velocity profiledevelopment, explore the Kelvin Helmholtz instability related to the flowand present the results resulting from including the difference of densitybetween layers.4.5.1 Entrance lengthEntrance length (Le) for internal flow is the distance from inlet where theflow no longer changes with axial direction and is said to be fully developed.Downstream of entrance length, the velocity profile and wall shear are bothconstant, and the pressure drops linearly with axial direction, White (2003).We set a criterion to calculate the entrance length as: the location atwhich the axial flow velocity reaches 98% of its fully developed value, figure684.5. Results4.4-c. The volume fraction distribution and axial velocity contours for thesame case are presented in parts a and b of figure 4.4, respectively.We predict the entrance length of the Newtonian case for the creepingand laminar flow regimes with Reynolds number in the range of Reh ∈(0.1 − 200) and predict the entrance length of the viscoplastic case for thecreeping and laminar flow regimes with the value of Reynolds number in therange of Reh ∈ (0.1− 10).Our study results are represented together with previous results of theliterature (Poole and Chhabra (2010) for pipe flow, McComas (1967) andMaia and Gasparetto (2003) for annular flow) in figure 4.4-d. A considerableamount of information can be gained from this figure. Firstly, the annularentrance length values are always considerably below the pipe values for thesame Reynolds number. There are noticeable differences between the valuesreported in the literature for the same cases (same radius ratio and Reynoldsnumber) for the entrance length in the laminar regime. The Newtonianfluid entrance length (Le/Dh) approaches to a unique value as the Reynoldsnumber goes to zero, the value is 0.28 in the case of radius ratio of 0.8.Visco-plastic fluids were found to have shorter entrance length for the sameReynolds number compared to Newtonian fluids in annuli.694.5. ResultsFigure 4.4: Entrance length. a) Volume fraction distribution: the red andblue colors denote the inner and outer fluids, respectively. b) axial velocitycontours with enlarged views of velocity vectors near the wall and in mid-gap, c) Centerline velocity normalized with analytical fully developed valueversus axial location, d) Entrance length for Newtonian and viscoplasticfluids in an annulus of radius ratio of 0.8.704.5. Results4.5.2 Tracking the interface locationTracking the interface location is an important task especially when it comesto fractionation efficiency. Figure 4.4-a shows the fluids volume fraction dis-tribution in the flow, the sharp interface can be easily distinguished in thefully developed region compared to entry region. This may be attributed tothe presence of transverse pressure in the inlet region. The fully developedlocation of the interface was at a radial location of 0.95. We have an ex-cellent agreement between the numerical and the analytical values with thepercentage difference falls below 1%.4.5.3 Axial velocity profile developmentFigure 4.5 shows the development of axial velocity profiles. The profiles atthe inlets are flats (uniform velocity) for both the inner and outer fluids. Itis clear that the flow is mostly yielded in the entrance region (z=0.05). Thisfact may be attributed to two factors. Firstly, while boundary layer approx-imation assumes constant lateral pressure (along each cross-section) in thefully developed region of the flow, the actual flow has a strong transversepressure especially at inlets.Only full elliptic equations can predict accurately the flow field charac-teristics in the entrance region (differences between parabolic and elliptic714.5. Resultsdifferential equations are detailed in Appendix C). This behaviour was re-ported by Vradis et al. (1993) who reported a yielded axial velocity profile inthe entrance region of a viscoplastic flow with a uniform inlet velocity profile.The numerical code used in this study (Ansys Fluent) has the capability tosolve the full elliptic equations.The second reason in our opinion for the flow to be yielded is somethingrelated to multi-layer viscoplastic flow nature, here, the inner and outerfluids compete in the entrance region exerting an additional pressure oneach other through the interface. Initially, the inner fluid fills only 1/5 of theannular volume and is restricted by the inner wall. The only way to satisfythe non-slip condition at the inner wall is to expand toward the outer wallexerting more pressure on the outer fluid through the interface and leavingless space for the outer fluid to occupy. Consequently, this condition leadsto more yielding in the outer fluid velocity profile.4.5.4 Kelvin Helmholtz instabilityWhen we consider parallel multi-layer flow with different layers velocities, itis important to check the Kelvin Helmholtz instability. The Kelvin Helmholtzinstability occurs when there is a difference in velocity across the interfacebetween two fluids of different densities. To check this instability, we set the724.5. ResultsFigure 4.5: Axial velocity profile development: a) Axial velocity profiles atinlets, entrance region and fully developed region. b) Axial velocity contours,Rez=437.17, Reθ=8.7, At=0, B(1)=0.058 and B(2)=0.035.734.5. ResultsRichardson number to 0.08 which is way below the critical number accordingto Richardson Number Criterion.The criterion states that the inequality Ri < 1/4 is a necessary condi-tion for linear instability of inviscid stratified parallel flows. However, thecriterion does not state that the flow is necessarily unstable if Ri < 1/4somewhere, or even everywhere, in the flow.Another requirement for this instability are Rayleighs inflection pointand Fjortofts stronger condition criteria. They states that a necessary (butnot sufficient) criterion for instability of an inviscid parallel flow is that thebasic velocity profile has a point of inflection and the flow should satisfy thecondition of (urr(u − ui) < 0) for all values of r, where urr is the secondderivative of velocity and ui is the velocity at the inflection points.Please refer to figure 4.6 that shows the velocity profile at the inlet andimmediately after the contact of the two fluids. Here due to the non-slipcondition at the walls, mass continuity and the continuous nature of thevelocity profile at the interface, the multi-layer velocity profile behaves asfollow:The slower stream (inner fluid) velocity profile adjusts instantaneouslyso that the minimum is zero at the wall and the maximum is at the in-terface. The faster stream (outer fluid) starts from zero value at the wall744.5. Resultsand increases to a maximum value (above the inlet velocity) somewhere be-tween the wall and the interface and then decreases to a lower value at theinterface. So essentially a zero velocity slope happens only once in the flow.Although the velocity profile has two inflection points at around r=0.85 and again around r = 0.9, the flow does not satisfy Fjortofts strongercritierion since (urr(u − ui) > 0). It is worth mentioning that the KelvinHelmholtz instability is based on the inviscid flow where bounding walls arenot in the picture. We can conclude that the presence of walls stabilizes theflow.Figure 4.6: Axial velocity profiles, Rez=5248, Reθ=7.25, At=0.05,B(1)=B(2)=0 and Ri=0.08.754.5. ResultsIn the following section we include the density difference effect on theflow as this is the general case for the real (industrial) fluids densities.4.5.5 Density difference effectsThe previous results for both analytical and numerical solutions were pro-duced based on the assumption of an iso-dense flow (which is most common)but there are some cases where the flow contains fluids of different densities.To address this, we set the Atwood number (At) to 0.2 to have the innerfluid as the heavier one (density unstable) or having the outer fluid as theheavier fluid (density stable). This value of Atwood number covers all rangesof Atwood number for the practical (industrial) combination of liquids.We study the effect of rotational Reynolds number Reθ and Binghamnumber, B to see the effect of the yield stress on rotating flow stability. Insingle fluid viscoplastic poisuielle flow in pipes and annuli, raising the yieldstress has a stabilizing effect on the flow by elevating the critical Reynoldsnumber, Rec where the laminar flow starts to change to the turbulent regime,Hanks and Dadia (1971).Figure 4.7 shows that Bingham number has a stabilizing effect on rotat-ing flow. Here we do not have any mixing in the Newtonian case (B(1) =B(2) = 0) for the flow without rotation (Reθ = 0) as shown in figure 4.7-a.764.5. ResultsWhen introducing the rotation (Reθ = 166) as shown in figure 4.7-b, anddue to the unstable density situation (heavier inner fluid), the Newtoniancase started to have some mixing. However, increasing Bingham numbersto (B(1) = B(2) = 0062) seems to prevent the mixing and the flow tendsto be more stable, figure 4.7-c. Increasing further the Bingham numbers to(B(1) = B(2) = 0064) as shown in figure 4.7-d leads to a fully stable flowwithout any mixing.For the density stable cases (with light inner fluid and heavy outer fluid)or for iso-dense cases (both fluids have the same density), the simulationsshowed that the flows are always stable.In this chapter, we solve the spiral multi-layer viscoplastic flow insidea full annulus including the entrance region. We report the hydrodynamicentrance length for the flow in the lower range of laminar regime and thecreeping regime for both Newtonian and viscoplastic fluids. We found thatthe entrance length of the viscoplastic flow is always smaller than the cor-responding Newtonian flow for the same range of Reynolds number.The Kelvin Helmholtz instability was not evident due to the presenceof the walls of this bounded flow. The flow was always found to be stablefor density stable or iso-dense cases. On the other hand, density unstablecases destabilize the flow but the flow can be stabilized by increasing the774.5. ResultsFigure 4.7: Density difference effect: Volume fraction distribution for At= 0.2 at: a) Rez = 403.26, Reθ = 0 and B(1) = B(2) = 0, b) Rez = 424.60,Reθ = 166 and B(1) = B(2) = 0, c) Rez = 20467.47, Reθ = 166 and B(1) =B(2) = 0.062, d) Rez = 30006.11, Reθ = 166 and B(1) = B(2) = 0.064. Thered and blue colors denote the heavy and light fluids, respectively.Bingham number (or increasing the yield stress).78Chapter 5The continuous fractionationdevice and fractionationefficiencyNow with all information we conclude from the studies in the previous chap-ters, i.e., the fully developed flow solution study, the full annulus solutionincluding the entrance region and the continuous fractionation study, we areable to design the full continuous fractionation device. We solve the flowinside the device and see the implications of the physical constraints of thedevice parts on the device operating parameters.Our aim is to obtain a numerical solution for the full system of Navierstokes equations of the pressure driven laminar multi-layer viscoplastic fluidflow problem inside the continuous fractionation device. The walls of thedevice are rotating walls at a constant angular velocity.795.1. Mathematical modeling and method of solution5.1 Mathematical modeling and method ofsolutionWe consider the fluid flow inside the continuous fractionation device geome-try, figure 5.1. The dimensions were chosen based on the constraints resultedfrom the consideration of the intended range of operation (fluids flow ratesand yield stresses values) and the available commercial pipe sizes.Figure 5.1: A schematic of the continuous fractionation device model ge-ometry.We follow the same approach that we use for the full annulus study(chapter 4) with the same assumptions and dimensionless groups. We fol-low also the same method of solution to solve numerically the flow underconsideration. Same boundary conditions are defined with the exception ofthe outlet and walls boundary conditions that are defined as follow:• At the outlet (z = 0.3952m):805.2. Grid independence testWe impose zero axial gradients at the (outflow condition in fluent).• At the walls:We impose the non-slip conditionuz = 0uθ = ωˆrˆWe conduct several experiments to study the behaviour of the flow atdifferent zones of the device. Our main parameter is the fluid flowrate ratios(Q1/Q2) that control the layers interface location, ri. The interface locationis a crucial factor that determines the efficiency of the fractionation as weshow shortly in the discussion of the results section. Details of the runsconditions are shown in table 5.1.5.2 Grid independence testTo achieve grid independence, the solution was computed for different gridsizes to assess the effects of mesh refinement. The grid sizes used were55,000 cells, 105,000 cells and 155,000 cells referred to as Grid-1, Grid-2,and Grid-3 respectively. Figures 5.2 shows the parameter used to checkthe grid independency. The parameter is the fully developed axial velocity.Percentage difference between all grids falls below 1%. Consequently, it is815.3. ValidationQ1/Q2 ri2.98 0.841.22 0.750.96 0.720.75 0.690.60 0.66Table 5.1: A summary of the simulation runs conditions at: Rez = 149.3,Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093.concluded that Grid-1 (105,000 cells) is satisfactory and gives an accuratesolution.5.3 ValidationTo check the adequacy of the present computer code, we compared the fullydeveloped axial velocity profile of the numerical solution with the corre-sponding analytical solution as seen in figure 5.3. The average percentagedifference between the numerical solution and the analytical solution lieswithin 2 percent difference reflecting a good agreement.825.4. ResultsFigure 5.2: Fractionation device model grid independence for Rez = 149.3,Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093 : (a) volume fraction distributionand (b) grid independence.5.4 ResultsIn this section we track the development of the axial velocity profiles indifferent regions of the fractionation device and estimate the entrance lengthof the mixing region, we investigate the behaviour of the flow mainly theinterface in the exit region of the device and predict the implications of thisbehaviour on the streams flowrates design and the size of the fractionationoperating window.835.4. ResultsFigure 5.3: Fractionation device model validation with analytical solutionat Rez = 149.3, Reθ = 7.88, B(1) = 0.155 and B(2) = 0.093. (a) volumefraction distribution. (b) validation.5.4.1 Axial velocity profile development and entrancelengthTracking the axial velocity profiles development is important to make surethat the flow reaches the fully developed condition early in the mixing zonein order to allow sufficient space for fractionation to take place. The mixingzone is the space where the two fluids meet between the inlet separator and845.4. Resultsthe outlet separator. Figure 5.4 shows the development of axial velocity pro-files at different locations before, inside and after the mixing zone. Althoughthe velocity profiles at the inlets are flat (uniform velocity), they developdownstream and adapt to fulfill the non-slip condition at walls (zero axialvelocity at walls).The velocity profiles reach the fully developed condition before and insidethe mixing zone in this case, however, in the exit region (after mixing zone)and specifically the inner fluid stream outlet, the velocity profile do not reachthe fully developed condition because the axial distance was not sufficientto cover the entrance length for that outlet.The entrance length (Le/Dh) was estimated to be in the order of 0.32.The value for the corresponding entrance length for a single viscoplastic fluidwith the same value of Reynolds number was found to be 0.20 as reportedin the previous chapter. This discrepancy is likely to be attributed to thefact that the current simulation is for the flow of two viscoplastic fluids flow.In either case the entrance length is way below the corresponding en-trance length of a viscoplastic pipe flow for the same range of Reynoldsnumber. It is preferable to have a short entrance length as what we havein this case. This is because the time scale for the particles residence timein the entry region is small and the fractionation occurs solely in the fully855.4. Resultsdeveloped region of the mixing zone.Figure 5.4: Axial velocity profile development: (a) Volume fraction dis-tribution, (b) axial velocity profile development (before, inside and aftermixing zone). Simulated at Rez = 149.3, Reθ = 7.88, B(1) = 0.155 andB(2) = 0.093.5.4.2 Outlet separator effects on interfaceThe presence of the outlet separator has an important effect on the flowcharacteristics especially the interface location as shown in figure 5.5. Inthe current simulation, the fluid interface in the fully developed region of865.4. Resultsthe mixing zone is located at a radial distance slightly higher than the radialdistance of the outlet separator centerline. This situation makes the interfacemove out towards the outer wall since the flow has to satisfy the no-slipcondition at the separator wall. Here, the interface tends to favor the areasof low shear stress in somewhere around the middle of the outer fluid outletgap compared to the areas of high shear stress near the separator walls.Figure 5.5: Fractionation device volume fraction distribution showing theeffect of outlet separator on interface location, Rez = 149.3, Reθ = 7.88,B(1) = 0.155 and B(2) = 0.093.This phenomenon may lead to poor fractionation as some of the non-targeted particles may leave with the targeted (fractionated) particles throughthe outer fluid outlet. The flow needs to be designed (choosing the twostreams flowrates) in a way that ensures that the interface of the fully devel-oped region of the mixing zone is located radially below the outlet separator.This selection of flowrates leads to a successful fractionation by making surethat only targeted particles penetrate to the outer stream before they reach875.4. Resultsaxially the outlet separator.5.4.3 Flow rates designAdditional simulations were carried out to investigate the behaviour of theflow at the exit region. Figure 5.6 shows the position of the interface forfour different cases (simulations) as a result of the presence of the outletseparator. The separator physical location is bounded between the radiallocations of (0.74 - 0.81). In the first two simulations when the values ofthe mixing zone fully developed interface were 0.84 and 0.75, the interfacediverged to the outer fluid outlet (a situation that potentially hinders thefractionation efficiency).Figure 5.6: Interface behaviour at the exit region, Rez = 149.3, Reθ = 7.88,B(1) = 0.155 and B(2) = 0.093 at: (a) ri = 0.84, Q1/Q2 = 2.98, (b)ri = 0.75, Q1/Q2 = 1.22, (c) ri = 0.72, Q1/Q2 = 0.96 and (d) ri = 0.69,Q1/Q2 = 0.75.885.4. ResultsWhen we designed the flowrates to have the fully developed interfacebe located at 0.72, we noticed that the interface diverged to the inner fluidoutlet (a favorable location in terms of fractionation efficiency) in close prox-imity to the separator. An interface of 0.69 would continue and exit freelythe device through the inner fluid outlet.We conclude that for this type of fractionation requirement and flu-ids rheology, the fluids flowrates need to be designed to have the maxi-mum bound of the mixing zone interface to be 0.72 which corresponds to aflowrates ratio (Q1/Q2) of 0.96.5.4.4 Exit region effects on fractionation diagramTo assess the device fractionation effectiveness especially with the distur-bance in the exit region caused by the outlet separator, we demonstrate aspecific fractionation case. We use the fractionation code that was devel-oped and presented in chapter three to generate the required fractionationdiagram.In order to construct a fractionation diagram, we need to solve the flowin thousands of cases (a flow case has a unique combination of streamsflowrates) in fluids of known properties. The analytical solution of chapterthree is a feasible way to accomplish this task.895.4. ResultsOn the other hand, the numerical solutions as presented in chapters fourand five would require several days to solve only one flow case. I generatea fractionation diagram for the case of fractionating a 2.5 mm sphere in aviscoplastic fluid suspension. The properties of the fractionation case aredetailed in table 5.2. As shown in the fractionation diagram, figure 5.8,the operating window occupies a noticeable area of the diagram giving theoperator flexibility with a wide range of flowrates to operate the device.The time ratio in figure 5.8 is defined as the ratio between the residence andthe settling times. The settling rate is the rate at which the particle travelsradially across the gap, figure 5.7. This rate and the gap size represent thesettling time in the chamber. The residence time on the other hand is thetime spent as the particle travels axially through the mixing zone of thedevice. For fractionation to take place, the time ration has to be greaterthan unity.When we include the effect of the outlet separator presence in the exitregion, figure 5.9, the operating window is affected drastically and its cover-age area shrunk considerably limiting the range of operation of the device.A relief of the operating window criterion from 50% (which is conservative inour opinion) to 40 or 30% would result in expanding the operating window.905.4. ResultsProperty Value UnitParticles diameter 2.5 mmParticles density 7800 kg/m3Suspension density 1000 kg/m3Annulus radius ratio 0.368 dimensionlessAnnulus outer radius 0.0476 mInner and outer fluids plastic viscosity 1.083 Pa.sInner fluid yield stress 10 PaOuter fluid yield stress 6 PaPressure drop/unit length range 750 to 5000 Pa/mRotation speed 40.22 rad/sTable 5.2: fractionation case properties.915.4. ResultsFigure 5.7: Settling and residence rates in the mixing zone of the device.This measure could be used if the operator elects to increase the range ofoperation. The fractionation diagram will look like figures 5.10 and 5.11 ifthe operating window criterion is reduced to 40% and 30%, respectively.In this chapter we design the continuous fractionation device based onthe information we have from the fully developed, full annulus and contin-uous fractionation model studies. We validate the current code with theanalytical solution. We track the velocity profile development to make surethat the flow reaches the fully developed condition early in the mixing zoneto allow sufficient space for fractionation to take place. We study the exit re-gion of the device and suggest certain guidelines to address the hurtful effectof the outlet separator on the layers interface and fractionation efficiency.925.4. ResultsFigure 5.8: Fractionation diagram for the case detailed in table 5.2. Timeratios in blue lines, interface locations in black dashed lines and dimension-less pressure drop rate in solid black lines.Figure 5.9: Fractionation diagram after including the effect of outlet sepa-rator. Case properties are per table 5.2. Time ratios in blue lines, interfacelocations in black dashed lines and dimensionless pressure drop rate in solidblack lines.935.4. ResultsFigure 5.10: Fractionation diagram for an operating window criterion of40%. Case properties are per table 5.2. Time ratios in blue lines, interfacelocations in black dashed lines and dimensionless pressure drop rate in solidblack lines.Figure 5.11: Fractionation diagram for an operating window criterion of30%. Case properties are per table 5.2. Time ratios in blue lines, interfacelocations in black dashed lines and dimensionless pressure drop rate in solidblack lines.94Chapter 6Summary and conclusionsIn this work we studied the application of the novel technology of particlesuspension fractionation that was proposed by Madani et al. (2010b) on thecontinuous operation process. The flow that is suggested to be employedin facilitating this new fractionation methodology is the spiral multi-layerpoisuille annular flow. Here we presented several studies of the flow indifferent circumstances/states namely: fully developed region, annular entryregion and full fractionation device. We focused on the stability of the flowin conjunction with the fractionation requirement.We found that spiral multi-layer Poisuille viscoplastic flow can be stablefor both annulus and continuous fractionation geometries. The range of sta-bility expands with increasing fluids yield stress as the yield stress supressesthe yielded interfacial instability.Gravity current instability was found to be evident in density unstablecases due to the centrifugal forces resulting from rotation and mixing can95Chapter 6. Summary and conclusionshappen in these cases. However, the flow can be stabilized by increasing thefluids yield stress (Bingham number). In the other hand, density current isnot a factor in density stable and iso-dense cases.Kelvin Helmholtz instability was not found on conditions tested. Thepresence of the walls in this internal flow has a stabilizing effect and preventsthe onset of Kelvin Helmholtz instability.All annular flows entrance lengths were found to be shorter than thecorresponding ones of pipe flows for the same values of Reynolds number.We reported the hydrodynamic entrance length for the creeping and lowrange laminar flow regimes for annuli of both Newtonian and viscoplasticfluids. The entrance length of the viscoplastic flow was found to be shorterthan the equivalent Newtonian flow for the same range of Reynolds number.The disturbance to the flow at the exit region of the continuous frac-tionation device hinders the fractionation process due to the presence ofthe outlet separator. This situation adds more restriction to the fractiona-tion operating window. In this regards, carful design of streams flowrates iscrucial to providing a successful fractionation process.96BibliographyMB Abbott. On the spreading of one fluid over another. La Houille Blanche,(6):827–856, 1961.U Ts Andres. Equilibrium and motion of spheres in a viscoplastic liquid. InSoviet Physics Doklady, volume 5, page 723, 1961.Byron E Anshus. Bingham plastic flow in annuli. Industrial & EngineeringChemistry Fundamentals, 13(2):162–164, 1974.DIH Barr. Densimetric exchange flow in rectangular channels. La HouilleBlanche, (6):619–632, 1967.GK Batchelor. Sedimentation in a dilute dispersion of spheres. Journal offluid mechanics, 52(02):245–268, 1972.T Brooke Benjamin. Gravity currents and related phenomena. Journal ofFluid Mechanics, 31(02):209–248, 1968.AN Beris, JA Tsamopoulos, RC Armstrong, and RA Brown. Creeping97Bibliographymotion of a sphere through a bingham plastic. Journal of Fluid Mechanics,158:219–244, 1985.SH Bittleston, J Ferguson, and IA Frigaard. Mud removal and cement place-ment during primary cementing of an oil well–laminar non-newtonian dis-placements in an eccentric annular hele-shaw cell. Journal of EngineeringMathematics, 43(2-4):229–253, 2002.Simon H Bittleston and Ole Hassager. Flow of viscoplastic fluids in a rotatingconcentric annulus. Journal of non-newtonian fluid mechanics, 42(1):19–36, 1992.Jason E Butler and Eric SG Shaqfeh. Dynamic simulations of the inhomo-geneous sedimentation of rigid fibres. Journal of Fluid Mechanics, 468:205–237, 2002.SS Chen, LT Fan, and CL Hwang. Entrance region flow of the binghamfluid in a circular pipe. AIChE Journal, 16(2):293–299, 1970.Raj P Chhabra. Bubbles, drops, and particles in non-Newtonian fluids. CRCpress, 2006.MR Chowdhury and Firat Yener Testik. Viscous propagation of two-98Bibliographydimensional non-newtonian gravity currents. Fluid Dynamics Research,44(4):045502, 2012.V Ciriello, S Longo, L Chiapponi, and V Di Federico. Porous gravity currentsof non-newtonian fluids within confining boundaries. Procedia Environ-mental Sciences, 25:58–65, 2015.JJ Derksen et al. Direct simulations of spherical particle motion in binghamliquids. Computers & Chemical Engineering, 35(7):1200–1214, 2011.V Di Federico, S Malavasi, and G Bizzarri. Viscous spreading of non-newtonian gravity currents in radial geometry. WIT Trans Eng Sci, 52:399–408, 2006a.V Di Federico, R Archetti, and S Longo. Spreading of axisymmetric non-newtonian power-law gravity currents in porous media. Journal of Non-Newtonian Fluid Mechanics, 189:31–39, 2012.Vittorio Di Federico, Stefano Malavasi, and Stefano Cintoli. Viscous spread-ing of non-newtonian gravity currents on a plane. Meccanica, 41(2):207–217, 2006b.F Durst, S Ray, B U¨nsal, and OA Bayoumi. The development lengths of99Bibliographylaminar pipe and channel flows. Journal of fluids engineering, 127(6):1154–1160, 2005.MP Escudier, PJ Oliveira, and FT Pinho. Fully developed laminar flow ofpurely viscous non-newtonian liquids through annuli, including the effectsof eccentricity and inner-cylinder rotation. International Journal of Heatand Fluid Flow, 23(1):52–73, 2002.EE Feldman, RW Hornbeck, et al. A numerical solution of laminar devel-oping flow in eccentric annular ducts. International journal of heat andmass transfer, 25(2):231–241, 1982.James Feng, Howard H Hu, and Daniel D Joseph. Direct simulation of initialvalue problems for the motion of solid bodies in a newtonian fluid part 1.sedimentation. Journal of Fluid Mechanics, 261:95–134, 1994.Michael McNeil Forbes, A. Apple, and B. Boat. Frequency of quality testingin syrup creation. Maple Science J., 255:139–144, 2010.Edmund J Fordham, Simon H Bittleston, and M Ahmadi Tehrani. Vis-coplastic flow in centered annuli, pipes, and slots. Industrial & Engineer-ing Chemistry Research, 30(3):517–524, 1991.100BibliographyArnold Fredrickson and R Byron Bird. Non-newtonian flow in annuli. In-dustrial & Engineering Chemistry, 50(3):347–352, 1958.IA Frigaard. Super-stable parallel flows of multiple visco-plastic fluids. Jour-nal of non-newtonian fluid mechanics, 100(1):49–75, 2001.Julio Gratton, Fernando Minotti, and Swadesh M Mahajan. Theory ofcreeping gravity currents of a non-newtonian liquid. Physical Review E,60(6):6960, 1999.IP Grinchik and A Kh Kim. Axial flow of a nonlinear viscoplastic fluidthrough cylindrical pipes. Journal of Engineering Physics and Thermo-physics, 23(2):1039–1041, 1972.Richard W Hanks. The axial laminar flow of yield-pseudoplastic fluids in aconcentric annulus. Industrial & Engineering Chemistry Process Designand Development, 18(3):488–493, 1979.Richard W Hanks and Bharat H Dadia. Theoretical analysis of the turbulentflow of non-newtonian slurries in pipes. AIChE Journal, 17(3):554–557,1971.J Happel and H Brenner. Low reynolds number hydrodynamics. 1965, 1965.Howard S Heaton, William C Reynolds, and William M Kays. Heat trans-101Bibliographyfer in annular passages. simultaneous development of velocity and tem-perature fields in laminar flow. International Journal of Heat and MassTransfer, 7(7):763–781, 1964.Benjamin Herzhaft and E´lisabeth Guazzelli. Experimental study of thesedimentation of dilute and semi-dilute suspensions of fibres. Journal ofFluid Mechanics, 384:133–158, 1999.Richard Holm, S Storey, Mark Martinez, and L Daniel So¨derberg. Visualiza-tion of streaming-like structures during settling of dilute and semi-diluterigid fibre suspensions. Physics of fluids, 2004.CK Huen, IA Frigaard, and DM Martinez. Experimental studies of multi-layer flows using a visco-plastic lubricant. Journal of non-newtonian fluidmechanics, 142(1):150–161, 2007.KOLF Jayaweera and BJ Mason. The behaviour of freely falling cylindersand cones in a viscous fluid. Journal of Fluid Mechanics, 22(04):709–720,1965.Lin Jianzhong, Shi Xing, and You Zhenjiang. Effects of the aspect ratioon the sedimentation of a fiber in newtonian fluids. Journal of AerosolScience, 34(7):909–921, 2003.102BibliographyLaurent Jossic and Albert Magnin. Drag and stability of objects in a yieldstress fluid. AIChE journal, 47(12):2666–2672, 2001.F Julien Saint Amand and B Perrin. Fundamentals of screening: Effectof rotor design and fibre properties. Proceedings of the Tappi PulpingConfernce, pages 941–955, 1999.Donald L Koch and Eric SG Shaqfeh. The instability of a dispersion ofsedimenting spheroids. Journal of Fluid Mechanics, 209:521–542, 1989.P Kumar and BV Ramarao. Enhancement of the sedimentation rates offibrous suspensions. Chemical Engineering Communications, 108(1):381–401, 1991.Peter B Laxton and John C Berg. Gel trapping of dense colloids. Journalof colloid and interface science, 285(1):152–157, 2005.Yu-Quan Liu and Ke-Qin Zhu. Axial couette–poiseuille flow of binghamfluids through concentric annuli. Journal of Non-Newtonian Fluid Me-chanics, 165(21):1494–1504, 2010.S Longo, V Di Federico, and L Chiapponi. Propagation of viscous grav-ity currents inside confining boundaries: the effects of fluid rheology andchannel geometry. In Proceedings of the Royal Society of London A: Math-103Bibliographyematical, Physical and Engineering Sciences, volume 471, page 20150070.The Royal Society, 2015.Sandro Longo, Vittorio Di Federico, Renata Archetti, Luca Chiapponi,Valentina Ciriello, and Marius Ungarish. On the axisymmetric spreadingof non-newtonian power-law gravity currents of time-dependent volume:an experimental and theoretical investigation focused on the inference ofrheological parameters. Journal of Non-Newtonian Fluid Mechanics, 201:69–79, 2013.Michael B Mackaplow and Eric SG Shaqfeh. A numerical study of thesedimentation of fibre suspensions. Journal of Fluid Mechanics, 376:149–182, 1998.A Madani, S Storey, JA Olson, IA Frigaard, J Salmela, and DM Martinez.Fractionation of non-brownian rod-like particle suspensions in a viscoplas-tic fluid. Chemical Engineering Science, 65(5):1762–1772, 2010a.A. Madani, D.M. Martinez, J.A. Olson, and I.A. Frigaard. The stability ofspiral poiseuille flows of newtonian and bingham fluids in an annular gap.Journal of Non-Newtonian Fluid Mechanics, 193:3 – 10, 2013.Ario Madani, James A Olson, D Mark Martinez, and Alan Fung. Novel104Bibliographymethods to characterize and fractionate papermaking fibres. Nordic Pulp& Paper Research Journal, 25(4):448–455, 2010b.Ario Madani, Harri Kiiskinen, James A Olson, and D Mark Martinez. Frac-tionation of microfibrillated cellulose and its effects on tensile index andelongation of paper. NORDIC PULP & PAPER RESEARCH JOURNAL,26(3):306–311, 2011.MCA Maia and CA Gasparetto. A numerical solution for the entranceregion of non-newtonian flow in annuli. Brazilian Journal of ChemicalEngineering, 20(2):201–211, 2003.R Marton and JD Robie. Characterization of mechanical pulps by a settlingtechnique. Tappi Tech Ass Pulp Pap Indus, 1969.ST McComas. Hydrodynamic entrance lengths for ducts of arbitrary crosssection. Journal of Basic Engineering, 89(4):847–850, 1967.Olivier FJ Meuric, Richard J Wakeman, Tak W Chiu, and Kenneth A Fisher.Numerical flow simulation of viscoplastic fluids in annuli. The CanadianJournal of Chemical Engineering, 76(1):27–40, 1998.Taegee Min, Hyoung Gwon Choi, Jung Yul Yoo, and Haecheon Choi. Lam-inar convective heat transfer of a bingham plastic in a circular pipeii. nu-105Bibliographymerical approach hydrodynamically developing flow and simultaneouslydeveloping flow. International Journal of Heat and Mass Transfer, 40(15):3689–3701, 1997.MA Moyers-Gonzalez, IA Frigaard, and C Nouar. Nonlinear stability of avisco-plastically lubricated viscous shear flow. Journal of Fluid Mechanics,506:117–146, 2004.C Nouar. Numerical and experimental investigation of thermal convectionfor thermodependent herschel-bulkley fluid in an annular duct with rotat-ing inner cylinder. In Progress and Trends in Rheology V, pages 297–298.Springer, 1998.Z Nowak and B. Gajdeczko. Laminar entrace region flow of the binghamfluid. Acta Mechanica, 49(3):191–200, 1983.Shinichi Ookawara, Kohei Ogawa, Norman Dombrowski, Esmail Amooie-Foumeny, and Ahmed Riza. Unified entry length correlation for newto-nian, power law and bingham fluids in laminar pipe flow at low reynoldsnumber. Journal of chemical engineering of Japan, 33(4):675–678, 2000.L Paavilainen. The possibility of fractionating softwood sulfate pulp accord-ing to cell wall thickness. Appita, 45(5):319–326, 1992.106BibliographyTasos C Papanastasiou. Flows of materials with yield. Journal of Rheology(1978-present), 31(5):385–404, 1987.JP Pascal. A two-layer model for a non-newtonian gravity current subjectedto wind shear. Acta mechanica, 162(1-4):83–98, 2003.Suhas Patankar. Numerical heat transfer and fluid flow. CRC press, 1980.S Pelipenko and IA Frigaard. Mud removal and cement placement dur-ing primary cementing of an oil well–part 2; steady-state displacements.Journal of Engineering Mathematics, 48(1):1–26, 2004.Jie Peng and Ke-Qin Zhu. Linear stability of bingham fluids in spiral couetteflow. Journal of Fluid Mechanics, 512:21–45, 2004.RJ Poole and RP Chhabra. Development length requirements for fully devel-oped laminar pipe flow of yield stress fluids. Journal of Fluids Engineering,132(3):034501, 2010.RJ Poole and BS Ridley. Development-length requirements for fully de-veloped laminar pipe flow of inelastic non-newtonian liquids. Journal ofFluids Engineering, 129(10):1281–1287, 2007.AMV Putz, TI Burghelea, IA Frigaard, and DM Martinez. Settling of an107Bibliographyisolated spherical particle in a yield stress shear thinning fluid. Physics ofFluids, 20(3):33102–33300, 2008.J Salmela, DM Martinez, and M Kataja. Settling of dilute and semidilutefiber suspensions at finite re. AIChE journal, 53(8):1916–1923, 2007.RK Shah. A correlation for laminar hydrodynamic entry length solutionsfor circular and noncircular ducts. Journal of Fluids Engineering, 100(2):177–179, 1978.VL Shah and RJ Soto. Entrance flow of a bingham fluid in a tube. AppliedScientific Research, 30(4):271–278, 1975.Lawrence F Shampine, Robert Ketzscher, and Shaun A Forth. Using adto solve bvps in matlab. ACM Transactions on Mathematical Software(TOMS), 31(1):79–94, 2005.ZP Shul’Man. Calculation of a laminar axial flow of a nonlinear viscoplasticmedium in an annular channel. Journal of engineering physics, 19(4):1283–1289, 1970.John E Simpson. Gravity currents in the laboratory, atmosphere, and ocean.Annual Review of Fluid Mechanics, 14(1):213–234, 1982.108BibliographyJohn E Simpson. Gravity currents: In the environment and the laboratory.Cambridge university press, 1999.C.M. Sloane. Kraft pulp processing - pressure screen fractionation . AppitaJ., 53(3):220226, 2000.Herve´ Tabuteau, Philippe Coussot, and John R de Bruyn. Drag force on asphere in steady motion through a yield-stress fluid. Journal of Rheology(1978-present), 51(1):125–137, 2007.D Vola, F Babik, and J-C Latche´. On a numerical strategy to compute grav-ity currents of non-newtonian fluids. Journal of computational physics,201(2):397–420, 2004.George C Vradis, John Dougher, and Sunil Kumar. Entrance pipe flow andheat transfer for a bingham plastic. International journal of heat andmass transfer, 36(3):543–552, 1993.Yeh Wang. Finite element analysis of the duct flow of bingham plasticfluids: an application of the variational inequality. International journalfor numerical methods in fluids, 25(9):1025–1042, 1997.YEH Wang. Axial flow of generalized viscoplastic fluids in non-circularducts. Chemical Engineering Communications, 168(1):13–43, 1998.109Frank M White. Fluid mechanics. 5th. Boston: McGraw-Hill Book Com-pany, 2003.110Appendix AFluid rheology and stabilityThe following figures A.1, A.2, A.3 and A.4 show the effect of the outer fluidyield stress, inner and outer fluids plastic viscosity and radius ratio on theflow stability.In addition, we present figure A.5 in which we specify the maximum flowrate Q that may be reached while maintaining a successful separation or inother words, having the residence to settling time ratio of unity.111Appendix A. Fluid rheology and stabilityFigure A.1: The solution of equations 3.22-3.23 for the conditions given asseries 2 in table 3.1. The colour map defines the stability of the flow state,as given in figure 3.3. The dashed line represents the interface position andthe solid line represents the pressure drop rate: (a) τ(2)y = 0, (b) τ(2)y = 0.5,(c) τ(2)y = 2.112Appendix A. Fluid rheology and stabilityFigure A.2: The solution of equations 3.22-3.23 for the conditions given asseries 3 in table 3.1. The colour map defines the stability of the flow state, asgiven in figure 3.3. The dashed line represents the interface position and thesolid line represents the pressure drop rate: (a) µ(1) = 0.375, (b) µ(1) = 0.75,(c) µ(1) = 1.5.113Appendix A. Fluid rheology and stabilityFigure A.3: The solution of equations 3.22-3.23 for the conditions given asseries 4 in table 3.1. The colour map defines the stability of the flow state, asgiven in figure 3.3. The dashed line represents the interface position and thesolid line represents the pressure drop rate: (a) µ(2) = 0.375, (b) µ(2) = 0.75,(c) µ(2) = 1.5.114Appendix A. Fluid rheology and stabilityFigure A.4: The solution of equations 3.22-3.23 for the conditions given asseries 5 in table 3.1. The colour map defines the stability of the flow state,as given in figure 3.3. The dashed line represents the interface position andthe solid line represents the pressure drop rate: (a) κ = 0.6, (b) κ = 0.8, (c)κ = 0.9.115Appendix A. Fluid rheology and stabilityFigure A.5: Maximum allowable flow rates Q in (l/s) for different stainlesssteel spherical particles in a fluid having a yield stress τy = 1 Pa. Theannulus axial length L = 50 cm, radius ratio, κ = 0.6 , outer radius, R =0.127 m, plastic viscosity, µ = 1.0833 Pa.s, spheres density, ρs = 7800kg/m3 and suspension density, ρf = 1000 kg/m3: a) fractionation curve, b)maximum allowable flow rates.116Appendix BFractionation computerprogram operationIn order to find the operating curves for a specific case i.e. (specific fluidsrheology and geometries and properties of targeted fractionated particles),we need to follow some steps to define these properties and parameters andfeed them to the program we developed in this concern:Inputs:• Set the annulus geometry parameters, radius ratio, κ (dimensionless)and outer radius, Rˆ (m).• Set the rotational speed, ωˆ (rad/s) required to fractionate specificparticles that have their specific size, Dˆ (m) and density, ρˆs .• Set the fluids physical and rheological properties, inner fluid yieldstress and plastic viscosity, τˆy1 (pa) , µˆ(1) (Pa.s), outer fluid yield117Appendix B. Fractionation computer program operationstress and plastic viscosity, τˆy2 (pa) , µˆ(2) (Pa.s) and fluids density ρˆf .• Set up ranges for independent variables, interface radius, ri (dimen-sionless) and pressure drop, Gˆ (Pa/m).• Run the program.Outputs:• Velocities profiles (axial and circumferential for both inner and outerfluids).• Flow rates.• Stability conditions.• The operating window.118Appendix CParabolic equations versuselliptic partial differentialequationsThe following general second order partial differential equation can be classi-fied as an elliptic or a parabolic partial differential equation based on the pa-rameter (b2−4ac). The partial differential equation is classified as parabolicif (b2 − 4ac) is zero and classified as elliptic if (b2 − 4ac) is negative.auxx + buxr + curr + dux + eur + fu = gWhere (uxx, uxr and urr) are second order partial derivatives and (ux andur) are first order partial derivatives. (a, b, c, d, e, f and g) are constants orfunctions of independent variables (x and r).The governing equations of the present study are elliptic partial differ-119Appendix C. Parabolic equations versus elliptic partial differential equationsential equations because both second order partial derivatives (urr and uxx)do exist. Elliptic partial differential equations govern the equilibrium prob-lems in which the solution of the partial differential equations is desired ina closed domain subject to a prescribed set of boundary conditions. In thiscase, the solution of the partial differential equations at every point in thedomain depends upon the prescribed boundary condition at every point onthe boundaries.Researchers mostly simplify the Navier stockes equations using the bound-ary layer approximation. They neglect the second order partial derivative(uxx) and the only second order partial derivative left in the equations was(urr). This led to changing the elliptic nature of the partial differentialequation to a parabolic nature.Parabolic partial differential equations govern the marching or propa-gation problems where the solution of the partial differential equations isrequired on an open domain subject to a set of initial conditions and a setof boundary conditions. The solution must be computed by marching out-ward from the initial data surface while satisfying the boundary conditions.120

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

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}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0305821/manifest

Comment

Related Items