MULTIPHASE FLOW THEORY APPLIED TO WATER TREATMENT SYSTEMS AND DEVELOPMENT OF A NEW WALL FILM MODEL FOR MEMBRANE FOULING By Sina Jahangiri Mamouri A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements f or the degree of Mechanical Engineering Docto r o f P hilosophy 2018 ABSTRACT MULTIPHASE FLOW THEORY APPLIED TO WATER TREATMENT SYSTEMS AND DEVELOPMENT OF A NEW WALL FILM MODEL FOR MEMBRANE FOULING By Sina Jahangiri Mamouri Pollution mitigation systems are rapidly evolving to address societal challenges and this work uses multiphase flow theory for studying and improving the performance of novel solid/liquid and liquid/liquid separat ion systems. The first part focuses on large (approx. 29,400 Mt) solid - in - liquid systems used for treating wastewater coming from Combined Sewer Overflows. Use of multiphase theory is needed to study the performance and design changes to a novel system int roduced in Dearborn (MI) and other locations around the USA due to the significant costs of these systems. A recently developed Combined Sewer Overflow detention system (called the Treatment Shaft system) is investigated and its performance for separation of suspended solids from water is evaluated using both Eulerian - Lagrangian and Eulerian - Eulerian multiphase flow approaches for various particle sizes and for a 10 - year, 1 - hour rainstorm event. In order to improve the efficiency of solid separation in th e Treatment Shaft system, the impact of a flocculating agent is evaluated using a Population Balance Model coupled with a multiphase flow solver to model particle aggregations and breakups. The effect of design changes on the performance of the system is e valuated. Multiphase flow theory is used again in the second part to investigate the performance of membranes for liquid - in - liquid separation. An emphasis is put on oil - in - water emulsion separation using crossflow filtration. Various membrane configuration s are first analytically studied, and new analytical solutions are introduced for estimating the flow field in such systems. In order to improve membrane efficiency, the impact of using charged membranes on oil droplet separation is then evaluated in vario us cylindrical membrane configurations. A significant issue in using membranes for oil/water separation is membrane fouling. To incorporate membrane fouling in a computational fluid dynamics simulation, a new model that employs a novel film model is devel oped and implemented in this work for the first time. The multiphase mixture model and the wall film model are coupled with a population balance model. The proposed approach is then used to model and predict the early stages of membrane fouling that accoun t for hydrodynamic interactions. A parametric study is carried out to demonstrate the features of the model and to evaluate the impact of different parameters such as transmembrane pressure, density ratio, viscosity ratio, and contact angle on the performa nce of a membrane separation system. The model is also validated using experimental data and found to correlate well. Copyright by SINA JAHANGIRI MAMOURI 2018 v This dissertation is dedicated to my parents, Azar and Abdollah, for their love, support, and their faith in me. It is also dedicated to my sister, Sahar, who has always been my inspiration, to my lovely niece, Selina , and to my grandmother, Batul, who passed away in 2015. vi ACKNOWLEDGEMENTS I would like to express my sincere gratitude to my advisor, Professor André Bénard , for his support, patience, and encouragement throughout my PhD studies. He has set an example of excellence as a researcher, mentor, and an instructor. In addition, I would like to thank Professor Charles A. Petty . He made me familiar with fundamentals of computational multiphase flow and the mixture model. I would like to thank Professor Farhad A. Jaberi for his guidance and comments as m y dissertation committee member. I would like to especially thank Professor Volodymyr V. Tarabara for his invaluable advice and guidance on membrane filtration technology. I would like to thank Institute for Cyber - Enabled Research (ICER) at MSU for providi ng High - performance computing service and their incredible support. I would like to express my gratitude to Professor Mehrdad T. Manzari and Professor Behshad Shafii my advisors during my M .Sc and B .Sc studies who played important roles in the research pat h I took . Financial support of this research from NSF PIRE project , MSU Sustainability, The James Dyson Foundation, and Dr. Saad Ghalib (from Water Engineering Technologies) is gratefully acknowledged. Last, but not least , I would like to thank my family members and my friends in East Lansing . Their support and encouragement were in the end what made this dissertation possible. vii TABLE OF CONTENTS LIST OF TABLES ................................ ................................ ................................ ................................ ..................... x LIST OF FIGURES ................................ ................................ ................................ ................................ .................. xi KEY TO SYMBOLS AND A BB REVIATIONS ................................ ................................ ................................ xvi Introduction ................................ ................................ ................................ ................................ .. 1 1.1. Overview ................................ ................................ ................................ ................................ ................. 1 1.2. Combined se wer overflows ................................ ................................ ................................ ............. 1 1.2.1. Problem statement ................................ ................................ ................................ ..................... 2 1.2.2. Objectives ................................ ................................ ................................ ................................ ....... 5 1.3. Crossflow filtration membranes ................................ ................................ ................................ .... 6 1.3.1. Produced Water ................................ ................................ ................................ ........................... 6 1.3.2. Problem statement ................................ ................................ ................................ .................. 10 1.3.3. Objectives ................................ ................................ ................................ ................................ .... 12 1.4. Methodology ................................ ................................ ................................ ................................ ....... 13 1.5. Outline of the thesis ................................ ................................ ................................ ......................... 13 Numerical approaches for multiphase flows with interacting dispersed phases ................................ ................................ ................................ ................................ ......................... 15 2.1. Introduction ................................ ................................ ................................ ................................ ........ 15 2.2. Population balance model ................................ ................................ ................................ ............. 17 2.2.1. Standard moment method ................................ ................................ ................................ ... 19 2.2.2. Quadrature method of moments ................................ ................................ ....................... 20 2.2.3. Direct quadrature method of moments ................................ ................................ .......... 22 2.2.4. Breakage models ................................ ................................ ................................ ...................... 24 2.2.5. Coalescence models ................................ ................................ ................................ ................ 27 2.3. Film drainage ................................ ................................ ................................ ................................ ..... 31 2.3.1. Film drainage models ................................ ................................ ................................ ............. 33 2.3.2. Drainage of thin liquid films with mobile interface subjected to surfactants .. 36 Numerical investigation on Treat ment Shaft system: separation efficiency and impact of flocculant polymer agents ................................ ................................ ................................ .......... 47 3.1. Introduction ................................ ................................ ................................ ................................ ........ 47 3.2. System features and geometry ................................ ................................ ................................ .... 49 3.3. Modeling approach ................................ ................................ ................................ .......................... 52 3.3.1. Residence Time Distribution study ................................ ................................ .................. 53 3.3.2. Natural sedimentation ................................ ................................ ................................ ........... 55 3.3.3. Impact of a flocculant on separation efficiency ................................ ............................ 57 3.3.4. The effect of various treatment shaft geometries ................................ ....................... 62 3.4. Numerical simulations details ................................ ................................ ................................ ..... 63 3.4.1. Boundary conditions ................................ ................................ ................................ .............. 64 viii 3.5. Results and Discussions ................................ ................................ ................................ ................. 66 3.5.1. Residence time distribution ................................ ................................ ................................ 66 3.5.2. Natural sedimentation ................................ ................................ ................................ ........... 68 3.5.3. Flocculation process ................................ ................................ ................................ ............... 73 3.5.4. The effect of various treatment shaft geometries ................................ ....................... 83 3.6. Summary and conclusions ................................ ................................ ................................ ............ 85 Analytical solutions for laminar flow in membrane channels with cylindrical symmetry: Single - and dual - membrane systems ................................ ................................ ................. 87 4.1. I ntroduction ................................ ................................ ................................ ................................ ........ 87 4.2. Approach ................................ ................................ ................................ ................................ .............. 89 4.2.1. Analytical solution approach ................................ ................................ ............................... 89 4.2.2. Numerical simulations and validation ................................ ................................ ........... 105 4.3. Results and Discussion ................................ ................................ ................................ ................. 107 4.3.1. Flow through a conventional membrane ................................ ................................ ..... 107 4.3.2. Tubular membrane with an impermeable cylinder at the center ...................... 109 4.3.3. Membrane at the center of a cylindr ical channel with solid (impermeable) wall ................................ ................................ ................................ ................................ ....................... 111 4.3.4. Two tubular co - axial membranes (dual membrane configuration) .................. 112 4.3.5. Effect of Reynolds number and permeability on the solutions ............................ 114 4.4. Conclusions ................................ ................................ ................................ ................................ ....... 117 4.5. Ackno wledgement ................................ ................................ ................................ .......................... 118 4.6. Nomenclature ................................ ................................ ................................ ................................ ... 118 Filtration of charged oil droplets in stationary and rotating tubular membranes ................................ ................................ ................................ ................................ ....................... 121 5.1. Abstract ................................ ................................ ................................ ................................ .............. 121 5.2. Introduct ion ................................ ................................ ................................ ................................ ...... 122 5.3. Modeling Approach ................................ ................................ ................................ ........................ 123 5.3.1. Continuous phase modeling ................................ ................................ .............................. 125 5.3.2. Dispersed phase modeling ................................ ................................ ................................ . 127 5.3.3. Boundary conditions ................................ ................................ ................................ ............ 128 5.4. Numerical Method ................................ ................................ ................................ .......................... 129 5.5. Results and discussions ................................ ................................ ................................ ............... 131 5.5.1. Flow pattern ................................ ................................ ................................ ............................. 131 5.5.2. Wall shear stress ................................ ................................ ................................ .................... 133 5.5.3. Grade efficiency ................................ ................................ ................................ ...................... 136 5.5.4. Effect of electric repelling forces on droplet separation ................................ ........ 136 5.5.5. Reduced Grade Efficiency ................................ ................................ ................................ ... 138 5.6. Summary and conclusions ................................ ................................ ................................ .......... 140 5.7. Nomenclature ................................ ................................ ................................ ................................ ... 140 5.8. Acknowledgments ................................ ................................ ................................ .......................... 142 Modeling membrane fouling using a wall - film model ................................ ............. 143 6.1. Introduction ................................ ................................ ................................ ................................ ...... 143 6.2. Modeling of fluid flow, droplet deposition, and membrane fouling ........................... 146 ix 6.2.1. Phenomenological description of membrane fouling mechanisms ................... 146 6.2.2. Multiphase flow modeling ................................ ................................ ................................ .. 148 6.2.3. Permeate flux o f a membrane ................................ ................................ ........................... 156 6.3. Parametric study ................................ ................................ ................................ ............................ 166 6.4. Results and discussion ................................ ................................ ................................ .................. 169 6.5. Summary and conclusions ................................ ................................ ................................ .......... 179 6.6. Appendix ................................ ................................ ................................ ................................ ............ 180 6.6.1. Droplet breakage ................................ ................................ ................................ ................... 180 6.7. Nomenclature ................................ ................................ ................................ ................................ ... 181 Summary and Conclusions ................................ ................................ ................................ . 185 7.1. The Treatment Shaft system as a CSO detention system ................................ ................ 185 7.1.1. Summary and conclusions ................................ ................................ ................................ .. 185 7.1.2. Limitations ................................ ................................ ................................ ............................... 187 7.1.3. Future work ................................ ................................ ................................ ............................. 188 7.2. Crossflow filtration membranes ................................ ................................ ............................... 189 7.2.1. Summary and conclusions ................................ ................................ ................................ .. 189 7.2.2. Limitations of membrane fouling model ................................ ................................ ...... 191 7.2.3. Future work ................................ ................................ ................................ ............................. 191 REFERENCES ................................ ................................ ................................ ................................ ..................... 194 x LIST OF TABLES Table 2 - 1: Examples of multiphase flows in natural, biological, and industrial systems [79] ................................ ................................ ................................ ................................ ................................ ..... 17 Table 3 - 1 : Residence time distribution key terms for characterization, from [ 152 , 153 ] .... 55 Table 3 - 2 : Residence time distribution key terms for characterization of present study ..... 67 Table 4 - 1 : Definition of dimensionless variables ................................ ................................ .................. 96 Table 4 - 2 : Deta ils of the different cases and their boundary conditions ................................ ... 105 Table 4 - 3 : Normalization of flow field variables ................................ ................................ ................. 107 Table 4 - 4: Effect of membrane permeability values on the accuracy of axial velocity profile at ................................ ................................ ................................ ................................ ......... 117 Table 5 - 1 : Additional accelerations on the particles ................................ ................................ .......... 1 28 Table 5 - 2 : Parameters for different simulations ................................ ................................ ................. 131 Table 6 - 1 : Algorithm to calculate the membrane permeate flux ................................ .................. 166 Table 6 - 2: details of variations used for parameter study. For each parameter, two variations of base case are introduced and compared with the base case. ................. 169 Table 6 - 3 : Steady state times and fluxes for different cases studie d for parametric study 173 xi LIST OF FIGURES Figure 1 - 1 : Distribution of Combined Sewer Systems (CSSs) in the United States [ 1 ] .............. 2 Figure 1 - 2: Map of produced water sample locations from conventional hydrocarbon wells as of December 08, 2017. Adapted from USGS, Energy Resource Program .................... 7 Figure 2 - 1: Number of publications (a), and citations (b) per year with keyword of ................................ ................................ ................................ ................. 16 Figure 2 - 2 : (a) nondeformable and (b) deformable droplets [ 135 ] ................................ ............... 32 Figure 2 - 3 : (a) Immobile, (b) Partially mobile, (c) fully mobile interfaces [ 135 ] ..................... 32 Figure 2 - 4 : Schematic of a liquid film between two deformable droplets ................................ .. 36 Figure 3 - 1: Schematic of Treatment Shaft of present study. a) Front view, b) isometric view from top corner b) top view, c) isometric view from bottom corner. Units are in meters. ................................ ................................ ................................ ................................ ..................... 51 Figure 3 - 2: Flow rate at the inlet of the treatment shaft in a 5 - hour timespan. This represents the varying flow rate of CSO for a 10 - year, 1 - hour rainstorm, and lastin g for about 5 hours. ................................ ................................ ................................ ................................ 52 Figure 3 - 3: Particle size distribution used at the inlet of the Treatment Shaft evaluating various geometries ................................ ................................ ................................ .............................. 63 Figure 3 - 4 : Particle size distribution of CSO samples; mean value adapted from [ 162 ] . ...... 65 Fig ure 3 - 5: Resident time distribution of Treatment Shaft. The numerical results of the present study is compared with the available data from [7]. ................................ ............. 67 Figure 3 - 6: Separation efficiency of the Treatment Shaft system during a captured rainstorm CSO for various particle sizes ranging from to . ................... 68 Figure 3 - 7: Separation efficiency (area - normalized) of the Treatment Shaft system during a captured rainstorm CSO for various particle sizes ranging from 1 µ m to 10320 µ m. 69 F igure 3 - 8: The average settling time and velocity in the Treatment Shaft system during a captured rainstorm CSO for various particle sizes ranging from 1 µ m to 10320 µ m. 70 Figure 3 - 9: The average particle time and velocity in the Treatment Shaft system at the end of a rainstorm CSO for various particle sizes ranging from 1 µ m to 10320 µ m. ........... 71 Figure 3 - 10: Contour of the separation efficiency for various steady - state operation conditions velocities, for particle sizes ranging from 1 µ m to 512 µ m ............................. 72 xii Figure 3 - 11: Particle size distribution of CSO samples used for flocculation analysis; mean value adapted from [162] . Highlighted area is based on the results of natural sedimentation. ................................ ................................ ................................ ................................ ...... 73 Figure 3 - 12: Velocity contour (m/s) and streamlines of cross - section of the Treatment Shaft at various timeframes during a rainstorm event. Sewer overflow enters the shaf t from top left of the domain. (a) Velocity contour at t = 500s, and (b) Streamlines at t = 500. (c) Velocity contour at t = 2000s, and (d) Streamlines at t = 2000s. (e) Velocity contour at t = 5000s. (f) Streamlines at t = 5000s. .......................... 74 Figure 3 - 13: Q - criterion of the Treatment Shaft at various timeframes during a rainstorm event. Sewer overflow enters the shaft from top left of the domain. (a) t = 500s, (b) t = 2000s, (c) t = 5000s, (d) t = 10000s. ................................ ................................ ......................... 76 Figure 3 - 14: Concentration of polymer flocculant (kg/kg) at the cross - section of the Treatment Shaft at various timefram es during a rainstorm event. Sewer overflow enters the shaft from top left of the domain. (a) t = 500s, (b) t = 2000s. ....................... 78 Figure 3 - 15: Sauter mean diameter (m) of polymer flocculant at the cross - section of the Treatment Shaft at various timeframes during a rainstorm event. Sewer overflow enters the shaft from top left of the domain. (a) t = 500s, (b) t = 2000s, (c) t = 5000s, (d) t = 10000s. ................................ ................................ ................................ ................................ ....... 79 Figure 3 - 16:Volume fraction of suspended solids at the cross - section of the Treatment Shaft at various timeframes during a rainstorm event. Sewer overflow enters the shaft from top l eft of the domain. (a) t = 500s, (b) t = 2000s, (c) t = 5000s, (d) t = 10000s. ................................ ................................ ................................ ................................ ................................ ..... 81 Figure 3 - 17: Distribution of Sauter mean diameter in the Treatment Shaft during a 10 - year, 1 - hour rainstorm. In various time frames. ................................ ................................ ................ 82 Figure 3 - 18: Separation efficiency comparison of original Treatment Shaft design and Treatment Shafts with smaller and larger baffle gap ................................ ............................ 83 Figure 3 - 19: Separation efficiency comparison of original Treatment Shaft design x = 0m), and Treatment Shafts with baffle - = +2.9m) ................................ ................................ ................................ ................................ .................. 84 Figure 3 - 20: Separation efficiency comparison of original Treatment Shaft design and Treatment Shaft wit h straight inlet and outlet channels ................................ ..................... 85 Figure 4 - 1: Geometry of the four membrane systems considered. Front (left) and cross sectional (right).views are given for : (a) single channel tubular membrane, (b) a tubular membrane with a co - axial solid cylinder at its center, (c) a tubular membrane located at the center of a channel with impervious walls, and (d) two tubular co - axial membranes. ................................ ................................ ................................ ........... 91 xiii Figure 4 - 2: Comparison of the present analytical method with numerical results and analytical approach of Karode [185] for case A: (a) Axial flow rate, (b) Axial velocity profile in diff erent cross sections, (c) Radial velocity profile at z/L=0.5, and (d) pressure changes along the axial axis in channel. Inlet flow rate is 6.67×10 - 6 m 3 /s. ................................ ................................ ................................ ................................ ................................ ... 109 Fig ure 4 - 3: Comparison of the analytical method (this work), numerical results (also this work), and analytical approach of Karode [185] for case B: (a) Axial flow rate, (b) Axial velocity profile in different cross sections, (c) Radial velocity profile at z /L=0.5, and (d) pressure changes along the axial axis in channel. Inlet flow rate is 6.67×10 - 6 m 3 /s. ................................ ................................ ................................ ................................ ........................ 110 Figure 4 - 4: Comparison of the analytical method (this work), numerical results (also this work), and analytical approach of Karode [185] for case C: (a) Axial flow rate, (b) Axial velocity p rofile in different cross sections, (c) Radial velocity profile at z/L=0.5, and (d) pressure changes along the axial axis in channel. Inlet flow rate is 6.67×10 - 6 m 3 /s. ................................ ................................ ................................ ................................ ........................ 112 Figure 4 - 5: Comparison of the analytical method (this work), numerical results (also this work), and analytical approach of Karode [185] for case D: (a) Axial flow rate, (b) Axial velocity profile in different cross sections, (c) Radial velocit y profile at z/L=0.5, and (d) pressure changes along the axial axis in channel. Inlet flow rate is 6.67×10 - 6 m 3 /s. ................................ ................................ ................................ ................................ ........................ 113 Figure 4 - 6: Comparing of normalized axial velocity at z/L=0.8 for different Reynolds numbers in (a) Case A, (b) Case B , (c) Case C, and (d) Case D ................................ .......... 114 Figure 4 - 7: Comparison of normalized axial velocitie s of present analytical approach with numerical results for different permeability values of 10 mD, 20mD, and 30mD - (a) Case A, (b) Case B , (c) Case C, and (d) Case D ................................ ................................ ........ 116 Figure 5 - 1 : Schematic of the computational domain (not to scale) ................................ ............. 124 Figure 5 - 2 : Schematic of t he computational domain of OCFF (not to scale) ............................ 125 Figure 5 - 3: Normalized Axial velocity at longitudinal position of z/L=0.5 in the membrane ................................ ................................ ................................ ................................ ................................ ... 132 Figure 5 - 4: Normalized radial velocity at longitudinal position of z/L=0.5 in the membrane ................................ ................................ ................................ ................................ ................................ ... 132 Figure 5 - 5: Normalized tangential velocity at longitudinal position of z/L=0.5 in the membrane ................................ ................................ ................................ ................................ ............. 133 Figure 5 - 6: Average wall shear stress components for SCFF and OCFF membranes, , ................................ ................................ ................................ ................................ ........ 135 xiv Figure 5 - 7 : Overall wall shear stress and its direction for SCFF and OCFF membranes ...... 135 Figure 5 - 8: Grade efficiency of membranes with and without charge on the droplets. The number on the legend shows the charge density. ................................ ................................ . 137 Figure 5 - 9: Grade efficiency of membranes with charge densities of 2 and 4 C/m 3 on the droplets. The number on the legend shows the charge density. ................................ ..... 137 Figure 5 - 10: Grade efficiency of a rotating OCFF membrane (1500 rpm) versus stokes number ................................ ................................ ................................ ................................ ................... 139 Figure 6 - 1 : Schematic of the operation of a crossflow filtration (CFF) membrane ............... 144 Figure 6 - 2: Droplet - membrane interaction scenarios in a computational cell; (a) There is no droplet in the cell, (b) there are droplets in the cell, but no droplets sitting on the membrane surface, (c) there are droplets both in the cell and sitting on the membrane su rface, but no droplet blocks the pores, (d) there are droplets in the cell, sitting on the membrane, and blocking and entering into the pores, (e) a film is partially formed on the membrane surface, and (f) a film is formed and entirely covered the cell, b locking and entering into the pores. ................................ ...................... 148 Figure 6 - 3 : Schematic of the top view and side view of droplets in a computational cell. .. 158 Figure 6 - 4: Schematic of the top view of droplets in a computational cell on the membrane surface (on the left); (a) the top view of a single droplet shows the diameter of the deposited dr oplet ( ), and the diameter of the wetting area ( ). (b) the side view of a droplet shows the radii of the sitting droplet and the wetting area and their relationship with angle , which is defined as , where is the contact angle. 160 Figure 6 - 5: Schematics of (a) physical representation of droplets on membrane surface; (b) computational representation of droplets on membrane surface. The grey checkerboard background mimi cs the computational faces on membrane surface. ................................ ................................ ................................ ................................ ................................ ... 161 Figure 6 - 6: Schematic of a thin two - dimensional film with the thickness of on a surface. Momentum changes du e to (I) Normal bulk flow pressure gradient, (II) Normal gravity - driven pressure gradient, (III) Tangential pressure gradient, (IV) Surface curvature, (V) Bulk flow shear stress, (VI) Wall shear stress, (VIII) Droplet deposition on surface. See Eq. ( 6 - 46 ) for more details. ................................ ....................... 163 Figure 6 - 7 : Droplet size distribution of sample emulsion used for simulations. .................... 167 Figure 6 - 8: Schematic of computational domain demonstrating the boundary conditions ................................ ................................ ................................ ................................ ................................ ... 167 Figure 6 - 9: Effect of various parameters on the permeate flux of a CFF membrane. (a) Transmembrane pressure difference, (b) Membrane resistance, (c) Feed volume xv fraction, (d) Feed velocity. See Table 6 - 2 for variations used to study each parameter. Percentages ind icate . ................................ ................................ ............... 170 Figure 6 - 10: Effect of various parameters on the permeate flux of a CFF membrane. (a) Contact angle, (b) Droplet size distribution, (c) Rela tive density, (d) Relative viscosity. See Table 6 - 2 for variations used to study each parameter. Percentages indicate . ................................ ................................ ................................ ................................ 172 Figure 6 - 11: Contour of mass per area of the deposited droplets on the membrane surface in various timeframes: 250s, 500s, 1000s, 2000s. ................................ ................................ 174 Figure 6 - 12: Contour of permeate flux at the mem brane surface in various timeframes: 250s, 500s, 1000s, 2000s. ................................ ................................ ................................ ............... 175 Figure 6 - 13: Contour of Sauter mean diameter near the membrane surface at steady state permeate flux ................................ ................................ ................................ ................................ ...... 176 Figure 6 - 14: Contour of the oil volume fraction near the membrane surface for three different inlet velocities. The plots are not in scale. The plots are not in - scale. The vertical axis shows 10% of the entire channel height (i.e. 0.3mm), and the horizontal axis is the entire channel length (i.e. 80mm). The lower edge represents the membrane surface. ................................ ................................ ................................ ............................ 177 Figure 6 - 15 : Sauter mean diameter at the outlet of the membrane channel. .......................... 178 Figure 6 - 16 : Sauter mean diamete r of droplets on the membrane surface. ............................. 179 xvi KEY TO SYMBOLS AND ABBREVIATIONS Shear rate (s - 1 ) Stress tensor (Pa) Film thickness Electric field vector (V/m) Force vector (N) Acceleration vector (m 2 /s) B ody force (N) N ormal vector V elocity vector (m/s) A rea (m 2 ) Hamaker constant (J) Capillary number Force (N) Grade Efficiency L ength of the membrane (m) xvii P ressure (Pa) V olumetric flow rate (m 3 /s) Membrane resistance (m - 1 ) R adius (m) Reynolds Number E igenvalue Volume (m 3 ) Location (m) Droplet diameter (m) Coalescence frequency Breakage kernel Permeate Flux (LMH) permeate velocity (m/s) Mass (kg) Number density Pressure (Pa) Charge (C) xviii V olumetric flow rate (m 3 /s) R adius (m) T hickness (m) Time (s) V elocity (m/s) Volume of a droplet (m 3 ) Coalescence Kernel P ermeabili ty (m 2 ) Volume fraction Probability density function M embrane thickness (m) Permittivity (F/m) A ngular coordinate Collision efficiency D ynamic viscosity (Pa.s) K inematic viscosity (m 2 /s) xix D ensity (kg/m 3 ) Surface tension (N/m) Contact angle Membrane blockage Electric potential (V) Stress tensor (Pa) CFF - Crossflow filtration CFV - Crossflow velocity DQMOM - Direct quadrature method of moments DSD - Droplet size distribution MPE - Mean percentage error OCFF - Optimized CrossFlow Filtration PBE - Population balance equation PBM - Population balance model PRESTO - Pressure staggering option QUICK - Quadratic Upstream interpolatio n for convective kinetics RGE - Reduced Grade Efficiency xx SCFF - Simple CrossFlow Filtration TGE - Total Grade Efficiency TMP - Transmembrane pressure UDF - User - defined function UDM - User - defined memory variable 1 Introduction Overview Clean Water is essential for population health, and daily operations. R eturn of thermally, organically, or chemically polluted water to the environment continues to be a major concern. Water pollution directly and indirectly harms public health and damages our limited sources of freshwater. To address this issue, p ollution mitigation systems are rapidly evolving. Water pollutants are in different forms such as solutes, suspended solids, and liquid droplets. Suspended solids a nd liquid droplets as two common sources of pollutions, despite their different characteristics , raise a common issue which is the separation of dispersed phase (i.e. pollutant) from a continuous phase (i.e. water). This study in two main parts utilizes th e multiphase flow theory and develops computational fluid dynamics models to simulate two different water treatment systems, predict their performance, and improve their design and operation conditions. The first part focuses on the Treatment Shaft technol ogy to separate suspended solids coming from combined sewer overflows. The second part investigates the crossflow membrane systems for separation of oil - in - water emulsions . The details are discussed below. Combined sewer overflows The first part of the th esis focuses on the use of multiphase flow theory and Computational Fluid Dynamics (CFD) for evaluating a water treatment technology for solid - in - liquid systems. 2 1.2.1. P roblem statement Combines sewer systems (CSS) and Separate sanitary systems (SSS) are two main types of public wastewater collection systems in the United States. CSSs were the earliest sewer systems constructed in the U.S. in 1900 - 1950. In the contrary to SSS, CSS was designed to collect both wastewater and stormwater in a single pipe. Rain an d snow events may cause excess flow through the system and can result in Combined Sewer Overflows (CSOs). Most CSS facilities discharge CSOs to surface waters such as rivers, coastal waters, etc. As of 2015 , there were 859 active permits for CSO discharges in the U.S. Figure 1 - 1 shows the distribution of CSSs in the United States [ 1 ] . The annual discharge of CSO in the U.S. is about 850 billion gallon [ 1 ] . The U.S. Environmental Protection Agency estimates that more than 770 cities in the United States suffer from water pollution resulting from CSOs. An approach to mitigate the effect of combined sewer overflows is to store rainfall runoffs temporarily through storage 3 and treatment of solids in flows that exceeds storage capac ity by sedimentation [ 2 ] . Content of wastewater solids are characterized as Total Suspended S olid (TSS), which assesses the quality of wastewater after and prior to treatment. The co ncentration of TSS is evaluated by standard laboratory tests [ 3 ] . The removal of these solid particles is important in order to reduce the pollution in receiving waters. This can be achieved by settli ng, filtration, coagulation and flocculation [ 4 ] . Settling tanks are one of the most practical ways of pollution management since they provide volumetric control of flows and treatment by settling [ 5 ] . The efficiency of a sedimentation tank to remove suspended solids depen ds on configuration of tank (solid concentration, water flow rate), and characteristics of the particles (shape, size, drag and buoyancy forces) [ 6 ] . The detention facility can contribute to removing solid particles from the overflow through sedimentation. Detention tanks can be designed in s everal ways [ 7 ] . One approach is to use in - system storage with real - time monitoring structures [ 8 ] , but apparentl y is vulnerable to malfunction in harsh environment. Another approach is to use offline detention basins which are typically constructed below - grade, and work based on sedimentation basins or clarifiers [ 9 ] . An issue with these systems is that populated urban environments limit the site availability for implementation of such systems. Another approach is based on deep storage t unnels built far below - grade where other utilities are not located. Such tunnels are subjected to surges and trapped air during the filling process. This may lead to formation of geysers in shafts and causes return of low - quality water to grade [ 10 ] . Another issue in a long tunnel is removal of solids when it is dewatered [ 11 ] . Inadequate tr eated combined sewer overflows (CSOs) and discharges of urban stormwater are among major causes of long - term persistence of poor water quality in receiving waters [ 12 ] . 4 A Treatment Shaft system [ 13 , 14 ] provides an alternative solution that can mitigate mos t of the aforementioned concerns. It was designed to achieve the least possible contact time under flow rates over 1500 million gallons per day (i.e. ~65 m 3 /s). A Treatment Shaft system is a detention basin, with cylindrical shape which has a relatively sm all footprint and is able to store a large volume of wastewater by extending underground. It employs baffles and partitions to ensure a very slow upward velocity flow, which enhances the settling efficiency. In addition, accumulated solids at the bottom of the shaft are discharged into the sanitary sewer system using chopper pumps. A flushing system uses high - pressure nozzles to clear the bottom of the Treatment Shaft during a final rinse cycle [ 7 ] . A Treatment Shaft system has been in operation in Dearborn, Michigan since 2010 and was designed to meet Michigan Department of Environmental Quality (MDEQ) presumptive criteria of 10 - minute detention of 10 - year, 1 - hour peak flow [ 7 ] . The issue in the Treatment Shaft systems is that their design and geometry details strongly depend on operation flow rates, particle size distribution of suspended solids, and size and scale of the Treatment Shaft. These parameters may significantly change from one site to another and can remarkably influence its settling performance. Moreover, the eff iciency of such systems varies depending on the rainstorm event and its subsequent combined sewer overflow. Unfortunately, due to high construction costs (more than $10M), experimental analysis of a Treatment Shaft of 1:1 scale is practically impossible. A s an alternative, performing experiments on prototype systems (e.g. 1:19 scale [ 7 ] ) is beneficial; however, although scaling based on dynamic similarities may preserve some flow features, it defectively captures solid - fluid interaction s , and cannot illustrate the impact of treatment additives such as flocculating agents. 5 In order to address the issues above, computer simulations of 1:1 scale are needed . Such simulations can incorporate variety of P article S ize D istributions (PSDs) , CSO flow rates, and possible geometry modifications, based on site details and rainstorm event. A better understanding of particle settling in a Treatment Shaft system could lead t o improved designs based on PDSs and common rainstorm events in a given location. In addition, computer simulations help to evaluate the impact of aggregation additives such as flocculating polymer agents. 1.2.2. Objectives The Treatment Shaft is a relatively new technology with a few systems constructed and in operation. A very limited information about its performance for particle settling is available , whereas its characteristics may change for various rainstorms and geometry modifications. On the other hand, although flocculating agents are recognized to improve particle settling, their impact on particle separation in the Treatment Shaft as a CSO detention system is unknown. The performance of the Treatment Shaft system is modeled using computational fluid dy namics (CFD) for different design alternatives such as various shaft sizes, and different baffle locations. The multiphase flows arising in these systems are modeled using a Lagrangian - Eulerian CFD approach to simulate the interaction of the continuous and dispersed phases for an actual rainstorm overflow situation. The simulations are carried out based on a 10 - year rainfall storm event. The efficiency of Treatment Shaft for particle separation is reported for the several design variations. Additionally, a Treatment Shaft is thoroughly evaluated to be use d as a CSO detention system with a flocculating agent. For this 6 purpose, an Eulerian - Eulerian CFD approach is employed, which is coupled with a Population Balance Model (PBM) to simulate particle - particle in teractions such as coalescence and breakage of flocs. Using CFD, the use of Treatment Shaft system with a flocculating agent during a 10 - year, 1 - hour rainstorm event is modeled, and the efficiency enhancement is evaluated. Crossflow filtration membranes Th e second part of the thesis focuses on the use of multiphase flow theory and Computational Fluid Dynamics (CFD) for evaluating a water treatment technolo gy for liquid - in - liquid systems, with an emphasis on oil - water mixtures. 1.3.1. Produced Water Water plays an important role in all production methods employed for fossil fuels extraction. It is either a byproduct to be cleaned up, or it enables the extraction of methane and natural gas through hydraulic fracturing, and maintain flow by sustaining reservoir pressu re [ 15 , 16 ] . Oil and gas industry in United States hav e to treat 7 to 8 barrels of water (called produced water) per each barrel of oil extracted [ 15 ] . The rapid industrial growth, such as in oil and gas, petrochemical, pharmaceutical, metallurgical and food industries, has led to the large production of oily wastewater. Necessity to treat the oily wastewater is an inevitable challenge. Produced water i s the largest waste stream generated in oil and gas industries . Because of the rapidly increasing volume of waste all over the world in the current decade, the outcome and effect of discharging produced water on the environment has lately become a signific ant environmental issue [ 17 ] . Produced water contains various 7 org anic and inorganic components and discharging it to the environment can pollute surface underground water and soil. The permit oil and grease (O&G) limit for treated produced water discharge offshore in Australia are 30 mg/L (on daily average), and 50 mg/L instantaneous [ 18 ] . The maximum daily limit in United States is 42 mg/L and the monthly average of 29 mg/ L. Figure 1 - 2 shows the map of produced water sample locations from conventional hydrocarbon wells as of December 2017 . Reynolds [ 19 ] presented different factors which can influen ce the amount of produced water production from a well: (i) Method of well drilling: a horizontal well produces a higher rate than a vertical well at similar drawdown. (ii) Location of well within homogenous or heterogeneous reservoirs (iii) Different type of completion 8 (iv) Single zone and commingled (v) Type of water separation technologies (vi) Water injection or water flooding for enhancing oil recovery (vii) Poor mechanical integrity (viii) Underground communications Produced water is a mixture of inorganic and organic components. The followin g factors affect the physical and chemical properties of produced water [ 20 ] : (1) geological location of the field ; (2) geological formation of the field; (3) lifetime of the reservoir; (4) type of produced hydrocarbons. The composition of produced water from different sources can differ quantitatively but is similar qualitatively [ 21 ] . The major compounds of produced water are [ 22 ] : Dissolved and dispersed oil compounds Dissolved formation minerals Production ch emical compounds Production solids (e.g. formation solids, bacteria, waxes, asphaltenes) Dissolved gases One of the components is dispersed oil, which consists of small droplets of oil suspended in the produced water. The amount of these droplets depends o n the density of oil, the shear 9 history of the droplet, the amount of oil precipitation, and interfacial tension between oil and water [ 23 ] . Polyaromatic hydrocarbons (PAHs) and some of the heavier alkyl phenols are less soluble in the water and appear as dispersed oil in produced water [ 24 ] . The concentration of PAHs and C6 C9 alkylated phenols is correlated to amount of di spersed oil droplets in produced water. Produced water is considered as oilfields waste and should be minimized, reused or recycled, or if neither is practical, disposal is the final option [ 25 ] . Some of the available options for produced water management are ass follow [ 26 ] : (1) injection (reusing the produced water); (2) discharge after meeting the environmental reg ulations; (3) reuse in other oil and gas operations; (4) consume in other beneficial uses after treatment, such as irrigation. The general objectives of produced water treatments are: De - oiling: removing dispersed phase of oil Soluble organics removal Suspended particles and sand removal Dissolved gas removal Desalination Softening Removing other components such as Naturally occurring radioactive materials (NORM) 10 De - oiling of produced water represents the largest scale routine separation of oil - water mi xtures. Centrifugation [ 27 ] , air flotation [ 28 , 29 ] and hydrocyclone separation [ 30 ] are the current methods of oil removal from produced water [ 31 - 34 ] . H oweve r , the efficiency of these methods decreases dramatically for droplets smaller than approximately 20 µm [ 35 ] . Consequently, to comply with environmental regulations multiple - step conventional treatment is applied, which results in high costs. More effective separation of oil - water mixtures into water and oil phases is a key technology to decrease the environmental footprint of the oil in dustry and to recover large amounts of oil that are now lost due to spilling or by disposing of oil as a waste stream. 1.3.2. Problem statement Membranes can potentially remove very small oil droplets from wat er [ 36 - 39 ] . Examples of promising membrane - based oil - water separation processes include reverse osmosis [ 36 ] , flocculation followed by microfiltration or ultrafiltration [ 37 , 40 ] , microfiltration [ 41 ] , membrane distillation [ 42 ] and ultrafiltration [ 43 , 44 ] . A common membrane configuration is when the feed passes tangentially over membrane surface (crossflow filtration) and particles larger than the membrane pores are rejected, while the continuous phase flows through the membrane pores. Crossflow filtration (CFF) membranes have applications in petroleum industries, wastewater treatment, etc. [ 45 - 49 ] . Studies on the filtration of oil droplets using membranes observed that the droplets rejected by a membrane form a deposit layer at the membrane surface [ 50 - 56 ] , which is usually referred fouling mechani sms can be observed on a membrane surface, which are thoroughly 11 discussed by Dickhout et al [ 57 ] . A lthough numerous studies investigated oil - water separation by membranes, the research was mostly empirical and fell short of expectations due to membrane fouling. While membranes show high flux an d rejection initially, the separation performance deteriorates rapidly with an onset of fouling that remains the single most significant barrier for the commercial success of membrane - based approaches to oil - water separation. A better understanding of memb rane fouling by oil droplets can provide routes to improve existing techniques and develop new technologies. Motin et al. [ 58 ] discussed the effect of membrane rotation on the performance of the membrane, and efficiency of the droplet separation. In addition to experimental efforts on the membrane fouling, other analytical modeling and numerical sim ulations were also presented to understand the fouling phenomenon and flux decline in membranes. Grenier et al. [ 59 ] proposed a methodology to quantify the fouling mechanism. Darvishzadeh et al. [ 60 ] studied the behavior of an oil droplet at a membrane pore entrance using Computational Fluid Dynamics. Won et al. [ 61 ] computationally investigated the effect of patterned membranes on particle depositions. Hou et al. [ 62 ] provided a combined blocking and cake filtration model. Zhang et al. [ 63 ] simulated the membrane fouling under a baffle - filled flow. Several other analytical or numerical studies were carried out in which the fluid dynamics o f the flow in a CFF membrane [ 64 - 66 ] and membrane fouling [ 67 - 71 ] were discussed . An extensive study of hydrodynamics of crossflow filtration membranes and membrane foulin g is essential to identify the enhanced operation conditions to mitigate fouling. A physical model of membrane fouling in oil - in - water emulsions could provide a better understanding of this phenomenon in membrane filtration and could lead to an improved CF F membrane system designs. 12 1.3.3. Objectives Despite the mentioned studies which mainly focus on modeling the fouling process using the mass transfer models, there is a lack in literature for a CFD approach to model membrane fouling of interacting droplets in a c ontinuous phase. Mass transfer equations are extensively used for filtration of solutes, and in some cases, they are empirically calibrated for separation of emulsions; however , they are not capable of effectively modeling the interaction of droplets within the continuous phase. In the second part of this thesis, the m ultiphase flow theory is used to investigate the performance of membranes for liquid - in - liquid separation; specifically, for oil - in - water emulsion separation using crossflow filtration. V arious membrane configurations are analytically studied, and new analytical solutions are introduced for estimating the flow field in such systems. To improve membrane efficiency, the impact of using charged membranes on oil droplet separation is then eval uated in various cylindrical membrane configurations. To incorporate membrane fouling in CFD simulation, a new model that employs a novel film model is developed and implemented in this work for the first time. The multiphase mixture model and the wall fil m model are coupled with a population balance model. The proposed approach is then used to model and predict the early stages of membrane fouling that account for hydrodynamic interactions. A parametric study is carried out to demonstrate the features of t he model. The proposed model is also validated using experimental dat a. Suggestions for future works are discussed at the end. 13 Methodology This research is based on analytical, theoretical, and numerical study of separation of suspended solids and liquid droplets from water. For provided analytical solutions, Wolfram Mathematica is employed. For both parts of this thesis , CFD simulations are carried out using ANSYS Fluent software. Pre - processing of the simulations are conducted using ANSYS DesignModeler, ANSYS ICEM CFD, and ANSYS Meshing tool s. In problems involving turbulence flow, the Reynolds Averaged Navier - Stokes (RANS) equations are modeled using Realizable model with scalable wall function. In Lagrangian - Eulerian approach, Discrete Phase Model (DPM) is employed to compute the force balance, and trajectories of dispersed phase. In Eulerian - Eulerian approach, the multiphase mixture model is utilized in combination with Population Balance Model (PBM) to account for droplets coalescence and breakag e. The User Defined Functions (UDFs) are used to introduce new models and modifications to ANSYS Fluent software. Post - processing of data is carried out using ANSYS CFD - Post, MATLAB, and in - house Python codes. Outline of the thesis This thesis is organize d in seven chapters. An overview of the discussed problems, including a brief background and problem statement were discussed in the present chapter. Chapter 2 provides a review of numerical approaches for multiphase flows with dispersed phases. Details of population balance model and its solution methods are also discussed in Chapter 2. In Chapter 3, a comprehensive numerical investigation on performance of the Treatment Shaft system is presented for a 10 - year, 1 - hour rainstorm event. Additionally, the imp act of polymer flocculant on the settling efficiency of the Treatment Shaft system is 14 studied. In Chapter 4, cylindrical crossflow filtration membranes are studied, and novel analytical solutions are derived for various configurations of cylindrical CFF me mbrane systems. THe introduced solutions are validated using numerical simulations. Chapter 5 addresses the filtration of charged oil droplets in stationary and rotating cylindrical CFF membranes, and the effect of electrical repelling forces on oil separa tion is discussed. In Chapter 6, a novel model coupled with CFD and population balance model is presented which can effectively predict the early stages of membrane fouling and permeate decline in CFF membrane systems. In Chapter 7, a summary of conclusion s is provided, and future research paths are suggested. 15 Numerical approaches for multiphase flows with interacting dispersed phases Introduction Particles and droplets are present in variety of the natural and industrial flows. They can significantly affect the behavior of such systems. There are different approaches to model multiphase flows. Eulerian - Lagrangian approach is useful for dilute flows and particles are tracked in the level of a single particle. Therefore, equation of motion for a particle is solved. In Eulerian - Eulerian approach, all phases are treated as continuous phases, and the dispersed phase equations are averaged in each cell. This approach is mostly suitable for dense flows. The other approach, which in fact belongs to Eulerian - Eul erian framework, is Volume of Fluid (VOF) approach. In contrast to Eulerian - Eulerian approach, in VOF, the phases are not interpenetrating. However, dispersed phase modeling is essential for all the aforementioned approaches. The simplest, and most common assumption is to describe dispersed phase as spherical particles. The study of these complex particulate flows requires an understanding of the behavior of the population of particles. The population is defined as the density of an extensive variable. Dens ity is the variable that usually employed for this purpose, however, sometimes variables such as mass or volume of the particles are incorporated [ 72 ] . Therefore, in multiphase flows involving a size distribution, in addition to momentum, mass and energy balance equations, a balance equation is essential to describe the changes in particle population. Populat ion balance model has many applications in engineering sciences such as: crystallization [ 73 ] , dissolution [ 74 ] , granulation [ 75 ] , drying [ 76 ] , mixing [ 77 ] , polymerization [ 78 ] , multiphase flows [ 79 ] , reacting flows [ 80 ] , 16 fermentation [ 81 ] , cell growth, division and death [ 82 ] . Figure 2 - 1 (a) and (b) shows the number of publications and ci [ 83 ] . The analysis of flows with gas bubbles, liquid drops and solid particles, and development of generic computational approaches are inevitable in multiphase flow problems. Some 17 examples of multiphase flows in natural, biological, and industrial systems are tabulated in Table 2 - 1 [ 79 ] . Gas - liquid Liquid - liquid Gas - particle Liquid - particle Three - phase Natura Rain droplets, mist formation Oil - water mixture Sand storms Soil erosion Oil - water - sand mixture Biological Aerosols Water and hydrophobic phases (lipids or membranes) Dust particles Nanoparticles in blood flow Plasma Industrial Desalination systems Fuel - cell systems Spray drying Sewage treatment plants Air lift pumps Population balance model The population balance equation (PBE) is a balance equation for the number density, , which is the sized droplet per unit volume of continuous phase. The general form of the population balance equation is given as ( 2 - 1 ) where the terms in the LHS provide the material derivative of number density. The terms in the RHS of Eq. ( 2 - 1 ) are birth due to coalescence, death due to coalescence, birth due to breaka ge, and death due to breakage of droplets, respectively. The birth due to coalescence is given by [ 84 ] 18 ( 2 - 2 ) where is the coalescence kernel. The death due to coalescence is given by ( 2 - 3 ) The birth due to breakage of droplets is given by ( 2 - 4 ) where is the breakage kernel, and is the daughter size distribution function. The death due to breakage is given by ( 2 - 5 ) in which is coalescence efficiency, is collision frequency, is droplet number densi ty, is number of daughter droplets , is daughter droplet size distribution function, and is breakage frequency. In Eq. ( 2 - 1 ) , the terms in RHS of the equation are droplets birth due to coalescence, death due to coalescence, birth due to breakage and death due to breakage of droplets , respectively. The population balance equation can be solved using different methods. There are two main category of methods: 1 - discrete (sectional) methods [ 85 - 87 ] , 2 - method of moments [ 84 , 88 - 90 ] . In discrete metho ds, the continuous particle size distribution (PSD) is treated as discretized bins or classes of particles. Discrete methods provide a PSD at every timestep in 19 each computational cell, however, to achieve an accurate PSD representation, a high number of bi ns is required. On the contrary, in a Method of Moments (MOM) approach, variables of the problem are reduced by introducing the moment transfer method , however, the number of variables may not be sufficient to provide an exact PSD in each computational cel l . Another classification of solution methods of PBE divides the solutions to homogenous methods, and inhomogeneous methods. In homogenous methods, all bins (or moments) move with the same velocity field, while in inhomogeneous methods, bins (or moments) move with various velocities . It should be noted that multiple velocity field s for transporting particle bins (or moments) is only important when particle size has a wide range, and where segregation is significant in the flow. For problems involving parti cles in a small size - range, a homogeneous method is computationally more efficient , and provides a reasonable accuracy. In current study, due to the presence of gravity force as a major player in particle/droplet motions, an inhomogeneous method of moments approach is employed. 2.2.1. Standard moment method In the Standard Moment Method ( SMM ), th e k th moment of the number density is defined as ( 2 - 6 ) According to the definition above, Eqs. ( 2 - 2 ) to ( 2 - 5 ) can be modified to reduce the number of variables [ 84 ] . The birth due to coalescence is given by [ 84 ] 20 ( 2 - 7 ) where superscript refers to the term associated for k th moment. The death due to coalescence is given by ( 2 - 8 ) The birth due to breakage of droplets is given by ( 2 - 9 ) The death due to breakage is given by ( 2 - 10 ) Solving PBE using SMM requires additional information about the number density function and is applicable only with further simplifications [ 90 ] or using a closure proposed by Diemer and Olson [ 91 ] . 2.2.2. Quadrature method of moments In the Quadrature Method of Moment s (QMOM), a quadrature approximation proposed by McGraw [ 92 ] is used as [ 89 ] ( 2 - 11 ) 21 Therefore, Eq. ( 2 - 6 ) can be modified to ( 2 - 12 ) where and are abscissas and weights, respectively. In this method, these values can be obtained from the lower - order moments. Therefore, in order to build a quadrature approximation of order , the first moments are sufficient. Therefore, for instance, to build a QMOM approximation of order 3, the fi rst six moments are needed (i.e. ). The algorithm to find and is the Product - Difference (PD) algorithm introduced by Gordon [ 93 ] , and is based on canonical moments theory [ 94 ] . A matrix is constructed with components from the moments . T he components of the first column are given by ( 2 - 13 ) where is the Kronecker delta. The component of the second column are given by ( 2 - 14 ) The remaining components (i.e. except columns 1 and 2) are ( 2 - 15 ) A Jacobi Matrix is created using and , such that 22 ( 2 - 16 ) ( 2 - 17 ) where is defined as ( 2 - 18 ) and and are the diagonal and the co - diagonal, respectively. The eigenvalues of this matrix are the abscissas. The weights are found as ( 2 - 19 ) where is the first component of the eigenvector, . 2.2.3. Direct quadrature method of moments In Direct Quadrature Method of Moments (DQMOM), in the contrary to QMOM, each node is treated as a distinct phase. The transport equation for the weights and weighted abscissas can be given as [ 88 ] ( 2 - 20 ) ( 2 - 21 ) where is the velocity of the dispersed phase , and . and can be calculated through a linear system of equations obtained from the first moment s as [ 88 ] 23 ( 2 - 22 ) where the matrix is a matrix as ( 2 - 23 ) and and are defined as ( 2 - 24 ) ( 2 - 25 ) and is defined as ( 2 - 26 ) and the RHS of Eq. ( 2 - 22 ) is calculated from the known source terms (due to coalescence and breakage events) as ( 2 - 27 ) where the source term of the moment in each computational cell is defined as 24 ( 2 - 28 ) The weights and abscissas are related to the volume fractions of phases, , and their effective length , according to [ 88 ] ( 2 - 29 ) ( 2 - 30 ) where is the volumetric shape factor. For spherical particles, . By replacing Eqs. ( 2 - 29 ) and ( 2 - 30 ) in Eqs. ( 2 - 20 ) and ( 2 - 21 ) , respectively, the transport equations can be written as ( 2 - 31 ) ( 2 - 32 ) 2.2.4. Breakage models The number of droplets in a c ontinuous phase may change due to droplet breakage and/or coalescence. Breakage of a droplet in a turbulent dispersion depends on the relative magnitude of two stresses: an external stress on the droplet due to turbulence in the continuous phase and a surf ace stress due to interfacial tension. The external stress causes deformation of the droplet, while the surface stress restores the spherical shape of the droplet and minimizes the interfacial area. If the external stress (due to turbulence) overcomes the restorative stress, the droplet breaks apart [ 95 ] . The balance between the 25 stress components leads to a prediction of the maximum stable droplet diameter. Breakage mechanisms can be categorized into four main categories [ 96 ] : (1) turbulent fluctuation and collision, (2) viscous shear stress, (3) shearing - off process, (4) interfacial stability. In the first case, the breakage is mainly due to turbulent pressure fluctuations along the surface, or by particle - eddy collisions . The particle surface becomes unstable when the amplitude of oscillation is close to a critical value; so deformation of the particle is started and stretching in one direction leads to a breakage of particle into two or more daughter particles. In this c ase, the dominant external force initiating the process is dynamic pressure difference around the particle. This breakage mechanism is a balance of dynamic pressure and its surface stress. In second scenario, viscous shear stress causes a velocity gradient around the interface and deforms the fluid particle, which eventually results in particle breakage. In this case, the breakage is a balance of external viscous stress and surface tension. In third scenario, the breakage is due to velocity difference acros s the interface. This process is a balance of viscous shear force and the surface tension. Forth scenario is different than the other three cases. In all those cases, the breakage occurs due to fluid dynamics characteristics in continuous phase. However, e xperiments show that in absence of such conditions, breakage can be caused by interfacial instabilities. Different models have been presented for breakage frequency, , due to turbulent fluctuation and collision [ 97 - 103 ] , viscous shear stress [ 104 , 105 ] , shearing - off [ 106 ] , and surface instability [ 87 , 107 ] . A mathematical model for the breakage frequency of droplets is given by Coulaloglou and Tavlarides [ 97 ] as 26 ( 2 - 33 ) where and are constants, is the dissipation rate, is the parent diameter, is the surface tension, and is the volume fraction of the dispersed phase, and is the density of the dispersed phase. and are empirically chosen as 0.00481 and 0.08, respectively. For solid flocs a breakage rate can be defined as a function of shear rate and particle diameter [ 108 , 109 ] given by ( 2 - 34 ) where is shear rate and is the volume of particle. and are empirical parameters which can be obtained by fitting the experimental data to the simulations. However, according to [ 109 - 111 ] and are chosen as 0.0047 and 1.6, respectively. There are different approaches to obtain a daughter droplet size distribution function, , including empirical models [ 112 ] , phenomenological models [ 113 , 114 ] , and statistical models [ 97 , 115 - 117 ] . A statistical model assumes that the size of daughter droplet is a random variable and its probability distribution satisfies either n ormal, beta , delta and uniform distribution functions [ 96 ] . In present study, for solid flocs a binary fragment distribution is used, in which a floc breakage result s in two flocs of identical volume. 27 2.2.5. Coalescence models Coalescence of droplets is the other phenomenon contributing in evolution of different droplet sizes in fluid flow. The coalescence of droplets is divided into thre e subprocesses [ 118 ] : (1) two bubbles collide and trap a small amount of liquid (continuous phase) between them; (2) bubbles keep in contact till the liquid film drains out to a critical thickness; (3) the rapture of the film results in coalescence. Several mechanisms promote collisions of droplets [ 118 ] : (1) Turbulent fluctuations of the continuous phase; (2) mean velocity gradients in the flow; (3) different droplet rise velocities, due to buoyancy or body forces; (4) droplet capture in an eddy; (5 ) wake interactions. In turbulence flows, one some or all the mechanisms above could be involved. But in laminar flows, shear - induced coalescence is the dominant mechanism. To implement coalescence in population balance equation, a coalescence kernel is re quired . The coalescence kernel , can be obtained as ( 2 - 35 ) in which two functions of collision frequency, , and collision efficiency, , and should be det ermined. 2.2.5.1. Droplet coalescence in laminar flow For shear - induced droplet collisions in a laminar flow, Friedlander [ 119 ] introduced a collision frequency given by ( 2 - 36 ) 28 where and are radii of colliding droplets, and is shear rate. The collision efficiency is obtained from a function originally introduced by Coulaloglou [ 97 ] as ( 2 - 37 ) where is the drainage time of droplets, and is the contact time. The drainage time is the time required for the intervening film between the droplets to thin to a critical thickness. The contact time is the interaction time of two colliding droplets. A model presented by Chesters [ 120 ] is used to obtain the drainage time for partially mobile interfaces as ( 2 - 38 ) where is interaction force, is equivalent radius obtained as , and and are initial and critical film thicknesses, respectively. The interaction force is given by [ 120 ] ( 2 - 39 ) The initial film thickness is considered as [ 114 ] . The critical (final) film thickness is adopted from Chesters [ 120 ] and presented by Venneker et al. [ 121 ] as ( 2 - 40 ) 29 where is Hamaker constant , which ranges between 10 - 20 J and 10 - 19 J. Hamaker constant is 10 - 20 J in this study. A model presented by Luo [ 122 ] is employed to calculate the contact time as follows: ( 2 - 41 ) where . The parameter is the added mass coefficient and is set to 0.65 in this study. is normally taken to be a constant between 0.5 and 0.8 [ 123 ] . Using Eqs. ( 2 - 36 ) to ( 2 - 41 ) , the coalescence kernel is calculated. 2.2.5.2. Suspended solid s coalescence in turbulence flow The overall collision frequency of two colliding solid particles is given by [ 79 , 124 ] ( 2 - 42 ) is collision frequency due to Brownian motion and is important for particles smaller than , and therefore for suspended solids is neglected in present study. is collision frequency due to shear for particles with diameters of and and is given by [ 125 ] ( 2 - 43 ) where is turbulent energy dissipation rate, and is kinematic viscosity. 30 in Eq. ( 2 - 42 ) is collision frequency due to gravity - induced sedimentation and becomes important for large particles when the density difference is significant. Since in this study the effect of flocculatio n on small particles is studies, this term is neglected. The collision efficiency, is the ratio of the successful encounters (leading to coalescence) to all the encounters of particles of sizes and . Only some of the interactions leads to coales cence, and therefore, maximum collision efficiency is 1.0 and . In solid - in - liquid suspensions, a polymer flocculant can enhance the collision efficiency by attaching to a floc. However, this process is a complex phenomenon and several factors can in fluence on it, such as charge density distribution on flocs, flocculant concentration in water, flow temperature, etc. [ 126 ] . Therefore, following assumptions are made to simplify the problem [ 127 ] : (1) The attachment of flocculants to the flocs are fast and adsorbed polymer is in equilibrium with the amount in water, (2) the adsorbed polymer is small compared to the total available polymer in the domain and it does not change the amount of f locculant in water. Assuming an exponential relation between the collision efficiency and the concentration of polymer flocculant, the collision efficiency , used in Eq. ( 2 - 35 ) , is given by ( 2 - 44 ) where is the concentration of polymer flocculant in water, and is a critical polymer concentration. A transport equation should be solved in order to obtain the concentration of flocculant 31 ( 2 - 45 ) where i s the concentration of flocculant , is the density of continuous phase , is velocity of continuous phase , is the molecular diffusivity of the polymer. refers to source or sink of the concentration. Film drainage Thin liquid films (TLF s ) are unstable structures formed in emulsions, foams, suspensions, and colloidal dispersed systems. Interaction between the phases in such systems is significant, since it affects the aggregation and breakage of the dispersed phase particles/droplets [ 128 ] . The kinetic of drainage of thin liquid film is widely studied in years [ 129 - 134 ] , in order to investigate the mechanical properties of the geometry of the film surface on the drainage [ 128 ] . As mentioned earlier, the drainage time, depends on rigidity of the surface of droplet, and mobility of the surface. The droplet can be either deformable or non - deformable, and its surfa ce can be immobile, partially mobile, or fully mobile (See Figure 2 - 2 and Figure 2 - 3 for details). 32 In most of the film drainage analysis in the literature, the lubrication theory is used to derive a thinning equation. In the next section, some of the existing models for film drainage of droplets are presented, according to rigidity and mobility of the droplet. 33 2.3.1. Film drainage models 2.3.1.1. Non - deformable rigid sphere In a case that the droplets viscosity is much higher than viscosity continuous phase, or the droplet diameter is very small, their interface at large distances behave such as rigid spherical droplets. Chesters [ 120 ] and Chester and Hofman [ 136 ] used Poiseuille relation to obtain the following drainage time: ( 2 - 46 ) where , is the continuous phase viscosity, and and are initial and critical film thicknesses. However, when the droplets are larger, deformation of droplet surface during the collision is not negligible. One of the ways to consider the deformation of form to parallel disks during the coalescence (see Figure 2 - 2 (b)). 2.3.1.2. Deformable droplets with immobile interfaces For immobile interfaces, the film drainage is contro lled by viscous thinning [ 118 ] . The velocity in the film is a parabolic profile with no slip boundary conditions on the interfaces (see Figure 2 - 3 (a)). For this case, the drainage time for different droplet sizes can be obtained as [ 120 , 137 ] ( 2 - 47 ) 34 Assuming the initial and critical thicknesses ( and ) as constants, one can get the simplified and commonly used model of Coulaloglou and Tavlarides [ 97 ] . 2.3.1.3. Deformable droplets with partially mobile interfaces In many problems, since the contribution of additional flow within the film is comparab ly smaller, the drainage is controlled by the motion of the film surface [ 118 ] . The drainage time by assuming a quasi - steady creeping flow for par tially mobile interfaces was calculated by Chesters [ 120 ] as ( 2 - 48 ) One of the commonly used models, introduced by Davis et al. [ 138 ] . They calculated the relation between the force during film drainage and drainage velocity ( 2 - 49 ) where is drainage velocity, and , and . 2.3.1.4. Deformable droplets with fully mobile interfaces This is most complicated drainage regime where the drainage is controlled by both inertia and viscous terms. Chesters [ 139 ] proposed the following model for this case 35 ( 2 - 50 ) in which, , and is the radius of the deformed droplet surface ( see Figure 2 - 2 ). Since there is no analytical solution for Eq. ( 2 - 50 ) , two limiting cases of high viscous liquids, and inertia - controlled drainage are often considered. As seen above, in drainage time relations, one need to embed the interaction force of collis ion, , which is usually not a constant. The force is assumed to be proportional to the mean square velocity difference at ends of the eddy size of the equivalent diameter [ 97 , 114 , 118 , 140 ] ( 2 - 51 ) For locally isotropic turbulence, the turbulent force is given by Narsimhan [ 141 ] as ( 2 - 52 ) for inertial subrange, and ( 2 - 53 ) for viscous subrange. In the equations above, is energy dissipation defined as , is dynamic viscosity of the mixture defined as , and is dynamic density of the mixture defined as . It should be noted 36 that is volume fraction of dispersed phase. The net force is sum of turbulent and colloidal forces. The colloidal force is given by Narsim han [ 141 ] as ( 2 - 54 ) where is initial film thickness, which is suggested [ 114 ] as . 2.3.2. Drainage of thin liquid films with mo bile interface subjected to sur factants In this section, a new model for drainag e of thin liquid films with mobile interface, subjected to surfactants are presented. Figure 2 - 4 shows a schematic of the liquid film thinning. The following assumptions are being made for the calculations: (i) Quasi - steady state process 37 (ii) Small film thickness (iii) Small Reynolds number (iv) Surface pro perties of the droplets are the same (v) Small Peclet number (vi) Diffusion - controlled adsorption (vii) Small deviations of concentration, and adsorption (viii) Droplet surfaces in spherical parts are approximated as paraboloids In the lubrication approximation [ 142 ] , the pressure in the continuous phase depends only on the radial coordinates, and the radial component of t he velocity is calculated as follow: ( 2 - 55 ) where is the dynamic viscosity, and and are the surface velocities on the lower and upper droplets, respectively. The integrated bulk continuity equation is given by [ 142 ] ( 2 - 56 ) where is a sum of Couette and Poiseuil le average velocities, given by 38 ( 2 - 57 ) According to lubrication theory, the momentum equation in a radial direction is given by By equating the flow rate of squeezed liquid film (due to relative velocity of droplets, i.e. rate of film thinning ) and the flow rate of the obtained velocity profile, the pressure gradient reads ( 2 - 58 ) where ( 2 - 59 ) in which . It should be noted that the spherical parts of the droplet surface is approximated as paraboloid [ 138 ] . In the lubrication approximation, the tangential stress boundary conditions at the film interfaces are simplified to [ 143 ] ( 2 - 60 ) ( 2 - 61 ) 39 By substituting Eq. ( 2 - 55 ) in Eqs. ( 2 - 60 ) and ( 2 - 61 ) , t he balance of bulk viscosity stress and surface tension gradient (Marangoni effect) on the inte rface of the droplet results in the following equation ( 2 - 62 ) where is the mean surface velocity, is surface tension, and is total surface viscosity defined as ( 2 - 63 ) where is interfacial shear viscosity , and is dilatational viscosity. For the thin film liquids the bulk surfactant concentration depends weekly on the vertical coordinate. Thus , where . Also, the Peclet number is small. Therefore, the diffusion equation in the film phase can be written as ( 2 - 64 ) where is the bulk diffusion coefficient. By integrating Eq. ( 2 - 64 ) over and from to : ( 2 - 65 ) By multiplying Eq . ( 2 - 56 ) by and adding it to Eq . ( 2 - 65 ) we get 40 ( 2 - 66 ) The equation of surfactant species balance at film interfaces in a diffusion - controlled adsorption is: ( 2 - 67 ) where is surfactant absorption on the film surface, is the tangential velocity of film interfaces during drainage, is the concentration of the surfactant on the film, and is surface diffusion coefficient. By adding the surfactant mass balance and surfactant species balance (i.e. Eq . ( 2 - 66 ) and ( 2 - 67 ) ), the total mass balance equation can be obtained as ( 2 - 68 ) We assume a quasi - steady state condition for Eq. ( 2 - 68 ) . Also small deviations from equilibrium conditions for both concentration and adsorption coefficient are assumed as follow [ 144 ] : ( 2 - 69 ) So the following equation can be obtained [ 142 , 144 ] 41 ( 2 - 70 ) To non - dimensionalize the equation above, the following variables are defined [ 143 ] ( 2 - 71 ) Therefore, the following non - dimensionalized equation is derived ( 2 - 72 ) where ( 2 - 73 ) and . By defining , two limiting cases can be considered. 2.3.2.1. Limiting case 1, In the case of large interface viscosity, large film thickness and low deformability ( ), the contribution of the interface viscosity on hydrodynamic pressure gradient is negligible. Therefore, the pres sure gradient is mainly developed by parabolic part of velocity profile. Thus, Eq. ( 2 - 72 ) can be simplified to 42 ( 2 - 74 ) By substituting Eq. ( 2 - 73 ) into Eq. ( 2 - 74 ) and applying the following boundary conditions ( 2 - 75 ) we get the following solution for (parallel disk region) as ( 2 - 76 ) And for (meniscus region), we get ( 2 - 77 ) In order to find , by equating the velocity gradients at ( 2 - 78 ) we get 43 ( 2 - 79 ) By substituting Eq. ( 2 - 79 ) into Eqs. ( 2 - 76 ) and ( 2 - 77 ) , the velocity profiles in both regions (parallel disk, and meniscus) are obtained. The hydrodynamic drag force, , is given as ( 2 - 80 ) By integrating by part, one can conclude that the following relation can also be used: ( 2 - 81 ) By substituting Eq. ( 2 - 58 ) into Eq. ( 2 - 81 ) and applying the variable changes of Eq. ( 2 - 71 ) , we get ( 2 - 82 ) It should be noted that in Eq. ( 2 - 82 ) , the first term is contribution of the interface mobility, and the second term accounts for immobile interfaces. By substituting Eqs. ( 2 - 76 ) and ( 2 - 77 ) into Eq. ( 2 - 82 ) the hydrodynamic force can be obtained for both parallel disk and meniscus region. Fo r parallel disk region ( ): 44 ( 2 - 83 ) For meniscus region ( : ( 2 - 84 ) 2.3.2.2. Limiting case 2, In the case of small interfacial viscosity and very thin film, we have . Thus, Eq. ( 2 - 72 ) can be simplified to ( 2 - 85 ) we get the following solution for (parallel disk region) as ( 2 - 86 ) And for (meniscus region), we get 45 ( 2 - 87 ) Similar to the pervious section, the hydrodynamic force can be obtained for both parallel disk and meniscus region. For parallel disk region ( ): ( 2 - 88 ) For meniscu s region, after some algebra, the hydrodynamic force is calculated as ( 2 - 89 ) The integral of Eq. ( 2 - 89 ) cannot be solved analytically. By assuming , we get ( 2 - 90 ) By applying Taylor series expansion to , we have ( 2 - 91 ) Therefore, the solution for the meniscus region ( ) is given a s 46 ( 2 - 92 ) 47 Numerical investigation on Treatment Shaft system: separation efficiency and impact of flocculant polymer agents Introduction Inadequate treated combined sewer overflows (CSOs) and discharges of urban stormwater are among major causes of long - term persistence of poor water quality in receiving waters [ 12 ] . The pollution is caused by different pollutants (such as heavy metals and solid particles) present in urban wastewater and washoff of urban surfaces by rain [ 145 - 147 ] . A combined sewer is a sewage collection facility that collects sewage and surface runoff in the same system. These facilities may cause serious water pollutions when a Combined Sewer Overflow (CSO) oc curs, i.e. the flows exceed the capacity of the treatment plant. An approach to reduce the number of combined sewer overflows is to store rainfall runoff temporarily [ 2 ] . After the e vent, the storage facility can be emptied by pumping the contents into the sanitary sewer system. The excess flow can be discharged to the receiving water body, but treatment (such as skimming, settling, screening and disinfection) is needed before release . Through sedimentation, solids are separated in these facilities. The detention facilities can also contribute to removing solid particles from the overflow through sedimentation. Solid particles are described as total suspended solid (TSS), which is a pa rameter to assess the quality of wastewater after treatment. The concentration of TSS is evaluated by standard laboratory tests [ 3 ] . The removal of these solid particles is important in order to reduc e the pollution in receiving water. This can be achieved by settling, filtration, coagulation and flocculation [ 4 ] . Settling tanks are one of the most practical and simplest ways of pollution management in water 48 since they provide volumetric control of flows and treatment by settling [ 5 ] . The efficiency of a sedimentation tank to remove suspended solids depends on details of rainstorm (solid concentration, water flow rate), and characteristics of the particles (shape, size, drag and buoyancy forces) [ 6 ] . Detention tanks can be designed in various ways [ 7 ] . One approach is to use in - system storage with real - time monitoring structures [ 8 ] , but it is vulnerable to malfunction in harsh environment. Another app roach is based on offline detention basins which are typically constructed below - grade, and work using sedimentation basins or clarifiers [ 9 , 148 ] . These systems however have large environmental footprints, and are difficult to constru ct in populated urban environments. Another approach is based on using deep storage tunnels far below - grade where other utilities are not located. Such tunnels however become filled (e.g. rainstorm events), and are prone to surges and trapped air during th e filling process, which may lead to formation of geysers in shafts and results in return of low - quality water to grade [ 10 , 149 ] . When it is dewatered, removal of solids i s another issue in a long tunnel [ 11 ] . A novel alternative approach, called the Treatment Shaft, appropriately addresses several of the aforementioned concerns. A Tre atment Shaft is a detention basin with a cylindrical shape which has relatively small footprint that is able to store a large volume by extending deep underground. The Treatment Shaft is a patented technology (U.S. Patent No. 6,503,404 [ 13 , 14 ] and other patents) that includes the necessary CSO control and treatment, with reduced footprints . A main baffle wall is located at the center to ensure a very slow upward velocity flow during filling , which enhances the settling efficiency. In addition, 49 accumulated solids at the bottom of the shaft are discharged into a sanitary sewer system. A flushing system uses high - pressure nozzles to clear the bottom of Treatment Shaft during in a final rinse cycle [ 7 ] . In this study, a Treatment Shaft is thoroughly evaluated as an innovative and a relatively new techn ology to be used as a CSO detention system with a flocculating agent. Using Computational Fluid Dynamics (CFD), the filling process of the Treatment Shaft system during a 10 - year, 1 - hour rainstorm event is modeled and the separation efficiency of the CSO s ystem is evaluated. The residence time distribution is obtained to ensure the contact time is sufficient to provide disinfection. A Lagrangian - Eulerian approach is used to estimate the natural sedimentation efficiency of different particle sizes. Then, a combination of an Eulerian - Eulerian approach and Population Balance Model (PBM) is presented to study the impact of a flocculant, and coalescence and emerging larger flocs on the separation efficiency of smaller particles. Lastly, different design alternat ives such as various shaft sizes, and different baffle locations is evaluated. System features and geometry A Treatment Shaft consists of a diverging inlet and a converging outlet (with flared channel angle of and length of ). The inlet channel is a r ectangular channel diverging into a cylinder with diameter of 28.96m. The cylinder is attached to a cone - shaped part at the bottom. A baffle of length 39.07m and thickness of 0.91m is located in the center of the shaft. The gap at the bottom of the shaft w ith the length of is the turning point in the shaft, 50 providing an upward flow to the outlet . In original design, dimensions are: , and . However, to reduce the environmental footprint, the following dimensions is considered as the main design in this study: , , and . It is shown later in this study that the efficiency difference is insignificant between these two different designs. The schematic of Treatment Shaft is illustrated in Figure 3 - 1 . Treatment Shaft should be able to handle a 10 - year, 1 - hour rainstorm event. This rainstorm event leads to a surface flow and eventually a CSO of 5 hours long, including an increasing flow (w hile storm is happening) followed by a decreasing flow (after storm). The CSO has a peak flow rate of 41.6 m 3 /s. 51 (a) (b) (c) (d) 52 Figure 3 - 2 presents the variation of flow rate into the shaft during a 5 - hour timespan. The data is obtained using Storm Water Management Model (SWMM) [ 150 ] for CSO 017 (located in Dearborn, Michigan). Particles are assumed as solid particles with density of 2650 kg/m 3 [ 151 ] suspended in water with density of 998.2 kg/m 3 and viscosity of 0.001 Pa.s. The size of s uspended particles ranges from 1 µ m to 10.32 mm. Modeling approach Depending on the particle size and flow pattern in a CSO detention system, larger particles may naturally sediment, while smaller particles tend to follow the fluid flow and unfavorably rea ch to the outlet of the detention system. By providing a slow upward flow 53 section in the Treatment Shaft, particles of different sizes sediment without any external force or moving parts. However, adding flocculants promotes coalescence of particles and th erefore, improves the separation efficiency of the system. In this section, first the Residence Time Distribution (RTD) is obtained for the treatment shaft. RTD is an important factor in detention systems which evaluates if the system provides sufficient time for disinfection. Then , natural sedimentation of particles (settling without any mechanical means, mixing, etc.) is studied using a Lagrangian - Eulerian approach. Thus, diameter in which the natural sedimentation occurs for at least 50% of particles ( ) is computationally obtained. A n Eulerian - Eulerian approach in combination wit h Population Balance Model (PBM) is then utilized to evaluate the impact of flocculant and coalescence and possible breakage of flocs, for particles of size and smaller. Lastly, the impact of various geometries for a sample particle size distribution is investigated. 3.3.1. Residence Time Distribution study Residence Time Distribution is a probability function of the amount of time a species remains inside a fluid domain. In present study, a tracer species analysis is used to obtain RTD. The transport of a ch emical species in a fluid flow is given by the following convection/diffusion equation [ 152 ] : ( 3 - 1 ) 54 where is any scalar variable, is species density, is turbulence viscosity, is velocity field, is the molecular diffusivity of the scalar , and is turbulent Schmidt number of the scalar . The term refers to source or sink of the scalar variable, which is equal to 0 for only - transport problems. There are two approaches to determine the residence time distribution: pulse approach and step approach. In pulse approach, the tracer is injected at for only one computational time - step, and then the species concentration at the outlet, , is me asured over time. In step approach, the injection will be continued over all time - steps, and the species concentration at the outlet is monitored. The dimensionless resident time, , for pulse injections is given by ( 3 - 2 ) The dimensionless normalized species concentration, at the outlet for a step injection of at the inlet is given by ( 3 - 3 ) It should be noted that step and pulse responses are related by 55 ( 3 - 4 ) In this study, the step approach is used. The mean residence time of the domain is calculated by ( 3 - 5 ) where is the concentration of the tracer at the outlet, and is time passed since injection of the tracer. The key terms to characterize the RTD curve are presented in Table 3 - 1 [ 152 , 153 ] . Table 3 - 1 : Residence time distribution key terms for characterization, from [ 152 , 153 ] Term Definition Theoretical hydraulic retention time (volume of the domain divided by flow rate) Time at which tracer first appear Mean residence time , , Time at which 10%, 50%, and 90% of tracer exits the domains at the outlet Index of short - circuiting (equals to 1.0 for plug flow) Index of mean retention time 3.3.2. Natural sedimentation In order to assess the natural settling of non - interacting solid particles, a Lagrangian - Eulerian approach is employed. In this approach the motion of particles is modeled using a Lagrangian formulation, in which provides a one - way coupling, and particles do not have any impact on fluid flow. The v olume fraction of particles in a CSO is less than 10%, and therefore, 56 this a pproach is an appropriate choice for studying the natural settling of particles [ 65 , 154 , 155 ] . To model the continuous phase (i.e. water), Navier - Stokes equation of isothermal, incompressible, single - phase turbulence flow is solved using an Eulerian approach. The Reynolds - Averaged Navier - Stokes (RANS) equation is given as [ 156 ] : ( 3 - 6 ) where is the Reynolds stress tensor, and and are density of viscosity of continuous phase, respectively. is the spatially - averaged velocity. The Reynolds stress te nsor is calculated using a Realizable model for turbulent viscosity. To model the dispersed phase (i.e. solid particles), the trajectory of solid particles, in a Lagrangian frame of reference can be given as ( 3 - 7 ) The velocity of a solid particle of RHS of Eq. ( 3 - 7 ) is obtained by integrating the following ( 3 - 8 ) 57 where and are instantaneous velocities of continuous phase and the solid particle, respectively. The first term in RHS of Eq. ( 3 - 8 ) is due to drag force. is the drag force per unit mass of particles given by ( 3 - 9 ) where is the drag coefficient, is the diameter of a solid particle, and is density of solid particles. The second and the third terms of RHS of Eq. ( 3 - 8 ) are acc eleration due to pressure gradient force, and acceleration due to virtual mass force , respectively. 3.3.3. Impact of a flocculant on separation efficiency To investigate the impact of flocculants on particle settling in Treatment Shaft, as mentioned earlier, a different approach is employed. An Eulerian - Eulerian approach is used, where continuous and dispersed phase are modeled using an Eulerian approach. Additionally, a Population Balance Model is utilized to account for particle coalescence and breakage. 3.3.3.1. M ulti phase flow modeling The Mixture model is used to simulate the flow of CSO in Treatment Shaft. W ater is the continuous (or primary) phase, and the dispersed (or secondary) phase consists of solid particles . The mixture flow variables are defined as 58 ( 3 - 10 ) ( 3 - 11 ) ( 3 - 12 ) where subsc ripts , , and refer to mixture, continuous, and dispersed phases, respectively. The mass conservation for the mixture is given by ( 3 - 13 ) The momentum equation for a mixture is given by ( 3 - 14 ) where is a body force, is turbulence shear stress and is the visco sity of the mixture. The terms and represent the drift velocity of continuous and dispersed phases, respectively. The drift velocity of any phase is defined as ( 3 - 15 ) 59 Relative (slip) velocity is modeled using a model proposed by Manninen et al. [ 157 ] , and drag force on particles is modeled using a model presented by and Schiller and Naumann [ 158 ] . 3.3.3.2. Flocculation model In order to study the particle interactions, a population balance equation is used as [ 90 ] ( 3 - 16 ) where is number density function, and is particle characteristic length. is birth of particles due to coalescence, is death due to coalescence, is birth of particles due to breakage, and is death of particles due to breakage. Birth of particles due to coalescence is given by [ 90 ] ( 3 - 17 ) where is coalescence kernel. The coalescence kernel is often defined as a product of two functions [ 118 ] ( 3 - 18 ) 60 in which is coalescence frequency, and is collision efficiency. The overall collision frequency of two particles is given by [ 79 , 124 ] ( 3 - 19 ) is collision frequency due to Brownian motion and is important for particles smaller than , and therefore is neglected in present study. is collision frequency due to shear for particles with diameters of and is given by [ 125 ] ( 3 - 20 ) where is turbulent energy dissipation rate, and is kinematic viscosity. in Eq. ( 3 - 19 ) is collision frequency due to gravity - induced sedimentation and becomes important for large particles when the density difference is significant. Since in this study the effect of flocculation on small particles is studies, this term is neglected. The collision effici ency, is the ratio of the successful encounters (leading to coalescence) to all the encounters of particles of sizes and . Only some of the interactions leads to coalescence, and therefore, maximum collision efficiency is 1.0 and . A polyme r flocculant enhances the collision efficiency by attaching to a floc. However, this process is a complex phenomenon and several factors can influence on it, such as charge density distribution on flocs, flocculant concentration in water, temperature, etc. [ 126 ] . Therefore, following assumptions are made to simplify the problem [ 127 ] : (1) The 61 attachment of flocculants to the flocs are fast and adsorbed polymer is in equilibrium with the amoun t in water, (2) the adsorbed polymer is small compared to the total available polymer in the domain and it does not change the amount in water. Assuming an exponential relation between the collision efficiency and the concentration of polymer flocculant, t he collision efficiency used in Eq. ( 3 - 18 ) is given by ( 3 - 21 ) where is the concentration of polymer flocculant in water, and is a critical polymer concentrati on. A species transport equation, Eq. ( 3 - 1 ) , should be solved in order to obtain the concentration of flocculant in Treatment Shaft. The death of particles due to co alescence is obtained by [ 90 ] ( 3 - 22 ) The birth of particles due to breakage is given by ( 3 - 23 ) where is breakage frequency , and is probability density function that contains information on the fragments produced by a breakage event. The death due to breakage is obtained by 62 ( 3 - 24 ) The breakage rate of flocs is a function of shear rate and particle diameter [ 108 , 109 ] and is give n as ( 3 - 25 ) where is shear rate and is the volume of particle. and are empirical parameters which can be obtained by fitting the experimental data to the simulations. Ho wever, according to [ 109 - 111 ] and are 0.0047 and 1.6, respectively. In present work, a binary fragment distribution is used, in which a floc breakage results in two flocs of identical volumes. 3.3.4. The effect of various treatment shaft geometries Various Treatment Shaft configurations are eval uated and compared. Specifically, the size of baffle, its location, and the channels geometry are discussed below. For this purpose, the particle size distribution (PSD) is a sample based on a report of New Jersey Department of Environmental Protection (NJ DEP), representing the various particles that would be associated with typical stormwater runoff [ 159 , 160 ] . Figure 3 - 3 shows the particle size distribution (PSD) used in the simulations . 63 In this work, following cases are investigated: The effect of the length of vertical baffle (or size of the gap, , at the bottom of the baffle size) on the separation efficiency (B=5.42m, B=9.32m, B=13.22m) . The effect of the location of the vertical baffle on the separation efficiency (2.9m toward inlet, center, 2.9m toward outlet) . The shape of the inlet/outlet channels (Flared: , , and Straight: , ). Numerical simulations detail s In present study, all the numerical simulations are carried out by ANSYS ® Fluent (Academic Research, Release 18.1). To obtain Residence Time Distribution (see Section 3.3.1 ), a single - phase turbulence flow solver in combination with species transport model is used. For studying the natural sedimentation (see Section 3.3.2 ) , and also to study the effect of various geometries (see Section 3.3.4 ) , a single - phase turbulence flow solver is employed 64 in combination with Discrete Phase Mo del (DPM) as a Lagrangian particle tracking method. To calculate the trajectory of particles, Equation ( 3 - 8 ) is integrated using an automated scheme selection by employing a Trapezoidal scheme as high order scheme, and implicit method as low order scheme [ 161 ] . The maximum number of time steps per iteration is 500, and to control the accuracy the tolerance of the solution is set to 10 - 5 . In order to track solid particle settling in the computational domain, active particle sampling is employed, which continuously samples the trapped particles in shaft surfaces. This active sampling method provides large data - files (including particle size, velocity , residence time, etc.) that are analyzed using in - house Python codes. To study the impact of flocculation (see Section 3.3.3 ), a multiphase mixture model is use d in combination with species transport model to obtain the concentration of the flocculant in the domain. Additionally, Direct Quadrature Method of Moments (DQMOM) is used to solve the population balance model which accounts for coalescence and breakage o f flocs. Flow variables of the continuous phase (such as pressure and velocity) are obtained by solving the Navier - Stokes equation. PISO pressure - velocity coupling is used for transient simulations using a bounded 2 nd order implicit scheme for time - marchin g. The Realizable k - turbulence model with scalable wall function is employed. Pressure interpolation is performed using Body Force Weighted scheme. Spatial discretization for pressure, momentum, and turbulence variables is given by 2 nd order upwind schem es. 3.4.1. Boundary conditions The boundary condition at the inlet of the Treatment Shaft is a varying uniform velocity - inlet (see Figure 3 - 2 ), which is introduced to the so ftware using a User - Defined Function 65 (UDF). The boundary condition at the outlet of the Treatment Shaft is a zero pressure at the outlet. The computational domain is extended at the outlet to ensure there is no reverse flow on the physical outlet. Top wall s of Treatment Shaft are slip - walls (zero - shear) to resemble an open - channel flow simulation without dealing with expensive computational costs of interface tracking. The walls of the cylinder, cone - shaped part, baffle wall, and inlet and outlet channel wa lls are no - slip boundary condition. In flocculation process, a polymer flocculant agent is injected at a rate of of water at the inlet. Due to the small concentration of polymer in flocculant solution, its density and viscosity is assumed same as water. The critical polymer concentration is 0.1%. Suspended solids are uniformly injected at the inlet surface. Figure 3 - 4 shows the Particle Size Distribution (PSD) of suspended solids, which is the mean PSD of CSO samples adapted from [ 162 ] . This PSD is used to assess the natural sedimentation, and impact of flocculation. 66 Results and Discussions In both natural sedimentation analysis and floccu lation process analysis, the separation efficiency is defined as the ratio of mass of the suspended solid s in the Treatment Shaft ( ) at the end of a rainstorm to mass of total solid particles entered the shaft ( ) during a rainstorm event (5 - hour timespan): ( 3 - 26 ) In this section, first the residence time distribution is obtained and compared with available experimental data. Then, natural sedimentation without the effect of flocculant agents is investigated. Later , the impa ct of flocculation is discussed. Lastly, t he influence of various geometries is assessed. 3.5.1. Residence time distribution In present study, the step approach is used to obtain the residence time distribution. A combined steady - state/transient algorithm is employed. For RTD evaluation, the inlet flow r ate is assumed the peak flow of 41.6 m 3 /s, which is corresponding to inlet velocity of 1.2 m/s. A steady - state solution is obtained for the fluid flow without any tracer at the domain boundaries or interior. The solution scheme is changed to the transient and a tracer (with the same material of the main fluid, i.e. water) with mass fraction of 1.0 is injected to the domain at inlet. Species transport equation is only being solved in this step. Therefore, the 67 tracer uses the existing velocity field of the fl uid in the domain. The tracer species concentration is monitored at the outlet boundary. The results are illustrated in Figure 3 - 5 . According to the experimental results from [ 7 ] , the mean breakthrough time, i.e. was measured as t=660s, which is in good agreement with the numerical results of present study, at t=728s. Table 3 - 2 tabulates the key values for characterization of RTD of Treatment Shaft. Table 3 - 2 : Residence time distribution key terms for characterization of present study Term value 745.5 s 335 s 890.7 s 424 s 728 s 1375 s 4299 s 0.449 0.976 68 3.5.2. Natural sedimentation In order to investigate the natural sedimentation of suspended solids, 41 groups of particles of various sizes are injected at the inlet of the Treatment Shaft system, such that , where is the diameter of a particle in i th group, and . Therefore, the largest solid particle is . Figure 3 - 6 shows the separation efficiency of the Treatment Shaft system during a rainstorm CSO for various particles sizes. It also illustrates the spots where the particles settle. It can be seen that most of solid particles settle on the large cone - shaped part of the shaft. Additionally, the particle diameter that has at least 50% of settling rate ( ) is shown in the graph. is indicated as 161.3 µ m. 69 Figure 3 - 6 also show s that separation efficiency is about 10% for particles smaller than 50 µ m. These particles are mainly remained in the shaft at the end of the CSO event and are easily separated during drainage of the shaft. On the other hand, separation efficiency is higher than 90% for particles larger than 500 µ m. Figure 3 - 7 provides the area - normalized separation efficiency values, which normalize the efficiency of each part (e .g. Small cone - shaped part) by area of that part. It is shown that flat bottom of the shaft has the highest accumulation of mass per area during the settling process. Figure 3 - 8 shows the Average S ettling Time (AST) and Average Settling Velocity (ASV) for various particle sizes during a 5 - hour CSO of a rainstorm. The average settling time is defined as the mean value of settling time for a specific particle size. The average settling velocity, on th e other hand, is defined as the mean value of velocity of particles settling in the Treatment Shaft. In an ideal condition, if particles reach the upward section of the Treatment 70 Shaft, they will have a higher AST and lower ASV compared to settling in the downward section. It is observed that in particles larger than 64 µ m, freefall due to gravity force overcomes the drag force and consequently, larger particles settle with higher velocity, and in a shorter time. On the contrary, particles smaller than 64 µ m tend to get trapped in eddies and follow fluid flow, due to less gravity force and a dominant drag force applied on them. Figure 3 - 9 shows the Average Particle Time (APT) and Average Particle Velocity (APV) for various particle sizes at the end of a 5 - hour CSO of a rainstorm. APT is defined as the average time a particle has been in the shaft, and APV is the average velocity of tha t particle at the end of the 5 - hour CSO rainstorm. It is seen that particles larger than 161.3 µ m have a small APT, which demonstrates that these particles were recently entered the shaft. A higher APV for these particles which is higher than fluid veloci ty at that timeframe shows that these particles are moving with a significantly higher velocity compared to smaller particles. 71 On the other hand, smaller particles remain in the shaft for a long time, while their APV is close to fluids velocity at that t imeframe. In order to investigate the efficiency of natural sedimentation in the shaft for PSD introduced in Figure 3 - 4 , the results provided in Figure 3 - 6 is employed. For this purpose, since particle - particle interactions and impact of particles on fluid flow are both neglected, a volume - based averaging approach is used to calculate the overall efficiency of the Treatment Shaft during a rainstorm event without flocculation. The overall efficiency of the Treatment Shaft for the given PSD ( Figure 3 - 4 ) is obtained as 33.9%. For the specified PSD, only particles smaller than 161.3 µ m only 20.8% of particles settle, while for particles larger than 161.3 µ m, the separation efficiency is 81.7%. This also highlights the importance of a remedy for smaller partic les. 72 In addition to the investigation of the Treatment Shaft system for a rainstorm event, it is performance under the steady - state condition is also evaluated. For this purpose, the aforementioned particle groups are assessed under steady - state operation of the shaft with various inlet velocities, ranging from 0.1 m/s to 2.0 m/s (i.e. flow rates of 3.42 m 3 /s to 68.38 m 3 /s). Figure 3 - 10 provides a contour map of the s eparation efficiency of various particle sizes, under various stead - state velocity fields. It is shown that for inlet velocity of 0.1 m/s, all the particles of size 80 µ m and larger are settled (i.e. ). By increasing the inlet velocity, the efficiency o f separation rapidly drops, since the force applied on particles due to the flow field is bigger compared to the gravity force. Consequently, for inlet velocity of 2.0 m/s, only particles of size 406 µ m or larger have a separation efficiency of 100% (i.e. ). It should be noted that for particles larger than 512 µ m, for all evaluated velocities, and therefore, omitted from the contour plot (see Figure 3 - 10 ). 73 3.5.3. Flocculation process Results from previous section indicates that . As mentioned earlier in this section, the flocculation process is only evaluated for particles smaller than , where improving the settling efficiency is required. Therefore, the PSD at the inlet of the Treatment Shaft is different that Figure 3 - 4 . Highlighted section in Figure 3 - 11 shows the PSD used to investigate the impact of flocculation. The flow pattern in the shaft is continuously changing because of turbulence and varying flow rate at inlet of the shaft. Figure 3 - 12 shows velocity contours and streamlines at a cross - section passing through the symmetry plane of the Treatment Shaft in different timeframes during a rainstorm event. It is seen that location of vortexes and velocity profile is continuously changing, whic h can affect the instantaneous settling efficiency. It is seen that the velocity of the downward flow is higher than the velocity of the upward flow in all the 74 timeframes, thus providing sufficient time for some solid particles to settle at the bottom of t he shaft. Additionally, surface streamlines at different timeframes are illustrated in Figure 3 - 12 . The baffle wall induces vertical vortexes near the baffle s surfac e, which can trap and significantly extend the time a solid particle would remain in the shaft. (a) (b) 75 (c) (d) (e) (f) Figure 3 - 13 illustrates Q - criterion of the flow field in various timeframes, colored with velocity magnitude. function is the second invariant of the velocity gradient tensor, defined as [ 163 ] 76 ( 3 - 27 ) where is the vorticity tensor, and is the strain rate tensor. represents the local balance between the strain rate, and vorticity magnitude, and is an indicator of turbulence structure. (a) (b) 77 (c) (d) It is seen that in addition to two main vortexes near the inlet channel, many smaller vortexes are formed in various spots in the channel, increasing the residence time of particles in the Treatment Shaft. Figure 3 - 14 illustrates the variations of the concentration of polymer flocculant, and Figure 3 - 15 shows the changes of Sauter mean diameter at a cross - section of the Treatment Shaft. Figure 3 - 14 shows the concentration of flocculant in various timeframes. It should be noted that after about 34 minutes, the flocculant is saturated in entire domain. 78 (a) (b) Figure 3 - 15 shows a higher Sauter mean diameter in regions with higher concentration of flocculant, demonstrating the i mpact of a flocculating agent on promoting particle aggregation. 79 (a) (b) (c) (d) 80 Figure 3 - 16 illustrates the volume fraction of suspended solids in a cross - section of the Treatment shaft during a rainstorm event in various timeframes. It is seen that after the peak flow (~3000s), the volume fraction on the upward section of the shaft is decre asing and it is increasing at the bottom of the shaft. It highlights the importance of slow upward flow in settling the suspended solids during after the peak flow. Additionally, a comparison of the results with Figure 3 - 15 shows that as expected, the region with accumulated particles at the bottom of the shaft is mainly occupied with the largest particles. 81 (a) (b) (c) (d) 82 Figure 3 - 17 shows the distribution of Sauter mean diameter in the shaft, and at the outlet, for various timeframes. Figure 3 - 17 illustrates that not only is larger than its value at the inlet channel, the distribution of Sauter mean diameter is changing toward larger droplets over time. Since only particles smaller than 161.3 µ m is investigated for the impa ct of flocculation, Figure 3 - 17 also illustrates the impact of flocculant in forming particles that are several times larger than the maximum particle size at the inlet. Results at the end of the rainstorm event shows that a flocculating agent for particles smaller than 161.3 µ m can incr ease the separation efficiency from 20.8% to 49 .3 %. 83 3.5.4. The effect of various treatment shaft geometries Various geometry comparison is carried out in this section. Note that the original design refers to the case with flared input/output channels, and , and baffle wall located at the center. 3.5.4.1. The effect of the length of vertical baffle A larger baffle wall, providing a gap distance of , and a smaller baffle wall, providing a gap distance of are compared with the original design. The comparison is presented in Figure 3 - 18 . It is seen that decreasing the gap distance lowers the separation efficiency, since it increases the flow ve locity near the bottom of the shaft. On the other hand, a larger baffle gap does not significantly change the separation efficiency. The separation efficiency for the given PSD (see Figure 3 - 3 ) is 33.2%, 34.8%, and 35.4% for B=5.42m, B=9.32m, and B=13.22m, respectively. 84 3.5.4.2. The effect of the location of vertical baffle A baffle wall moved from center toward inlet ( ), and a baffle wall moved toward outlet ( ) are compared with the original design ( ). The comparison is presented in Figure 3 - 19 . It is shown that moving the baffle wall toward inlet slightly increases the separation efficiency for particles larger than . However, the difference is insignificant. The separation e fficiency for the given PSD (see Figure 3 - 3 ) is 34.8 %, 34.8%, and 35 .6 % for = - 2.9 m, , and , respectively. 3.5.4.3. The shape of the inlet and outlet channels Flared inlet channel enables a lower velocity of fluid at the entrance of the shaft, while flared outlet channel provides a longer fluid residence time in the outlet channel. However, flared inlet/outlet channels require more construction efforts, they are costly, and they have larger environmental footprint. The effect of shape of the inlet and outlet channels is evaluated in this section. Figure 3 - 20 compares the efficiency ratio of the original Treatment 85 S haft and the shaft with straight inlet and outlet channels . It can be seen that changing the shape of the inlet and outlet channel has no remarkable influence on the separation efficiency. Therefore , this design variation is recommended, due to significant construction costs reduction, and reduced environmental footprint. The separation efficiency for the given PSD (see Figure 3 - 3 ) is 34.8%, and 35.2% for the shaft with flared channels, and the shaft with straight channels, respectively. Summary and conclusions In this study, the Treatment Shaft technology as a CSO det ention system was evaluate. For this purpose, the residence time distribution was measured to assess the contact time of the system . Then, the natural sedimentation was investigated, and it was seen that 50% of particles with diameter of 161.3 µ m settle at the bottom of the shaft. Using the result from natural sedimentation, it was concluded that for a CSO event with a common PSD, for particles larger than 161.3 µ m, the separation efficiency is 81.7%, while this value is only about 20.8% for particles smaller than 161.3 µ m. Therefore, the impact of using a flocculating 86 agent was estimated using a population balance model coupled with the multiphase mixture model. Results showed that the incorporated flocculant can significantly improve the particle aggregation and increase the Sauter mean diameter across the domain. That leads to an increase in separation efficiency for small particles (smaller than 161.3 µ m ) from 20.8% to 4 9.3 %. Lastly , an investigation on design alternatives showed that the Treatment Shaft with a straight inlet and outlet channels, gap distance of B=9.32m, and baff le located at center provides the highest efficiency with reduced cost and environmental footprint compared to the original design. 87 Analytical solutions for laminar flow in membrane channels with cylindrical symmetry: Single - and dual - membrane systems Introduction Crossflow membrane filtration (CFF) has a broad range of applications in food [ 49 ] , biotechnology [ 48 ] , and petroleum [ 45 ] industries , as well as in desalination [ 164 ] and waste water treatment [ 165 ] . To ev aluate the performance of a given CFF system, the flow field must be known. Although computational fluid dynamics (CFD) offers reliable techniques and tools for finding velocity and pressure fields, analytical solution can be very useful for rapidly estim ating the performance of such systems under various operational conditions. Approximate analytical solutions can also be used to initialize the flow field variables in order to enhance the convergence of numerical simulations. Available mathematical models can predict the permeate flux but not the velocity and pressure profiles [ 166 , 167 ] . In 1953, Berman [ 168 ] developed analytical solutions for the Navier Stokes equations of fluid flow in a rectangular slit channel with two porous walls. He used perturbation method and assumed a steady state laminar flow, in which the velocity of fluid leaving the porous walls (permeate velocity) did not vary with the position along the channel. Afterwards, Berman [ 169 ] as well as Terrill and Shrestha [ 170 ] developed solutions for flow through cylindrical channels with porous walls. Yuan and Finkelstein [ 171 ] employed a similar approach to that of Berman [ 168 ] for a cylindrical channel with uniform permeate velocity. Kozinski et al. [ 172 ] proposed an arbitrary exponential function of wa ll velocity to calculate the pressure drop along a porous tube. Moussy and Snider [ 173 ] proposed an analytical expression for 88 low - Reynolds number laminar flow over pipes with suction and injection through the porous wall. Erd ogan and Imrak [ 174 ] developed a solution based on Bessel function of the first kind for velocity field in a uniformly porous pipe. Mellis et al. [ 175 ] presented analytical solutions for flow in a porous tube with two differe nt boundary conditions on the porous surface, but the solutions were not valid for wall Reynolds numbers larger than 0.5; moreover, they concluded that the models were useful qualitatively, but could not predict the performance accurately. Many authors use mathematical models [ 176 - 184 ] not depend on the pressure drop across the pipe. For higher recovery systems this assumption may be too strong. To address this issue, Karode [ 185 ] presented an analytical approach to derive the pressure drop in a cylindrical membrane channel. However which results in errors for radial systems in which the thickness of the membrane is comparable to the diameter of the tube. Analytical approximations, as d iscussed above, have drawbacks that should be addressed. The currently available solutions have at least one of the following issues: 1) They are derived for flows with very low Reynolds numbers (and not for all laminar flow range); 2) They are not valid f or wall Reynolds numbers of higher than 0.5; 3) They do not represent the pressure drop across the channels; 4) They show high errors for radial systems; 5) They involved perturbation or similarity solutions with complex algebraic calculations. 89 In this stu dy, approximate analytical solutions for the laminar flow are presented that further accuracy in radial flows. Moreover, solutions are presented for tubular membrane system s of several configurations, using different Reynolds numbers and different values of membrane permeability. The accuracy of the derived solutions is validated by comparing them with numerical simulation results. Approach 4.2.1. Analytical solution approach (d)) is a dual - membrane configuration. The membranes are modeled as porous media. 90 (a) (b) 91 (d) 4.2.1.1. Flow through a conventional membrane for crossflow filtration The fully developed velocity field for a cylindrical channel with imper meable walls and no - slip boundary condition (Poiseuille flow) is given by ( 4 - 1 ) 92 where is the axial velocity component, is dynamic viscosity, and is pressure of the fluid. The flow rate can be obtained by integrating over the cross - sectional area of the flow domain: ( 4 - 2 ) The derivative of Eq. ( 4 - 2 ) with respect to the axial coordinate is ( 4 - 3 ) which should be equal to zero for a cylinder with impermeable wall. However, for a cylinder with permeable wall, with very small fluid velocity through the membrane ( ), Eq. ( 4 - 3 ) can be modified to yield ( 4 - 4 ) where is the velocity of the fluid leaving the permeable wall of the cylinder (i.e. , and the pressure is ass umed to vary with , i.e. . ( 4 - 5 ) where is permeate flow rate, which is the radial flow through the membrane. is permeability, and is the area of the membrane available for permeation. is the unit vector 93 normal to the inner surface of the membrane. For the case of a one - dimensional flow in the membrane and constant wall permeability, the velocity of fluid t hrough the membrane due to radial pressure difference in the membrane (i.e. permeate flux), is given by ( 4 - 6 ) where is thickness of the membrane and is the pressure on the permeate side of the membrane. For radial flow through a membrane, Eq. ( 4 - 5 ) can be rewritten as and simplified for a cyl indrical membrane as ( 4 - 7 ) By rearranging Eq. ( 4 - 7 ) and integrating from to : ( 4 - 8 ) where ter radius, and is the pressure at the distance from the centerline. Hence, the permeate flux the inner surface of the membrane ( ) can be calculated as . So Eq. ( 4 - 8 ) can be rewritten as ( 4 - 9 ) It should be noted that a key assumption is that the permeate velocity, , is constant along the channel. The flow inside the membrane is due to the pressure difference of 94 permeate side and wetted side (exposed to the axial flow). The pressure at the permeate side is constant, and the pressure changes at the wetted side is very small. Therefore, this assumption is valid. Substituting and rearranging Eq. ( 4 - 4 ) to obtain gives ( 4 - 10 ) By using the change of variables and , we integrate . ( 4 - 11 ) The solution of the Eq. ( 4 - 11 ) can be easily found as ( 4 - 12 ) where and are integration constants. Imposing the inlet flow rate at gives ( 4 - 13 ) and assuming the inlet pressure at , one can find . The following solution i s obtained ( 4 - 14 ) 95 The pressure is then given by . ( 4 - 15 ) The pressure distribution is required in order to find the axial and radial velocity fields. Using Eq. ( 4 - 9 ) , the radial velocity of the fluid on the inner surface of membrane can be found as , ( 4 - 16 ) and the radial flow rate at each axial point can be obtained by ( 4 - 17 ) With the radial flow rate known, the axial flow rate can be simply calculated as ( 4 - 18 ) So a corrected velocity field, that accounts for the weak flow through the wall, at each section in the fluid field, is given by ( 4 - 19 ) To find the radial velocity profile in the tubular membrane with no rotati onal velocity, we use the continuity equation 96 ( 4 - 20 ) By integrating the equation above, the radial velocity in the cylindrical channel can be obtained as ( 4 - 21 ) Eq. ( 4 - 21 ) needs one boundary condition to calculate the integration constant. To obtain dimensionles s equations, the dimensionless variables are defined in Table 4 - 1 . Dimensionless variable Definition Aspect ratio of the membrane Radial location Axial location Outer radius of the inner membrane Outer radius of the outer membrane Inner radius of the inner membrane Eigenvalue Flow rate Axial velocity Radial velocity Darcy number (outer membrane) Darcy number (inner membrane) Reynolds number Pressure Using values from Table 4 - 1 , the dimensionless forms of Eqs. ( 4 - 15 ) , ( 4 - 16 ) , and ( 4 - 19 ) are ( 4 - 22 ) 97 ( 4 - 23 ) ( 4 - 24 ) where . Furthermore, the ratio of the flow through permeable membrane to the ratio of the inlet flow, , is given by ( 4 - 25 ) 4.2.1.2. Flow through an annulus channel with porous walls In these configurations (Cases B, C, and D) a cylindrical porous/impervious channel of radius is located at the center of a tubular porous/impervious channel of radius and length . The benefit of these configurations is an increased shear at the wall. A laminar flow with rate of enters at and leaves the tube at . The fully developed velocity field for the mentioned flow is given by [ 186 ] : ( 4 - 26 ) By integrating the velocity field, the flow rate in the channel can be calculated as ( 4 - 27 ) 98 4.2.1.2.1. Tubular membrane with an impermeable cylinder at center In the case of an impermeable cylinder at the center of a tubular membrane, the r adial (permeate) flow only passes through the outer membrane. Thus, assuming the constant permeate velocity , the derivative of flow rate is given by ( 4 - 28 ) Substituting Eq. ( 4 - 27 ) in Eq. ( 4 - 28 ) leads to ( 4 - 29 ) Eq. ( 4 - 29 ) can be simplified to Eq. ( 4 - 11 ) , by using and defining . By imposing the boundary conditions of inlet flow rate at , and outlet pressure at , and replacing with , the following solution can be obtained for pressure distribution along the membrane channel: ( 4 - 30 ) 99 The radial velocity in the membrane, radial flow rate, and axial flow rate at each point along the membrane axis can be obtained by Eqs. ( 4 - 16 ) to ( 4 - 18 ) . The corrected velocity profile is given by ( 4 - 31 ) The radial velocity profile in the channel can be found using Eq. ( 4 - 21 ) . Using Table 4 - 1 , the dimensionless forms of Eqs. ( 4 - 30 ) and ( 4 - 31 ) are ( 4 - 32 ) ( 4 - 33 ) where . Furthermore, the ratio of the flow through the membrane to the ratio of the inlet flow, , is given by 100 ( 4 - 34 ) 4.2.1.2.2. Membrane at the center of a cylindrical channel with solid (impermeable) wall In the case of a tubular membrane at the center of a cylindrical channel with impermeable wall, the radial (permeate) flow only passes through the inner membrane. Thus, assuming the constant permeate velocity, , the derivative of flow rate is given by ( 4 - 35 ) Substituting Eq. ( 4 - 27 ) in Eq. ( 4 - 35 ) leads to ( 4 - 36 ) Eq. ( 4 - 36 ) can be simplified to Eq. ( 4 - 11 ) , by using and defining . By imposing the boundary conditions of inlet flow rate at , and outlet pressure at , and replacing with , the pressure along the membrane can be calculated. Eq. 101 ( 4 - 15 ) can be used to obtain the pressure by replacing by . Moreover, the radial velocity in the membrane, permeate flow rate, and axial flow rate at each point along the me mbrane centerline can be obtained by Eqs. ( 4 - 16 ) to ( 4 - 18 ) . In addition, the axial and radial velocity profiles can be calculated using Eqs. ( 4 - 19 ) and ( 4 - 21 ) , respectively. Using Table 4 - 1 , the dimensionless forms of pressure profile is ( 4 - 37 ) where . The dimensionless velocity profile can be obtained using Eq. ( 4 - 33 ) . Moreover, the ratio of the flow through permeable membrane to the ratio of the inlet flow, , is given by ( 4 - 38 ) 4.2.1.2.3. Two tubular co - axial membranes (dual membrane configuration) In this configuration, the radial flow rate passes through both the inner and the outer membranes. Assuming constant permeability of the both outer ( ), and inner cylindrical 102 membranes ( ), and the permeate pressures of and for outer membrane and inner membrane, respectively, the derivative of flow rate is given by ( 4 - 39 ) in which, , , and are permeability, outer radius, and inner radius of the inner me mbrane, respectively. Substituting Eq. ( 4 - 27 ) in Eq. ( 4 - 39 ) leads to ( 4 - 40 ) By replacing with , Eq. ( 4 - 40 ) can be converted to ( 4 - 41 ) Eq. ( 4 - 41 ) can be simplified to ( 4 - 42 ) By defining , and as 103 , ( 4 - 43 ) Eq. ( 4 - 42 ) can be easily solved as ( 4 - 44 ) By imposing the boundary conditions of inlet flow rate at , and inlet pressure at , and replacing with the following solution can be obtained for pressure along the membrane: ( 4 - 45 ) The radial velocity on the porous surface of the membranes can be found as ( 4 - 46 ) and ( 4 - 47 ) 104 Thus, the radial (permeate) flow rate at each axial point can be obtained by ( 4 - 48 ) Having the radial flow rate, one can find the axial flow rate using Eq. ( 4 - 18 ) and the corrected velocity using Eq. ( 4 - 31 ) . Moreover, the radial velocity profile in the channel can be calculated by Eq. ( 4 - 21 ) . The dimensionless form of pressure equation is ( Table 4 - 1 ): ( 4 - 49 ) where and . The velocity profile can be obtained using Eq. ( 4 - 33 ) . Moreover, the ratio of the flow through permeable surfaces to the ratio of the inlet flow, , is given by 105 ( 4 - 50 ) 4.2.2. Numerical simulations and validation T able 4 - 2 : Details of the different cases and their boundary conditions Setup Parameter A B C D Length (mm) 250 250 250 250 Membrane radius (mm) 3 3 3 3 Inline membrane radius (mm) - 0.5 0.5 0.5 Outer membrane - Permeable Permeable - Permeable - Permeate pressure (kPa) 0 0 - 0 - Radius (mm) 5 5 - 5 - Permeability (m 2 ) 10 - 14 10 - 14 - 10 - 14 - Thickness (mm) 2 2 2 Inner membrane - - Impermeable solid Permeable Permeable - Radius (mm) - - 0.25 0.25 - Permeate pressure (kPa) - - 0 0.5 - Permeability (m 2 ) - - 10 - 14 10 - 14 - Thickness (mm) - 0.25 0.25 0.25 Flow rate (m 3 /s) 6.67×10 - 6 6.67×10 - 6 6.67×10 - 6 6.67×10 - 6 Mesh size - ~2.3×10 6 ~1.5×10 6 ~1.2×10 6 ~1.8×10 6 Mesh quality - 0.79 0.98 0.98 0.99 106 4.2.2.1. Boundary conditions The inlet velocity of the CFF is assumed to correspond to a fully developed flow in all the cases. Thus, Eqs. ( 4 - 1 ) and ( 4 - 26 ) are used for the inlet velocity. The transmembrane pressure is set to 50 kPa for all the simulations. The gauge pressure for the outer permeate surface, , is assumed to be zero, and for inner permeate surface, is assumed zero or 0.5 kPa in different cases (see Table 4 - 2 for details). No - slip wall boundary condition is imposed on all the walls. 4.2.2.2. Numerical solution procedure The Eulerian approach is employed for modeling the fluid flow. The velocity and pressure profiles in the membranes are obtained by solving the Navier - Stokes equation for lam inar flow. The computational domain is a hexahedral mesh, generated by ANSYS ICEM CFD 16.2 [ 187 ] . Mesh sizes and qualities are also ta bulated in Table 4 - 2 . For The solutions for coarser and finer grids are performed and mesh independency is found for all solutions. The Navier - Stokes equation for la minar flow is solved using a pressure - based algorithm of ANSYS® Fluent 16.2. Since the simulation is steady - state and single phase flow, the coupled algorithm is chosen due to its robustness [ 174 ] . The simulations are carried out for three dimensional fully devel oped flows in steady state and isothermal conditions. Water is used as the continuum phase, with density of 998.2 kg/m 3 , dynamic viscosity of 0.001003 Pa s. The momentum equation is discretized using the Quadratic Upstream Interpolation for Convective Kinetics (QUICK) scheme [ 175 , 188 ] , which is more suitable for hexahedral grid and provides higher accuracy [ 174 , 188 ] . Pressure interpolation is performed using Pressure 107 Staggering Option (P RESTO) scheme [ 189 ] , which is a more appropriate estimation for flows involving filter media, flows in curved domains and high swirling flows [ 174 , 188 ] . To ensure the stable solution, the C ourant Friedrichs Lewy (CFL) number is set to 200, and the explicit relaxation factor is set to 0.75 for both pressure and velocity. To achieve accurate solutions, the scaled residuals of the simulations for continuity and momentum equations are set to be less than 10 - 8 . Results and Discussion The accuracy of the analytical solutions is investigated to validate the presented analytical approach. Different cases (A, B, C and D) are studied, and normalized velocities, pressures and flow rates are compared. Th e mentioned variables are normalized using the parameters introduced in Table 4 - 3 . Parameters , , and are 0 mm, 3 mm, and 250 mm for case A; and 0.5 mm, 3 mm and 250 mm for cases B, C, and D, respectively Moreover, normalized radial velocities are compared to the solution approach of Karode [ 185 ] . In Table 4 - 1 , re fers to dimensionless pressure of exterior membrane ( ). In case of no exterior membrane, the pressure on interior membrane is used ( ). Table 4 - 3 : Normalization of flow field variables Para meter Flow rate Axial velocity Radial velocity Pressure Radius Axial position Normalized 4.3.1. Flow through a conventional membrane In this section the analytical solution for the flow through a conventional membrane is validated. A solution by Karode [ 185 ] was presented for this problem. Figure 4 - 2 shows axial 108 flow rate along the membrane channel, axial velocity profile in different cross sections, radial velocity profile, and pressure along the axial axis, respectively. Figure 4 - 2 (a) and (d) demonstrate a good agreement of the analytical solution with the numerical results. However, Figure 4 - 2 (b) uncovers the fact that by moving along the longitudinal axis of the CFF, the normalized axial velocity has comparatively high deviation from the numerical results, which is due to the approximation assumed in Eq. ( 4 - 19 ) , which is a fully developed parabolic approximation for the velocity profile. Figure 4 - 2 (c) illustrates the comparison of the present [ 185 ] and nu merical results. Both analytical results show good agreement with the numerical results when the normalized radius is smaller than 0.4. For bigger values of normalized radius, the presented solution both qualitatively and quantitatively. However, noting the velocity profile slopes near the wall, a good agreement approach. . Moreover, the results are accurate for wa ll Reynolds number of Wall Reynolds number is given by [ 175 ] : ( 4 - 51 ) Mellis et al. [ 175 ] concluded that their presented method provides a distorted parabolic velocity pr ofile for . 109 4.3.2. Tubular membrane with an impermeable cylinder at the center In this section the analytical solution for the flow through a membrane with a solid impermeable cylinder at the center of the tube is validated. This case is applicable for the filtration systems in which electric fields due to cylindrical electrodes at t he center improves the filtration [ 190 ] . Figure 4 - 3 shows axial flow rate along the membrane channel, axial 110 velocity profile in diffe rent cross sections, and radial velocity profile, respectively. Figure 4 - 3 (a), (b) and (c) show the high accuracy of the presented analytical solution. The small deviations in axial velocity are due to the approximation in Eq. ( 4 - 36 ) , whic h is the approximation of fully developed flow in the channel. Moreover, Figure 4 - 3 (c) demonstrate the high accuracy of the presented approach compared to that of t [ 185 ] . Wall Reynolds number is 0 .97 for this case. 111 4.3.3. Membrane at the center of a cylindrical channel with solid (impermeable) wall In this section the analytical solution for the flow through a tubular membrane positioned inside a cylindrical channel with impermeable walls is validated. Figure 4 - 4 s hows axial flow rate along the membrane channel, axial velocity profile in different cross sections, and radial velocity profile, respectively. Figure 4 - 4 (a), (b) and (c) show the high accuracy of the presented analytical solution. The small deviations in axial velocity compared to numerical results are due to approximation in Eq. ( 4 - 36 ) . It should be noted that the deviations in normalized axial velocity of this case are more noticeable compared to the other cases. For the same boundary conditions, the radial velocity (through the membrane) in the current case is almo st 4 times bigger than the other case, which results in a higher radial pressure gradient near the membrane surface. The high radial pressure gradient is in contrast with the assumption of one - dimensional (i.e. only axial) pressure gradients in the flow. Figure 4 - 4 (c) demonstrate the accuracy of the presented approach compared to that [ 185 ] . Wall Reynolds number is 0.70 for this case. 112 4.3.4. Two tubular co - axial membranes (dual membrane configuration) In this section the analytical solution for the flow through two tubular co - axial membranes is validated. Figure 4 - 5 illustrates axial flow rate along the membrane channel, axial velocity profile in different cross sections, and radial velocity profile, respectively. Figure 4 - 5 (a), (b) and (c) demonstrate the high accuracy of the presented analytical solution. 113 There are deviations in axial velocity when the normalized radius is smaller than 0.5 (i.e. in the vicinity of the inner cylindrical membrane). For values larger than 0.5, the normalized axial velocity is relatively accurate. Figure 4 - 5 (c) demonstrate the accuracy of the presented [ 185 ] . It is evident that the presented [ 185 ] has an error of up to 40%. Wall Reynolds number is 0.97 for this case. 114 4.3.5. Effect of Reynolds number and permeability on the solutions To investigate the ef fect of Reynolds number on the accuracy of the results, the cases A,B, C and D are solved for three different feed flow rates: 1.67×10 - 6 m 3 /s, 3.33×10 - 6 m 3 /s, 6.67×10 - 6 m 3 /s, and 1.0×10 - 5 m 3 /s. The axial velocities of each case are plotted and compared to the numerical results at ( Figure 4 - 6 ). 115 It is evident that by increasing the Reynolds number from 353 to 2112 for case A, and from 302 to 1810 for cases B, C and D, the deviation of the normalized axial velocity from the numerical solution increases. The reaso n is the fact that by increasing the Reynolds number, approximations in Eqs. ( 4 - 19 ) and ( 4 - 36 ) are no longer justifiable. In order to study the effect of permeability on the accuracy of the solutions, cases A and B are solved for three different permeability coefficients of 1 × 10 - 14 m 2 , 2 × 10 - 14 m 2 , and 3 × 10 - 14 m 2 (i.e. 10 mD, 20 mD, and 30 mD, respectively). Accordingly, the ratio of the permeate flow to inlet flow ( ) for different cases, varies from 0.16 to 0.69. It should be mentioned that the definition of Darcy number was presented in Table 4 - 1 . Darcy number is a dimensionless number representing the permeability of the membrane. However, the commonly used unit, L/(m 2 .h.bar), is also provided in the results (see Table 4 - 4 ). This unit represents the permeability divided by viscosity. The normalized axial velocities of the mentioned problems at are illustrated in Figure 4 - 7 . 116 Increasing the permeability of the membrane results in a higher deviation between the numerical and analytical solution s (i.e. the analytical solution has lower accuracy). The reason is that by increasing the permeability, radial flow rate increases and the effect of radial pressure gradient adversely affects the accuracy of derivative of axial flow rate and therefore the axial velocity profile. Table 4 - 4 presents the quantitative results, and mean 117 percentage error (MPE) for different cases and different permeability v alues. In cases where , the suction flow occurrs at outlet. Table 4 - 4 : Effect of membrane permeability values on the accuracy of axial velocity profile at Case Permeability Specific permeate flux (L/(m2.h.bar)) Mean percentage error of axial velocity (%) A 10 mD 0.23 0.97 1758 2.59 20 mD 0.46 1.95 3516 7.99 30 mD 0.69 2.93 5274 19.57 B 10 mD 0.23 0.97 1758 3.53 20 mD 0.46 1.95 3516 6.60 30 mD 0.69 2.93 5274 10.50 C 10 mD 0.16 0.70 13044 8.95 20 mD 0.33 1.42 26904 14.86 30 mD 0.50 2.13 40764 20.23 D 10 mD 0.39 0.97 2725 6.18 20 mD 0.79 1.95 5521 12.19 30 mD 1.18 2.93 8246 23.03 In cases A and B for , the maximum MPEs are 19.57% and 10.50%, respectively. The maximum MPEs for case B and C, are 20.23% and 23.03%, respectively. In the latter case, the radial flow is the dominant flow regime. Conclusions The flow field in different membrane configurations was analytica lly solved using an approximate method. The present analytical solution is applicable for micro - filtration (MF), Ultra - filtration (UF), and Nano - filtration (NF), where the flow rates are in laminar flow range, and as long as no - slip boundary conditions are valid. It should be mentioned that in the present solution, salute and foulant is not taken into consideration. The resulting expressions are simple and compare favorably to numerical simulations for different Reynolds numbers and membrane permeabilities . The introduced approximate solutions are 118 accurate for all laminar flows even when the radial flow is relatively large. The solution has small deviations for higher Reynolds numbers (up to for single tubular membrane, and up to for othe r configurations), and provides an approximation with accuracy superior to that of the solutions reported earlier. The analysis of the effect of membrane permeability on the solution accuracy shows that higher permeability translates into increased deviati on of the approximate solutions from the well - converged numerical solutions. Higher permeability of the porous medium increases the radial flow rate and radial pressure gradient, which undermines the assumption of one - dimensional (i.e. axial) pressure grad ient in the membranes. However, it is shown that for high permeability values in which the radial flow is significantly higher than the axial flow, the present analytical solution deviates from the numerical solutions by 25% only. We conclude that the prop osed analytical approach gives an acceptable solution for velocity and pressure profiles in laminar fluid flows in cylindrical membranes and dual - membrane systems. Acknowledgement This material is based upon work supported in part by the National Science Foundation Partnerships for International Research and Education program under Grant IIA - 1243433 . Nomenclature body force (N) normal vector velocity vector (m/s) 119 area (m 2 ) length of the membrane (m) pressure (Pa) volumetric flow rate (m 3 /s) radius (m) eigenvalue permeate velocity (m/s) volumetric flow rate (m 3 /s) radius (m) thickness (m) velocity (m/s) Greek symbols permeability (m 2 ) membrane thickness (m) angular coordinate 120 dynamic viscosity (Pa.s) kinematic viscosity (m 2 /s) density (kg/m 3 ) Subscripts entrance exterior membrane inner interior membrane outer radial coordinate axial coordinate Abbreviations CFF Crossflow filtration MPE Mean percentage error 121 Filtration of charged oil droplets in stationary and rotating tubular membranes Abstract Cross flow filtration (CFF) is a common membrane separation process with applications in food, biochemical and petroleum industries. In particular, m embranes can be used for liquid - liquid separation processes such as needed in oil - water separation. A major challenge in cross flow filtration is membrane fouling. It can decrease significantly the permeate flux g can be mitigated by inducing shear on the possible approach to improve membrane efficiency consists of repelling droplets/particles from the porous surface towar d the centerline using a repulsive electric force. For this purpose, the surface of the membrane can be exposed to electric potential and droplets/particles are also induced to have the same electric charge. In this work, numerical simulations of charged n on - deformable droplets moving within an axially rotating charged tubular membrane are performed. The results show that by increasing the electric potential on the membrane surface, the repelling force increases which obviously improves the grade efficiency of the membrane. However, the electric field gradients found in the flow field require large potentials on the membrane surface to observe a noticeable effect. Hence, a smaller solid cylinder is located in the centerline of the flow channel with zero pote ntial. This solid cylinder enhances the electric field gradient in the domain which results in higher repelling forces and larger grade efficiency of the membrane 122 at small potentials. The addition of a small cylinder in the flow field also improves the gra de efficiency increases due to the higher shear stress near the membrane surface. Introduction Cross - flow filtration (CFF) is a common membrane separation process which has bro ad applications in food industry [ 191 ] , biotechnology [ 192 ] , desalination [ 193 ] and wastewater treatment [ 47 ] . Membranes can be used for liquid - liquid separation processes and potentially for oil - water separation [ 17 ] . A major challenge in oil - water separation is fouling of the membrane by the oil. It can dramatically decrease the permeate flux and therefore sho rten the membrane life and its efficiency [ 194 ] . Depending on the application, CFF is possible under steady - state or time - dependent operations. For example, for solid - liquid separations, the deposited layer of particles thickens and results in time - varying permeate flux [ 195 ] . However, the shear stress generated by the crossflow can decrease the particle accumulation on the membrane surface [ 196 ] In this regard, the induction of swirl in a stationary membrane was shown to be significantly ben eficial [ 58 ] . The flow swirl also improves the separation efficiency of the mixture, when the d ensity of droplets are less than that of continuous phase (e.g. oil droplets in water). In this situation, droplets move toward the core of a vortex due to centrifugal force. The effect of rotational flow in rotating disk membrane separation was studied by [ 194 , 197 - 199 ] and the effect of rotational flow through rotating porous annuli was investigated by [ 200 - 203 ] . 123 Dispersed liquid droplets are often charged and surrounded by ions of opposite sign from the continuous phase. The arrangement of the electric charges on the oil droplet surfaces with the balancing charge in the continuous phase is often called the electrical double layer (EDL). Numerous studies to infer the electric charge by electrophoresis techniques and measure the zeta potential of oil drop lets in different aqueous emulsions or in ionic surfactant solutions [ 204 ] . It was observed that the electrical characteristics of the droplets in an emulsion can be easily modified by adsorbing surfac tant onto the liquid - liquid interfaces [ 205 , 206 ] . The addition of surfactant and/or alkali can increase the zeta potential of the oil droplets and make the particles negatively or po sitively charged [ 207 ] . This leads to the possibility of improving membrane efficiency by u sing a repulsive elect ric force [ 208 ] . In this study, the effect of an electric repelling force on the efficiency of the membrane separation is investigated. Moreover, the effectiveness of swirling flow and electric repelling force are compared. The efficiency of electric repelling force on dif ferent droplet sizes, and in two different membrane geometries are also studied. Modeling Approach The studied cross flow filtration system is a rotating tubular membrane. The geometry of the system is shown in Figure 5 - 1 . The CFF system has a membrane rotating around the longitudinal z - axis with angular velocity of . The membrane has length of 250 mm, is modeled as a filter with a thickness of 2 mm, inner diamete r of 6 mm, and permeability of 1×10 - 14 which is based on the ceramic membranes manufactured by TAMI industries. 124 The Eulerian - Lagrangian approach is employed for modeling the continuous phase (water) and dispersed phase (oil droplets). After a three - phase separator used in produced water treatment, the concentration of oil droplets in water are typically much lower than 10% [ 17 ] . Hence, the droplets are assumed to not interact wi th each other and have no effect on fluid domain [ 154 , 155 ] . The velocity and pressure profiles in the membrane are obtained by solving the Navier - Stokes equation for laminar flow and the dispersed phase is solved by tracking oil droplets through the calcula ted flow field using a Lagrangian approach. The electric field in the fluid flow must be computed to calculate the electric repelling forces acting on the droplets. Two different CFF systems are considered in this regard. In the S imple CFF (SCFF) , the elec tric potential is applied on the membrane surface and fluid inlet and outlet in the tube have no potential. End - effects change the electric field in the domain. In the O ptimized CFF (OCFF) , a small solid cylinder is located inside the channel to increase t he electric field gradients in the fluid flow. Figure 5 - 2 illustrates the schematics of an OCFF. 125 The computational domain is a hexahedral mesh, generated by ANSYS ICEM CFD 15 [ 187 ] . For the SCFF, the total number of cells is 2.13×10 6 and the minimum orthogonal quality of the mesh is 0.79. For the OCFF, the number of cells is 1.50×10 6 and the minimum orthogonal quality of the mesh is 0.98. 5.3.1. Continuous phase modeling The fluid system studied is assumed to be incompressible, isothermal and Newtonian. It is described by continuity and Navier - Stokes equations . The continuity equation for an incompressible fluid can be simplified to: ( 5 - 1 ) where is the velocity vector. The Navier - Stokes equation has a source term, included accounting for the electrohydrodynamic forces: ( 5 - 2 ) 126 In Eq. ( 5 - 1 ) , , and The electrohydrodynamic force is given by [ 209 ] : ( 5 - 3 ) where is the electric field, is the charge density, is the permittivity, and is the temperature. The first and second terms represent electrophoreti c and dielectrophoretic forces, respectively. The last term is called electrostriction force. The electrophoretic term is due to the action of the electric field on the free volume charge and can be neglected when the net free charge in bulk fluid is zero. The dielectrophoretic force (2 nd term) is due to effect of field on the polarization volume charge, which is due to nonunifor mity of the permittivity. Since the permittivity of fluid and droplets are assumed to be constant, this force is only active on droplet interfaces on which there is a sharp change in permittivity. The electrostriction force (3 rd term) is due to variation of electric permittivity with fluid density [ 209 , 210 ] and since we assume an incompressible fluid this term can be neglected [ 209 ] . To calculate the electric field in the domain, it can be shown that the electric field is both irro tational ( ) and divergence free ( ). This leads to the definition of electric potential, , that satisfies the Laplace equation. ( 5 - 4 ) Eqs. ( 5 - 2 ) and ( 5 - 4 ) need to be solved to obtain the velocity field and electric field in the domain. The follow ing assumptions are made: 127 The free charges in the domain is zero. The permittivity of the both fluid and droplets are constant and uniform. Since the size of the droplets are much smaller compared to the characteristic length of the domain (~1/100), the ef fect of droplets on the continuous phase is neglected. Therefore, the change of permittivity on the interface of fluid and droplets are neglected as well. Thus Eqs. ( 5 - 2 ) and ( 5 - 4 ) are decoupled and can be solved separately during the solution procedure. 5.3.2. Dispersed phase modeling The motion of droplets is obtained by tracking their trajectories in the flow field. The trajectory of the dispersed phase, can be predicted in a Lagrangian frame of reference as ( 5 - 5 ) The dispersed phase velocity in the right hand side of Eq. ( 5 - 5 ) is obtained by ( 5 - 6 ) where and are instantaneous velocities of continuous and dispersed phases, respectively. Also, is the drag force per unit mass of droplet, 128 is the drag coefficient and is the diameter of droplet. and are densities of continuous and dispersed phases, respectively. The term is overall acceleration due to high pressure gradient and virtual mass force. Table 5 - 1 shows the detailed definition of each of these accelerations. Presence of electric field, can absorb or repel a charged particle. The force acted on a particle due to an electric field of is given by ( 5 - 7 ) in which, is the charge of the droplet in Coulomb. 5.3.3. Boundary conditions The inlet velocity of the CFF is assumed to be fully developed flow. For the axial flow in a SCFF the fully developed velocity profile is ( 5 - 8 ) where is the inlet volumetric flow rate, and is the radius of the cylinder. On the other hand, for the axial flow between two concentric cylinders the fully developed velocity profile can be calculated using the following equation [ 186 ] 129 ( 5 - 9 ) where and are the radii of inner and outer cylinders, respectively. The coefficient is given by ( 5 - 10 ) The transmembrane pressure is set to 50kPa for all the simulations. The gauge pressure for the permeate surface is assumed to be zero. The feed flow rate of 6.67×10 - 6 m 3 /s at the inlet boundary is used for the simulations. This flow rate corresponds to Reynolds numbers of 1415 for SCFF and 1207 for OCFF. In addition, to induce the swirl to the flow, the membrane rotates about the longitudinal axis at the angular velocity of . The following velocities are used for different cases: 0 (no rotation), 500 and 1500 rpm. In order to impose the electric field in the fluid flow, an electric potential of is applied to the membrane surface. For both SCFF and OCFF, the potentia l at the inlet and outlet boundaries are assumed to be zero. However, for the OCFF, the internal solid cylinder is assumed to have zero potential. Numerical Method The Navier - Stokes equation for laminar flow is solved using a pressure - based algorithm of AN SYS® Fluent 1 5 .0 . Since the simulation is steady - state and single phase flow, the coupled algorithm is chosen due to its robustness [ 211 ] . The si mulations are carried out for 130 a three dimensional fully developed swirling flow in a steady state and isothermal conditions. The momentum equation is discretized using the Quadratic Upstream Interpolation for Convective Kinetics (QUICK) scheme [ 211 , 212 ] , which is more suitable for hexahedral grid and provides higher accuracy [ 211 , 213 ] . Pressure interpolation is performed using Pressure Staggering Option (PRESTO) scheme [ 189 ] , which is more appropriate estimati on for flows involving filter media, flows in curved domains and high swirling flows [ 211 , 213 ] . To ensure the stable solution, the Courant Friedrichs Lewy (CFL) number is set to 200, and the explicit relaxation factor is set to 0.75 for both pressure and velocity. To achieve the accurate solutions, the scaled residuals of the simulations for continuity, momentum and electric potential equations are set to be less than 10 - 8 . To simulate the motion of dispersed droplets, groups of droplets with specific gravity of 0.73 (diesel oil) and different diameters are injected at the inlet plane. To obtain the trajectory, Eq. ( 5 - 6 ) is integrated using a 5 th order Runge Kutta scheme [ 214 ] . The maximum number of time steps is set from 5×10 6 to 5×10 7 , depending on the boundary condit ions and droplet sizes. Also the tolerance of the solution is set to 10 - 6 to control the accuracy. In this work, three different types of problems are solved: Type I: the rotating membrane has no electric potential (SCFF). Type II: the rotating membrane ha s electric potential on its surface (SCFF). Type III: the rotating membrane has electric potential on its surface, and a smaller conductive solid cylinder with zero potential is located at the centerline of the flow channel (OCFF). 131 The different types of s olved cases are tabulated in Table 5 - 2 . Case ID Type Membrane Type Angular Velocity (rpm) flow rate (m3/s) Electric Potential on the membrane Inner cylinder SCFF - I - 0 I SCFF 0 N/A N/A SCFF - I - 500 I SCFF 500 N/A N/A SCFF - I - 1500 I SCFF 1500 N/A N/A SCFF - II - 0 II SCFF 0 + N/A SCFF - II - 500 II SCFF 500 + N/A SCFF - II - 1500 II SCFF 1500 + N/A OCFF - 0 III OCFF 0 + + OCFF - 500 III OCFF 500 + + OCFF - 1500 III OCFF 1500 + + Results and discussions The flow pattern and wall shear stress components are compared for all the mentioned types of membranes. In addition, the grade efficiency is investigated and compared for different cases and differ ent particle diameters. 5.5.1. Flow pattern The axial, tangential and radial velocity components are presented for two membrane types of SCFF and OCFF for different angular speeds. Figure 5 - 3 shows the normalized axial velocity at the position of z/L=0.5 in the membrane for different angular velocities. The normalized axial velocity is given by . It is evident that the maximum velocity occurs at the different radius for SCFF and OCFF membranes. So the maximum axial velocity occurs at for SCFF membranes and at for the presented OCFF membranes. 132 Figure 5 - 4 demonstrates the normalized radial velocity of SCFF and OCFF membranes for different angular velocities. The normalized radial velocity is given by , so that . It should be noted that is the overall permeate flow rate. Although the radial velocity is relative ly larger for SCFF membranes, the value is almost identical for both types of membranes. 133 Figure 5 - 5 demonstrates the normalized tangential velocity of SCFF and OCFF membranes for different angular velocities. The normalized tangential velocity is given by . It is observed that for lower angular velocity, OCFF membrane has the lower tangential velocity compared to a SCFF membrane. However, by increasing the angular velocity, the tangential velocity of the both SCFF and OCFF membranes tends to be identical. 5.5.2. Wall shear stress The shear force acting on the membrane surface is calculated from ( 5 - 11 ) is the normal component and it can be ignored since it is normal to the membrane surface. The terms and are the shear components in longitudinal and tangential directions, respectively. These terms can be obtained using the following expressions 134 ( 5 - 12 ) ( 5 - 13 ) In addition, the average longitudinal and tangential wall shear stresses are numerically calculated by the following equations ( 5 - 14 ) ( 5 - 15 ) The overall wall shear stress and its angle are calculated as ( 5 - 16 ) ( 5 - 17 ) 135 Figure 5 - 6 shows the average shear stress components for both SCFF and OCFF membranes for different angular velocities. It is evident that increasing the angular velocity, increases both axial and radia l components of wall shear stress. In addition, both axial and radial components of shear stress for OCFF are higher than the ones for SCFF. Figure 5 - 7 also demonstrates the overall wall shear stress ( ) and its direction ( ) for both SCFF and OCFF membranes for different angular velocities. 136 5.5.3. Grade efficiency The total grade efficiency is the ratio of the total mass of droplets that leave the membrane channel to the retentate stream. Since the droplets are assumed to have no interaction with fluid domain and also to have no break up or coalescence, the total grade efficiency can be simplif ied to: ( 5 - 18 ) in which and are the initial number of droplets injected to the flow at inlet, and the number of droplets at outlet with retentate stream, respectively. The grade efficiency presents the effect of both membrane filtration and the vortex separation. 5.5.4. Effect of electric repelling forces on droplet separation To investigate the effect of electric repelling force on the droplet separation, electric potential is app lied on the membrane surface on problems type II and III (see Table 5 - 2 ) . For the constant electric potential of on the membrane surface, different droplet c harge densities of 0, 2, 4 and 8 C/m 3 are used. In addition, the effect of droplet size ( ), on the effectiveness of droplets separation is studied. 137 Figure 5 - 8 shows that grade efficiency of the non - rotating membrane for both of SCFF and OCFF membranes are almost the same. For a r otating membrane (500 rpm), Although the grade efficiency of both systems are similar in droplets smaller than 50 micron, OCFF shows a slight advantage compared to SCFF for droplets bigger than 50 micron. However, for charged droplets, because of higher el ectric gradient in the flow domain, OCFF is more efficient and all the droplets of will leave the membrane with retentate stream. While for the similar system of SCFF, the grade efficiency is 1 for droplets with . 138 Figure 5 - 9 , compares the grade efficiency of the rotating membranes (500 rpm) with charge densities of 4 and 8 C/m 3 . The results show that by increasing the charge density for OCFF membranes, the cut - o ff diameter decreases and for charge densities of 4 and 8 C/m 3 , cut - off diameters are and , respectively. In addition, both Figure 5 - 8 and Figure 5 - 9 demonstrate that the grade efficiency of all cases improves by increasing the droplet size, which highlights the effectiveness of electric repell ing for on separation of smaller droplets. 5.5.5. Reduced Grade Efficiency To examine the effect of vortex induced forces (in absence of membrane filtration) or electric repelling forces, reduced grade efficiency can be defined. The reduced grade efficiency el iminates the effect of main flow stream on particle separation and quantifies the performance of the separation due to only swirling flow and/or electric repelling force. The reduced grade efficiency is given by ( 5 - 19 ) where . The Stokes number is defined as the ratio of the characteristic time of a droplet ( to a characteristic time of the flow ( . For high Stokes numbers, a droplet will continue along its initial trajectory, while for low Stokes number, a droplet will be affected by flow inertia and follow fluid streamlines. The Stokes number is defined as 139 ( 5 - 20 ) The parameter is the characteristic length of the fluid flow. For a SCFF, , and for an OCFF, . The effect of electric repelling force on the reduced grade efficiency is studied below. Figure 5 - 10 illustrates that reduced grade efficiency of an OCFF membrane with angular velocity of 1500 rpm and different charge density of particles. It indicates that the reduced grade efficiency rapid ly increases by increasing the charge density on the droplet. In addition, according to reduced grade efficiency for droplets with no charge, it is evident that rotating membrane and induced swirl is only effective for . For the present te st case, electric repelling forces can enhance the performance of the membrane and expand the range droplet size of existing crossflow filtration system. 140 Summary and conclusions Numerical simulations of charged non - deformable droplets moving within an axia lly rotating charged tubular membrane are performed. It is observed that by increasing the angular velocity, due to the increase in wall shear stress, the grade efficiency of the separation improves rapidly. However, this increase is more significant for t he larger droplets (corresponding to higher Stokes numbers). Hence, the electric repelling force is employed to enhance the separation efficiency. In this regard, for the constant flow rate, it is shown that the electric repelling force is more effective f or a membrane system with an electrode, compared to a membrane without an electrode in its center. The results show that, the repelling force improves significantly the grade efficiency of a membrane when an electrode is introduced at the center of the mem brane. Electric repelling forces can enhance the performance of the membrane and expand the range droplet size of existing crossflow filtration system especially for droplets with very small Stokes number (< 10 - 4 ) which are notoriously difficult to separat e and will foul membranes. Nomenclature Acceleration vector (m 2 /s) Electric field vector (V/m) Force vector (N) Grade Efficiency Pressure (Pa) 141 Charge (C) Time (s) Velocity vector (m/s) Location (m) Permittivity (F/m) Electric potential (V) Dynamic viscosity (m 2 /s) Density (kg/m 3 ) Stress tensor (Pa) Subscript C Continuous phase d Droplet D Dispersed phase el Electrohydrodynamic f Flow i Inlet o Outlet 142 Abbreviations CFF CrossFlow Filtration OCFF Optimized CrossFlow Filtration RGE Reduced Grade Efficiency SCFF Simple CrossFlow Filtration TGE Total Grade Efficiency Acknowledgments S upport for this work from the National Science Foundation Grant IIA - 1243433 is gratefully acknowledged. 143 Modeling membrane fou ling using a wall - film model Introduction Separation of dispersed droplets from a continuous phase is an important stage in applications found in biology [ 215 ] , chemistry [ 216 ] , food sciences [ 217 ] , and petroleum engineer ing [ 54 , 55 ] . In particular, in the oil and gas industry, separation of oil from an oil - in - water emulsion can become a challenge depending on droplet sizes, flow regime, and mixture components, which can be widely different for various oil fields [ 57 , 207 ] . For separation of oil droplets suspended in water, common methods of droplet removal includes centrifugation, air flotation, and hyd rocyclones [ 58 ] . Although these methods are effective for a wide range of flows, small droplets can escape from these systems, and their efficiency declines for droplets smaller than approximately 15 µm [ 54 , 218 ] . Membrane separation is a technology that selectively separates a secondary material from a primary phas e via small pores and gaps fabricated in its continuous structure. Membrane technology has shown a remarkable performance in separation of droplets [ 219 ] , particles [ 220 ] , and solutes [ 2 21 , 222 ] . In recent years, the research on the use of membrane separation techniques for oil - water emulsions has been growing [ 17 , 57 ] . Microfiltration [ 35 , 223 , 224 ] , and Ultrafiltration membranes [ 50 , 225 , 226 ] were studied to evalua te the ir capability for separation of droplets of various sizes. Different membrane configurations such as flat - sheet membranes [ 51 ] , tubular membranes [ 227 , 228 ] , slotted pore filters [ 223 ] , and hollow fiber membranes [ 229 ] were evaluated. Crossflow filtration (CFF) membranes have many applications in various industries such as biotechnology [ 48 ] , food [ 49 ] , petroleum industries [ 45 ] , 144 desalination [ 164 , 230 ] , and wastewater treatment [ 46 , 47 ] . Figure 6 - 1 shows a schematic of the operation of a CFF membrane. Studies on the removal of oil droplets from water using membranes reported that the droplets rejected by a membrane form a layer at the membrane surface [ 50 - 56 ] , which is Membrane fouling, is the single largest issue in membrane - base treatment [ 54 ] ; it can cause a significant decrease in the permeate flux and reduces its efficiency. The permeate flux and the rate of rejection is correlated with operating parameters of the CFF membrane system such as feed concentration, temperature, crossflow velocity (CFV), and transmembrane pressure (TMP) [ 226 , 231 , 232 ] . Dickhout et al. [ 5 7 ] discussed membrane fouling by oil droplets and categorized them into five main mechanisms: (1) complete pore blocking (pore blocked by a large droplet); (2) standard (or partial) pore blocking (pore blocked by coats of small droplets); (3) intermedia te blocking (a layer of droplets narrow 145 (5) oil layer (impermeable oil layer) blocking. One of the strategies to reduce fouling is to decrease droplet deposition on the membrane surface using hydrodynamic methods. Tummons et al. [ 54 , 55 ] presented a direct visualization of the real - time images of a membrane surface fouling by oil droplets in crossflow filtration. Motin et al. [ 58 ] discussed the effect of membrane rotation on the performance of the membrane, and efficiency of the droplet separation. In addition to experim ental efforts on the membrane fouling, other analytical modeling and numerical simulations were also presented to understand the fouling phenomenon and flux decline in membranes. Grenier et al. [ 59 ] proposed a methodology to quantify the fouling mechanism of solid particles in dead - end filtration. They identified that surface blocking may reduce f iltering surface up to 80% . Darvishzadeh et al. [ 60 ] studied the behavior of an oil droplet at a membrane pore entrance using direct numerical simulations. Won et al. [ 61 ] computationally investigated the effect of patterned membranes on particle depositions. Hou et al. [ 62 ] provided a combined blocking and cake filtration model. Zhang et al. [ 63 ] simulated the membrane fouling under a baffle - filled flow. Several other analytical or numerical studies were carried out in which the fluid dynamics of the flow in a CFF membrane [ 64 - 66 ] and membrane fouling [ 67 - 71 ] are discussed. There is currently no model that can provide a comprehensive view of system design, flow patterns and droplet concentration near the membrane surface, droplet - droplet and droplet - fluid interactions in CFF membranes, and droplet transport to membrane surface. Such a model can offer valuable information for membrane operation and performance under various conditions. In this study, a novel model is introduced that can simulate early 146 stag es of membrane fouling and predict the permeate flux and the time of its occurrence in a CFF membrane system. The model is the first model based on an Eulerian - Eulerian multiphase simulation approach that uses the concept of Eulerian wall film model for tr ansport of droplets from bulk flow to the membrane surface. Additionally, a Population Balance Model (PBM) is utilized to account for droplet size distribution and droplet coalescence in bulk flow, providing more information about the droplet size distribu tion and its impact on membrane fouling . Modeling of fluid flow, droplet deposition, and membrane fouling The wall film model was originally developed to model droplet - surface interactions of diesel droplets with cylinder walls of internal combustion engin es [ 233 - 235 ] . Subsequently, it has been extensively used for modeling droplet - surface interactions in various applications [ 236 , 237 ] . The wall film model benefits from lubrication theory approximation [ 238 ] which results in a two - dimensional velocity profile in the film. In the current study, the wall film model is presented and adapted to simulate droplets transport from bulk flow to a membrane surface. 6.2.1. Phenomenological description of m embrane fouling mechanisms Dickhout e t al. [ 57 ] discussed membrane fouling b y oil droplets. Figure 6 - 2 shows various scenarios of droplet - membrane interaction in a computational cell (dashed box). In Figure 6 - 2 (a), there is no droplets in near - membrane cell. In Figure 6 - 2 (b) there a re some droplets in the cell, however, there is not any droplet deposited on the membrane surface. In Figure 6 - 2 (c), in addition to droplets in the cell, some droplets are deposited on the 147 membrane surface, but no pore blockage occurred. In Figure 6 - 2 (d), in addition to dropl ets in the cell, some droplets are deposited on the membrane surface. However, in contrast to the case of Figure 6 - 2 (c), due to high transmembrane pressure differenc e, droplets enter the pores and there is no permeate flux of continuous phase through the membrane. However, dispersed phase gradually enters the pores. In Figure 6 - 2 (e), besides the deposited droplets, due to coalescence of the droplets on the membrane surface, a thin film is formed on the membrane surface; the cell is partially covered and there is a permeate flux of continuous phase. In Figure 6 - 2 (f), the film formed in Figure 6 - 2 (e) is expanded and covers the entire cell. Moreover, the high pressure - difference leads to entrance of oil into pores. In this study, a model for droplet deposition during stages of (a) through (c) is introduced. Developing of a criterion for ne 148 6.2.2. Multiphase flow modeling Multiphase flow models are nee ded to estimate the transport of droplets in the continuous phase, droplet - droplet interactions, and transport of droplets to the membrane surface. The averaging procedure of the flow field equations used in study is discussed below. Averaging of droplets in the bulk flow is first reviewed. Then, droplet - droplet 149 interactions are modeled using a population balance model. Lastly, the averaging method on a surface is proposed. 6.2.2.1. Fluid averaging in bulk flow A mixture of - phase flow is considered as one continu ous phase (subscript ), and dispersed phases. The void fraction is defined as the ratio of the volume occupied by phase to the total volume as ( 6 - 1 ) where is a control volume, and takes a value of 1.0 where phase is present, and is zero otherwise. It is found that . Any other property, can be averaged as ( 6 - 2 ) For instance, the effective averaged velocity of phase is given by ( 6 - 3 ) The mass conservation equation of phase can be given as ( 6 - 4 ) 150 where is the density, and is the mass transfer to and from phase . Since the mass that leaves one phase must add to another phase, , where is the mixture mass source term. Due to permeate flux at membrane surface, this term is nonzero for a near - membrane computational cells and is zero everywhere else in the domain. The continuity equation of the mixture is then obtained as ( 6 - 5 ) where and are averaged mixture density and averaged mixture velocity, respectively. and are given by ( 6 - 6 ) ( 6 - 7 ) The momentum conservation equation for phase is given by ( 6 - 8 ) where , , and are viscous stress, turbulence stress, and the effect of surface tension force, respectively. By summing the momentum equation of all phases, the conservation of momentum equation for the mixture can be obtained as 151 ( 6 - 9 ) where is the effect of surface tension force on the mixture and is defined as . The last term in the RHS of Eq. ( 6 - 9 ) is the mixture momentum source term. In the RHS of Eq. ( 6 - 9 ) , is average viscous stress defined as ( 6 - 10 ) is the average turbulence stress defined as ( 6 - 11 ) where is the drift velocity of phase and is defined as ( 6 - 12 ) is the average diffusion stress due to the phase slip given by ( 6 - 13 ) Manninen et al. [ 157 ] proposed the slip velocity as 152 ( 6 - 14 ) where is the droplet relaxation time which represents how quickly a droplet adjusts its velocity to the surrounding fluid flow. The droplet relaxation time is given by ( 6 - 15 ) where is the diameter of droplets. In Eq. ( 6 - 14 ) , is the acceleration of droplets given by ( 6 - 16 ) and is drag function of a droplet in free stream and is modeled using [ 158 ] ( 6 - 17 ) where is particle Reynolds number. 6.2.2.2. Population Balance Model A Population Balance Model (PBM) is employed to account for droplet - droplet interactions in the bulk flow and track changes in the distribution of droplet size . The Population Balance Equation (PBE) is a balance equation for the number density, , which is the number of sized droplet per unit volume of continuous phase. The general form of the po pulation balance equation for a flow without particle growth (nucleation) is given as [ 239 ] 153 ( 6 - 18 ) where is particle volume. The first term on LHS is transient and convection terms combined, the terms on RHS are birth due to coalescence, death due to coalescence, birth due to breakage, and death due to breakage, respectively. The birth due to coalescence is ( 6 - 19 ) where is coalescence kernel , which is defined as a product of coalescence frequency and collision efficiency . The equation above represents coalescence of droplets with volume of and droplets of volume , and forming of droplets of volume . The death due to coalescence is ( 6 - 20 ) The birth due to breakage is ( 6 - 21 ) where is the breakage kernel or the frequency of disruption of droplets of volume , and is probability density function that contains information on the fragments produced by a breakage event. The death due to breakage is ( 6 - 22 ) 154 In this study, a Method of Moment (MOM) approach is employed to solve the population balance equation. In MOM approaches, Droplet Size Distribution (DSD) is modeled by solving a balance equation for few moment terms which reduces the computational cost. Although the number of mo ments may not be sufficient to provide an exact DSD in each computational cell, it can provide Sauter mean diameter ( ) for each cell, defined as ( 6 - 23 ) Due to the importance of gravity force on droplet deposition, an inhomogeneous velocity approach is utilized to obtain a more accurate droplet transport in bulk flow. Therefore, Direct Quadrature Method of Moments (DQMOM) is used to solve population balance equation . DQMOM provides an inhomogeneous velocity field for various droplet groups [ 84 , 89 , 90 ] and has been extensively used for multiphase flow simulations [ 239 ] . It should be noted that in current study only droplet coalescence in bulk flow is considered, since droplet breakage is negligible because of flow conditions. See Appendix A for more details. 6.2.2.2.1. Droplet coalescence in bulk flow Different mechanisms promote collision of droplets and may cause droplet coalescence as a result. These collisions can be due to turbulent fluctuat ions, velocity gradients, eddy - capture, buoyancy and wake effect [ 118 ] . In this study, a laminar flow is investigated and therefore, only the coal escence of droplets due to velocity - gradients and buoyancy (shear - induced coalescence) is considered. Coalescence kernel is defined as [ 118 ] 155 ( 6 - 24 ) in which is coalescence frequency, and is collision efficiency. The model introduced by Friedlander [ 119 ] is used for coalescence. It was developed for a shear - induced collision in laminar flow given by ( 6 - 25 ) where and are diameters of colliding droplets, an d is shear rate. The collision efficiency is obtained from a function originally introduced by Coulaloglou [ 97 ] as ( 6 - 26 ) where is the drainage time, and is the contact time. The drainage time is the time required for the intervening film between the droplets to thin to a critical thickness. The contact time is the interaction time of two colliding droplets. A model presented by Chesters [ 120 ] is used to obtain the drainage time for partially mobile interfaces: ( 6 - 27 ) where is interaction force, is equivalent radius obtained as , and and are initial and critical film thicknesses, respectively. The interaction force is given by [ 120 ] ( 6 - 28 ) 156 The initial film thickness is considered as [ 114 ] . The critical (final) film thickness is adopted from Chesters [ 120 ] and presented by Venneker et al. [ 121 ] as ( 6 - 29 ) where is Hamaker, which ranges between 10 - 20 J and 10 - 1 9 J. Hamaker constant is set to 10 - 20 J in this study. A model presented by Luo [ 122 ] is employed to calculate the contact time as follows: ( 6 - 30 ) where . The parameter is the added mass coefficient and is set to 0.65 in this study. is normally taken to be a constant between 0.5 and 0.8 [ 123 ] . Using Eqs. ( 6 - 25 ) to ( 6 - 30 ) , the coalescence kernel is calculated. The droplet coalescence on the membrane surface is neglected for the current study. 6.2.3. Permeate flux of a membrane Permeate flux of a membrane may decrease by deposition of droplets and blockage of [ 64 ] permeate flux can be obtained as ( 6 - 31 ) 157 where is the transmembrane pressure difference (i.e. ), is the overall resistance , and is the surface blockage ratio defined as the surface area with blocked pores to the total surface area of the porous medium. The overall resistance term can be given as ( 6 - 32 ) where is the membrane resistance, and is the cake resistance. Membrane resistance ( ) is normally available through membrane manufacturer or can be measured based on the initial flux of a clean membrane. Cake resistanc e ( ) is due to the presence of oil droplets in a near - membrane cell (see Figure 6 - 2 ). Although these droplets are not deposited on the membrane, and they do not enter pores, they can decrease the flow rate of continuous phase reaching to the membrane surface. For solid parti cles in dilute suspensions, the cake resistance, is given by ( 6 - 33 ) where is specific dispersed phase resistance estimated by Carmen - Konezy equation [ 240 ] as . i s the void fraction, and is the particle surface area per unit volume of particles. For a dispersed phase layer of incompressible spheres of diameter , these parameters are given as , , and [ 63 , 241 ] . The equivalent thickness of dispersed phase layer, in Eq. ( 6 - 33 ) , is calculated based on the volume fraction of the dispersed phase in the near - membrane cell as ( 6 - 34 ) 158 It should be noted that oil droplets in an oil - in - water emulsion behave differently than solid particles in a suspension. However, due to low capillary numbers (i.e. insignificant droplet deformation in bulk flow), and low volume fractions of oil in such e mulsions, it can be assumed that the model above provides a reasonable estimation of cake resistance in early stages of membrane fouling. However, it should be mentioned that the numerical results show that this term is small compared to the membrane resis tance in early stages of fouling. 6.2.3.1. A veraging o f the droplets on the membrane surface A similar approach of averaging in bulk flow is employed for multiphase averaging on the membrane surface. In early stages of membrane fouling, a layer of droplets deposit on the membrane surface. Figure 6 - 3 illustrates a schematic of deposited droplets on the membrane surface. 159 The surface void fraction is defined as the ratio of the surface occupied by droplets of phase on the membrane surface to the total membrane surface as ( 6 - 35 ) where is a control surface, takes a value of 1.0 where phase is present, and is zero everywhere else. When a group of droplets of diameter deposit on a surface, Eq. ( 6 - 35 ) can be approximated as ( 6 - 36 ) where is the number of droplets on the surface, and is the diameter of the wetted area on the surface by the droplet. Figure 6 - 4 shows the schematic of a droplet on the membrane surface. Assuming the spherical shape of droplet on the surface, the diameter of the deposited droplet ( ) and the diameter of droplet before deposition ( ) is related by ( 6 - 37 ) is defined as , where is the contact angle of droplet on the surface. Receding and advancing contact angles are assumed equal. 160 Since , the diameter of wetted area and the diameter of the droplet befo re deposition are related by ( 6 - 38 ) Therefore, Eq. ( 6 - 36 ) is rewritten as ( 6 - 39 ) In order to calculate using Eq. ( 6 - 39 ) , the number of droplets ( ) must be known. Hence, it is required to obtain the mass of deposited droplets on the membrane surface. 161 By fluid averaging method on the surface, an averaged repre sentation of droplets deposited on a surface is provided. Therefore, the computationally intense and practically impossible problem of Figure 6 - 5 (a) can be simplifi ed to the problem of Figure 6 - 5 (b) with a reduced computational cost. In order to obtain the mass of deposited droplets, the equivalency of droplets and film repres entations (see Figure 6 - 5 ) is needed. The mass of droplets deposited on the membrane surface ( ) is given by ( 6 - 40 ) 162 Figure 6 - 5 (b) illustrates that the deposited droplets are averaged as a film layer on the membrane surface, i.e. , where is the mass of film layer. Therefore, Eq. ( 6 - 40 ) can be rewritten as ( 6 - 41 ) where is film thickness. Therefore, in order to find and calculate using Eq. ( 6 - 39 ) , the film thickness ( ) is needed in the averaged representation. In early stages of membrane fouling, is smaller than size of a single droplet, and therefore, the film thickness is muc h smaller than the characteristic length of the channel, i.e. . Additionally, the film thickness is much smaller compared to the curvature of the membrane surface (if any). Hence, the flowing assumptions are made for the film: - The film is a spatially - averaged continuum of droplets that are deposited on the membrane surface. - The film is approximated as a two - dimensional thin film of thickness , since its thickness is small compared to characteristic length of the channel . - The velocity field in the film is small. However, to model film flow on the membrane surface, a quadratic velocity profile is assumed for the film , which is adapted from lubrication theory approximation [ 238 ] and is based on experiments by Mudawar [ 242 ] . In order to calculate the film thickness, mass and momentum balance equations are required. Figure 6 - 6 illustrates a schematic of forces and mass fluxes on a control volume of a film formed on the membrane surface. 163 The finite volume mass conservation of liquid film at the wall is given by [ 234 ] ( 6 - 42 ) where is the area of computational cell on the wall, is mean film velocity , is the surface normal, is the cell size, is the number of cell sides, is film thickness, and is source term (due to transport of droplets to the film). The finite volume momentum conservation of liquid film is obtained as ( 6 - 43 ) 164 where is the pressure, is gravity acceleration tangential to the film, is the number of cell edges, and is its corresponding shear stress. is velocity vector of the secondary phase. The differential form of Equations ( 6 - 42 ) and ( 6 - 43 ) is also obtained. The mass conservation is given by ( 6 - 44 ) where is surface gradient operator defined as ( 6 - 45 ) In Eq. ( 6 - 44 ) , is mass source per unit wall area . The momentum conservation for a thin liquid film is obtained as ( 6 - 46 ) The first and second terms in the LHS are transient and convection terms, respectively. The tensor indicates the differential advection term. Term (I) is due to the pressure gradient of continuous phase; term (II) is the effect of gravity component normal to the wall; term (III) is due to gravity component tangential to the film; term (IV) is due to th e surface curvature; term (V) and (VI) are shear stress at the interface of continuous phase and film, and wall shear stress, respectively. These terms are given by 165 ( 6 - 47 ) ( 6 - 48 ) where , and . The aforementioned laminar quadratic velocity profile is used to obtain the wall shear stress. Term (VII) in Eq. ( 6 - 46 ) is the source term, which is due to droplet transport from bulk flow to the membrane surface. The mass source term of Eq. ( 6 - 44 ) is given by ( 6 - 49 ) where is the volume fraction of dispersed phase. is the dispersed p hase velocity normal to the wall surface , and is w all surface area . Accordingly, a m omentum source term for wall film can be given as ( 6 - 50 ) where is velocity vector of the secondary phase. Using the source terms obtained from Eqs. ( 6 - 49 ) and ( 6 - 50 ) , the balance equations for film thickness, i.e. Eqs. ( 6 - 44 ) and ( 6 - 46 ) , can be solved to find . The algorithm to calculate the membrane permeate flux is summarized Table 6 - 1 . 166 For number of timesteps Calculate the permeate flux of the membrane surface using Eq. ( 6 - 31 ) Calculate the source terms of Eqs. ( 6 - 5 ) and ( 6 - 9 ) Solve mass and momentum equations of the mixture model, Eqs. ( 6 - 5 ) and ( 6 - 9 ) Calculate Sauter mean diameter using PBM, Eqs. ( 2 - 1 ) and ( 6 - 23 ) For all computational cells located next to the membrane surface Calculate mass and momentum source terms of the mixture using Eqs. ( 6 - 49 ) and ( 6 - 50 ) Calculate using Eqs. ( 6 - 44 ) and ( 6 - 46 ) Calculate the number of deposited droplets using Eq. ( 6 - 41 ) Calculate the blockage ratio using Eq. ( 6 - 39 ) Update the blockage ratio End For End For Parametric study A parametric study is carried out to investigate the effect of various parameters on membrane fouling, permeate flux, and droplet s ize distribution on the membrane . A base case is introduced that consists of a rectangular channel of dimensions of 8 0mm× 3 0mm×3mm. A flat - sheet membrane with initial resistance of is embedded at the lower surface of the channel. The average crossflow velocity is 0.1 m/s at the inlet , while the volume fraction of droplets is 0.1% (1000 ppm) . The droplet size distribution (DSD) at the inlet is shown in Figure 6 - 7 . The Sauter mean diameter of the DSD is . The contact angle and the interfacial tension of the membrane - droplet combination is 130°, and 0.0 05 N/m, respectively. 167 The fluid flow is a multiphase flow of water as continuous phase, and oil as dispersed phase. The densities of water and oil are 998.2 kg/m 3 and 8 00 kg/m 3 , respectively. The gravity acceleration is 9.81 m/s 2 . The viscosities of water and oil are 0.001 Pa.s and 0.003 Pa.s, respectively. The transmembrane pressure ( ) is zero, while the outlet pressure is adjusted to provide the required transmem brane pressure difference ( ), which is 20 kPa for the base case. The lower surface is a permeable wall (membrane surface), while other boundaries are impermeable. Figure 6 - 8 shows a schematic of computational domain. A three - dimensional hexahedral mesh of 74 0,000 cells is generated using ANSYS Meshing module. The Navier - Stokes equation for laminar flow is solved using a pressure - 168 based algorithm by ANSYS Fluent 1 8.1. Mixture model with implicit body force implementation is empl oyed to perform the multiphase flow simulation. Eulerian wall film model is used to model the transport of droplets from bulk flow to the membrane surface. In - house User - defined functions (UDFs) and User - defined memory variables (UDMs) are used to enforce the permeate flux and membrane blocking process. First Order Explicit discretization scheme is used for transient term, and Second Order Upwind scheme is utilized for both continuity and momentum equations of Eulerian wall film. The pressure - velocity coupl ing is carried out using PISO scheme with Green - Gauss Node Based scheme for gradient terms, Pressure Staggering Option (PRESTO) for pressure, and Quadratic Upstream Interpolation for Convective Kinetics (QUICK) scheme for momentum and volume fraction equations. A Bounded Second Order Implicit formulation is used for transient time - marching of Navier - Stokes equation, providing an improved accuracy and better stability compared to the first order formulatio n. To achieve accurate transient solutions, the scaled residuals of the simulations are set to be less than for continuity, momentum, and volume fraction equations. The simulations are carried out for 90 minutes of flow. T he effect of variation of fo llowing parameters are investigated: Transmembrane pressure difference ( ), resistance of membrane ( ), relative density ( ), relative viscosity ( ), contact angle ( ), Feed oil volume fraction ( ), feed velocity, and droplet s ize distribution. Therefore, a set of simulations for various cases are carried out and the results are compared with the base case introduced above. The cases investigated in parametric study are tabulated in Table 6 - 2 . 169 Parameter Case 1 (Base) Case 2 Case 3 Operational conditions Transmembrane pressure (Pa) 20,000 10,000 30,000 Membrane resistance (m - 1 ) 5.0×10 11 2.5×10 11 1.0×10 12 Feed volume fraction (%) 0.1 0.01 1.0 Average Feed velocity (m/s) 0.1 0.05 0.15 Fluid Properties Contact angle (deg.) 130 110 150 Droplet size (m) DSD ( Figure 6 - 7 ) 0.5 × DSD 1.5 × DSD Relative density 0.8 0.7 0.9 Relative viscosity 3.0 1.0 5.0 Results and discussion The effect of different parameters on permeate flux is investigated. Four constants are defined to characterize the permeate flux of early stages of membrane fouling: (1) Steady - state permeate flux, , (2) time of steady - state flux, (3) median steady - state flux, (4) time of median - steady - state flux, . Median steady - state flux is defined as ( 6 - 51 ) where refers to initial permeate flux at . Figure 6 - 9 (a) shows that a higher transmembrane pressure difference leads to a higher permeate flux. However, the flux drop (i.e. ratio of steady - state flux to initi al flux , ) is larger for that case, since it creates more suction flow on membrane surface, resulting in more droplet transport to the surface . 170 (a) (b) (c) (d) Figure 6 - 9 (b) shows that a higher membrane resistance is associated with less permeate flux; however, the flux drop (i.e. the ratio of steady - state flux to initial flux) is higher for smaller membrane resistance, when all other parameters are constant. Figure 6 - 9 (c) shows that a higher volume fraction of droplets at the inlet increases the rate of fouling. Figure 6 - 9 (d) demonstrates that a higher inlet velocity does not have significant impact on 171 the steady - state permeate flux; however, for higher inlet velocity (in laminar regime), the steady - state condition occurs in a shorter period of time. Figure 6 - 10 (a) s hows that a lower contact angle increases the fouling and decreases the stea dy - state permeate flux. A droplet with lower contact angle covers more area on membrane surface and consequently blocks more pores. Figure 6 - 10 (b) demonstrates that compared to smaller droplets, larger droplets lead to smaller fouling rate and an increase in permeate flux. Buoyancy force on t hese droplets pushes them away from the membrane and they are less likely to reach to the surface and block the pores. Figure 6 - 10 (c) shows the effect of relative density ( ) on the permeate flux. The higher relative density (i.e. the density difference between oil and water is smaller) mitigate the effect of buoyancy force, and therefore, droplets are more likely to follow the streamlines, move toward gravity and approach the membrane surface. Hence, more fouling and less permeate flux is expected. Figure 6 - 10 (d) illustrates the effect of relative viscosity ( ) on the fouling. It is shown that a higher viscosity of oil droplets results in an increase in membrane fouling and dec rease in permeate flux. A higher viscosity lowers the droplet coalescence by decreasing the collision frequency (see Eqs. ( 6 - 26 ) and ( 6 - 27 ) ). Therefore, the mean droplet size will remain smaller for the case of higher viscosity, resulting in more fouling on the membrane (see Figure 6 - 10 (b) for the effect of droplet size distribution on the membrane fouling). 172 (a) (b) (c) (d) Table 6 - 3 presents the mid - steady - state flux, its occurrence time, steady - state flux, and its occurrence time for the cases discussed above. 173 Parameter (s) Base case 144 133 123 477 2390 Transmembrane pressure Lower 72 69 67 1356 5296 Higher 215 195 176 314 1973 Contact angle Lower 144 128 112 468 3542 Higher 144 139 134 484 2599 Droplet size Smaller 144 120 97 226 2013 Larger 144 139 135 1161 4715 Relative density Lower 144 135 127 768 3898 Higher 144 130 117 281 1669 Membrane resistance Lower 287 256 226 244 1651 Higher 72 69 67 1367 5300 Relative viscosity Lower 144 139 134 539 2064 Higher 144 127 110 408 3000 Feed volume fraction Lower 144 136 128 930 3359 Higher 144 100 56 60 313 Feed velocity Lower 144 132 121 1248 3814 Higher 144 133 122 290 1727 Figure 6 - 11 shows mass per area of the deposited droplets on the membrane surface for the base case in various timeframes. It is seen that the proposed model can predict the deposition of droplets on the membrane surface, as well as its advance over time. 174 t = 2 5 0s t = 500s t = 1 0 00s t = 2 0 00s Mass per area (kg/m 2 ) Figure 6 - 11 also indicates that the droplet deposition mainly occurs in the center of the membrane, and it has a lower rate near inlet of the crossflow channel. In addition, its advancing is lower near t he walls of the crossflow channel. Figure 6 - 12 shows the permeate velocity at the membrane surface in various timeframes. It is seen that the permeate flux is decreas ing over time, and its reduction is in agreement with the droplet mass deposited on the membrane surface (see Figure 6 - 11 ). Flow direction 175 t = 2 5 0s t = 50 0s t = 1 0 00s t = 2 0 00s Permeate velocity (m/s) Figure 6 - 13 illustrates the Sauter mean diameter near the membrane surface at steady state permeate flux. It is seen that the has its maximum value near the center of membrane surface. Flow direction 176 Sauter Mean Diameter (m/s) Figure 6 - 14 shows the volume fraction of oil near the membrane surface at steady - state permeate flux. It compares the base case, the case with a lower inlet velocity, and the case with the higher inlet v elocity. The volume fraction contour is plotted for the symmetry plane of the crossflow channel . Since the impact of permeate flux on the volume fraction is insignificant far from the membrane surface, Figure 6 - 14 illustrates only a portion of the channel height. The vertical axis shows only 10% of the channel height (0.3mm), and the horizontal axis shows the entire channel length ( 8 0mm). The lower edge represents the membrane surface, where the volume fraction is expected to be higher. Flow direction 177 Volume fraction Figure 6 - 14 illustrates that a thin layer is form ed near the membrane surface, which has a relatively high volume fraction of oil, and at some spots the volume fraction is about 10 times of its value at the inlet of crossflow channel. It is seen that in a CFF with a higher inlet velocity the thickness of this thin layer is smaller compared to the one with a lower inlet velocity. A higher inlet velocity tends to push the oil layer toward the membrane surface. Figure 6 - 15 shows a histogram of Sauter mean diameter ( ) at the outlet of the membrane channel for cases with different DSDs at the inlet (see Table 6 - 2 for details). I t can Flow direction 178 be seen that of only a small percentage of each DSD has changed. Figure 6 - 15 indicates that there is a small fraction of droplets with of higher than the initial value at the outlet. That specifies that as expected, coalescence has only affected a small fraction of droplets due to a very low volume fraction of droplets at the inlet. However, this could be different for droplets near the surface, where the volume fraction of droplets is higher (see Figure 6 - 14 ). Figure 6 - 16 shows Sauter mean diameter of the droplets of the base case near the membrane for various timeframes. It can be seen that has significantly changed compared to the initial value of for this case. Additionally, Figure 6 - 16 indicates that over time, a higher number of large droplets tend to be formed near the membrane. 179 Summary and conclusions In the present chapter , a new model for membrane fouling is proposed. The new model for the first time incorporates a wall film model, the mixture model, and an inhomogeneous population balance model that accounts for droplet coalescence in bulk flow. The results demonstrated that the proposed model is capable of predicting the varying permeate flux of the crossflow membrane system and provides a detailed information about the deposited droplets on the membrane surface. Since the model is fully coupled with computational fluid dynamics, it can be utilized to provid e system design and operation information. A parametric study was carried out to exhibit the ability of the proposed model to be used for several operation conditions, including various crossflow rates, droplet size distributions, various droplet densities and viscosities, and various transmembrane pressure differences. Experimental validations of the proposed model are proceeding and will be published 180 shortly. The current study lacks a model for droplet coalescence and film formation on the membrane surfac e. Additionally, a future study to model the droplet separation from the membrane surface and droplet entrance into the pores are recommended to expand the features of the proposed model in current study. Appendix 6.6.1. Droplet breakage A droplet tends to deform and eventually breakup due to both crossflow and permeate flow in a crossflow filtration membrane. The tendency of a droplet to deform can be estimated by capillary number, , which is the ratio of viscous and interfacial tension forces. The capillary number due to a crossflow is given by [ 55 ] ( 6 - 52 ) The capillary number due to the permeate flux is obtained as [ 55 ] ( 6 - 53 ) where is the permeate flux. Tummons et al. discussed that even though the viscous and body forces are able to change droplet shape, they are incapable of breaking droplets up, especially far from the channel walls, where shear rate is smaller compared to near - w all regions. Hence, the droplet breakage is neglected in regions far from the wall. Another possibility is droplet breakage on the pore entry [ 55 ] , in which a part of droplet permeates 181 through the pore and causes droplet breakup. A critical capillary number as a breakup criterion was introduced by Darvishzadeh et al. [ 60 ] ( 6 - 54 ) where . For the cases used in this study, both and are much smaller than , and therefore, the droplet breakage on the surface is neglected as well. Nomenclature Hamaker constant (J) Capillary number Droplet diameter (m) Force (N) Coalescence frequency Breakage kernel Film thickness Permeate Flux (LMH) Mass (kg) Number density Pressure (Pa) Membrane resistance (m - 1 ) 182 Droplet radius (m) Reynolds Number time (s) velocity vector (m/s) Volume (m 3 ) Volume of a droplet (m 3 ) Greek symbols Volume fraction Probability density function Contact angle Shear rate (s - 1 ) Coalescence Kernel Collision efficiency Viscosity (Pa.s) Membrane blockage Density (kg/m 3 ) Surface tension (N/m) Stress tensor (Pa) 183 Subscripts Continuous phase Drift Dispersed phase arbitrary secondary phase Mixture Transmembrane Abbreviations CFF Crossflow filtration CFV Crossflow velocity DQMOM Direct quadrature method of moments DSD Droplet size distribution TMP Transmembrane pressure PBE Population balance equation PBM Population balance model PRESTO Pressure s taggering o ption QUICK Quadratic Upstream interpolation for convective kinetics UDF User - defined function 184 UDM User - defined memory variable 185 Summary and Conclusions The multiphase flow theory was applied to two industrial multiphase flows. The first part tries to investiga te and characterize the Treatment Shaft system to identify its performance for separation of suspended solids from water. The second part studies various crossflow filtration membranes and aims to understand their flow patterns and improve the separation e fficiency. For this purpose , t he multiphase flow theory for interacting droplets and particles were reviewed and population balance model was explained as the base for interaction of dispersed phase in current study. Moreover , appropriate coalescence and breakage models were provided for interactions of liquid droplets, and interactions of solid particles in a continuous phase. For each part, summary of findings, conclusions, limitations of introduced models, and suggestions for future works are discussed below. The Treatment Shaft system as a CSO detention system 7.1.1. Summary and conclusions The Treatment Shaft was initially designed to be employed as a CSO detention system and needed to provide a short contact time during a 1 - year, 1 - hour rainstorm event. The main feature of the shaft is the baffle wall, which creates a slow directed downward and upward velocity fields to promote settling of suspended solids. R esults show that the length of baffle wall affects its efficiency. A shorter baffle wall, providing a larger gap at the bottom of the shaft can reduce the velocity at the bottom of the shaft, and increases the efficiency of settling. However, the change is less than 1% in the efficiency. On the other hand, a longer 186 baffle wall, not only increases the const ruction costs, but also decrease the efficiency due to inducing a higher velocity at the bottom of the shaft. Therefore, considering the construction costs, and efficiency, the gap length of 13.22m is recommended. Further study on the location of the shaf t (±10% from the center) showed an insignificant change in the efficiency . The shape of the inlet and outlet channels were also studied. Although the original shape of the Treatment Shaft provides a slower flow at its inlet and outlet, the environmental fo otprint is large, and its construction is more expensive. Therefore, using CFD simulations a new design alternative was proposed in which the channels are straight (instead of flared channels of original design). The results showed that the change in the o verall efficiency of the shaft during a 10 - year, 1 - hour rainstorm is insignificant, and the proposed design can provide the same performance with reduced environmental trace. Further study was conducted on the contact time of the Treatment Shaft. A residen ce time distribution analysis demonstrated that the retention time is 745.5s, while the mean residence time is 890.7s. The settling efficiency of the shaft was studied with and without the use of flocculating agents. Without flocculant, at least 50% of par sedimentation of smaller particles is not substantial. An analysis on a sample particle size distribution of CSOs, showed that this results in settling of only 20.8% of particles in the sample. The fact that most of the p the importance of using flocculants in such systems. Hence, a population balance model coupled with CFD was utilized to study the impact of a polymer flocculant on the separation efficiency of the Tre atment Shaft. It was shown by use a flocculating agent, the separation efficiency of particles smaller than 1 4 9.3 %, compared to 20.8% (without using flocculant). The simulations illustrat e that at 34 187 minutes of rainstorm event, flocculant is vastly mixed with the fluid, and increases the rate of coalescence and consequently, the Sauter mean diameter of particles in the Treatment Shaft. Additionally, it was illustrated the Sauter mean diameter in the shaft inc reases over time, which eventually enhances the settling efficiency of the system. 7.1.2. Limitations The current study provides information for design and operation purpose for the high flow rate ), the amount of available computational resources, and characteristics of used models entailed some limitations on the accuracy. details are omitted to grant better convergence for simulations without remarkable accuracy loss. In addition. t his study focused on 2 - phase simulations (i.e. sewage particles, water) , where the Treatment Shaft was assumed to be initially filled with water. This delivers a significant computational benefit, compared to a three - phase model, where interact ions of air with the 2 - phase flow is considered. Although 3 - phase flow (i.e. air, sewage particles, and water) lead s to capture the free surface , due to high ratio between the kinetic viscosity of air and other materials, a very small time - step is required , which can increase the simulation time for up to 50 times. Additionally, although the turbulence model used in present study is effective for most industrial applications, it performs poorly for complex flows and where strong streamline curvatures are pr esent. The choice of current model was based on available computational resources. 188 The impact of flocculating agent on the particle coalescence were evaluated based on characteristic data and constants gathered from other studies. These constants could be different in various situations and might need calibration for a given problem. 7.1.3. Future work The performance of the Treatment Shaft technology was demonstrated , but there are many opportunities for extending the current stud y: Design alternatives : The Treatment Shaft system is based on a simple design due to its inexpensive construction costs. Although some design enhancements were introduced in current study, more investigation is needed to improve the efficiency. Some suggestions are: Including all design details, adding inclined baffles attached on main baffle, adding circular baffles on the main cylinder, etc. Modeling free surface : by modeling free surface during the filling process of the Treatment Shaft, a better understanding of settling du ring first 20 minutes of rainstorm event will be provided. It can also help to identify the spots with higher erosion rate, which is valuable from aspect of construction and maintenance of the shaft. Turbulence model : A turbulence model with higher accurac y can provide a better picture of vortexes, flow curvatures, and consequently a more accurate settling efficiency. The use of Reynolds Stress Model by implementing anisotropic Reynolds stress tensor is suggested for future studies. 189 Improving coalescence and breakage models : coalescence and breakage models can be improved by calibrating their constants using experimental studies. A combination of experimental measurements and CFD simulations is recommended for future study of the Treatment Shaft as a CSO d etention system. Crossflow filtration membranes 7.2.1. Summary and conclusions Crossflow filtration membranes of various configurations were studied. In order to present a better understanding of flow field in CFF membranes, various configurations of tubular CFF membranes were investigated. Existing analytical solutions were reviewed, and a lack in the literature for accurate solutions for membranes with relatively high permeate flux was noted . In order to resolve this problem, a set of new analytical solutions fo r four different configurations of tubular systems with permeate flux were introduced in this study. The provided approximate solutions are applicable for different filtration methods of laminar flows and are able to predict the pressure drop across the cr ossflow channel. C ontrary to other studies, the proposed analytical solutions could accurately obtain the flow field variables even with high permeate flux, where radial flux is as high as 50% of inlet flux. To improve the separation efficiency of tubular CFF membrane system s, specifically for oil - in - water emulsions, rotating membrane and charging of membrane surface were studied using CFD simulations . Two different configurations were analyzed: (1) a regular tubular membrane channel, and (2) a tubular memb rane channel with a circular rod at the centerline. Results demonstrated that for nonrotating systems, the second configuration 190 slightly improves the efficiency. By rotating membrane systems, a noticeable increase in separation efficiency is seen. For drop membrane has separation efficiency of 94%, while a rotating membrane (500rpm) can reach 100% of separation efficiency. By applying charge on membranes, the increase is more significant, and it can decrea se the cut - membrane configuration. These studies were not pursued further but could lead to devices with high separation performance. Membrane fouling is a major issue in membrane filtration of oil - in - wat er emulsions . For the first time, a model was introduced which employs the wall film model, and the multiphase mixture model , fully coupled with population balance model to include droplet coalescence in bulk flow. The proposed model can identify the perme ate flux of a CFF membrane system, blockage ratio, blockage advancing over time, and distribution of Sauter mean diameter across the membrane, and in bulk flow. A parametric study was carried out to demonstrate features of the proposed model . The effect of various operational conditions including transmembrane pressure, membrane resistance, feed volume fraction, and crossflow velocity were evaluated. Moreover, the impact of oil properties including contact angle, droplet size distribution, density, and visc osity were investigated. Simulations show that the advancing front of the blockage form over time until the permeate flux reaches its steady - state value. In addition, it is seen that the model can effectively predict the high concentration region near the membrane surface, as well as it changes over time, and for various crossflow velocities. 191 7.2.2. Limitations of membrane fouling model The current study provides novel analytical solutions, and the first CFD model for membrane fouling model due to oil droplets. Ho wever, several assumptions are made that limit the use of proposed solutions and models and highlights the need for future study in this field. Although the proposed analytical solutions were derived for CFF membrane systems, the impact of foulant, or solu te is neglected due to its low volume fraction in the flow. Therefore, the permeate flux remains unchanged in proposed solutions, which limits the use of these analytical solutions for transient systems. In addition , the approximate solutions are limited t o laminar flows and are not suitable for turbulence flows in tubular CFF membrane systems. Moreover, the velocity profile at the inlet of these tubular systems are assumed as fully - developed, while this assumption may not be valid for some industrial cases . Membrane fouling is a very complicated phenomenon. In the fouling model proposed in this study, droplet coalescence on the membrane surface was not considered. Droplets on the surface were represented as an averaged film. Therefore, drag and lift forces on droplets, and droplet detachment from the surface is not particularly modeled. The proposed fouling model represents the early stages of membrane fouling wher e only a single layer of droplets are deposited. 7.2.3. Future work While this thesis provided new solutions and models for CFF membrane systems, numerous opportunities are available for further research in this field. 192 Analytical solutions for mass transfer : Curr ent study was conducted to obtain analytical solutions for flow field of a single - phase flow in tubular CFF membrane systems. These solutions can be used for to derive mass transfer approximate solutions in tubular CFF systems. Such solutions could be bene ficial for understanding the mass transfer phenomena during filtration of solutes from a continuous phase using a CFF membrane. Impact of charged membranes on oil separation : this study demonstrated the benefits of charged membranes for oil - in - water separ ation. The same approach is suggested to use for various membrane configurations, including flat sheet membranes, dead - end setups, etc. The model can be highly improved by including the effect of surfactant concentration in the mode l , which directly impact s the charge on oil droplets. Droplet coalescence on membrane surface : In present study , droplet coalescence on the membrane surface is not modeled. Droplet coalescence can create larger droplets, and eventually form spots with a combination of films and d roplets. Experimental efforts in combination with physical models are suggested to develop a model for coalescence on the surface. Droplet re - entrainment from the membrane surface : due to changes in force balance of a droplet on the membrane surface, some droplets may detach from the surface and get washed away by retentate flow. Including these features in proposed fouling model will enhance its accuracy and is suggested for future works. 193 Droplet entrance in the pores : in current study, entrance of droplet s into the pores are not implemented. This phenomenon occurs due to high transmembrane pressure difference ( ) and can completely block membrane pores. Adding this feature to the proposed model is suggested for future studies. Mul tilayer membrane fouling : The present study only models the early stages of membrane fouling, in which a single layer of droplets is deposited on the membrane surface. However, in CFF membrane systems, permeate flux continues to decline by addition of more layers of droplets on the first deposited layer. This can significantly improve the proposed fouling model and is recommended for future research. 194 REFERENCES 195 REFERE NCES [ 1] U. EPA, "Report to Congress: Impacts and Control of CSOs and SSOs," Office of Water, 2004. [2] [3] APHA, AWWA, and WEF, "Standard methods for the examination of water and wastewa ter," in American Public Health Association (APHA): Washington, DC, USA , ed, 2005. [4] I. Oke, K. Oladepo, A. Olajumoke, and E. Ajayi, "Settlement properties of solids in a domestic - institutional wastewater," Journal of Applied Sciences Research, vol. 2, p p. 385 - 90, 2006. [5] P. Piro, M. Carbone, N. Penna, and J. Marsalek, "Characterization of the settling process for wastewater from a combined sewer system," Water research, vol. 45, pp. 6615 - 6624, 2011. [6] A. M. Goula, M. Kostoglou, T. D. Karapantsios, an d A. I. Zouboulis, "The effect of influent temperature variations in a sedimentation tank for potable water treatment A computational fluid dynamics study," Water research, vol. 42, pp. 3405 - 3414, 2008. [7] S. J. Wright, S. Ghalib, and A. Eloubaidy, "Treat ment Shaft for Combined Sewer Overflow Detention," Water Environment Research, vol. 82, pp. 434 - 439, 2010. [8] D. M. Hudson, "Inflatable gates improve storage, reduce CSO," Public Works, vol. 129, 1998. [9] M. Abu - Orf, G. Tchobanoglous, H. D. Stensel, R. T suchihashi, F. Burton, G. Bowden , et al. , Wastewater engineering: treatment and Resource recovery : McGraw Hill Education, 2014. [10] S. J. Wright, J. Lewis, and J. G. Vasconcelos, "Mechanisms for stormwater surges in vertical shafts," Contemporary modeling of urban water systems, pp. 109 - 132, 2007. [11] J. Dettmar and P. Staufer, "Modelling of flushing waves for optimising cleaning operations," Water science and technology, vol. 52, pp. 233 - 240, 2005. [12] J. Suarez and J. Puertas, "Determination of COD, BO D, and suspended solids loads during combined sewer overflow (CSO) events in some combined catchments in Spain," Ecological Engineering, vol. 24, pp. 199 - 217, 2005. [13] S. A. Ghalib, "Wastewater treatment system and method," ed: Google Patents, 2003. 196 [14] S. A. Ghalib, "Method for treating wastewater," ed: Google Patents, 2008. [15] J. Veil and C. Clark, "Produced water volume estimates and management practices," SPE Production & Operations, vol. 26, pp. 234 - 239, 2011. [16] J. A. Veil, "Produced water mana gement options and technologies," in Produced Water , ed: Springer, 2011, pp. 537 - 571. [17] - Razi, A. Pendashteh, L. C. Abdullah, D. R. A. Biak, S. S. Madaeni, and Z. Z. Abidin, "Review of technologies for oil and gas produced water treatment," J ournal of Hazardous Materials, vol. 170, pp. 530 - 551, 2009. [18] J. M. Neff, Bioaccumulation in marine organisms: effect of contaminants from oil well produced water: Elsevier, 2002. [19] R. R. Reynolds, "Produced Water and Associated Issues," Oklahoma Geo logical Survey, 2003. [20] J. A. Veil, M. G. Puder, D. Elcock, and R. J. Redweik Jr, "A white paper describing produced water from production of crude oil, natural gas, and coal bed methane," Argonne National Laboratory, Technical Report, 2004. [21] J. Fil lo, S. Koraido, and J. Evans, "Sources, characteristics, and management of produced waters from natural gas production and storage operations," in Produced Water , ed: Springer, 1992, pp. 151 - 161. [22] B. Hansen and S. Davies, "Review of potential technologies for the removal of dissolved components from produced water," Chemical engineering research & design, vol. 72, pp. 176 - 188, 1994. [23] M. Stephenson, "A survey of produced water studies," in Produc ed water , ed: Springer, 1992, pp. 1 - 11. [24] P. Ekins, R. Vanner, and J. Firebrace, "Zero emissions of oil in water from offshore oil and gas installations: economic and environmental implications," Journal of Cleaner Production, vol. 15, pp. 1302 - 1315, 20 07. [25] J. Veil, "Research to Improve Water - use Efficiency and Conservation: Technologies and Practice," ed, 2007. [26] J. D. Arthur, B. G. Langhus, and C. Patel, "Technical summary of oil & gas produced water treatment technologies," All Consulting, LLC, Tulsa, OK, 2005. [27] A. Cambiella, J. Benito, C. Pazos, and J. Coca, "Centrifugal separation efficiency in the treatment of waste emulsified oils," Chemical Engineering Research and Design, vol. 84, pp. 69 - 76, 2006. 197 [28] A. Al - Shamrani, A. James, and H. Xiao, "Separation of oil from water by dissolved air flotation," Colloids and Surfaces A: Physicochemical and Engineering Aspects, vol. 209, pp. 15 - 26, 2002. [29] K. Bensadok, M. Belkacem, and G. Nezzal, "Treatment of cutting oil/water emulsion by coupling coagulation and dissolved air flotation," Desalination, vol. 206, pp. 440 - 448, 2007. [30] K. A. Hashmi, H. A. Hamza, and J. C. Wilson, "CANMET hydrocyclone: an emerging alternative for the treatment of oily waste streams," Minerals engineering, vol. 17, p p. 643 - 649, 2004. [31] J. Hargreaves and R. Silvester, "Computational fluid dynamics applied to the analysis of deoiling hydrocyclone performance," Chemical engineering research & design, vol. 68, pp. 365 - 383, 1990. [32] D. Colman and M. Thew, "Correlation of separation results from light dispersion hydrocyclones," Chemical Engineering Research and Design, vol. 61, pp. 233 - 240, 1983. [33] M. Thew, "Hydrocyclone redesign for liquid - liquid separation," Chemical engineer, pp. 17 - 23, 1986. [34] D. Wolbert, B. F Hydrocyclones using trajectory analysis," AIChE Journal, vol. 41, pp. 1395 - 1402, 1995. [35] I. Cumming, R. Holdich, and I. Smith, "The rejection of oil by microfiltration of a stabil ised kerosene/water emulsion," Journal of Membrane science, vol. 169, pp. 147 - 155, 2000. [36] T. Mohammadi, M. Kazemimoghadam, and M. Saadabadi, "Modeling of membrane fouling and flux decline in reverse osmosis during separation of oil in water emulsions," Desalination, vol. 157, pp. 369 - 375, 2003. [37] J. Zhong, X. Sun, and C. Wang, "Treatment of oily wastewater produced from refinery processes using flocculation and ceramic membrane filtration," Separation and Purification Technology, vol. 32, pp. 93 - 98, 2003. [38] R. O. Dunn Jr, J. F. Scamehorn, and S. D. Christian, "Use of micellar - enhanced ultrafiltration to remove dissolved organics from aqueous streams," Separation science and technology, vol. 20, pp. 257 - 284, 1985. [39] A. Chen, J. Flynn, R. Cook, an d A. Casaday, "Removal of oil, grease, and suspended solids from produced water with ceramic crossflow microfiltration," SPE production engineering, vol. 6, pp. 131 - 136, 1991. 198 [40] J. Benito, G. Ríos, E. Ortea, E. Fernández, A. Cambiella, C. Pazos , et al. , "Design and construction of a modular pilot plant for the treatment of oil - containing wastewaters," Desalination, vol. 147, pp. 5 - 10, 2002. [41] H. Ohya, J. Kim, A. Chinen, M. Aihara, S. Semenova, Y. Negishi , et al. , "Effects of pore size on separation m echanisms of microfiltration of oily water, using porous glass tubular membrane," Journal of membrane science, vol. 145, pp. 1 - 14, 1998. [42] M. Gryta and K. Karakulski, "The application of membrane distillation for the concentration of oil - water emulsions ," Desalination, vol. 121, pp. 23 - 29, 1999. [43] M. Bodzek and K. Konieczny, "The use of ultrafiltration membranes made of various polymers in the treatment of oil - emulsion wastewaters," Waste management, vol. 12, pp. 75 - 84, 1992. [44] K. Karakulski, A. Ko zlowski, and A. Morawski, "Purification of oily wastewater by ultrafiltration," Separations Technology, vol. 5, pp. 197 - 205, 1995. [45] - Razi, A. Pendashteh, L. C. Abdullah, D. R. A. Biak, S. S. Madaeni, and Z. Z. Abidin, "Review of technologies for oil and gas produced water treatment," J. Haz. Mat., vol. 170, pp. 530 - 551, 2009. [46] F. Carstensen, A. Apel, and M. Wessling, "In situ product recovery: submerged membranes vs. external loop membranes," Journal of membrane science, vol. 394, pp. 1 - 3 6, 2012. [47] P. Le - Clech, V. Chen, and T. A. Fane, "Fouling in membrane bioreactors used in wastewater treatment," Journal of membrane science, vol. 284, pp. 17 - 53, 2006. [48] C. Charcosset, "Membrane processes in biotechnology: an overview," Biotechnol. Adv., vol. 24, pp. 482 - 492, 2006. [49] G. Daufin, J. - P. Escudier, H. Carrere, S. Berot, L. Fillaudeau, and M. Decloux, "Recent and emerging applications of membrane processes in the food and dairy industry," Food Bioproduct Process., vol. 79, pp. 89 - 102, 2 001. [50] S. Lee, Y. Aurelle, and H. Roques, "Concentration polarization, membrane fouling and cleaning in ultrafiltration of soluble oil," Journal of Membrane Science, vol. 19, pp. 23 - 38, 1984. [51] P. Lipp, C. Lee, A. Fane, and C. Fell, "A fundamental study of the ultrafiltration of oil - water emulsions," Journal of Membrane Science, vol. 36, pp. 161 - 177, 1988. [52] Y. Matsumoto, T. Kawakatsu, M. Nakajima, and Y. Kikuchi, "Visualization of filtration phenomena of a suspended solution including O/W emulsion or solid particle and membrane separation properties of the solution," Water Research, vol. 33, pp. 929 - 936, 1999. 199 [53] D. - Q. Cao, E. Iritani, and N. Katagiri, "Properties of filter cake formed durin g dead - end microfiltration of O/W emulsion," Journal of Chemical Engineering of Japan, vol. 46, pp. 593 - 600, 2013. [54] E. N. Tummons, J. W. Chew, A. G. Fane, and V. V. Tarabara, "Ultrafiltration of saline oil - in - water emulsions stabilized by an anionic su rfactant: Effect of surfactant concentration and divalent counterions," Journal of Membrane Science, vol. 537, pp. 384 - 395, 2017. [55] E. N. Tummons, V. V. Tarabara, J. W. Chew, and A. G. Fane, "Behavior of oil droplets at the membrane surface during cross flow microfiltration of oil water emulsions," Journal of Membrane Science, vol. 500, pp. 211 - 224, 2016. [56] H. J. Tanudjaja, V. V. Tarabara, A. G. Fane, and J. W. Chew, "Effect of cross - flow velocity, oil concentration and salinity on the critical flux of an oil - in - water emulsion in microfiltration," Journal of Membrane Science, vol. 530, pp. 11 - 19, 2017. [57] J. Dickhout, J. Moreno, P. Biesheuvel, L. Boels, R. Lammertink, and W. de Vos, "Produced water treatment by membranes: A review from a colloidal per spective," Journal of colloid and interface science, vol. 487, pp. 523 - 534, 2017. [58] A. Motin, V. V. Tarabara, and A. Bénard, "Numerical investigation of the performance and hydrodynamics of a rotating tubular membrane used for liquid liquid separation," Journal of Membrane Science, vol. 473, pp. 245 - 255, 2015. [59] A. Grenier, M. Meireles, P. Aimar, and P. Carvin, "Analysing flux decline in dead - end filtration," Chemical Engineering Research and Design, vol. 86, pp. 1281 - 1293, 2008. [60] T. Darvishzadeh, V. V. Tarabara, and N. V. Priezjev, "Oil droplet behavior at a pore entrance in the presence of crossflow: Implications for microfiltration of oil water dispersions," Journal of membrane science, vol. 447, pp. 442 - 451, 2013. [61] Y. - J. Won, S. - Y. Jung, J. - H. Jang, J. - W. Lee, H. - R. Chae, D. - C. Choi , et al. , "Correlation of membrane fouling with topography of patterned membranes for water treatment," Journal of Membrane Science, vol. 498, pp. 14 - 19, 2016. [62] L. Hou, Z. Wang, and P. Song, "A precise combine d complete blocking and cake filtration model for describing the flux variation in membrane filtration process with BSA solution," Journal of Membrane Science, vol. 542, pp. 186 - 194, 2017. [63] W. Zhang, X. Ruan, Y. Ma, X. Jiang, W. Zheng, Y. Liu , et al. , "Modeling and simulation of mitigating membrane fouling under a baffle - filled turbulent flow with permeate boundary," Separation and Purification Technology, vol. 179, pp. 13 - 24, 2017. [64] S. J. Mamouri, V. V. Tarabara, and A. Bénard, "Analytical solution s for laminar flow in membrane channels with cylindrical symmetry: Single - and dual - membrane systems," Journal of Membrane Science, vol. 523, pp. 373 - 384, 2017. 200 [65] S. J. Mamouri, V. V. Tarabara, and A. Bénard, "Numerical Simulation of Filtration of Charge d Oil Particles in Stationary and Rotating Tubular Membranes," in ASME 2015 International Mechanical Engineering Congress and Exposition , 2015, pp. V07AT09A041 - V07AT09A041. [66] A. Pak, T. Mohammadi, S. Hosseinalipour, and V. Allahdini, "CFD modeling of po rous membranes," Desalination, vol. 222, pp. 482 - 488, 2008. [67] M. Yang, D. Yu, M. Liu, L. Zheng, X. Zheng, Y. Wei , et al. , "Optimization of MBR hydrodynamics for cake layer fouling control through CFD simulation and RSM design," Bioresource technology, v ol. 227, pp. 102 - 111, 2017. [68] A. Abrahamse, A. Van der Padt, R. Boom, and W. De Heij, "Process fundamentals of membrane emulsification: simulation with CFD," AIChE journal, vol. 47, pp. 1285 - 1291, 2001. [69] M. Rahimi, S. Madaeni, M. Abolhasani, and A. A. Alsairafi, "CFD and experimental studies of fouling of a microfiltration membrane," Chemical Engineering and Processing: Process Intensification, vol. 48, pp. 1405 - 1413, 2009. [70] A. Parvareh, M. Rahimi, S. Madaeni, and A. Alsairafi, "Experimental and CFD study on the role of fluid flow pattern on membrane permeate flux," Chinese Journal of Chemical Engineering, vol. 19, pp. 18 - 25, 2011. [71] B. Marcos, C. Moresoli, J. Skorepova, and B. Vaughan, "CFD modeling of a transient hollow fiber ultrafiltration system for protein concentration," Journal of Membrane Science, vol. 337, pp. 136 - 144, 2009. [72] D. Ramkrishna, "Population balances: theory and applications to particulate systems in engineering. 2000," Academic, New York . [73] K. A. Oucherif, S. Raina, L. S. Taylor, and J. D. Litster, "Quantitative analysis of the inhibitory effect of HPMC on felodipine crystallization kinetics using population balance modeling," CrystEngComm, vol. 15, pp. 2197 - 2205, 2013. [74] S. Zhan, M. Li, J. Zhou, J. Yang, and Y. Zh ou, "CFD simulation of dissolution process of alumina in an aluminum reduction cell with two - particle phase population balance model," Applied Thermal Engineering, vol. 73, pp. 805 - 818, 2014. [75] A. Chaudhury, H. Wu, M. Khan, and R. Ramachandran, "A mecha nistic population balance model for granulation processes: effect of process and formulation parameters," Chemical Engineering Science, vol. 107, pp. 76 - 92, 2014. [76] S. T. F. Mortier, K. V. Gernaey, T. De Beer, and I. Nopens, "Development of a population balance model of a pharmaceutical drying process and testing of solution methods," Computers & Chemical Engineering, vol. 50, pp. 39 - 53, 2013. 201 [77] M. Sen and R. Ramachandran, "A multi - dimensional population balance model approach to continuous powder mix ing processes," Advanced Powder Technology, vol. 24, pp. 51 - 59, 2013. [78] Á. Bárkányi, S. Németh, and B. G. Lakatos, "Modelling and simulation of suspension polymerization of vinyl chloride via population balance model," Computers & Chemical Engineering, vol. 59, pp. 211 - 218, 2013. [79] G. H. Yeoh, C. P. Cheung, and J. Tu, Multiphase flow analysis using population balance modeling: bubbles, drops and particles: Butterworth - Heinemann, 2013. [80] E. K. Yapp, C. G. Wells, J. Akroyd, S. Mosbach, and M. Kraft, "Modelling PAH curvature in a laminar premixed ethylene flame using a detailed population balance model," in Proc. Combust. Inst. , 2015. [81] Y. K. Ho, P. Doshi, H. K. Yeoh, and G. C. Ngoh, "Int erlinked population balance and cybernetic models for the simultaneous saccharification and fermentation of natural polymers," Biotechnology and bioengineering, vol. 112, pp. 2084 - 2105, 2015. [82] A. Bertucco, E. Sforza, V. Fiorenzato, and M. Strumendo, "P opulation balance modeling of a microalgal culture in photobioreactors: Comparison between experiments and simulations," AIChE Journal, vol. 61, pp. 2702 - 2710, 2015. [83] D. Ramkrishna and M. R. Singh, "Population balance modeling: current status and futur e prospects," Annual review of chemical and biomolecular engineering, vol. 5, pp. 123 - 146, 2014. [84] D. L. Marchisio, R. D. Vigil, and R. O. Fox, "Implementation of the quadrature method of moments in CFD codes for aggregation breakage problems," Chemical Engineering Science, vol. 58, pp. 3337 - 3351, 2003. [85] M. Hounslow, R. Ryall, and V. Marshall, "A discretized population balance for nucleation, growth, and aggregation," AIChE Journal, vol. 34, pp. 1821 - 1832, 1988. [86] M. M. Attarakih, C. Drumm , and H. - J. Bart, "Solution of the population balance equation using the sectional quadrature method of moments (SQMOM)," Chemical Engineering Science, vol. 64, pp. 742 - 752, 2009. [87] T. Wang, J. Wang, and Y. Jin, "Population balance model for gas - liquid flows: influence of bubble coalescence and breakup models," Industrial & engineering chemistry research, vol. 44, pp. 7540 - 7549, 2005. [88] R. Fan, D. L. Marchisio, and R. O. Fox, "Application of the direct quadrature method of moments to polydisperse gas solid fluidized beds," Powder technology, vol. 139, pp. 7 - 20, 2004. 202 [89] D. L. Marchisio, J. T. Pikturna, R. O. Fox, R. D. Vigil, and A. A. Barresi, "Quadrature AIChE Journal, vol. 49, pp. 1266 - 1276, 200 3. [90] D. L. Marchisio, A. A. Barresi, and M. Garbero, "Nucleation, growth, and agglomeration in barium sulfate turbulent precipitation," AIChE Journal, vol. 48, pp. 2039 - 2050, 2002. [91] R. Diemer and J. Olson, "A moment methodology for coagulation and b reakage problems: Part 2 moment models and distribution reconstruction," Chemical Engineering Science, vol. 57, pp. 2211 - 2228, 2002. [92] R. McGraw, "Description of aerosol dynamics by the quadrature method of moments," Aerosol Science and Technology, vol. 27, pp. 255 - 265, 1997. [93] R. G. Gordon, "Error bounds in equilibrium statistical mechanics," Journal of Mathematical Physics, vol. 9, pp. 655 - 663, 1968. [94] H. Dette and W. J. Studden, The theory of canonical moments with applications in statistics, pr obability, and analysis vol. 338: John Wiley & Sons, 1997. [95] A. Motin, J. M. Walsh, and A. Bénard, "Modeling Droplets Shearing and Coalescence Using a Population Balance Method in Produced Water Treatment Systems," in ASME 2015 International Mechanical Engineering Congress and Exposition , 2015, pp. V07BT09A023 - V07BT09A023. [96] Y. Liao and D. Lucas, "A literature review of theoretical models for drop and bubble breakup in turbulent dispersions," Chemical Engineering Science, vol. 64, pp. 3389 - 3406, 2009. [97] C. Coulaloglou and L. Tavlarides, "Description of interaction processes in agitated liquid - liquid dispersions," Chemical Engineering Science, vol. 32, pp. 1289 - 1297, 1977. [98] E. G. Chatzi, A. D. Gavrielides, and C. Kiparissides, "Generalized model for prediction of the steady - state drop size distributions in batch stirred vessels," Industrial & engineering chemistry research, vol. 28, pp. 1704 - 1711, 1989. [99] G. Narsimhan, J. Gupta, and D. Ramkrishna, "A model for transitional breakage probability of droplets in agitated lean liquid - liquid dispersions," Chemical Engineering Science, vol. 34, pp. 257 - 265, 1979. [100] C. - H. Lee, L. Erickson, and L. Glasgow, "Bubble breakup and coalescence in turbulent gas - liquid dispersions," Chemical Engineering Comm unications, vol. 59, pp. 65 - 84, 1987. 203 [101] C. M. - B. An, "On the breakup of an air bubble injected into a fully developed turbulent flow. Part 1. Breakup frequency," J. Fluid Mech, vol. 401, pp. 157 - 182, 1999. [102] F. Lehr and D. Mewes, "A transport equat ion for the interfacial area density in two - phase flow," in Second European Congress of Chemical Engineering, Montpellier, France , 1999. [103] columns," AIChE Journal, vol. 48, pp. 2426 - 2443, 2002. [104] H. P. Grace, "Dispersion phenomena in high viscosity immiscible fluid systems and application of static mixers as dispersion devices in such systems," Chemical Engineering Communications, vol. 14, pp. 225 - 277, 1982. [105 ] P. Elemans, H. Bos, J. Janssen, and H. Meijer, "Transient phenomena in dispersive mixing," Chemical engineering science, vol. 48, pp. 267 - 276, 1993. [106] X. Fu and M. Ishii, "Two - group interfacial area transport in vertical air water flow: I. Mechanisti c model," Nuclear Engineering and Design, vol. 219, pp. 143 - 168, 2003. [107] X. Sun, S. Kim, M. Ishii, and S. G. Beus, "Modeling of bubble coalescence and disintegration in confined upward two - phase flow," Nuclear Engineering and Design, vol. 230, pp. 3 - 26 , 2004. [108] J. Pandya and L. Spielman, "Floc breakage in agitated suspensions: theory and data processing strategy," Journal of Colloid and Interface Science, vol. 90, pp. 517 - 531, 1982. [109] B. Han, S. Akeprathumchai, S. Wickramasinghe, and X. Qian, "Flocculation of biological cells: Experiment vs. theory," AIChE journal, vol. 49, pp. 1687 - 1701, 2003. [110] S. Wickramasinghe, B. Han, S. Akeprathumchai, A. Jaganjac, and X. Qian, "Modeling flocc ulation of biological cells," Powder technology, vol. 156, pp. 146 - 153, 2005. [111] AIChE journal, vol. 42, pp. 1612 - 1620, 1996. [112] R. P. Hesketh, A. W. Etchells, and T. F. Russell, "Experimental observations of bubble breakage in turbulent flow," Industrial & Engineering Chemistry Research, vol. 30, pp. 835 - 841, 1991. [113] H. Zhao and W. Ge, "A theoretical bubble breakup model for slurr y beds or three - phase fluidized beds under high pressure," Chemical engineering science, vol. 62, pp. 109 - 115, 2007. [114] C. Tsouris and L. Tavlarides, "Breakage and coalescence models for drops in turbulent dispersions," AIChE Journal, vol. 40, pp. 395 - 4 06, 1994. 204 [115] K. J. Valentas, O. Bilous, and N. R. Amundson, "Analysis of breakage in dispersed phase systems," Industrial & engineering chemistry fundamentals, vol. 5, pp. 271 - 279, 1966. [116] M. Konno, M. Aoki, and S. Saito, "Scale effect on breakup pr ocess in liquid - liquid agitated tanks," Journal of Chemical Engineering of Japan, vol. 16, pp. 312 - 319, 1983. [117] columns," AIChE Journal, vol. 36, pp. 1485 - 1499, 1990. [118] Y. Liao and D. Lucas, "A literature review on mechanisms and models for the coalescence process of fluid particles," Chemical Engineering Science, vol. 65, pp. 2851 - 2864, 2010. [119] S. K. Friedlander, "Smoke, dust and haze: Fundamentals of aerosol behavior," New York, Wiley - Interscience, 1977. 333 p., 1977. [120] A. K. Chesters, "The modelling of coalescence processes in fluid - liquid dispersions: a review of current understanding," Chemical engineering research & design, vol. 69, pp. 259 - 270, 1991. [121] B. C. Venneker, J. J. Derksen, and H. E. Van den Akker, "Population balance modeling of aerated stirred vessels based on CFD," AIChE Journal, vol. 48, pp. 673 - 685, 2002. [122] H. Luo, "Coalescence, breakup and liquid circulation in bubble column reac tors," 1995. [123] S. Jeelani and S. Hartland, "Effect of approach velocity on binary and interfacial coalescence," Chemical engineering research & design, vol. 69, pp. 271 - 281, 1991. [124] R. I. Jeldres, F. Concha, and P. G. Toledo, "Population balance mo delling of particle flocculation with attention to aggregate restructuring and permeability," Advances in colloid and interface science, vol. 224, pp. 62 - 71, 2015. [125] P. Saffman and J. Turner, "On the collision of drops in turbulent clouds," Journal of Fluid Mechanics, vol. 1, pp. 16 - 30, 1956. [126] V. Runkana, P. Somasundaran, and P. Kapur, "A population balance model for flocculation of colloidal suspensions by polymer bridging," Chemical Engineering Science, vol. 61, pp. 182 - 191, 2006. [127] K. Samara s, A. Zouboulis, T. Karapantsios, and M. Kostoglou, "A CFD - based simulation study of a large scale flocculation tank for potable water treatment," Chemical Engineering Journal, vol. 162, pp. 208 - 216, 2010. [128] S. I. Karakashev and D. S. Ivanova, "Thin li quid film drainage: Ionic vs. non - ionic surfactants," Journal of colloid and interface science, vol. 343, pp. 584 - 593, 2010. 205 [129] D. Ekserova, I. Ivanov, and A. Sheludko, "Issled. v Obl. Poverkhn. Sil, Akad. Nauk SSSR, Inst," in Fiz. Khim., Sb. Dokl. na V toroi Konf., Moscow , 1964, pp. 158 - 163. [130] B. Radoev, D. S. Dimitrov, and I. Ivanov, "Hydrodynamics of thin liquid films effect of the surfactant on the rate of thinning," Colloid and Polymer Science, vol. 252, pp. 50 - 55, 1974. [131] T. Traykov, E. Mane v, and I. Ivanov, "Hydrodynamics of thin liquid films. Experimental investigation of the effect of surfactants on the drainage of emulsion films," International Journal of Multiphase Flow, vol. 3, pp. 485 - 494, 1977. [132] E. Ruckenstein and A. Sharma, "A n ew mechanism of film thinning: enhancement of Reynolds' velocity by surface waves," Journal of colloid and interface science, vol. 119, pp. 1 - 13, 1987. [133] C. Stubenrauch and R. Von Klitzing, "Disjoining pressure in thin liquid foam and emulsion films ne w concepts and perspectives," Journal of Physics: condensed matter, vol. 15, p. R1197, 2003. [134] S. I. Karakashev and A. V. Nguyen, "Effect of sodium dodecyl sulphate and dodecanol mixtures on foam film drainage: Examining influence of surface rheology and intermolecular forces," Colloids and Surfaces A: Physicochemical and Engineering Aspects, vol. 2 93, pp. 229 - 240, 2007. [135] M. Simon, "Koaleszenz von Tropfen und Tropfenschwärmen," Dissertation, TU Kaiserslautern, Germany, 2004. [136] A. Chesters and G. Hofman, "Bubble coalescence in pure liquids," in Mechanics and Physics of Bubbles in Liquids , ed: Springer, 1982, pp. 353 - 361. [137] D. Chappelear, "Models of a liquid drop approaching an interface," Journal of Colloid Science, vol. 16, pp. 186 - 190, 1961. [138] R. H. Davis, J. A. Schonberg, and J. M. Rallison, "The lubrication force between two viscous drops," Physics of Fluids A: Fluid Dynamics (1989 - 1993), vol. 1, pp. 77 - 81, 1989. [139] A. K. Chesters, "The applicability of dynamic - similarity criteria to isothe rmal, liquid - gas, two - phase flows without mass transfer," International Journal of Multiphase Flow, vol. 2, pp. 191 - 212, 1975. [140] A. Hasseine, A. H. Meniai, M. B. Lehocine, and H. J. Bart, "Assessment of drop coalescence and breakup for stirred extracti on columns," Chemical engineering & technology, vol. 28, pp. 552 - 560, 2005. [141] G. Narsimhan, "Model for drop coalescence in a locally isotropic turbulent flow field," Journal of colloid and interface science, vol. 272, pp. 197 - 209, 2004. 206 [142] I. Ivanov , D. Dimitrov, P. Somasundaran, and R. Jain, "Thinning of films with deformable surfaces: diffusion - controlled surfactant transfer," Chemical engineering science, vol. 40, pp. 137 - 150, 1985. [143] K. D. Danov, D. S. Valkovska, and I. B. Ivanov, "Effect of surfactants on the film drainage," Journal of colloid and interface science, vol. 211, pp. 291 - 303, 1999. [144] I. Ivanov and D. Dimitrov, Thin film drainage : Marcel Dekker, New York, 1988. [145] T. J. Pettersson, "Characteristics of suspended particles in a small stormwater pond," in Global Solutions for Urban Drainage , ed, 2002, pp. 1 - 12. [146] J. Vaze and F. H. Chiew, "Nutrient loads associated with different sediment sizes in urban stormwater and surface pollutants," Journal of Environmental Engineering , vol. 130, pp. 391 - 396, 2004. [147] G. W. Characklis, M. J. Dilts, O. D. Simmons, C. A. Likirdopulos, L. - A. H. Krometis, and M. D. Sobsey, "Microbial partitioning to settleable particles in stormwater," Water Research, vol. 39, pp. 1773 - 1782, 2005. [148] J. G. Li, H. Horneck, D. Averill, J. A. McCorquodale, and N. Biswas, "High - rate retention treatment basins for CSO control in Windsor, Ontario," Water quality research journal of Canada, vol. 39, pp. 449 - 456, 2004. [149] Q. Guo and C. C. Song, "Surging in urban storm drainage systems," Journal of hydraulic engineering, vol. 116, pp. 1523 - 1537, 1990. [150] United States Environment Protection Agency, 2004. [151] Q. Guo, Correlation of total suspended solids (TSS) and suspended sediment concentration (SSC) test methods, 2006. [152] J. Bridgeman, B. Jefferson, and S. A. Parsons, "Computational fluid dynamics modelling of flocculation in water treatment: a review," Engineering Applicati ons of Computational Fluid Mechanics, vol. 3, pp. 220 - 241, 2009. [153] J. C. Crittenden and H. Montgomery Watson, Water treatment principles and design . Hoboken, N.J.: J. Wiley, 2005. [154] Model for Eulerian Multiphase Flows," The Canadian Journal of Chemical Engineering, vol. 83, pp. 829 - 834, 2005. [155] P. Zhang, R. M. Roberts, and A. Bénard, "Computational guidelines and an empirical model for particle deposition in curved p ipes using an Eulerian - Lagrangian approach," Journal of Aerosol Science, vol. 53, pp. 1 - 20, 2012. 207 [156] J. Hinze, "Turbulence McGraw - Hill," New York, vol. 218, 1975. [157] M. Manninen, V. Taivassalo, and S. Kallio, "On the mixture model for multiphase flow ," ed: Technical Research Centre of Finland Finland, 1996. [158] L. Schiller and Z. V. Naumann, "Deutsch," Ing. 1935. pp. 77 - 318, 1935. [159] Q. Guo, "Correlation of total suspended solids (TSS) and suspended sediment concentration (SSC) test methods," Rutgers University, Department of Civil and Environmental Engineeering, 2006. [160] NJDEP, "Total Suspended Solids Laboratory Testing Procedures, Division of Science, Research and Technology," New Jersey Department of Environmental Protection (NJDEP), Tren ton, New Jersey2003. [161] A. Fluent, "Theory Guide 17.2," Ansys Inc. USA, 2016. [162] W. Zhang and T. Li, "Particle size distributions in combined sewer overflows in a high - intensity urban catchment in Shanghai, China," Desalination and Water Treatment, v ol. 56, pp. 655 - 664, 2015. [163] J. Jeong and F. Hussain, "On the identification of a vortex," Journal of fluid mechanics, vol. 285, pp. 69 - 94, 1995. [164] flow holl ow fiber membrane distillation devices integrated with a heat exchanger," AIChE J., vol. 57, pp. 1780 - 1795, 2011. [165] treatment in a membrane bioreactor: a review," Environ. Prog., vol. 23, pp. 59 - 68, 2004. [166] H. Yeh, H. Wu, and J. Dong, "Effects of design and operating parameters on the declination of permeate flux for membrane ultrafiltration along hollow - fiber modules," J. Membr. Sci., vol. 213, pp. 33 - 44, 2003. [167] Y. S. Polyakov, "Approximate method for nonlinear differential and integrodifferential equations," AIChE J., vol. 52, pp. 3813 - 3824, 2006. [168] A. S. Berman, "Laminar flow in channels with porous walls," J. Appl. Phys., vol. 24, pp. 1232 - 1235, 1953. [169] A. S. Berman, "Laminar flow in an annulus with porous walls," J. Appl. Phys., vol. 29, pp. 71 - 75, 1958. [170] R. M. Terrill and G. M. Shrestha, "Laminar flow through parallel and uniformly porous walls of different permeability," Zeitschrift für angewandt e Mathematik und Physik ZAMP, vol. 16, pp. 470 - 482, 1965. 208 [171] S. Yuan, "Laminar pipe flow with injection and suction through a porous wall," DTIC Document1955. [172] A. Kozinski, F. Schmidt, and E. Lightfoot, "Velocity profiles in porous - walled ducts," I nd. Eng. Chem. Fundam., vol. 9, pp. 502 - 505, 1970. [173] Y. Moussy and A. Snider, "Laminar flow over pipes with injection and suction through the porous wall at low Reynolds number," J. Membr. Sci., vol. 327, pp. 104 - 107, 2009. [174] Int, J. Non - Linear Mech., vol. 43, pp. 292 - 301, 2008. [175] R. MELLIS, W. N. Gill, and G. BELFORT, "Fluid dynamics in a tubular membrane: theory and experiment," Chem. Eng. Comm., vol. 122, pp. 103 - 125, 19 93. [176] D. Bhattacharyya, S. Back, R. Kermode, and M. Roco, "Prediction of concentration polarization and flux behavior in reverse osmosis by numerical analysis," J. Membr. Sci., vol. 48, pp. 231 - 262, 1990. [177] A. Da Costa, A. Fane, and D. Wiley, "Spac er characterization and pressure drop modelling in spacer - filled channels for ultrafiltration," J. Membr. Sci., vol. 87, pp. 79 - 98, 1994. [178] E. Pellerin, E. Michelitsch, K. Darcovich, S. Lin, and C. Tam, "Turbulent transport in membrane modules by CFD s imulation in two dimensions," J. Membr. Sci., vol. 100, pp. 139 - 153, 1995. [179] K. Madireddi, R. Babcock, B. Levine, J. Kim, and M. Stenstrom, "An unsteady - state model to predict concentration polarization in commercial spiral wound membranes," J. Membr. Sci., vol. 157, pp. 13 - 34, 1999. [180] V. Geraldes, V. Semião, and M. Norberta Pinho, "Numerical modelling of mass transfer in slits with semi - permeable membrane walls," Eng. Comp., vol. 17, pp. 192 - 218, 2000. [181] J. Pharoah, N. Djilali, and G. Vickers, "Fluid mechanics and mass transport in centrifugal membrane separation," J. Membr. Sci., vol. 176, pp. 277 - 289, 2000. [182] F. Chedevergne, G. Casalis, and J. Majdalani, "Direct Numerical Simulation and biglobal stability investigations of the gaseous moti on in solid rocket motors," J. Fluid Mech., vol. 706, pp. 190 - 218, 2012. [183] N. Tilton, D. Martinand, E. Serre, and R. M. Lueptow, "Incorporating Darcy's law for pure solvent flow through porous tubes: Asymptotic solution and numerical simulations," AICh E J., vol. 58, pp. 2030 - 2044, 2012. 209 [184] O. Makinde, S. Khamis, M. Tshehla, and O. Franks, "Analysis of heat transfer in Berman flow of nanofluids with Navier slip, viscous dissipation, and convective cooling," Adv. Math. Phys., vol. 2014, 2014. [185] S. Karode, "Laminar flow in channels with porous walls, revisited," J. Membr. Sci., vol. 191, pp. 237 - 241, 2001. [186] J. O. Wilkes, Fluid Mechanics for Chemical Engineers with Microfluidics and CFD: Pearson Education, 2006. [187] ANSYS, ANSYS ICEM CFD Tutori al Manual: ANSYS Inc., 2011. [188] M. Gupta and M. Goyal, "Viscous incompressible flow between two coaxial rotating circular cylinders with small uniform injection at the inner cylinder," Ind. J. Pure Appl. Math., vol. 2, pp. 657 - 72, 1972. [189] S. Patankar, Numerical heat transfer and fluid flow: CRC press, 1980. [190] R. J. Wakeman and E. Tarleton, "Membrane fouling prevention in crossflow microfiltration by the use of electric fields," Chem. Eng. Sci., vol. 42, pp. 829 - 842, 1987. [191] G. Daufin, J. - P. Escudier, H. Carrere, S. Berot, L. Fillaudeau, and M. Decloux, "Recent and emerging applications of membrane processes in the food and dairy industry," Food and Bioproducts Processing, vol. 79, pp. 89 - 102, 2001. [192] C. Charcosset, "Membrane process es in biotechnology: an overview," Biotechnology Advances, vol. 24, pp. 482 - 492, 2006. [193] flow hollow fiber membrane distillation devices integrated with a heat e xchanger," AIChE Journal, vol. 57, pp. 1780 - 1795, 2011. [194] S. Ri, X. Zhenliang, Z. Ying, and Y. Kim, "Experimental study on revolving cross - flow microfiltration of highly viscous liquids," Chinese Journal of Chemical Engineering, vol. 16, pp. 961 - 964, 2 008. [195] adhesion," AIChE journal, vol. 39, pp. 1292 - 1302, 1993. [196] N. S. Hanspal, A. Waghode, V. Nassehi, and R. J. Wakeman, "Development of a predictive mathematical model for coupled stokes/Darcy flows in cross - flow membrane filtration," Chemical Engineering Journal, vol. 149, pp. 132 - 142, 2009. [197] I. Hamedan, "A new technique for solving steady flow and heat transfer from a rotating disk in high permeability media," Int. J. of Appl. Math and Mech, vol. 8, pp. 1 - 17, 2012. 210 [198] A. Brou, L. Ding, P. Boulnois, and M. Y. Jaffrin, "Dynamic microfiltration of yeast suspens ions using rotating disks equipped with vanes," Journal of Membrane Science, vol. 197, pp. 269 - 282, 2002. [199] F. Liebermann, "Dynamic cross flow filtration with Novoflow's single shaft disk filters," Desalination, vol. 250, pp. 1087 - 1090, 2010. [200] O. Anwar Bég, O. Makinde, and J. Zueco, "Hydromagnetic Viscous Flow in a Rotating Annular High - Porosity Medium with Nonlinear Forchheimer Drag Effects: Numerical Study," World Journal of Modelling and Simulation, vol. 8, pp. 83 - 95, 2012. [201] N. A. Gumerov a nd R. Duraiswami, "Modeling of particle motion in viscous swirl flow between two porous cylinders," in Proc. ASME Fluids Engineering Division Summer Meeting , 1998, pp. 21 - 25. [202] J. Morel, Z. Lavan, and B. Bernstein, "Flow through rotating porous annuli, " Physics of Fluids (1958 - 1988), vol. 20, pp. 726 - 733, 1977. [203] E. C. Johnson and R. M. Lueptow, "Hydrodynamic stability of flow between rotating porous cylinders with radial and axial flow," Physics of Fluids (1994 - present), vol. 9, pp. 3687 - 3696, 1997 . [204] Y. Gu and D. Li, "Electric charge on small silicone oil droplets dispersed in ionic surfactant solutions," Colloids and Surfaces A: Physicochemical and Engineering Aspects, vol. 139, pp. 213 - 225, 1998. [205] M. Berthiaume and J. Jachowicz, "The eff ect of emulsifiers and oil viscosity on deposition of nonionic silicone oils from oil - in - water emulsions onto keratin fibers," Journal of colloid and interface science, vol. 141, pp. 299 - 315, 1991. [206] G. Stalidis, A. Avranas, and D. Jannakoudakis, "Inte rfacial properties and stability of oil - in - water emulsions stabilized with binary mixtures of surfactants," Journal of Colloid and Interface Science, vol. 135, pp. 313 - 324, 1990. [207] J. Sheng, Modern chemical enhanced oil recovery: theory and practice: G ulf Professional Publishing, 2010. [208] M. - m. Kim and A. L. Zydney, "Effect of electrostatic, hydrodynamic, and Brownian forces on particle trajectories and sieving in normal flow filtration," Journal of Colloid and Interface Science, vol. 269, pp. 425 - 43 1, 2004. [209] N. - T. Nguyen, Micromixers: fundamentals, design and fabrication: William Andrew, 2011. [210] A. Esmaeeli and M. N. Reddy, "The electrohydrodynamics of superimposed fluids subjected to a nonuniform transverse electric field," International Jo urnal of Multiphase Flow, vol. 37, pp. 1331 - 1347, 2011. 211 [211] A. F. Ansys, "14.0 Theory Guide," ANSYS inc, pp. 390 - 1, 2011. [212] B. P. Leonard, "A stable and accurate convective modelling procedure based on quadratic upstream interpolation," Computer meth ods in applied mechanics and engineering, vol. 19, pp. 59 - 98, 1979. [213] B. Andersson, R. Andersson, L. Håkansson, M. Mortensen, R. Sudiyo, and B. Van Wachem, Computational fluid dynamics for engineers : Cambridge University Press, 2011. [214] J. R. Cash a nd A. H. Karp, "A variable order Runge - Kutta method for initial value problems with rapidly varying right - hand sides," ACM Transactions on Mathematical Software (TOMS), vol. 16, pp. 201 - 222, 1990. [215] H. N. Joensson, M. Uhlén, and H. A. Svahn, "Droplet s ize based separation by deterministic lateral displacement separating droplets by cell - induced shrinking," Lab on a Chip, vol. 11, pp. 1305 - 1310, 2011. [216] A. B. Theberge, F. Courtois, Y. Schaerli, M. Fischlechner, C. Abell, F. Hollfelder , et al. , "Micro droplets in microfluidics: an evolving platform for discoveries in chemistry and biology," Angewandte Chemie International Edition, vol. 49, pp. 5846 - 5868, 2010. [217] A. Sarkar, K. K. Goh, and H. Singh, "Colloidal stability and interactions of milk - protei n - stabilized emulsions in an artificial saliva," Food Hydrocolloids, vol. 23, pp. 1270 - 1278, 2009. [218] T. Frankiewicz, "The Dirty Dozen," in The Produced Water Seminar, Houston, TX , 2001. [219] J. Cho, G. Amy, and J. Pellegrino, "Membrane filtration of n atural organic matter: factors and mechanisms affecting rejection and flux decline with charged ultrafiltration (UF) membrane," Journal of Membrane Science, vol. 164, pp. 89 - 110, 2000. [220] K. - J. Kim, V. Chen, and A. G. Fane, "Ultrafiltration of colloidal silver particles: flux, rejection, and fouling," Journal of colloid and interface science, vol. 155, pp. 347 - 359, 1993. [221] E. Gibbins, M. D'Antonio, D. Nair, L. S. White, L. M. F. dos Santos, I. F. Vankelecom , et al. , "Observations on solvent flux and solute rejection across solvent resistant nanofiltration membranes," Desalination, vol. 147, pp. 307 - 313, 2002. [222] drinking water by nanofiltration membranes," Separat ion and Purification Technology, vol. 42, pp. 137 - 144, 2005. 212 [223] A. Ullah, R. G. Holdich, M. Naeem, and V. M. Starov, "Stability and deformation of oil droplets during microfiltration on a slotted pore membrane," Journal of membrane science, vol. 401, pp . 118 - 124, 2012. [224] M. Cheryan and N. Rajagopalan, "Membrane processing of oily streams. Wastewater treatment and waste reduction," Journal of membrane science, vol. 151, pp. 13 - 28, 1998. [225] D. Lu, T. Zhang, and J. Ma, "Ceramic membrane fouling during ultrafiltration of oil/water emulsions: roles played by stabilization surfactants of oil droplets," Environmental science & technology, vol. 49, pp. 4235 - 4244, 2015. [226] B. Chakrabarty, A. Ghos hal, and M. Purkait, "Ultrafiltration of stable oil - in - water emulsion by polysulfone membrane," Journal of Membrane Science, vol. 325, pp. 427 - 437, 2008. [227] N. Nabi, P. Aimar, and M. Meireles, "Ultrafiltration of an olive oil emulsion stabilized by an a nionic surfactant," Journal of Membrane Science, vol. 166, pp. 177 - 188, 2000. [228] J. Pope, S. Yao, and A. Fane, "Quantitative measurements of the concentration polarisation layer thickness in membrane filtration of oil - water emulsions using NMR micro - ima ging," Journal of Membrane Science, vol. 118, pp. 247 - 257, 1996. [229] S. Yao, M. Costello, A. Fane, and J. Pope, "Non - invasive observation of flow profiles and polarisation layers in hollow fibre membrane filtration modules using NMR micro - imaging," Journ al of Membrane Science, vol. 99, pp. 207 - 216, 1995. [230] N. Shamsuddin, D. B. Das, and V. M. Starov, "Filtration of natural organic matter using ultrafiltration membranes for drinking water purposes: circular cross - flow compared with stirred dead end flow ," Chemical Engineering Journal, vol. 276, pp. 331 - 339, 2015. [231] A. Salahi, M. Abbasi, and T. Mohammadi, "Permeate flux decline during UF of oily wastewater: Experimental and modeling," Desalination, vol. 251, pp. 153 - 160, 2010. [232] B. Hu and K. Scott , "Microfiltration of water in oil emulsions and evaluation of fouling mechanism," Chemical Engineering Journal, vol. 136, pp. 210 - 220, 2008. [233] D. W. Stanton and C. J. Rutland, "Multi - dimensional modeling of thin liquid films and spray - wall interaction s resulting from impinging sprays," International Journal of Heat and Mass Transfer, vol. 41, pp. 3037 - 3054, 1998. [234] D. W. Stanton and C. J. Rutland, "Modeling fuel film formation and wall interaction in diesel engines," SAE Technical Paper 0148 - 7191, 1996. [235] P. J. O'Rourke and A. Amsden, "A spray/wall interaction submodel for the KIVA - 3 wall film model," SAE Technical Paper 0148 - 7191, 2000. 213 [236] T. Dhanasekaran and T. Wang, "Numerical model validation and prediction of mist/steam cooling in a 180 - degree bend tube," International Journal of Heat and Mass Transfer, vol. 55, pp. 3818 - 3828, 2012. [237] R. Ingle, R. Yadav, H. Punekar, and J. Cao, "Modeling of particle wall interaction and film transport using Eulerian wall film model," in ASME 2014 Gas Turbine India Conference , 2014, pp. V001T03A006 - V001T03A006. [238] I. Ivanov, Thin liquid films vol. 29: CRC Press, 1988. [239] H. M. Omar and S. Rohani, "Crystal population balance formulation and solution methods: A review," Crystal Growth & Design, 2017. [240] P. Carman, "Fundamental principles of industrial filtration," Trans. Inst. Chem. Eng, vol. 16, pp. 168 - 188, 1938. [241] G. Belfort, R. H. Davis, and A. L. Zydney, "The behavior of suspensions and macromolecular solutions in crossflow microfiltr ation," Journal of Membrane Science, vol. 96, pp. 1 - 58, 1994. [242] I. Mudawar and R. A. Houpt, "Measurement of mass and momentum transport in wavy - laminar falling liquid films," International Journal of Heat and Mass Transfer, vol. 36, pp. 4151 - 4162, 1993 .