University of Southern California Dissertations and Theses

Flow of fluids with yield stress in porous media

(USC Thesis Other)

Flow of fluids with yield stress in porous media

PDF

Download

Open document

Flip pages

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content FLOW OF FLUIDS WITH YIELD STRESS IN POROUS Copyright 2004 MEDIA by Min Chen A Dissertation Presented to the FACULTY OF THE GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulfillment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (CHEMICAL ENGINEERING) December 2004 Min Chen R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. UMI Number: 3155392 INFORMATION TO USERS The quality of this reproduction is dependent upon the quality of the copy submitted. Broken or indistinct print, colored or poor quality illustrations and photographs, print bleed-through, substandard margins, and improper alignment can adversely affect reproduction. In the unlikely event that the author did not send a complete manuscript and there are missing pages, these will be noted. Also, if unauthorized copyright material had to be removed, a note will indicate the deletion. ® UMI UMI Microform 3155392 Copyright 2005 by ProQuest Information and Learning Company. All rights reserved. This microform edition is protected against unauthorized copying under Title 17, United States Code. ProQuest Information and Learning Company 300 North Zeeb Road P.O. Box 1346 Ann Arbor, Ml 48106-1346 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. DEDICATION This dissertation is dedicated to my parents, grandparents, and my husband for their encouragement and support. ii R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. ACKNOWLEDGEMENTS I would like to express my sincere gratitude and appreciation to my advisor Professor Yannis C. Yortsos, who helped me to formulate and tackle the problem and shared his knowledge, experience and insight. His professionalism inspired me throughout the work. Special thanks also go to my other committee members, Prof. Charlie Sammis and Prof. Iraj Ershaghi for their assistance, and Prof. Ching-An Peng and Prof. Katherine Shing for serving on my Qualifying committee. Helpful suggestions from Prof. William R. Rossen are deeply appreciated. I am also very grateful for support and assistance from Prof. Muhammad Sahimi, and staff members in chemical engineering department, Karen Woo, Brendan Char, and Frank and Maria Valenzuela. R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. TABLE OF CONTENTS DEDICATION ii ACKNOWLEDGEMENTS iii LIST OF FIGURES vi ABSTRACT xvii 1. INTRODUCTION AND LITERATURE REVIEW 1 1.1. Introduction........................................................................................................1 1.2. Literature Review..............................................................................................8 1.2.1. Rheo logical Background..................................................................... 9 1.2.2. Foams in Porous M edia......................................................................13 1.2.2.1 Foam Texture........................................................................13 1.2.2.2 Foam Stability.......................................................................16 1.2.2.3 Foam Generation...................................................................19 1.2.2.4 Foam Mobility...................................................................... 22 2. FLOW AND DISPLACEMENT OF FLUIDS WITH YIELD STRESS 26 2.1. Introduction..................................................................................................... 26 2.2. Algorithm for Flow of a Fluid with Yield Stress....................................... 34 2.2.1. Invasion Percolation with Memory (IPM )......................................37 2.2.2. Extension of IPM to Account for Dynamic Effects...................... 41 2.3. Single Phase Flow .......................................................................................... 47 2.3.1. The Full Occupancy Case.................................................................. 48 2.3.2. The Partial-Occupancy C ase.............................................................64 2.4. Immiscible Displacement of a Bingham fluid by a Newtonian fluid 68 2.5. Summary and Conclusions............................................................................ 84 iv R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3. FOAM GENERATION AND MOBILIZATION IN POROUS MEDIA 86 3.1. Introduction.....................................................................................................86 3.2. Lamella-Creation Mechanisms..................................................................... 88 3.3. Models for Foam Generation....................................................................... 92 3.4. Pore Network Model Development..............................................................98 3.5. Lamella Generation during Displacement.................................................103 3.6. Lamella Generation Due to Capillary-Pressure Fluctuations................. 115 3.7. Repetitive Snap-off at Germination Sites..................................................129 3.7.1. State of Incipient Mobilization (Minimum Pressure Gradient) . 130 3.7.2. Sequentially Increasing Pressure Gradient....................................151 3.8. Summary and Conclusions..........................................................................154 4. FLOW OF A BINGHAM PLASTIC AROUND A CYLINDER IN A POROUS MEDIUM USING A CONTINUUM MODEL 157 4.1. Introduction...................................................................................................157 4.2. Problem Formulation....................................................................................159 4.3. Results............................................................................................................161 4.4. Conclusions...................................................................................................166 5. CONCLUDING REMARKS 167 5.1. Summary and Conclusions..........................................................................167 5.2. Recommendations for Future W ork...........................................................172 BIBLIOGRAPHY 174 APPENDIX A FLOW OF FLUID WITH YIELD STRESS IN PARALLEL CAPILLARIES 186 APPENDIX B SIMULTANEOUS TWO-PHASE FLOW IN A PORE NETWORK 188 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. LIST OF FIGURES 1.1. Shear stress and shearing rate relationship for (a) Bingham body; the slope /? gives the apparent viscosity, (b) Pseudo-plastic or shear-thinning fluid, (c) Dilatant or shear-thickening fluid. Dashed lines show Newtonian behavior (After Ref. [122])....................................................................................................10 1.2. Photomicrograph of foams in a transparent glass micro-model (After Ref. [71]).......................................................................................................................... 15 1.3. The ideal relationship of disjoining pressure 77 and the thickness of the film h (After Ref. [3])........................................................................................................ 17 1.4. Experimental results (After Ref. [7]).................................................................... 17 2.1. The current-voltage relation of an individual resistor in the network problem considered by Roux and Herrmann [108]. An analogous relationship for the flow rate-pressure drop dependence of a Bingham fluid in a pore throat will be used in the pore-network study in this chapter..............................................31 2.2. Description of the invasion rules in IPM, before (a), and after (b) a growth step. The growth site FG is that for which the sum S = F1 (F )+ rF ris minimum. F denotes a site currently on the front, F ' a site which is a nearest neighbor o f F. Function, F,(F)is recursively determined for all invaded sites, as described in the text. Invasion is from right to left (After Ref. [62])......... 40 2.3. The process of the opening of a new path after the applied pressure difference (P) is incrementally increased, so that only one new path opens. Curve ACB denotes an already open path. For an incrementally higher pressure, an additional path (curve ACD) opens. Here it is assumed to be connected to the MTP. Flow is from right to left............................................................................ 43 2.4. The PDF of the yield-stress threshold for a pore throat corresponding to the Rayleigh pore throat-size distribution in equation (2.11)................................. 45 vi R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 2.5. The patterns of the open paths for the dynamic case, corresponding from left- to-right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three-fourths of the sites belonging to the open paths. Full-occupancy case, uniformly distributed thresholds (n -1) and viscous model (2.5). A more layered structure of the mobilized paths is apparent in the dynamic case 46 2.6. The patterns of the open paths for the static case, corresponding from left-to- right, top-to-bottom, to the MTP, one-fourth o f the sites, half of the sites and three-fourths of the sites belonging to the open paths. Full-occupancy case, uniformly distributed thresholds (n=1) and viscous model (2.5). A more layered structure of the mobilized paths is apparent in the dynamic case 47 2.7. The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient, for the dynamic (solid line) and the static (dashed line) cases (A) and the dimensionless flow velocity-pressure gradient relationship. Full occupancy case, uniformly distributed thresholds («=T) and viscous model (5). Note the quadratic dependence near the onset of mobilization (region A) and the approach to a linear relationship at higher flow rates (point B)................................................51 2.8. Plot of the dimensionless permeability of the flowing fluid as a function of the fraction of sites belonging to the open paths. Solid line: full-occupancy case, uniformly distributed thresholds (n=l) and viscous model (2.5). Dashed line: full-occupancy case, strongly-disordered system (n=5) and viscous model (2.5). Dot-dashed line: full-occupancy case, Rayleigh throat-size distribution (2.11) and viscous model (2.6). The results appear to be well fitted with a cubic relationship (y=0.2994x3, y=0.5641x3, and y=0.9866v3 for the respective cases). For the Rayleigh case, the cubic fit pertains only to fraction of sites less than 0.9................................................................................. 54 2.9. The patterns of the open paths at different stages for the dynamic case, and for «=0.1 corresponding from left-to-right, top-to-bottom, to the MTP, one- fourth of the sites, half o f the sites and three-fourths of the sites belonging on open paths. Full-occupancy case and viscous model (2.5)............................... 55 2.10. The patterns of the open paths at different stages for the dynamic case, and for n— 5, corresponding from left-to-right, top-to-bottom, to the MTP, one- fourth of the sites, half of the sites and three-fourths of the sites belonging on open paths. Full-occupancy case and viscous model (2.5)............................... 56 2.11. The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient for the dynamic case (panel A) and the dimensionless flow velocity-pressure gradient relationship (panel B). Full-occupancy case, strongly-disordered system (n=5) and viscous vii R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. model (2.5). Note again the quadratic dependence near the onset of mobilization and the approach to a linear relationship at higher flow rates... 58 2.12. The patterns of the open paths at different stages for the dynamic case corresponding from left-to-right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three-fourths of the sites belonging to open paths. Full-occupancy case, uniformly distributed thresholds (n— 1, range (0.1, 1)) and viscous model (2.6).........................................................................................60 2.13. The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient for the dynamic case (panel A) and the dimensionless flow velocity-pressure gradient relationship. Full-occupancy case, uniformly distributed thresholds (n=1, range (0.1, 1)) and viscous model (2.6).........................................................................................61 2.14. The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient for the dynamic case (panel A) and the dimensionless flow velocity-pressure gradient relationship (panel B). Full-occupancy case, Rayleigh throat-size distribution (2.11) and viscous model (2.6). The large conductances dominate the flow behavior and the flow rate-pressure drop relation rapidly approaches a straight line 63 2.15. The open paths at two different stages of mobilization (Panel A) for the partial-occupancy case (S=0.5), Rayleigh throat-size distribution (2.11), and viscous model (2.6). Gray denotes the fluid, the occupancy of which was determined by Invasion Percolation without trapping, black denotes the open paths. Note the existence of a number o f dead-end pores.................................66 2.16. The fraction of sites (pores) belonging to the open paths as a function of the applied pressure gradient (Panel A) and the flow rate-pressure drop relation for the partial-occupancy case (<S=0.5), Rayleigh throat-size distribution (2.11), and viscous model (2.6)............................................................................ 67 2.17. The fraction of sites (pores) belonging to the open paths as a function of the applied pressure gradient (Panel A) and the flow rate-pressure drop relation for the partial-occupancy case (<S-0.5), Rayleigh throat-size distribution (2.11), and viscous model (2.6)............................................................................ 68 2.18. Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime A o f Figure 7B (onset of mobilization) for M=oo. Gray denotes mobilized Bingham fluid, white denotes non mobilized Bingham fluid and black is the displacing Newtonian fluid 73 2.19. Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime A of Figure 7B (onset of mobilization) viii R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. for M= 1. Gray denotes mobilized Bingham fluid, white denotes non mobilized Bingham fluid and black is the displacing Newtonian fluid 74 2.20. Patterns o f the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime A o f Figure 7B (onset of mobilization) for M=0.05. Gray denotes mobilized Bingham fluid, white denotes non mobilized Bingham fluid and black is the displacing Newtonian fluid 75 2.21. The variation of flow rate with time during the displacement of a Bingham fluid by a Newtonian fluid, corresponding to point A of Figure 7B, and for different values of the viscosity ratio. Results from pore-network simulations (panel A), simple analytical model (panel B)..................................................... 81 2.22. Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime B of Figure 7B (after the mobilization of most of flow paths) for M = < x> . Gray denotes mobilized Bingham fluid, white denotes non-mobilized Bingham fluid and black is the displacing Newtonian fluid. The relationship to DLA patterns is apparent...................... 82 2.23. Phase diagram o f the displacement of a Bingham fluid by a Newtonian fluid. The yield-stress dominant regime corresponds to the needle-like structure of the MTP pattern and may be affected also by the threshold distribution. The viscous dominant regime is of the DLA or the PL type depending on the values of the viscosity ratio (large or small, respectively). In the figure, the viscosity ratio is assumed large (hence the DLA-type regime shown). When the Bingham fluid displaces the Newtonian fluid, the yield-stress dominant regime is PL, as in [65]..........................................................................................83 3.1. Schematic of lamella creation by the leave-behind mechanism. Gray diamonds represent sand grains; gaps between them represent pore bodies and throats (After Ref. [113]).......................................................................................89 3.2. Schematic of a single lamella-division event. (After Ref. [103])..................... 89 3.3. Schematic of lamella division on the pore-network scale. In the first panel, one moving lamella produces 14 lamellae in place by division. In the second panel, a second moving lamella (taking a different path) produces 9 more in place, plus displaces an extra 2 moving ahead of it. In the third panel, a third moving lamella (taking a different path) produces 2 more in place, plus displaces an extra 5 moving ahead of it (After Ref. [113]).............................. 90 3.4. Schematic of snap-off in a pore throat. Black denotes pore-throat wall, gray the water and white the gas (After Ref. [113])................................................... 90 IX R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3.5. Schematic of Roof snap-off as gas invades a pore body. Black is pore wall, gray is water and white is gas............................................................................... 92 3.6. Schematic of the results of the percolation model for foam generation by Rossen and Gauglitz. / is the fraction of pore throats through which gas can flow (not blocked by water-filled pores, lenses or lamellae). Here 7=0.25 represents the percolation threshold, where the pressure gradient for foam mobilization approaches zero. The diagonal line on right side of figure represents foam generation during continuous injection of gas and liquid (starting from no-foam state of continuous gas flow). The diagonal line on the left side of the figure represents mobilization o f discontinuous-gas foam (After Ref. [106])....................................................................................................95 3.7. Schematic of results of the percolation model o f Rossen et al. [41] modified from that of Rossen and Gauglitz [104]. The two halves of the model diverge from each other at the percolation threshold (/= 0.25)....................................... 97 3.8. Sequence of patterns illustrating liquid displacement in a 2-D network and the creation of lamellae. The top panel shows the conditions at breakthrough, the middle panel conditions at a higher value of the applied pressure drop. The bottom panel is an enlargement of a region in the middle panel showing a number of lamellae generated by snap-off and leave-behind mechanisms. The liquid flows at a constant rate corresponding to Ca=10' . Gas is white, liquid and rock grains are black. Extensive lamella formation by leave-behind occurs as the liquid is being displaced...............................................................106 3.9. Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability 7s=0.1, M=100 and Ca=10'3. Relative permeabilities (kr) as a function o f the applied pressure drop. The liquid saturation Sw falls as p rises....................................................107 3.10. Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability f s=0.1, M - 100 and Ca= 10‘3. The fraction F of pore throats containing gas without lamellae, leave-behind lamellae, snap-off lamellae, and lamellae created by division, and the flow velocity-pressure drop relation, as a function of the applied pressure drop. The liquid saturation Sw falls as Ap rises................................. 108 3.11. Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability / s =0.7, M=100 and Ca=10'3. Relative permeabilities (kr), the fraction F of pore throats containing gas without lamellae, leave-behind lamellae, snap-off lamellae, and lamellae created by division, and the flow velocity-pressure drop relation, as a function of the applied pressure drop. Note the existence of a critical gas saturation x R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. (and a critical pressure gradient), below which the gas relative permeability vanishes..................................................................................................................110 3.12. Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability f s=0.7, M=100 and Ca=10'3. The fraction F of pore throats containing gas without lamellae, leave-behind lamellae, snap-off lamellae, and lamellae created by division, and the flow velocity-pressure drop relation, as a function of the applied pressure drop. Note the existence of a critical gas saturation (and a critical pressure gradient), below which the gas relative permeability vanishes I l l 3.13. Steady-state relative permeability results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability f s=0.1 (top panel) and f s=0.7 (bottom panel), M= 100 and Ca=10' . The increase in the liquid rate (capillary number) leads to increased relative permeabilities in both fluids.............................................................................................................. 114 3.14. The minimum pressure gradient corresponding to the state of incipient mobilization, in which lamellae are generated by capillary pressure fluctuations. (fs=0.3, />c=5.1677, pw=7.6677, 0 2 =7 .6 7 7 x 10"4, Z=40). The value fluctuates between two extreme states, one in which the gas has a continuous path (corresponding to a zero value of the minimum pressure gradient) and a state in which the gas is completely discontinuous. In-between states correspond to the displacement of lamellae from the blocked paths and the creation of a continuous path. The number of steps corresponds to the number of lamellae moved.................................................................................. 117 3.15. Lamella patterns corresponding to Figure 3.14. The system fluctuates between two extreme states, one in which the gas has a continuous path (shown above as a white line) and a state in which the gas is completely discontinuous. When the first state is reached, the capillary pressure drops and new lamellae are generated in the continuous path. These are advected downstream until a continuous gas path opens. For emphasis, in the top figure gas in the continuous path is white and the remaining gas gray; otherwise, gas is white, liquid and rock grains black................................................................................ 118 3.16. The minimum pressure gradient corresponding to the state o f incipient mobilization, in which lamellae are generated by capillary pressure fluctuations. (£=0.7,/>c=5.1677,/?/0 =7.6677, Ca=7.677xl0'4, L=40). The value fluctuates between two extreme states, one in which the gas has a continuous path (corresponding to a zero value of the minimum pressure gradient) and a state in which the gas is completely discontinuous. In-between states correspond to the displacement of lamellae from the blocked paths and the creation of a continuous path. The number of steps corresponds to the number of lamellae moved.................................................................................................119 xi R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3.17. Lamella patterns corresponding to Figure 3.16. The system fluctuates between two extreme states, one in which the gas has a continuous path (shown above as a white line) and a state in which the gas is completely discontinuous. When the first state is reached, the capillary pressure drops and new lamellae are generated in the continuous path. These are advected downstream until a continuous gas path opens. For emphasis, in the top figure gas in the continuous path is white and the remaining gas gray; otherwise, gas is white, liquid and rock grains black................................................................................ 120 3.18. The minimum pressure gradient at incipient mobilization corresponding to the maximum value of the fluctuating cycle between the states, as a function of the parameter p c-pio (for f s~0-3), at the top panel, and of the initial snap-off probability f , (for p c-pw=-2.5 corresponding to point A in the top panel) at the bottom panel.......................................................................................................... 121 3.19. Lamellae patterns corresponding to two different states under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations (fs=0.7, p c=5A677, pio= 7.6677, |Yp|=0.17, G z = 6 .8 x 10'4). The system fluctuates between two extreme states, one in which the gas has a continuous path (shown in the bottom figure as a white line) and a state in which the gas is completely discontinuous (top figure, corresponding to the initial state, but which is also similar to one of the two extremes of the asymptotic state). Gas is gray or white, liquid and rock grains black...................................................................................................125 3.20. The gas flow rate as a function of time for the process of Figure 3.19, namely under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations. One complete cycle is shown. At the end of the cycle, a continuous gas path opens, new lamellae are generated by capillary pressure and the process is repeated. The increasing flow rate as a function of time results from the mobilization and the advection downstream of lamellae in a pore network in which the gas phase is disconnected (fs=0.7, p c=5A677, pi0 =7.6677, |Vp|=0.17, Ca=6.8xl0’4 )...................................................................................... 126 3.21. Lamellae patterns corresponding to two different states under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations (fs=0.7, p c=5A677, pw=7.6677, |Vp|=0.275, C a = l.lx l0 '3 ). The system fluctuates between two extreme states, one in which the gas has a continuous path (shown in the bottom figure as a white line) and a state in which the gas is completely discontinuous (top figure, corresponding to the initial state, but which is also similar to one of the two extremes of the asymptotic state). Gas is gray or white, liquid and rock grains black...................................................................................................127 xii R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3.22. The gas flow rate as a function of time for the process of Figure 3.21, namely under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations. One complete cycle is shown. At the end of the cycle, a continuous gas path opens, new lamellae are generated by capillary pressure and the process is repeated. The increasing flow rate as a function o f time results from the mobilization and the advection downstream of lamellae in a pore network in which the gas phase is disconnected (fs=0.7, p c— 5.\677, pio=7.6677, |Vp|=0.275, C a = l.lx l0 '3 )....................................................................................128 3.23. Foam generation by repetitive snap-off at germination sites (fs=0.2, uniform pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization in panel (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is continuous in all states......................................................................................... 134 3.24. Foam generation by repetitive snap-off at germination sites (fs- 0.2, uniform pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states of the Figure 3.23a. The gas is continuous in all states....................................................................................................................... 135 3.25. Foam generation by repetitive snap-off at germination sites (fs- 0.2, uniform pore-size distribution): Simplified schematic illustrating the concept of a lamella streamtube. The arrows in the latter indicate the mobilized lamellae and the direction in which they are mobilized. In the first panel, the moving lamella divides into two unoccupied downstream throats. In the second panel, the mobilized lamella moves to the throat on the left. In the third and consecutive panels (not shown here for simplicity), the lamella moves in the direction of the arrow and occupies consecutive pores without any new lamella generation. The pressure gradient fluctuates as a result of the fluctuating pore size. In the fourth and last panel, the lamella exits the system, thus recreating the top panel. This cycle is then repeated...............................136 3.26. Foam generation by repetitive snap-off at germination sites (fs=0.7, uniform pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is discontinuous in all states.................................................................................... 137 3.27. Foam generation by repetitive snap-off at germination sites (fs- 0.7, uniform pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.26a. Gas is white, liquid and rock grains are black..................................................................................................... 138 x iii R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3.28. The minimum pressure gradient as a function of the final fraction of the pore throats occupied by a lamella for the case o f foam generation by repetitive snap-off at germination sites and uniform pore size distribution. Compare to the plot in Figures 3.6 and 3.7 (noting that F = \-f and that here the results correspond to a network with different coordination number than that in Figures 3.6 and 3.7).............................................................................................. 139 3.29. Foam generation by repetitive snap-off at germination sites (£=0.2, uniform pore-size distribution, lamellae moved out flow back in from inlet): The minimum pressure gradient corresponding to the state o f incipient mobilization (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is discontinuous in all states.. 142 3.30. Foam generation by repetitive snap-off at germination sites (£=0.2, uniform pore-size distribution, lamellae moved out flow back from inlet): Two lamella patterns corresponding to the initial state and one o f the final states of Figure 3.29a. Gas is white, liquid and rock grains are black.......................................143 3.31. Foam generation by repetitive snap-off at germination sites (£=0.7, uniform pore-size distribution, lamellae moved out flow back in from inlet): The minimum pressure gradient corresponding to the state of incipient mobilization (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is discontinuous in all states.. 144 3.32. Foam generation by repetitive snap-off at germination sites (£=0.7, uniform pore-size distribution, lamellae moved out flow back from inlet): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.31a. Gas is white, liquid and rock grains are black.......................................145 3.33. Foam generation by repetitive snap-off at germination sites (£=0.2, unimodal pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization and the fraction of pore throats occupied by continuous gas or lamellae of the various types. The gas is continuous in all states....................................................................................................................... 146 3.34. Foam generation by repetitive snap-off at germination sites (£=0.2, unimodal pore-size distribution): Two lamella patterns corresponding to the initial state and one o f the final states of Figure 3.33a. The gas is continuous in all states. Gas is white, liquid and rock grains are black.................................................. 147 3.35. Foam generation by repetitive snap-off at germination sites (£=0.7, unimodal pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization and the fraction of pore throats occupied by continuous gas or lamellae of the various types. The gas is discontinuous in all states..................................................................................................................148 xiv R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3.36. Foam generation by repetitive snap-off at germination sites unimodal pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.35a. Gas is white, liquid and rock grains are black. The gas is discontinuous in all states....................................149 3.37. Foam generation by initial snap-off at specific germination sites, corresponding to f s=0.2 (top) and f s=0.7 (bottom). The minimum pressure gradient corresponding to the state of incipient mobilization is shown. The gas is continuous in all states, and very weak foam is observed....................150 3.38. Foam generation by repetitive snap-off at germination sites corresponding to f s=0.3 and a uniform pore-size distribution. The top panel shows the minimum pressure gradient at conditions of incipient mobilization for the specific f s value. The bottom panel is a plot of the flow velocity-pressure drop relationship. Because the gas phase is continuous, the fluid behaves initially as Newtonian, with a shear-thickening behavior later, as lamellae are mobilized and additional lamellae are created by division............................. 153 3.39. Foam generation by repetitive snap-off at germination sites corresponding to f s- 0.7 and a uniform pore-size distribution. Because of the high snap-off probability, the system rapidly becomes fully occupied by lamellae (through snap-off and lamella division); hence it becomes identical to a Bingham plastic, as discussed in chapter 2. The corresponding flow velocity-pressure drop relationship is identical to that of the flow o f a fluid with yield stress. 154 4.1. The change o f the yield surface with different s values. The arrow in the figure shows the direction of increasing e......................................................... 165 B. 1 Snapshots of invasion of the non-wetting phase, before (top panel) and after a step (bottom panel). Black part is non-wetting phase, white part is wetting phase....................................................................................................................... 198 B.2 Relative permeabilities of the two phases as functions of the wetting phase saturation, for drainage case, M= 1 and different capillary numbers. Solid line, IP; dotted line, Ca=10"4; short dash, Ca=10'2; and dash-dot line, Ca=0.5... 199 B.3 The necessary pressure gradient of the system as a function of the wetting phase saturation, for drainage case, M= 1 and different capillary numbers. Top panel, Ca=10'4; middle panel, Gz=10‘2; and bottom panel, Ca=0.5..............200 B.4 Change o f fraction of wetting phase flow rate with wetting phase saturation, for drainage case, M= 1 and different capillary numbers. Solid line, C<a=10"4, dash line, Ca= 10'2, and dotted line, Ca=0.5..................................................... 201 x v R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. B.5 Change o f fraction of wetting phase flow rate with wetting phase saturation, for drainage case, Ca=10"4 and different viscosity ratio M values. Solid line, M= 0.5, medium dash, M= 1, short dash, M= 2 and dotted line, M=5............. 202 B.6 The necessary pressure gradient of the system as a function of wetting phase saturation, for drainage case, Ca=10'4 and different viscosity ratio M values. Solid line, M= 0.5; medium dash, M= 1; short dash, M= 2; and dotted line M= 5. ................................................................................................................................ 203 B.7 Relative permeabilities of the two phases as functions of the wetting phase saturation for imbibition case, M= 1 and different capillary numbers. Solid line, Cb=T0'4; short dash, C<a=10"2; and dotted line, Ca— 0.5..........................204 B.8 The necessary pressure gradient of the system as a function of wetting phase saturation, for imibibition case, M= 1 and different capillary numbers. Top panel, Ca=10'4; middle panel, Ca=10'2; and bottom panel, Ca=0.5..............205 B.9 Fraction of wetting phase flow rate as a function of wetting phase saturation, for imbibition case, M= 1 and different capillary numbers. Solid line, Ca=10’4; short dash, Ca=10‘2; and dotted line, Ca=0.5................................................... 206 B.10 Fraction of wetting phase flow rate vs. wetting phase saturation, for imbibition case, Ca=10'4 and different viscosity ratio M values. Solid line, M -20; long dash, M= 5; short dash, M= 1; dash-dot, M= 0.2; and dotted line, M=0.05...................................................................................................................207 B.ll The necessary pressure gradient vs. wetting phase saturation, for imbibition case, Ca=10'4 and different viscosity ratio M values. Solid line, M= 20; long dash, M= 5; short dash, M= 1; dash-dot, M= 0.2; and dotted line, M=0.05. ...208 B.12 Different flow directions of the moving interface, downwards in (a) and upwards in (b)....................................................................................................... 209 B.13 koj/r*2 as a function of vapor saturation at different viscosity and density ratios....................................................................................................................... 210 x v i R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. ABSTRACT This dissertation deals with the flow of fluids with yield stress in porous media using numerical simulations in pore network model. A dynamic algorithm based on Invasion Percolation with Memory (IPM) is proposed. The main idea of this new algorithm is to take into account viscous flow effects and to update the thresholds correspondingly. It is found that the relation between the macroscopic pressure gradient to mobilize a fluid with yield stress and the flow rate of the fluid is consistent with the results obtained by Roux and Herrmann, which is also observed for different distributions of thresholds and conductances. The fraction o f the sites which belong to the open paths increases more slowly with the applied pressure gradient, compared to some other cases studied. In all cases, the qualitative features of the quadratic behavior of the flow rate near the onset of flow can be captured by a model of parallel capillaries. R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Based on the patterns obtained by applying the above algorithm to displacements, a phase diagram is constructed involving the capillary number and the yield stress number. The previous algorithm, appropriately modified is used to study different scenarios of foam generation and mobilization in porous media. Drainage displacement, lamella generation due to capillary pressure fluctuations and repetitive Roof snap-off are investigated where foam formation mechanisms are included correspondingly. In each process, different cases are considered, such as the incipient mobilization, and the application of constant or increasing pressure gradient. The relative permeabilities of the two phases at steady states are computed. Foam properties are found to be fluctuating in some cases with a well-defined period and in others in a random-chaotic manner. A new foam generation mechanism is presented and the result shows that it is possible to generate strong foams without repetitive Roof snap-off. Minimum pressure gradient calculations using the latter lead to results qualitatively similar to the theory of Rossen and Gauglitz. A continuum model for the flow of Bingham plastic around a circular cylinder is solved asymptotically in the last part of this dissertation. The result obtained is qualitatively consistent with observations by other researchers. xviii R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Chapter 1 INTRODUCTION AND LITERATURE REVIEW 1.1 Introduction Many problems of importance to industry and science involve flows in porous media, including, to name a few, filtration, ground water movement, production of gas and oil from reservoirs, chemical reactor design, polymer engineering and separation processes. Because of the obvious complexity of porous media, a number of studies are through direct experimentations. The fluids involved in the flows may be Newtonian or non-Newtonian. For Newtonian fluids, after many years of extensive studies, we have now a good understanding of characteristics and flow mechanisms. On the other hand, non-Newtonian fluid flow in porous media, which is encountered typically in the fields of chemical, environmental and petroleum engineering [27, 118], still requires significant additional investigation for its full understanding. In Enhanced Oil Recovery, extemal-drive mechanisms, such as water flood, polymer flooding, C 02 flood, steam injection etc, are employed to improve 1 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. recovery performance. Many fluids in these processes involve complex physical and chemical properties, such as shear-dependent viscosity, thixotropy and elasticity. Interactions happen between fluids, or even with the porous formation, when they come into contact. With length scales ranging from 10'6 to 103 m, the underground structure itself has significant heterogeneities [92]. These complexities impose many difficulties on investigations. A particular class o f complex fluids, of particular interest to this thesis, is fluids with yield stress. Even though the very subject of yield stress has been under debate in the last decade or so [5], it has been used quite often in different areas, including colloid science, food science, chemical and biochemical engineering, and blood flow. The simplest definition of such a fluid is one which exhibits a drastic change from a stationary solid state to a flowing liquid state, i.e. continual deformation, when a specific yield stress threshold has been exceeded. Typical examples of fluids with yield stress are Bingham plastics and foams. Bingham plastics are named after E. C. Bingham [11], Some examples of Bingham plastics are pain, slurries, pastes, and food substances like margarine, mayonnaise and ketchup [90]. Other fluids are approximated as Bingham plastics as well, including drilling mud [131] and heavy oils [17, 131]. Bird et al. [13] gave a review of systems which exhibit yield-stress behavior. Bingham plastics are special in that the fluid shows a yield stress followed by nearly Newtonian behavior. In EOR processes, one of the intentions of the recovery method is to reduce the mobility ratio between injected and in-place fluids to improve sweep efficiency 2 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. [19]. Foams can be used as mobility control agents in miscible floods and well treatment [77]. By blocking the pore throats due to the gas-liquid films, foams flowing in permeable media can drastically reduce the mobility o f a gas phase, which thus can be viewed as an enhancement o f effective viscosity of a single phase flow or as a decrease in the gas phase permeability. The mobility reduction can alleviate undesirable features, such as gravity override, gas channeling and viscous fingering caused by low viscosity and density o f gas. With potential in increasing oil recovery efficiency, foams are under wide investigation to better understand their properties, and the mechanisms by which they flow and interact with porous media and each other. Past models on transport processes involving Non-Newtonian fluids have been based on continuum representation [4, 6, 92]. The common feature of most of the applications of continuum models is the use of Darcy’s law or its straightforward extension, by introducing an effective viscosity in single-phase flow and a relative permeability in multi-phase flow. A similar approach is followed in the studies on Bingham plastics, with the fluid assigned an effective yield stress, the subsequent flow being approximated as that in a single capillary. Wu [130, 131] used such a method to study single-phase flow and fluid displacement. But there exist extensive debate on whether such continuum methods can give an effective account of the effect of microstructure, hence predict the behavior of the complex fluids when they flow through porous media. An alternative method is the use of micro-models which represent a simplified geometry but they allow for a direct conversion from 3 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. micro-model quantities to yield values of continuum variable. To extend the laws of motions for purely Newtonian fluids to Theologically complex ones, one has to have apparent fluid properties derived from a representation of the porous media. A common approach is to represent the medium by a bundle of parallel capillary tubes. In these simplified models, the macroscopic properties o f the fluids can be expressed in analytical form, e.g. as an apparent shear rate [21, 52, 121, 124, 128]. One model which has been used widely and which will also be used in this thesis is the pore-bond network model, in which the porous media is simplified to pores interconnected by throats with different coordination numbers. Pore-bond network models are ideally suited to computational techniques and give considerable insights to physical phenomena as well. These kinds of models have been extensively used to analyze multiphase transport phenomena [16, 18, 87, 88, 91]. Although originally developed to study purely Newtonian fluids, they have been extended to investigate the single-phase flow of non-Newtonian fluids, including power-law fluids [80, 92, 117], Xanthan gum [121], and Bingham plastics [62, 64], Pearson and Tardy [92] concluded that the tortuosity of the medium significantly affected transport properties. Lopez et al. [80] used a detailed description of the pore space derived from a sand pack and Berea sandstone to make predictions of properties in natural porous materials. Pore-network models can naturally account for such effects. Although there is a wealth of literature on Non-Newtonian fluid flows in the bulk, the same does not apply to the study of flows involving fluids with yield 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. stress in porous media. One topic which has attracted some attention is how to identify the minimum threshold path (MTP) [62, 64, 105, 109], Using EMA, Sahimi [109] showed that the flow rate-pressure drop relationship is linear when the flow rate is large enough and non linear when the flow rate is small. His results are consistent with Roux and Herrmann’s study on a network consisted with resistors [108]. K harabaf s work [62, 64] created a pioneering algorithm, termed Invasion Percolation with Memory (IPM) to identify the Minimum Threshold Path (MTP). Even though IPM is a novel and effective approach to studying MTP, its application to mobilization and displacement of Bingham plastics was limited by the fact that the authors did not consider effects introduced by viscous flow. For studies on foams, one focus is the foam generation processes, with the subsequent foam properties determined by the process of generation and flow in porous media [2, 37, 43, 49, 66, 73, 74, 89, 96, 99, 119, 123]. Independently, research work in relatively simple geometries has been carried out to understand the properties of foams without the complexity imposed by porous media [45, 54, 67, 85, 98, 102, 103, 132]. Population-balance model together with pore-network models, especially percolation theory and related model have been used for the modeling of foam generation and mobilization [10, 28, 29, 35, 39, 41, 58, 59, 65, 68, 69-71, 76, 82, 104]. One important feature of the flow o f foams is that it involves not only yield stress, but also the change of the number o f lamellae and their properties. This dynamic effect introduces additional difficulty. 5 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. This study focuses on pore-network modeling of flow and displacement of fluids with yield stress which leads to a more comprehensive understanding of the process in porous media. Motivated by the importance o f testing the various theories for foam generation at the pore-network level, we also modeled different foam generation mechanisms under various scenarios flow conditions. The dissertation is organized as follows: In this chapter, after the introduction, a brief review of the literature pertinent to flow and mobilization of fluids with yield stress and foam generation and mobilization is presented. In Chapter 2, an extended IPM algorithm is introduced to incorporate the effects introduced by flow of fluids with yield stress, especially Bingham plastics, in porous media. The algorithm involves updating the yield stress of each bond belonging to flowing paths and the iterative application of the IPM method. A macroscopic relation between applied pressure gradient and flow rate is obtained, which is found to be quadratic at small pressure gradient and linear at a high pressure gradient. The microstructure of the porous media is found to affect significantly the flow patterns of the mobilization. Wide distributions of the radii of the bonds result in tortuous flowing paths and narrow distribution giving clearly layered patterns. The “relative” permeability is found to have a cubic relationship with the fraction of sites belonging to open paths. A study o f the partial occupancy case suggests a simple upscaling method describing the two-phase flow involving fluids with yield stress. Application of the algorithm to displacement processes leads to interesting flow patterns corresponding to different system parameters. 2-D Patterns are obtained to identify 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. the flow regimes for two-phase displacement involving a yield-stress fluid, which is different from the Newtonian case proposed by Lenormand [79]. In Chapter 3, the method developed in Chapter 2 is applied to the study of the generation and mobilization of foams in porous media. Foams show yield stress when they are mobilized, which is caused by the temporary increase o f surface area necessary to push bubbles past each other in shear flow and is inversely related to bubble size [101]. To model foam formation and propagation in porous media, we should consider two parts, the flow of two immiscible fluids, where capillary and viscous forces are important, and the flow of a fluid exhibiting a yield stress, for example of a Bingham plastic. The dynamic change of the properties o f foams in the process is a very important factor to be considered. Two-phase flow has been extensively studied as reviewed recently by Blunt and Jackson [16]. We will introduce a simple method to simulate simultaneous two-phase flow, with more details shown in Appendix B. To account for the effect of yield stress, the dynamic algorithm is applied when necessary. In experimental works, a commonly used method is to inject the gas phase into a liquid-filled medium (e.g. Ref. [96]). Then, a foam phase will be generated when the gas phase invades and propagates further inside. The first case studied in this chapter is lamellae mobilization and generation in drainage displacement. Lamellae may be generated by snap off, leave-behind and lamella division mechanisms. Another case is lamella creation due to capillary-pressure fluctuations. The creation of lamellae is affected by liquid pressure gradients, the capillary pressure and the 7 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. snap off probability. The last one is that Roof snap-off occurs repeatedly at specific, randomly distributed germination sites. Foam properties were found to be fluctuating, in some cases with a well-defined period and in others in an apparently chaotic manner. The results show that strong foam generation by lamella mobilization and division and capillary fluctuations is possible without the necessity of repetitive Roof snap-off. Minimum pressure gradient calculations lead to results qualitatively similar to the theory of Rossen and Gauglitz [104], without their assumption o f a uniform pressure gradient in the pore network. In the last chapter, a continuum model for the flow of Bingham plastics is used and an asymptotic solution is given to clarify a special case where Bingham fluid is flowing around a circle. The yield surface is determined. 1.2 Literature Review First, a brief description of rheological aspects of fluids with yield stress and foams is presented. Then, the special properties of foams in porous media will be introduced. They include relevant research works on foam generation, foam texture, foam stability, and the flow and mobilization o f foams. Because of the large amount of literature available and also the limited space in this thesis, we do not intend to give a complete review of these topics. Only subjects that are closely related to the work presented will be discussed. 8 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 1.2.1 Rheological Background For an incompressible Newtonian fluid, the shear stress is proportional to the shear rate where the constant is the viscosity of the fluid. This relationship is a straight line through the origin in a plot of shear stress and shear rate. One simple definition of a non-Newtonian fluid is that for which the curve of shear stress versus shear rate is non-linear and/or does not pass through the origin. Then, the apparent viscosity of the fluid depends on the flow conditions. Generally, non-Newtonian fluids can be grouped into three main categories. One is time independent: for this kind of fluids the rate of shear at any point is only determined by the shear stress at that point and depends on nothing else. For more complex systems, the rate of shear can change with the duration of shearing and their kinematic history. And the last group is visco-elastic fluids which have both the characteristics of fluids and solids and show partial elastic recovery after deformation. In this thesis, we will not consider that last two kinds of complex fluids. Figure 1.1 shows the qualitative flow behavior of different time-independent fluids. Here we follow the categorizing method by Tanner [122]. For this kind of fluids, the behavior is described by a relation of the form Y = /(r ) or its reverse t = W ) 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Depending on the format of the function/ or/ / the fluids are categorized as three types: (a) Bingham plastic (and related); (b) Pseudo-plastic fluids; (c) Dilatant fluids. T / T T Figure 1.1: Shear stress and shearing rate relationship for (a) Bingham body; the slope P gives the apparent viscosity, (b) Pseudo-plastic or shear-thinning fluid, (c) Dilatant or shear-thickening fluid. Dashed lines show Newtonian behavior (After Ref. [122]). 10 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Pseudo-plastic or shear-thinning fluids are characterized by an apparent viscosity which decreases with increasing shear rate. Most pseudo-plastic polymers and melts exhibit Newtonian behavior at very high and very low shear rates. Dilatant or shear-thickening fluids show no yield stress like shear-thinning fluids, but their apparent viscosity increases with increasing shear rate. Both (b) and (c) can be described by power-law equations, with the power less than unity for shear- thinning fluids and the power greater than unity for shear-thickening fluids. Viscoplastic fluid behavior is characterized by a yield stress value. The fluid deforms or flows only after the stress becomes larger than the yield stress, and the flow curve may be linear or non-linear but will not pass the origin. One explanation to this phenomenon is that there is a rigid 3-D structure in the substance which will be broken down after the stress is high enough. There are different models to describe these behaviors. One simple model is the Bingham model, which describes the fluid behavior with a linear flow curve after the applied stress exceeds the yield stress. r = MBf for H > y = 0 for Id < Other popular equations to describe liquids with yield stress include the Casson model [81], which has often been used for describing the steady shear stress of many food stuff and biological materials 11 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. and the Herschel-Bulkley [51] model, sometimes also called the generalized Bingham model r = Tq + k(y)n for Id > \ t * I y - 0 for |r| < where the dimensions of k depend on the value o f n. The Herschel-Bulkley model is a combination of a Bingham model and a power-law model and has a better fit for some experimental data. Foam flow is a complex process especially in porous media. Because the properties o f foams are not steady, with repeated break and reform events, the size and shape of the foams change dynamically, as a result of the interaction of foams with the pore structure. At low Reynolds numbers, there are at least four microstructures of foams identified [71]. When the individual bubble size is smaller that the characteristic size of the pore space, the foam is called bulk foam. Bulk foam can be further divided into ball foams and polyhedral foams depending on the gas fraction of the foams. In the opposite case, where the characteristic size of the space is o f the same order or smaller than the bubble size, foams are called confined foams, and can be in two different quality states depending on the fraction of gas flow. At high fraction of gas flow, the bubbles are close to each other and in direct contact. This lamellae regime was identified by Hirasaki et al. [54], It is 12 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. different from the state at low fraction of gas flow, where the bubbles are far away from each other. The rheological behavior of foams in porous media depends on the liquid-gas ratio [94], Usually, foam displacement applications in porous media have low liquid-gas ratio and foams are in the form of polyhedral foams and individual lamellae. These foams are often considered as fluids with yield stress (see, for example, Ref. [104]). Different models [20] have been used to describe foam flow, including Bingham plastic models and Herschel-Bulkley models. 1.2.2 Foams in Porous Media Foams appear widely in petroleum industry, in both undesirable and desirable states. Significant understanding of foams in porous media has been achieved in the past several decades. 1.2.2.1 Foam Texture Foams in porous media are the same as ordinary foams in composition. They are a complex system involving gas phase and a liquid phase with surfactant. Surfactants are needed to stabilize the film of the foams. Their properties and concentration can largely affect the structure and stability of foams. The foaming agent tends to decrease the interfacial tension and to increase the interfacial viscosity which 13 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. provides a mechanical resistance to film thinning and rupturing. The choice of surfactant or the solubility of the surfactant and the concentration of the surfactant are important factors in foam-related processes. Apaydin et al. [2] found that decreasing surfactant concentration can lead to weaker foams, lower displacement efficiency and higher gas mobility. They stated that the concentration of surfactant affected the coalescence forces, which determine the bubble size, foam lifetime and foam texture. A similar observation was made by Owete [89] in heterogeneous porous media, where foams of various sizes were generated: small ones at high surfactant concentration and large ones at lower concentrations. Gauglitz et al. [45] found that the surfactants increased the time needed for bubbles to breakup. Ball foams and polyhedral foams are bulk foams, which are usually encountered in everyday life. The foam structure in porous media can be different from bulk foams [35]. There are observations [71] that the size of foams is of the same order as the pore size, which means that foams in porous media may not be bulk foams and be more close to the direct lamellae described by Hirasaki et al. [54]. A direct observation of foams spanning the porous media is shown in Figure 1.2. The size of foams in porous media is affected by more factors than traditional bulk foams, which are determined primarily by surfactant type and its concentration. Variables that may affect the texture of foams in porous media include the characteristics of the porous media, such as pore size and shape, permeability, and capillary pressure, and the velocities of gas and liquid [68], besides surfactant and its concentration. Foam texture in porous media is also a direct result of the competition of lamella 14 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. generation, trapping, mobilization and destruction [43]. Ettinger et al. [34] compared the effluent foam size distribution for pre-generated foams with that for foams generated in situ, under the same flow conditions. They found that the two are essentially the same and reached the conclusion that the porous medium shaped the foams to its “liking” no matter what type of foam was supplied. The experiments by Owete et al. [89] also showed that in heterogeneous porous medium, pore geometry determined the shape of large foams. A 2-D network model was used to simulate the coarsening of foam in porous media [29]. 100 [im Figure 1.2: Photomicrograph of foams in a transparent glass micro-model (After Ref. [71]). 15 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 1.2.2.2 Foam Stability For bulk foams, the stability of foams is determined both by bulk solution and interfacial properties, including gravity drainage, capillary suction, surface elasticity, viscosity (bulk and surface), electric double-layer repulsion, dispersion force attraction and steric repulsion [114]. Similarly, for foams in porous media, there are various factors affecting foam stability and foam texture. One variable studied extensively is the disjoining pressure (see, for example, Ref. [3]). The disjoining pressure is a function of the film thickness. Its value reflects the property of the film forces, with negative values indicating net attractive forces and positive values net repulsive forces. Figure 1.3 shows an ideal disjoining pressure and film thickness, while Figure 1.4 shows a corresponding experimental result. The capillary pressure is connected to the disjoining pressure through the augmented Young-Laplace equation, thus it governs the thickness of a lamella at equilibrium and at rest [57]. When capillary pressure rises, the disjoining pressure increases and the thickness o f the film decreases. Figure 1.4 shows that after the disjoining pressure reaches IJm ax, the film thickness is 4nm. An even higher capillary pressure will lead to the rupture of that film. For part of the film without the inner branch, the foam film will collapse when the disjoining pressure reaches T 7m ax. When the thickness of the film is close to hm a x, the foam is not stable to mechanical disturbance and easy to break. 16 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. + NBF r c* 'm i ci CBF 0 Film Thickness, h Figure 1.3: The ideal relationship of disjoining pressure 77 and the thickness of the film h (After Ref. [3]). 120 ------------3r - r ---------------- 1 -------------------1 --------' ---------- * A / ! / ! « | | AM A 100 A i A w 0* 00 1 M S u S m o[ 0.18 M NaCI ~ 8 0 I i o t & 60 ts 1I . 40 A 20 - ¥ 0 i i 10 h(nm) Figure 1.4: Experimental results (After Ref. [7]). 15 20 17 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. In a flowing system, the coalescence of foams is complex. Khatib et al. [66] studied the effect of capillary pressure on coalescence of foams flowing through porous media. They stated that the because of flow of foams, even when a limiting capillary pressure is reached, not all of the lamellae can break at the same time. The obvious effect is that the quality of foams becomes coarser. The limiting pressure of foams in porous media should change with gas velocity and the permeability of the medium, in addition to the surfactant and its concentration. Through experiments, a limiting capillary pressure is found to exist during foam displacement, with increasing gas fractional flow and a constant gas flow rate. Significant foam coalescence was marked by an abrupt drop of the capillary pressure and the change of foam texture to coarse foams. Khatib et al. [66] also compared the limiting capillary pressure with that of a bulk foam system. Their results showed that the limiting capillary pressure also depends on the electrolyte concentration. Aronson et al. [3] studied the influence o f disjoining pressure on foam stability. They measured the disjoining pressure isotherm for a single foam film for two different SDS concentrations and they found that the effect of addition of salt on rupture pressure depended on surfactant concentration. One interesting observation is that the limiting pressure for rapid foam coalescence in porous media was close to the rupture pressure obtained from the isotherm they measured. They reached the conclusion that high repulsive disjoining pressure in single foam films resulted in strong foam in porous media which demonstrated high flow resistance. Ettinger et al. [34] identified a regime, called the limiting-capillary-pressure 18 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. coalescence, where the coalescence and generation rates were almost the same. Another mechanism for foam coalescence is gas diffusion. For stationary foams, smaller bubbles have higher pressure than larger ones and gas diffuses through the lamellae into the large ones around it. The result of diffusion is that the small bubbles ultimately disappear [101]. 1.2.2.3 Foam Generation For foam applications in EOR processes, foams are generated before or during the process. Pre-generated foams are also used in some experimental research works [36], Foams generated by a foam generator are bulk foams, and they are usually in ball or polyhedral forms. After injection into experimental apparatus, like a bead pack, the shape of the foams may change [34], Foams generated in situ are created mainly by three mechanisms, snap-off, leave-behind and lamella division. The description of the three mechanisms can be found in Ref. [23]. These mechanisms were identified by Ransohoff et al. [96]. Roof [97] proposed the snap-off mechanism in a study of oil droplets in water-wet porous media. Ransohoff et al. [96] recorded visually the primary mechanisms of foam generation, i.e. snap-off, leave-behind and lamella division, in bead packs. They found that snap-off and lamella division could be dominant mechanisms depending on the capillary numbers. They also introduced the concept of germination sites, namely those which meet the static, geometric criterion for snap 19 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. off. They concluded that the importance of the snap-off mechanism depends on the number of germination sites. Chambers et al. [21] stated that for snap-off or germination sites in porous media, the size of the pore should be at least two times large as the throat size, and the slope should not be very sharp. One necessary condition is that there should be sufficient amount of liquid around the throat, so that it can flow back and snap off. In the experiments by Owete et al. [89] and Gauglitz et al. [45], snap-off mechanisms were confirmed as well. There are several studies [26, 72, 112, 133] on the calculation of times related to the snap-off process, including the times required for liquid flow and lamella moving. Kovscek and Radke [72] used a pore-level model to study foam generation by snap-off in a constricted, cornered pore. They analyzed the time to accumulate liquid and form an aqueous lens and the time to displace that lens so that snap-off can happen again. They found that the liquid continued to accumulate even when the capillary pressure was high, but they also stated that there was a static criterion that the pore- throat to pore-body ratio must be small enough. Kovscek et al. [73] did calculations in cornered, constricted pores, and focused on the initial volume of liquid phase in the comers of a pore and then on the possibility o f a minimum liquid content below which snap-off is inhibited. They found that smoothly constricted pores with relatively small throats connected to large pore bodies exhibit excellent ability to snap-off regardless of their initial liquid content. And the presence of this kind of stmcture can lead to strong foam. Rossen [102], on the other hand, argued that conventional theoretical methods based on single capillary experiments is not good 20 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. for the complex geometric conditions in porous media. He also stated that if repeated snap-off did happen in porous media, its frequency depended on some factors which were not accounted in previous studies. Leave-behind mechanism is another important mechanism for foam generation. Because of the position of the lamellae generated by leave-behind, this kind of lamellae cannot block the flow of gas. No separate bubbles are formed by leave- behind [96] and these lamellae create more dead ends and block gas flow paths. Lamella-division has the special property that a lamella should exist before this division can happen. This kind of lamellae cannot be generated from scratch. Another important condition is that the pressure should be high enough so that the lamella can be mobilized. The existence of foams can alter the flow of gas phase, thus two different cases are identified. One is that foams block pathways and change the flow of gas phase, but gas phase is still continuous. The other is that all of the pathways are blocked by foams, in which case the gas phase becomes discontinuous. In these two cases, the liquid phase is always occupying small pores and comers of the pores. The first case is possibly corresponding to the situation with small liquid saturation, so that foam generation decreases and continuous gas channels develop. There are arguments on the relative importance of these mechanisms, especially snap-off and lamella-division. In leave-behind mechanism, the lamellae are parallel to the flow direction of the gas. Even though it can happen very frequently, this generates relatively weak foam. Snap-off can occur repeatedly at one site and 21 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. influences the flow properties of the gas phase by creating lamellae and blocking gas flow paths. Hence, it generates strong foam [96]. Kovscek et al. [71] asserted that snap-off is the dominant foam generation mechanism. Based on the shear- thinning behavior of foam observed in experiments, they believed that lamellae division is not significant. Rossen [103] critically reviewed the evidences of Roof snap-off mechanism in experimental observations and concluded that there was no compelling evidence strongly supporting repeated Roof snap-off as the basis of steady-state foam generation. Tanzil et al. [123] also stated that lamella division could also lead to strong foam in porous media. 1.2.2.4 Foam Mobility Gas mobility is affected by the presence of foams, and is directly related to foam texture. The effect of foams in porous media on gas mobility can be conceptually divided into two parts. One is on the change of permeability, the other is viscosity. One observation [42, 56, 110] is that the presence of foams will not affect the relative permeability of liquid phase. The reason is that the liquid phase always chooses the smallest pores and occupies almost the same pores it would without foam. Thus, the relative-permeability and viscosity of liquid are not affected by foam. Flumerfelt and Prieditis [40] conducted gas-only injection experiments and found that the permeability of the bead pack to gas was 2 orders o f magnitude less than in the case without foam. They also concluded that there were 100 times less 22 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. channels present in the foam case than in the foam-free case. The reduction of permeability of gas can be related to gas trapping. Gas trapping is found in experimental studies [89], with tracer studies revealing that 80% of gas may be trapped at any time during flow [47], The gas trapped may be mobilized if the pressure imposed is large enough and there are channels with flowing lamellae. But flow in these channels is intermittent, because they can be plugged again by lamellae flowing [29], Hazlett and Furr [50] invoked a percolation model for the permeability reduction in porous media by continuous-gas foams. They used a foam efficiency parameter to model the complex process o f in situ foam generation and propagation. The flow resistance from the flowing bubble trains can be viewed as a change in gas viscosity. Hirasaki and Lawson [54] studied the apparent viscosity in smooth capillaries. They found that foam texture is very important on the mobility and that the apparent viscosity is a function of the equivalent bubble radius. They also concluded that the apparent viscosity is affected by three factors. One is from the slugs of liquid between bubbles, the second is the resistance to deforming the interface when a bubble passes through a capillary and the last is the surface tension gradient caused by the surfactant concentration difference between the front and the back o f the bubble. Falls et al. [36] obtained apparent viscosity from the measured mobility of gas and extended Hirasaki and Lawson’s theory to account for the capillary pressure imposed by the porous medium and constricted flow paths. Their results show that the apparent gas viscosity depends strongly on foam 23 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. size and foams of uniform texture are pseudoplastic. Xu and Rossen [132] predicted effective viscosity of foams in periodically constricted tubes using a dynamics model. Gas mobility is closely related to foam generation. As more foam is generated, the more flow paths are blocked and the apparent viscosity o f foams from moving bubble trains will change too. So, one observation of foam generation in laboratory usually is the sudden reduction of gas mobility [59]. According to the theory of Ransohoff et al. [96], there is a critical capillary number above which strong foams can form. Through experiments, Gauglitz et al. [43] found that there may be a minimum pressure gradient for foam generation in porous media, which varies with permeability. The onset of foam generation usually corresponds to the onset of mobilization [123]. Kovscek and Bertin [68] employed conventional and percolation network scaling ideas to describe the mobility in heterogeneous porous media. These factors, including permeability, porosity, capillary pressure, limiting capillary for foam collapse and liquid and gas velocity, can be combined to describe gas mobility upon assumption of local equilibrium of foam generation and coalescence rates. Their experiments showed that in-situ foam generation is effective in controlling gas mobility and diverting injected fluid to low permeability areas in heterogeneous porous media, and foamed gas can be more mobile in lower permeability porous media. Falls et al. [35] developed a mechanistic simulator in which the conservation equations are coupled with balances on the densities of 2 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. flowing and stationary bubbles in the foam. Foam is generated in situ by capillary snap-off. Another kind of model used is based on percolation [106]. Flow of foams in porous media was also studied in single capillaries (for example, Ref. [85] and Ref. [54]). Nguyen et al. [85] investigated the motion of foam films in axisymmetric diverging-converging channels, taking into account surface viscosity and elasticity. They postulated that because of the wide presence of diverging-converging channels in oil reservoirs, the large resistance observed may be largely due to the surface elasticity and viscosity o f the films. Population balances [39, 58, 59, 69] are widely used in modeling flow of foams, in addition to various experimental studies [39, 74, 75]. Foam in bulk has a yield stress [75, 93] and experiments [74] demonstrated that the foam in fracture could form bulk foam too. When the foam is strong in a porous media, the gas flowing paths are blocked and gas becomes discontinuous. For flow of this kind o f foam, there is a minimum pressure gradient below which gas is blocked, as though the gas phase were a homogeneous fluid with a yield stress [36, 99, 101, 123]. Different methods were used to find the minimum pressure gradient. Kharabaf et al. [65] determined it with Invasion Percolation with Memory. Rossen et al. [104] used percolation model and proposed the model for foam generation on Bethe tree network. Rossen [99] derived the minimum pressure gradient along a bubble train by analyzing in detail the passage of a single lamella across a single pore. 25 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Chapter 2 FLOW AND DISPLACEMENT OF FLUIDS WITH YIELD STRESS 2.1 Introduction The flow and displacement of non-Newtonian fluids with yield-stress rheology Bingham fluids) is an important subject, encountered in a variety of industrial applications. The particular rheology describes solid-like behavior, if the applied stress is below a yield value, and liquid-like behavior, when the yield value is exceeded. Typical Bingham fluids include paint, slurries, pastes, and food substances, as well as fluids, such as heavy oils and foams in porous media. The latter will be o f particular interest to this chapter. Heavy oil rheology has often been approximated as that o f a Bingham fluid [4, 6]. Likewise, foams in porous media effectively display yield-stress rheology. Foams are commonly used in underground oil reservoirs to improve reservoir displacement (sweep) efficiency, 2 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. block swept channels, and in gas storage and acidizing operations. The state of the art in the field has been summarized by Rossen [101]. The majority of studies in the flow and displacement o f Bingham fluids in porous media are phenomenological. Typically [4, 6], a continuum description is used, with the Bingham fluid assigned an effective yield stress, the subsequent flow being approximated as that in an equivalent capillary. Wu et al. [130, 131] used such an approach to investigate single-phase flow of a Bingham fluid in porous media, and to study by analytical means the 1-D displacement of a Bingham by a Newtonian fluid. While useful in extending our phenomenological understanding, however, these approaches cannot provide an account o f the effect of the microstructure, generally relying on ad hoc effective parameters. As a result, the state of the art in the area is incomplete. Important issues in need of clarification include how to upscale single-capillary expressions to a representative element of the porous medium, understanding the dependence of the effective yield stress and the flow parameters on the microstructure, and identifying and upscaling the displacement patterns in the case of two-phase displacements. Steps in this direction have already been taken by previous authors, mostly in the context of foam flow in porous media. Rossen and Gauglitz [104] developed an approximate model for foam mobilization, separately addressing continuous- and discontinuous-gas regimes. They employed percolation theory concepts to demonstrate that a minimum pressure gradient (equivalently, an effective yield stress or macroscopic threshold) exists, because it is necessary to displace foam 27 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. lamellae (liquid films) from pore throats. Rossen [98, 99] and Cox et al. [31] addressed the relationship between pore geometry and the macroscopic threshold for foam flow on a given pathway, without addressing its specific determination on the pore-network scale, however. Clearly, the macroscopic threshold (minimum pressure gradient) directly follows from the geometry of the path, along which mobilization first occurs, and the two are intertwined. In the following, for consistency in terminology, we will denote this path as the MTP (minimum threshold path). In a network with randomly distributed thresholds, the MTP is the connected (through nearest neighbors) path between two given boundaries (or points), along which the sum of thresholds is the minimum possible. Fundamental to this concept is the notion that specific, local thresholds must be exceeded across a given pore throat, for the fluid therein to be mobilized, and that these thresholds are distributed in the network. Rossen and Mamun [105] defined a simplified model, they called the “simple percolation model”, in which the MTP is approximated as the subset of the percolation cluster [38] that randomly samples those bonds (throats) with the smallest individual thresholds. Consider a 2D square network, in which the percolation cluster consists of bonds with threshold values in the lower 50% of the distribution. Then, the Rossen and Mamun model postulates that the minimum pressure gradient is the arithmetic average of the lower 50% of the distribution. The authors went on to show that the actual macroscopic threshold in a Bethe (or Cayley) tree (which is a lattice without reconnection) is, in fact, lower than this 28 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. estimate. The MTP occasionally samples bonds with a higher yield stress threshold than for bonds in the percolation cluster, when by doing so it can access a long cluster of low-threshold bonds. Both the “simple percolation model” and the numerical analysis on the Bethe tree exclude the concept o f geometric tortuosity, however - the fact that avoiding some high-threshold bonds could lengthen the path and raise the overall path threshold. Thus, the true threshold in a 2D or 3D network could be higher than estimated by Rossen and Mamun. One expects the geometric structure of the MTP to differ fundamentally from the backbone cluster in a conventional percolation process [38], because in minimizing overall thresholds the path has an incentive to minimize its length. Also, there are no dangling “deadend” clusters on the MTP, contrary to the percolation cluster in conventional percolation. Recognizing the need for an alternative approach, Kharabaf [62] and Kharabaf and Yortsos [64, 65] proposed a different variant of the percolation process, termed by the authors Invasion Percolation with Memory (IPM), to identify the MTP. Based on this approach, Kharabaf and Yortsos addressed static properties of various problems with yield stress, namely, the onset of the mobilization of a single-phase Bingham fluid in a porous medium [64], or the foam formation and propagation in porous media [65] in the absence o f flow effects. In addition, Kharabaf [62] studied displacement problems involving Bingham plastics, again in the quasi-static limit. Through IPM, these authors were able to elucidate various properties of the MTP, including its relation to the percolation cluster. However, their approach was limited, in that it 29 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. did not account for dynamic (viscous flow) effects, which will arise subsequent to mobilization. Viscous effects during the flow and displacement in porous media have been extensively analyzed as a function of the pore-network microstructure, when the fluids are Newtonian (e.g., see recent review by Blunt et al. [16]). Extensions to power-law fluids were considered by Sahimi [109], Shah and Yortsos [117], Tian and Yao [126] and Shah et al. [116]. Falls et al. [37] and Xu and Rossen [132] presented flow-rate pressure-gradient relations for foam flow in bundles of capillaries, but not in interconnected networks. Sahimi [109] derived an Effective Medium Approximation for Bingham plastic flow in a network, to find that in a certain range the macroscopic flow rate depends quadratically on the departure of the applied pressure difference from its minimum value. Such dependence was in fact first proposed by Roux and Herrmann [108]. The latter considered a square 2- D network of resistors, having the voltage-current relation o f Figure 2.1. The critical thresholds, below which flow does not take place, were randomly distributed on the lattice from a uniform distribution in the range between 0 and 1. Roux and Herrmann demonstrated the existence of an effective (macroscopic) threshold, F=(0.22±0.02)Z, linearly scaling with the lattice size L. They also found that the relation between current and applied voltage is generally non-linear, Q~(r-rJ (2.1) where, in the absence of finite-size effects, < 5= 2 for small and intermediate values of the current Q, and < 5 = 1 for large values of Q. A self-consistency argument was used 30 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. to explain the quadratic dependence near the onset of flow: Incremental increases in the driving force V-Vc have a combined effect, namely a proportional increase in the number of open paths, as well as a proportional increase of the flow rate in the paths already open. How the value of Vc might depend on the distribution of microscopic thresholds was not extensively discussed, however. V c Figure 2.1: The current-voltage relation of an individual resistor in the network problem considered by Roux and Herrmann [108], An analogous relationship for the flow rate- pressure drop dependence of a Bingham fluid in a pore throat will be used in the pore- network study in this chapter. 31 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. In these studies, the overall flow rate was determined without separately focusing on the effects of the trapping (non-mobilization) of a potentially large fraction o f the fluid and of the rheology of the mobilized fraction. Studies with gas- phase tracers in foam flow in porous media [42, 85, 95] have estimated that from 70% to 99% of the gas phase remains immobile (trapped), even at significant pressure gradients. Thus, attempts to model foam flow mechanistically [10, 24, 35, 71, 82] typically tend to separate the effect of gas trapping from the rheology of the gas (foam) that flows [10, 37], Various relative-permeability algorithms are then implemented to account for the effect of gas trapping [10, 24, 35, 71, 82]. At this point, it is not clear whether this division of gas mobility into effects of gas trapping and rheology is feasible or meaningful, however. It is one of the motivations of this work to assess the validity of such an approach. Another motivation is to understand the immiscible displacement of Bingham fluids in a porous medium. In addition to their generic interest, such displacements find important applications, for example in the recovery o f heavy oils with yield stress rheology [4], or in foam acid diversion [24, 84, 107]. In the latter, how the displacement goes forward defines the mobility of the acid and the extent of the diversion. Immiscible displacements involving Newtonian fluids have been investigated at great lengths in the past (e.g. see Lenormand [79], Yortsos et al. [134], and Blunt et al. [16] among a large number o f references). Lenormand [79] pioneered an approach that delineated and characterized three possible regimes that may develop during drainage (the displacement of a wetting by a non-wetting 32 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. fluid), namely capillary controlled, and stable and unstable viscous-controlled regimes. The state o f the art in non-Newtonian fluids is not nearly as advanced as with Newtonian fluids, however. Displacements involving power-law fluids at the pore-network level were addressed by Shah [115], who developed a phase diagram analogous to Lenormand’s (see also Shah et al. [116]). Displacements involving Bingham plastics were addressed by Kharabaf [62] who also proposed a phase diagram. However, these models did not include viscous effects, which arise following the onset of mobilization (but see also Kharabaf et al. [63]). Other studies in the literature involving non-Newtonian displacements were restricted to the solution of continuum equations, based on plausible empirical models (e.g. see [6, 130, 131]). The general objective of this chapter is to provide a comprehensive approach to the modeling of flow and displacement of Bingham plastics, or other yield-stress fluids, in porous media. For this purpose, we will extend the IPM approach of Kharabaf and Yortsos [64] to account for viscous flow effects, which were neglected in previous applications (e.g. in [65]). The new algorithm will then be used to explore, as a function of the pore microstructure, the macroscopic relations between the applied pressure drop and the resulting flow rate in single-phase flow, and patterns in the displacement of a Bingham by a Newtonian fluid in two-phase flow. The chapter is organized as follows: First, we present a brief review of the IPM algorithm and then extend it to incorporate dynamic flow effects. The new algorithm is subsequently applied to determine the macroscopic flow rate-pressure 33 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. drop relationship for single-phase flow in a pore network. The effect of various geometrical parameters on the macroscopic flow rate-pressure drop relation is discussed. Subsequently, we use the algorithm to study the displacement of a Bingham plastic by a Newtonian fluid, classify the patterns obtained and discuss their upscaling to macroscopic flows. 2.2 Algorithm for Flow of a Fluid with Yield Stress Consider the flow of a Bingham fluid in a porous medium, represented as a lattice or network of pores and throats. A pore-network representation has proved to be very effective in many applications of flow and transport in porous media [16]. The flow of a Bingham fluid in the bulk or in simple geometries is well-known. The bulk rheology, expressed as [12,14, 86] relates the deviatoric stress x to the rate-of-strain y through the yield stress xe and the plastic viscosity fxp. Its upscaling to flow in an infinitely long cylindrical capillary (subscript i) was obtained by Bird et al. [14], who derived the following relation t = Te - n py when x > xt y = 0 when x < x t (2.2) AP when r, > x. Qt = 0 when x t < r (2.3) 3 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. between the volumetric flow rate of the fluid Qi and the pressure gradient APi/l APR across the cylinder. Here, we defined the wall shear stress t, = — — where R; is 2/ the radius of the capillary. Using dimensionless expressions for the volumetric flow ^ A R t A P i A- R i — ~— , the pressure drop = ------- and the radius r. = —- j t R 2It , R is a reference radius, equation (2.3) further becomes < 7 , 1 = r} V - i ^ 3\Api \ri 3(Apiri) y Apt | when [ A p t |> — r qt = 0 when | Ap t |< — (2.4) ri In the above notation, we have alluded to the obvious fact that when flow occurs it is in the direction o f decreasing pressure. When the fluid is in a porous medium, the upscaling of rheology (2.2) to the smallest elements under consideration, namely the pore throats, is necessary. This problem is not trivial, and it is, in general, non-local. For further progress, a key approximation will be made, that the flow in any pore throat between two lattice nodes (pores) can be represented by expressions of the type (2.4). This approach is typical of lattice models. The underlying assumption is that the pressure drop-flow rate relation in a pore throat in the network is the same as if the geometric element is isolated, and subject to the same pressure drop. Furthermore, for the sake of simplicity, as well as the recognition that they were really derived for infinitely long capillaries, in the following we will use simpler versions o f equations (2.3) 35 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. and (2.4), in two different forms: In one model, we will assume that the conductivity o f each pore throat is the same for all throats, namely we will take | q t |=| Ap t | —— when | A.pt |> — rt r, qt - 0 when | Apt |< — (2.5) rt which is similar to the Roux and Herrmann [108] model (see Figure 2.1 and replace current with flow rate and voltage with pressure drop). In another, we will allow for a variable conductivity and take the model V i 1 = r f AP, I r, , 1 Apt | when | Ap t |> — r, qt = 0 when | A p i |< — (2.6) ri In both models, the dependence of the flow rate on the pressure difference is linear, upon the onset of flow (but remains overall non-linear). In addition, in postulating equation (2.6) we tacitly assumed a geometric similarity between the bond radius and its length, so that the effective conductance scales with the third, rather than the fourth, power of the radius. In all cases, radii will be distributed randomly from a given size distribution, a(r), therefore thresholds and conductances are also randomly distributed in the porous medium. 3 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. We will consider the mobilization and flow across two opposite boundaries of a lattice, subject to periodic lateral boundary conditions. Clearly, for flow of the Bingham fluid to occur, a minimum overall threshold must first be exceeded, (2-7) i * t where the summation is over all throats over a connected path. As noted, the latter problem was analyzed in detail by Kharabaf and Yortsos [64], who showed that there exists a specific path, the Minimum Threshold Path (MTP), over which the corresponding sum of thresholds is minimum. The authors developed the method of Invasion Percolation with Memory (IPM) to determine the MTP and its properties. Because we will make repeated use of this technique in the subsequent, we will first briefly review its main features. 2.2.1 Invasion Percolation with Memory (IPM) Consider a lattice consisting of sites and bonds (pores and throats). A dimensionless threshold r,y from a given distribution is randomly assigned on the bond ij connecting adjacent sites i and j. The first objective is to find the MTP, connecting any point o f one boundary (e.g. on the right in Figure 2.2) to any point of the opposite boundary (e.g. on the left in Figure 2.2), over which the sum of the thresholds is minimum. Invasion Percolation with Memory is a technique that allows the identification of this, as well as of higher-energy paths [64]. The 37 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. algorithm consists of performing a virtual invasion process, starting from the right boundary. The invading front invades one site at a time. Denote an arbitrary site currently on the front by F , its nearest neighbor by F ', the site from which the front will invade next by FG, and the site which will actually be invaded by FG. Each site that has been invaded is recursively assigned a value Vx (F) , using the following rule: At any time, and for any site F on the front, the sum S = Vx (F) + xFp. is formed. Then, the minimum such sum over all front sites F is found. The site corresponding to that minimum is identified as the growth site F0 , its neighbor as FG, and the assignment Vx (F0 ) = Vx (FG) + tfg f G is next made. Using the boundary condition Vx (R) = 0, for all sites R on the right boundary, determines recursively Vx (F) for all invaded sites. Then, the process is repeated (as shown schematically in steps (a) and (b) in Figure 2.2). The MTP is identified when the front first reaches the opposite (left) boundary for the first time. It is the path connecting adjacent sites whose difference in the values of 1 V is identical to the threshold of the bond separating the two sites (see Ref. [64] for a proof). Various properties of IPM were detailed in Ref. [64], In particular, the authors discussed geometrical properties of the MTP, showing that it is not a subset of the percolation cluster, or the minimum path of percolation, in general. They also explored the minimum value, Vl< n ia - Vc, 3 8 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. as a function of the threshold distribution. Its ratio with the lattice size denotes the minimum applied gradient for the onset of mobilization. Kharabaf and Yortsos [64] proceeded next to identify higher-energy paths, which will open at a higher value of the applied gradient Vx! L . This was done under static conditions, however, namely in the absence of viscous effects, which should occur following the onset of flow. Because in their calculations, flow in an open path did not affect the distribution of pressure, the identification of higher- energy paths was strictly a static (quasi-thermodynamic) process. In the case of Bingham fluids, this would correspond to a vanishing plastic viscosity jup. When viscous flow effects are important the algorithm for opening new paths must be modified. Now, the applied gradient to open a new path will not only be determined by the pre-assigned throat threshold — , but also by the velocity field, ri which affects the pressure distribution along the mobilized paths. Thus, even though the first path to be mobilized is indeed the MTP, identified using IPM, finding the subsequent fraction of the mobilized fluid requires a more elaborate approach. 3 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. FF AB (a) (b) Figure 2.2: Description of the invasion rules in IPM, before (a), and after (b) a growth step. The growth site FG is that for which the sum S = Fj(.F)+ tff, is minimum. F denotes a site currently on the front, F ' a site which is a nearest neighbor of F. Function, Vl (F)is recursively determined for all invaded sites, as described in the text. Invasion is from right to left (After Ref. [62]). 40 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 2.2.2 Extension o f IPM to Account for Dynamic Effects In general, after the onset of mobilization, an incrementally higher pressure is needed to open a new path, with only one path opened at a sufficiently small incremental pressure. Such path may emanate from an already mobilized path (being a branch of or a loop connected to it, as for example, path CD in Figure 2.3) or may be a completely new path. Because of the small increment in pressure, the flow rates in its pore throats will be infinitesimally small, thus not contributing any viscous flow effects in that path. For example, the pressure drop across a bond j in path CD of Figure 2.3 will still be equal to the corresponding yield stress (taken above as — ). However, over any bond k in the already mobilized path (e.g. path ri AC of Figure 2.3), the pressure drop, hence the effective threshold, would be different and equal to the sum of the yield stress and of the pressure drop due to viscous flow. Therefore, in order to find the next path to be open, one must apply IPM, except that the threshold along the open, flowing part must be updated. For example, the updated dimensionless threshold in the flowing path is now * k = 9 k + — (2-8) rk Thus, the MTP at the new applied incremental pressure, and the new paths to open, are determined using IPM, in which all thresholds in the flowing part are updated following (2.8), with all others having their pre-assigned values — . This process rj 41 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. works incrementally, by determining minimum-energy (threshold) paths successively, at small increments of the applied pressure drop. In doing so, paths determined at the prior step always remain minimum paths at the next increment. The smallest increment required to open only one new path may require an iterative process, which converges to a unique solution, however. The approach works well whether or not the new path is connected to a previously open path. To determine the pressure distribution due to viscous flow, we assume that the fluid is incompressible, equations (2.5) or (2.6) (or other appropriate local expression) apply for the momentum transfer, and write a mass balance for any site i along the flow path i > » = ° < 2-9> j where q{ j is the volumetric flow rate from site i to its neighbor site j , and Z is the coordination number o f the lattice (e.g. Z = 4 for a square or Z = 6 for a cubic lattice). The resulting equations for the pressure were solved using a combination of Conjugate Gradient and Successive Over-Relaxation methods. Because of the inherent non-linearities of the problem and the need for iterations, the algorithm is time-consuming. This limited the networks studied to relatively small sizes, e.g. 50x50. For this reason, some finite-size effects may exist in the results obtained below. 4 2 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0 B D Figure 2.3: The process of the opening of a new path after the applied pressure difference (P) is incrementally increased, so that only one new path opens. Curve ACB denotes an already open path. For an incrementally higher pressure, an additional path (curve ACD) opens. Flere it is assumed to be connected to the MTP. Flow is from right to left. We remark that in a site belonging to an open path and adjacent to bonds which may not be open (not yielding), viscous forces may make additional contributions to the pressure drop from that we assumed. Such were neglected in the present algorithm. We do not anticipate that this error would be significant, however. The next sections will present simulations of both single-phase and two-phase flow. In one set of simulations, we used model (2.5) and distributed the yield-stress 43 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. thresholds following the distribution of Kharabaf and Yortsos [64], for the sake of consistency. Namely, first, we assumed that r = /?", where /? is a random number, uniformly distributed in (0, 1), and where n> 0. The corresponding threshold probability density function (pdf), mean and standard deviation are then given by f(r) = I r (V»H 5 (r ) =r - J _ , = -----” (2.10) n 1 + n [l + n y jl + 2n respectively. With this model, the system is homogeneous with respect to flow conductances, but heterogeneous with respect to thresholds, the exponent n parametrizing the heterogeneity: For n=\, the distribution is uniform; for n -> 0, it approaches a delta function, peaked at 1, and corresponds to a homogeneous case); for n — > oo, it is peaked at 0 and corresponds to the most heterogeneous case. The results obtained were then normalized such that they correspond to the case where the mean threshold is always equal to one (see below). In another set of simulations, we also considered a Rayleigh-type pore-size distribution with a mean equal to one, a (r ) = ~ r exp f 2 \ nr hence / ( r ) = ~ j exP I t r \ n 4 r 2 (2.11) 4 4 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where relationship r tj = — was used. A plot of the corresponding threshold distribution is shown in Figure 2.4. The results were also normalized to correspond to a mean threshold of one. .0 0.8 0.6 0.4 0.2 0.0 2 4 6 8 10 12 14 0 Figure 2.4: The PDF of the yield-stress threshold for a pore throat corresponding to the Rayleigh pore throat-size distribution in equation (2.11). 4 5 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.5: The patterns of the open paths for the dynamic case, corresponding from left- to-right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three- fourths of the sites belonging to the open paths. Full-occupancy case, uniformly distributed thresholds (n=1) and viscous model (2.5). A more layered structure of the mobilized paths is apparent in the dynamic case. 4 6 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.6: The patterns of the open paths for the static case, corresponding from left-to- right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three-fourths of the sites belonging to the open paths. Full-occupancy case, uniformly distributed thresholds (n= 1) and viscous model (2.5). A more layered structure of the mobilized paths is apparent in the dynamic case. 2.3 Single-Phase Flow In this section we present results for the onset of mobilization and flow of a single phase, for two qualitatively different cases. The first (referred to in the following as the full-occupancy case) corresponds to the lattice being fully occupied by the 47 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Bingham fluid, the other (the partial occupancy case) to the fluid occupying only a fraction of the lattice. First, we will take all thresholds distributed according to (2.10) and use the viscous flow models of the type (2.5) to examine the effect of the threshold distribution at constant flow conductance. A comparison with the “thermodynamic” approach of [64] will be made. Then, the effect of a distributed flow conductance (model (2.6)), as well as of the Rayleigh model (2.11), will be considered. In addition to probing the microscale details, approximate macroscopic expressions will be sought, to upscale the results found to useful macroscopic expressions. 2.3.1 The Full-Occupancy Case Consider, first, the case of uniformly distributed thresholds (n=1) using model (2.5), which, as noted above, also corresponds to Roux and Herrmann [108]. Unless indicated otherwise, the results presented will correspond to a 50x50 lattice. Figure 2.5 shows mobilization patterns at four different stages, corresponding to the MTP at the onset of mobilization, and the cases where about l A, V 2 or 3 /4 of the lattice sites belong to mobilized paths. For comparison, the static patterns corresponding to [64] are also shown in Figure 2.6. The first open path (the MTP) is identical in both the static and the dynamic cases, as expected. Following the onset o f mobilization, the patterns o f open paths are quite different, however, with a more-layered structure apparent in the dynamic flow case. Because of the viscous friction associated with 48 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. flow in the open paths, the effective pressure drop across a throat transverse to the main flow direction is reduced, thus making more difficult its opening in the dynamic case compared to the static case. Figure 2.7A shows the fraction of the sites belonging to the open paths, X site, as a function of the dimensionless applied pressure gradient. In foam flow in porous media, this plot would correspond to the flowing gas fraction as a function of the applied pressure gradient. As noted above, for comparison purposes, the axes were rescaled to correspond to a mean value of the dimensionless threshold equal to one (in the simulations with n=1, the mean threshold is 0.5, e.g. see (2.10)). Clearly, the fraction o f sites on the mobilized paths is substantially smaller in the dynamic compared to the static case. In the former case, the additional viscous friction hinders the opening of new paths and fluid mobilization. In fact, mobilizing the fluid in all sites requires the doubling of the threshold value, compared to only a slight increase of the applied gradient beyond the minimum, in the static case. Even more difficult is the mobilization of fluid in pore throats (not shown here for lack of space). The open fraction of the latter increases much more slowly with the pressure gradient, as a result of the anisotropy in the mobilization alluded to in Figure 2.6. In the dynamic case, a bond between two rows transverse to the main flow direction can open only as a result of viscous pressure drop along the main direction. For a fairly uniform threshold distribution, such pressure difference will be provided only at very high flow rates. This is not so in the static case, where viscous effects are not relevant, and thresholds open sequentially. As a result, the 49 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. fraction of throats mobilized at a given pressure difference is substantially smaller in the dynamic case. In fact, in the latter case, the approach to the mobilization of the entire lattice (sites and bonds) is very slowly reached at very large values of the applied pressure gradients. It is instructive to compare Figure 2.7A with a typical drainage capillary- pressure curve in a porous medium [33]. In such a comparison, the applied pressure gradient would stand for the capillary pressure and the fraction of sites on open paths for the saturation (volume fraction) of the non-wetting phase. Viewed from such a perspective, the two curves are quite similar, both displaying macroscopic thresholds and the effect of heterogeneity (see also below). O f course, there exist obvious differences, for instance, in the scaling of the properties near the threshold. Nonetheless, the analogy is useful for qualitative purposes. Another, simple, way to understand the dependence in Figure 2.7A is to consider flow of a Bingham plastic in a bundle of parallel capillaries. As shown in the Appendix, in such a case the curve coincides with the cumulative density function of the threshold distribution, which for the specific problem here is a straight line, quite different than the curve on Figure 2.1 A. As with the capillary-pressure curve, the bundle of parallel capillaries model is at best a qualitative approximation. Nonetheless, it serves to show the effect of heterogeneity in the threshold distribution, with more heterogeneous distributions resulting into more spreading curves (see also below). 5 0 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.8 0 5 '55 X 0.4 0.2 — Dynamic IPM - Static IPM 0.0 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 R*|VP|/2t6 (A) 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0.0 0.0 0.4 0.8 1.2 1.6 2.0 R*|VP|/2xe (B) Figure 2.7: The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient, for the dynamic (solid line) and the static (dashed line) cases (A) and the dimensionless flow velocity-pressure gradient relationship. Full occupancy case, uniformly distributed thresholds («= 1) and viscous model (5). Note the quadratic dependence near the onset of mobilization (region A) and the approach to a linear relationship at higher flow rates (point B). 51 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. The relation between the dimensionless velocity and the dimensionless applied pressure gradient is shown in Figure 2.7B, where $ denotes the porosity of the porous medium. The minimum pressure gradient for the onset of flow is as predicted by the IPM algorithm [64], For the specific example, we find a minimum value of 2x0.296, compared to 2x0.22 in Roux and Herrmann [108] (recall the factor of two rescaling to normalize the mean thresholds). However, in the latter study, the lattice was rotated by 45 degrees, which may explain the difference observed (0.22 • V2 w 0.31). As the pressure gradient increases, fluid in more sites is mobilized (Figure 2.7A) and the flow rate increases as well. The relation has approximately a quadratic dependence in a region above the threshold, as predicted by Roux and Hermann [108], and linearizes at larger values of the pressure gradient. The self-consistency argument of Roux and Herrmann for the quadratic dependence is in agreement with our findings: in the neighborhood o f the threshold, the increase in the site fraction is proportional to the difference of the applied pressure from the minimum (Figure 2.7A). At the same time, the flow rate in the open paths increases proportionally. It is rather interesting that this regime essentially corresponds to flow of capillaries in parallel (see Appendix A). Note that the four snapshots in Figure 2.5 are all within the quadratic region (region A) of Figure 2.7A; that is, 75% of the sites are mobilized for values of (R*\ VP\/2re ) less than about 1. The linear behavior at large rates can be predicted by a simple equation from the observation that in this regime most sites are open and flow is 5 2 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. almost in series (Appendix A). Because of the uniform flow conductances, it is not difficult to derive the following relationship by using model (2.5), Equation (2.12) shows that the linear part of the curve should intercept the horizontal axis at the mean value of the dimensionless thresholds, equal here to 1. Rewritten in the more conventional form this equation can be used, under the restriction of homogeneous conductances, as a continuum approximation of the macroscopic velocity-pressure drop, far from the onset of mobilization. A different way to represent the relation of Figure 2.7B is to consider the sites, X sile, which are not “trapped.” The results are plotted in Figure 2.8, together with a cubic equation fit to the model results. Some versions of the “population balance” model for foam flow in porous media [71], [82] postulate a cubic relationship (with a pre-exponential factor of 1) between “relative permeability” of gas in foam and the flowing gas volume fraction, which corresponds to the fraction of open sites in our network. Figure 2.8 shows that the network behavior is fitted well with the cubic fit. The pre-exponential factor in the fit to Figure 2.8 would be (2.13) “relative Permeability” of the fluid, — ------- ( /> R 2 | V P 8 npu , as a function of the fraction of 53 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. accounted for in a population-balance model by an elevated effective viscosity of gas in the presence of foam, equal to the inverse of this factor. The effective viscosity here is a constant, unlike the (2/3)-power dependence assumed for foam [35, 71, 82]; this is expected given the linear flow dependence assumed at the pore level in equation (2.5). 3.0 ------------------------------------------------------------------------------------------ 2.5 - qT > 2.0 - 0.0 0.2 0.4 0.6 0.8 1.0 Xsite Figure 2.8: Plot of the dimensionless permeability of the flowing fluid as a function of the fraction of sites belonging to the open paths. Solid line: full-occupancy case, uniformly distributed thresholds («= 1) and viscous model (2.5). Dashed line: full-occupancy case, strongly-disordered system (n-5) and viscous model (2.5). Dot-dashed line: full-occupancy case, Rayleigh throat-size distribution (2.11) and viscous model (2.6). The results appear to be well fitted with a cubic relationship (y=0.2994x\ y=0.5641x3, and y=0.9866X3 for the respective cases). For the Rayleigh case, the cubic fit pertains only to fraction of sites less than 0.9. 54 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.9: The patterns of the open paths at different stages for the dynamic case, and for n=0.1 corresponding from left-to-right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three-fourths of the sites belonging on open paths. Full-occupancy case and viscous model (2.5). 55 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.10: The patterns of the open paths at different stages for the dynamic case, and for n=5, corresponding from left-to-right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three-fourths of the sites belonging on open paths. Full-occupancy case and viscous model (2.5). For different threshold distributions the patterns are somewhat different. Figures 2.9-2.10 shows results analogous to Figure 2.5, for two different values of the exponent n, always assuming equation (2.5). At small values o f n, the threshold distribution is almost unimodal (near 1), all paths are therefore straight, and the flow occurs as if in a bundle of parallel capillaries of very similar thresholds. The 56 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. almost-parallel flow results in very few transverse bonds becoming open, even at high flow rates. Mobilization of the sites and the approach o f the curve in Figure 2.7B to its linear limit (and to equation (2.13)) occurs very fast, as also suggested by equation (A.3) in Appendix A. At large n, on the other hand, the relative disorder is much greater, and the paths are more tortuous. It was shown in Kharabaf and Yortsos [64] that this case theoretically approaches the process of directed Invasion Percolation, where the minimum-energy paths belong to the percolation cluster. Now, open paths contain a larger fraction o f transverse bonds. The fraction of sites belonging to the mobilized fluid depends on the applied pressure gradient over a much larger range than in the previous case (compare Figure 2.11 A to Figure 2.7A), which is also consistent with the spread o f the cumulative density function, as predicted by equation (A.3). The corresponding flow rate-pressure drop relation is shown in Figure 2.1 IB. The quadratic dependence regime extends over a larger region of the pressure gradient, consistent with Figure 2.11 A. Ultimately, though, the curve becomes linear and is fitted well with an equation of the type (2.13). For completeness, the fit of the “relative” permeability to the empirical cubic relationship is shown in Figure 2.8. As with the more homogeneous cases, the empirical law appears to be quite good, although not as good as in the more homogeneous case (see also below). 5 7 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 1.0 0.8 0.6 0.4 0.2 0.0 2 3 4 5 0 1 R *| VP|/2xe (A) 4 3 2 1 0 2 3 4 5 0 1 R*|VP|/2re (B) Figure 2.11: The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient for the dynamic case (panel A) and the dimensionless flow velocity-pressure gradient relationship (panel B). Full-occupancy case, strongly-disordered system (n=5) and viscous model (2.5). Note again the quadratic dependence near the onset of mobilization and the approach to a linear relationship at higher flow rates. 58 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Next, we used model (2.6) to explore the effect o f a variable flow conductance, for the case o f uniformly distributed thresholds («=1). To obtain finite values for the mean and the variance of the conductance, we truncated from the threshold distribution the very low values, and considered instead a uniform distribution in the interval (0.1, 1). The corresponding conductance distribution is then /? « •)= J ^ r 4'3 wherel < (T < 103. (2.14) The resulting stretched distribution is typical of flows in porous media, it applies equally well to Newtonian fluids and has important implications to flow ermeability, as pointed out in earlier studies, e.g. by Katz and Thompson [61]. Simulations were conducted in a smaller lattice (30x30). Flow patterns, the fraction of sites mobilized and the flow-rate pressure gradient relationships are shown in Figures 2.12 and 2.13, respectively. Because very large flow conductances are associated with the open paths at the early stages, where smaller thresholds (large pores) are open, the viscous flow contribution to the pressure is much smaller than previously. One expects, therefore, that the pattern of open paths should be similar to the static case of Kharabaf and Yortsos [64]. Figures 2.12-2.13 show that this is the case during most of the range of the quadratic behavior. The fraction of open sites rises much faster than in the case of model (2.5), as the large conductances dominate the early behavior (Figure 2.13A), while because of the increased conductance, it takes relatively moderate pressure gradients to open all of the bonds in the lattice. The linear relationship between flow rate and pressure drop is also 59 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. reached faster (Figure 2.13B). Even though new paths are opened, the large- conductance clusters, once formed, dominate the flow behavior, much like in the Newtonian fluid case [61]. The approach to a linear model of the type (2.13), therefore, is likely to be a realistic description for typical porous media, where the variability in conductances is large [61]. Figure 2.12: The patterns of the open paths at different stages for the dynamic case corresponding from left-to-right, top-to-bottom, to the MTP, one-fourth of the sites, half of the sites and three-fourths of the sites belonging to open paths. Full-occupancy case, uniformly distributed thresholds (»=1, range (0.1, 1)) and viscous model (2.6). 6 0 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.8 0 ) (/) x 0.4 0.2 0.0 0 1 2 3 4 5 R*|VP|/2xe (A) 0 ) % - °'9 Q . ^ 0.6 0.3 0.0 0.0 0.4 0.8 1.2 1.6 2.0 R*|VP|/2t6 (B) Figure 2.13: The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient for the dynamic case (panel A) and the dimensionless flow velocity-pressure gradient relationship. Full-occupancy case, uniformly distributed thresholds (n=1, range (0.1, 1)) and viscous model (2.6). 61 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. For completeness, we also conducted simulations using the Rayleigh distribution (2.11), which results in the stretched-exponential conductance distribution /?(0 = 7 r 1 /3 e x p 6 f < 4/3A (2.15) 4 Corresponding results are shown in Figure 2.14. The findings are similar to the previous: the wide conductance distribution causes paths to be opened at a faster rate and at smaller values of the applied pressure gradient, at least near the onset of flow; the large conductances associated with the first few open paths dominate the flow behavior even at the later stages, so that the linear flow regime between flow rate and pressure drop is reached faster; and even though more paths are opened at later stages (note the long tail in the site distribution in Figure 2.14A), their conductance does not significantly add to the overall flow behavior. The “relative permeability” for this distribution is plotted in Figure 2.8. We note that the part of the curve corresponding to values of X s ite less than about 0.9, is fitted very well with the cubic fit. At higher values, near full site occupancy (high rates), the fit does not work very well, however. In that range, the overall flow conductance increases due to the fact that additional bonds open without the need for the opening o f new sites. 6 2 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.8 < 1 ) - 4 — » 0.4 0.2 0.0 0 1 2 3 4 R*|VP|/2xe (A) 4 3 * a: - © - 2 3 '4 ' 1 0 0.0 0.4 0.8 1.2 1.6 2.0 R*|VP|/2xe (B) Figure 2.14: The fraction of sites (pores) of the lattice belonging to the open paths as a function of the dimensionless applied pressure gradient for the dynamic case (panel A) and the dimensionless flow velocity-pressure gradient relationship (panel B). Full-occupancy case, Rayleigh throat-size distribution (2.11) and viscous model (2.6). The large conductances dominate the flow behavior and the flow rate-pressure drop relation rapidly approaches a straight line. 63 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 2.3.2 The Partial-Occupancy Case Consider, next, the flow behavior of a Bingham fluid, when it occupies and flows in a restricted fraction of the lattice. This problem essentially corresponds to a reduced lattice coordination number, so we intuitively expect a behavior similar to the previous, except that the pressure gradient requirements would be higher, since the phase connectivity is diminished. In the results to be shown, the phase occupancy was determined using Invasion Percolation (IP) without trapping [63], in which the fluid assumes the pattern of the invading phase after the breakthrough. We will use the “saturation”, S, to denote the fraction of the lattice occupied by the Bingham fluid. S= 1 corresponds to the full-occupancy case, discussed previously, and it is the pattern obtained after the completion of the invasion process. Smaller saturation values correspond to less-connected patterns, with fractal patterns arising near percolation ( S = 0). In the results shown below, the fluid occupancy remained fixed and identical to the initial (thus, mobilization and flow did not give rise to a displacement, which is a process further discussed in a subsequent section). Figures 2.15-2.16 show flow patterns, the fraction of sites mobilized and the flow-rate pressure gradient relationships for a partial-occupancy case corresponding to 5=0.5, the Rayleigh model of equation (2.11) and the viscous model (2.6). In the figures, the fluid-in-place occupied the largest pores, effectively simulating a drainage process. The results are quite instructive. For the specific saturation value, 64 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the fraction of sites on the mobilized paths increases near the onset of mobilization, but relatively rapidly reaches a plateau, which for the range examined represents only about 20% of the occupied sites. Correspondingly, the flow rate-pressure drop relation becomes linear almost right after the onset of flow. This behavior is partly due to the large conductances, which dominate the process, as discussed in the previous case. In addition, here the fluid occupies a large fraction of sites which are dangling and belong to dead-ends. The fluid in these sites will not be mobilized, regardless of the applied pressure gradient, just as it would not flow if the fluid were Newtonian. As a result, they do not contribute to the flow process. Figure 2.17 shows the variation of the minimum pressure gradient with saturation for two different cases, one in which the fluid configuration was determined by invading the largest pores (the drainage case) and another by invading the smallest pores (the imbibition case). In both cases, the minimum pressure gradient increases as the occupancy is reduced, and sharply diverges as the percolation threshold is approached, due primarily to the increased tortuosity of the pattern. The effect is similar to the increase in the gradient, in the full-occupancy case, as exponent n increases. The hysteresis between the two branches is due to the largest thresholds needed to mobilize the fluid in the imbibition case, as the fluid occupies smaller throats with correspondingly larger thresholds. As percolation is approached, the minimum pressure gradient diverges, following a critical scaling behavior, which is expected to have the same dependence as the chemical path of percolation theory. Note that because of finite size effects, the lowest saturation 65 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. value is not equal to zero, as theoretically expected. The latter is approached algebraically, as the lattice size increases, and in a rather slow dependence for the case of 2-D lattice, e.g., Sc ~ L~°'u [32], For simplicity, we will not proceed with further discussion o f these issues, however. Figure 2.15: The open paths at two different stages of mobilization (Panel A) for the partial-occupancy case (.5=0.5), Rayleigh throat-size distribution (2.11), and viscous model (2.6). Gray denotes the fluid, the occupancy of which was determined by Invasion Percolation without trapping, black denotes the open paths. Note the existence of a number of dead-end pores. 66 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.12 0.10 0.08 0.06 0.04 0.02 0 2 4 6 8 10 12 14 R*|VP|/2xe (A) 0.4 0.2 0.0 0 2 4 6 8 10 12 14 R*|VP|/2xe (B) Figure 2.16: The fraction of sites (pores) belonging to the open paths as a function of the applied pressure gradient (Panel A) and the flow rate-pressure drop relation for the partial- occupancy case (S=0.5), Rayleigh throat-size distribution (2.11), and viscous model (2.6). 6 7 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 1.8 1.6 <D L_» 1.4 C N c E 1.2 C L > 1.0 hi 0.8 0.6 0.4 — Drainage - Imbibition s v 0.3 0.4 0.5 0.6 0.7 s 0.8 0.9 1.0 Figure 2.17: The minimum pressure gradient for the onset of mobilization as a function of the volumetric fraction, S, of the Bingham fluid, put in place using a drainage (solid line) and an imbibition (dashed line) process. Partial-occupancy case, Rayleigh throat-size distribution (2.11), and viscous model (2.6). 2.4 Immiscible Displacement of a Bingham Fluid by a Newtonian Fluid Consider, next, the immiscible displacement of a Bingham fluid by a Newtonian fluid. As we noted before, the majority o f studies in this area have focused on Newtonian fluids. We will refer to two specific studies, a pioneering one by Lenormand [79], who classified displacement regimes at the pore-network level, 68 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. and its extension to fully-developed displacements by Yortsos et al. [134]. Lenormand [79] delineated three regimes that may occur during the displacement of one immiscible fluid by another: A percolation (IP) regime, at low values of the capillary number, Ca, where capillary forces predominate, and two viscous regimes, at high values of Ca , one corresponding to piston-like (PL) displacement, at small values of the viscosity ratio M, and another to viscous fingering (VF), at large values of M, respectively. Here, we defined the capillary number Ca = , r where v is the fluid velocity and y the interfacial tension between the fluids, and the ratio M - between the viscosities of the displaced fluid, denoted by subscript 2, M x and the displacing fluid, denoted by subscript 1, respectively. Yortsos et al. [134] extended Lenormand’s approach to fully-developed drainage, which they characterized in terms of two global regimes, a stabilized and a viscous-fingering regime. A number of related studies have addressed many other properties. The state o f the art in non-Newtonian fluids is not nearly as advanced as with Newtonian fluids. Consider the displacement of a Bingham plastic by a Newtonian fluid. In such displacements, in addition to the capillary number and the viscosity ratio, defined above, and where for the viscosity of the displaced fluid we will now use [ip, a third dimensionless number arises, representing the ratio of the yield stress to the capillary pressure. This number was defined in Kharabaf [62] and Kharabaf et al. [63] as the yield stress number 69 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (actually, in [62] and [63] the square root of the permeability, k , was used instead of yfk ). We also add that there is one more parameter of importance, namely the variance of the pore-size distribution, which will affect the distribution of the thresholds in the pore throats. For instance, the variance o f the pore-size distribution enters directly in affecting the tortuosity of the MTP. The existence of these dimensionless parameters implies that the phase diagram of the displacement must be portrayed in, at least, a three-dimensional space. Kharabaf [62] studied only the competition between capillary forces and the yield stress (see also the foam invasion model, Ref. [65]). Kharabaf et al. [63] discussed viscous flow aspects, but their findings were incomplete. Here, we will paraphrase their findings for the case when viscous forces are negligible, as follows: For a fixed lattice size, L, there exist two limiting regimes, depending on the relative magnitude of the yield stress compared to the capillary forces. The pattern will be of the IP type, if the Bingham fluid is readily mobilized, and capillary forces dominate, namely if condition tLNx « 1 (2.17) applies, where t is the tortuosity o f the MTP, a function of the threshold distribution [64], Clearly, as the front advances in a finite lattice, the IP pattern becomes more easily satisfied. In the opposite case, 70 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. tLNt » 1 (2.18) at least initially, the pattern coincides with the MTP, thus resulting into a very poor displacement. Depending on the conditions, however, as the displacement proceeds, it will become easier to mobilize the Bingham fluid, and a mixed pattern with some IP features, will emerge. Clearly, condition (2.18) is progressively better satisfied as the lattice size increases, thus, under these conditions and in the absence of viscous effects, the MTP will be the dominant pattern. Thus, at conditions of small capillary numbers, we expect a transition from an E P pattern to an MTP pattern, as the yield stress number increases. We note that this transition is opposite to that for foam invasion [65], in which it is the Bingham fluid that displaces a Newtonian fluid. In the latter case, the transition is from IP to a compact (PL) pattern, as the yield stress number increases. Consider, now, the limit of large capillary numbers, in which viscous and yield stress effects predominate and capillary forces are negligible. The objective is to identify the displacement patterns as a function of the flow rate and the viscosity ratio. In our pore-network simulations of this displacement, an interface separating the two fluids exists at all times, even though capillary forces were neglected. The relevant dynamic parameter now is the ratio - (see also [63]). For large Ca v//, values of the latter, the displacement will be more like the MTP, while in the opposite case, more like to that for Newtonian fluids. Effectively, the former regime would correspond to the onset of flow (e.g. point A in Figure 2.6B), while 71 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the latter to the regime where the fluid is more like Newtonian (point B in Figure 2.6B). We used simulations to test these expectations. For the rheology of the Bingham fluid we took model (2.6) but with a fourth-power conductance dependence Displacements were conducted over a fixed-size lattice, across which a constant pressure drop is applied. Periodic boundary conditions where applied on the lateral sides of the lattice. In the simulations, we first consider the case in which the pressure drop is maintained equal to that required for the onset of mobilization of N the initial fluid configuration (corresponding to the case —L » 1 ) . Clearly, this is Ca sufficient to sustain the flow of the mobilized fluid, the flow rate increasing constantly as the displacement proceeds. It follows from the above equations that the only relevant parameter for this problem is the viscosity ratio M, the sensitivity to which was examined in the following. q( = 0 when | Apt |< — (2.19) while the Newtonian fluid obeys the usual Poiseuille law q , = Mr* APi (2.20) 72 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.18: Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime A of Figure 7B (onset of mobilization) for M=00. Gray denotes mobilized Bingham fluid, white denotes non-mobilized Bingham fluid and black is the displacing Newtonian fluid. 73 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.19: Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime A of Figure 7B (onset of mobilization) for M=\. Gray denotes mobilized Bingham fluid, white denotes non-mobilized Bingham fluid and black is the displacing Newtonian fluid. Given an interface configuration, the pressure field was calculated in both the Newtonian and the Bingham plastic, and that interface in the particular pore throat that moves at the highest rate was advanced in the next time step. During the same time interval, all other interfaces were accordingly advanced. The algorithm is similar to the IPM discussed before, except that in the Newtonian phase, calculation of the flow rate does not require the opening o f a path (effectively 74 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. corresponding to zero thresholds), and only the effective thresholds of the Bingham fluid are updated. Patterns at different stages of the displacement and at different mobility ratios from these simulations are shown in Figures 2.18-2.20. In the figures, white denotes the Bingham fluid, gray the part of that fluid which is mobilized, and black the Newtonian fluid. Figure 2.20: Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime A of Figure 7B (onset of mobilization) for M=0.05. Gray denotes mobilized Bingham fluid, white denotes non-mobilized Bingham fluid and black is the displacing Newtonian fluid. 75 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.18 corresponds to M = 00, thus the pressure drop in the Newtonian fluid is zero. As a result, the pressure applied on the Bingham plastic remains constant during the process, and an increasingly larger fraction of the fluid becomes progressively mobilized, as the displacement proceeds. Correspondingly, the flow rate rises, as time increases (see below and also Figure 2.18), and when breakthrough is approached, it shows an almost exponential rise. The displacement efficiency is quite poor, however. The displacement pattern follows a relatively narrow path, close, although not identical, to the MTP. This pattern is certainly very different from DLA [38], which characterizes the same unfavorable mobility (M = 00) displacement of two Newtonian fluids. As the viscosity of the displacing fluid increases, the pressure in that phase drops, resulting into a smaller fraction of open paths (e.g. see Figures 2.19 and 2.20 for M = 1 andM = 0.05, respectively). Correspondingly smaller is the flow rate (Figure 2.18). Even though there is a slight improvement in the displacement efficiency as the viscosity ratio decreases, the pattern remains almost quasi 1-D, and far from the compact Piston-Like displacement obtained under favorable mobility conditions when both fluids are Newtonian. An interesting observation is that when M is small the flow rate appears to roughly stabilize as a function of time (Figure 2.18). This is due to the fact that as the interface advances, the pressure drops proportionally, but at the same time the force required to open a mobilization path ahead also decreases proportionally. For a suitable flow rate, the two 7 6 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. competing tendencies are roughly equal. This can be readily demonstrated with the following 1-D example. Consider the displacement of a Bingham fluid by a Newtonian fluid in a 1-D geometry. We will use equations (2.19)-(2.20) with radius r0 , replace pressure differences by pressure gradients and drop subscript i. At time t, the front is at location x f (t) , while the pressure distribution is where P0 > L denotes the dimensionless applied pressure. Applying (2.22) at x = L , where P = 0, gives the following differential equation for the front evolution, P = Pq ^~rx in 0 < x < x f (t) Mr* (2.21) and Mr0(rQ P0 - ( L - x f )) x f = — T77T----------;------------- M (L - X f ) + Xj- (2.23) where we took into account the identity x / = q! r0 2 . Integration o f (2.23) gives the implicit relationship / X f + P0r0 —L + X f r0Mt Po^o-L (2.24) v 11 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Based on (2.23) and (2.24), the flow rate can be calculated as a function of time. Results are shown in Figure 2.2 IB for different values of the viscosity ratio M. The features are very similar to those from the pore-network simulation in Figure 2.21A. As liquid displaces foam in foam-acid diversion [24, 84, 107], the yield stress number Nr in equation (2.16) is of order 1 or somewhat less; the resistance to movement per lamella is of the order of the capillary entry pressure [98, 99] and lamellae are spaced on the order of one per one or a few pores [101]. The capillary number Ca, as defined above based on water viscosity, is much less than 1, and M is large. These displacements would be best represented by Figure 2.20. We conclude that for large values of the ratio - , yield stresses predominate Ca vjux and the patterns are quite different from those for two Newtonian fluids, being more like 1-D, MTP-like, regardless of the viscosity ratio. As a result, the displacement has very low efficiency. The MTP like patterns are in fact similar to those observed in CT images by Nguyen [83, 84], where liquid displaces foam. As N the ratio decreases, yield-stress effects become weaker, and at sufficiently small values, the flow regime corresponds to point B in Figure 2.7B, in which case the patterns approach those of the corresponding Newtonian displacement. An illustration o f this effect is shown in the displacement patterns o f Figure 2.22, for the case M - o o , and where the applied pressure is sufficiently large to mobilize a 78 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. large fraction of the Bingham fluid (corresponding to regime B in Figure 2.7B). The features of this displacement are quite similar to DLA, which characterizes displacements involving Newtonian fluids at infinite viscosity ratio. Similar is the analogy with other values of the viscosity ratio. The above can be summarized in the phase diagram o f Figure 2.23. The yield- stress dominant regime corresponds to the needle-like structure o f the MTP pattern and will also be affected by the threshold distribution (e.g. see Figures 2.5, 2.9 and 2.10). The viscous dominant regime is of the DLA or the PL type depending on the values of the viscosity ratio (large or small, respectively). The boundaries separating the yield-stress dominant regime from the viscous regime at large values of the capillary number are suggestive. Further analysis is needed to define precisely these boundaries. We also point out that when the Bingham fluid displaces the Newtonian fluid, the yield-stress dominant regime is piston-like, rather than MTP-like (see also [65]). We close this section by noting that the findings o f the previous section for the partial occupancy case (e.g. see Figures 2.15-2.17) can be used to obtain a conventional upscaling o f the two-phase flow problem. These findings suggest that in a limiting regime, the flow rate-pressure drop relation can be approximated as that of an equivalent Bingham fluid, but with relative permeability and critical thresholds being functions of the saturation of the fluid, namely g j(S )\= krj(S) A,pj | when | A p . |> r AS) 7 9 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. q. = 0 when | Api |< x} (S ) (2.25) Here, we used dimensionless notation for the flow of phase j, and denoted by krj(S) and by t -(S) the relative permeability and the effective threshold of the Bingham phase, as a function of its saturation S. Figures 2.15-2.17 show some of the features of the above two functions, obtained from limited investigations. Determining the detailed properties of these functions will require additional study using the above-introduced algorithms. O f considerable interest would be the analysis of displacements in which the yield-stress fluid saturation approaches zero, in which the case the minimum pressure gradient for its mobilization diverges. The corresponding mathematical problem should be of certain interest and will be worth exploring. 8 0 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 1e+2 1e+1 1e+0 1e-1 0 1e-2 1e-3 1e-4 1e-5 0 500 1000 1500 2000 2500 Tim e (A) 1e+2 1e+1 1e+0 1e-1 0 1e-2 1e-3 1e-4 1e-5 0 500 1000 1500 2000 2500 T im e (B) Figure 2.21: The variation of flow rate with time during the displacement of a Bingham fluid by a Newtonian fluid, corresponding to point A of Figure 7B, and for different values of the viscosity ratio. Results from pore-network simulations (panel A), simple analytical model (panel B). 81 W = O Q M =1 M=0.2 M=0.1 M=0.05 M = » M =1 M=0.2 M = 0 .1 MM D.05 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 2.22: Patterns of the displacement of a Bingham fluid by a Newtonian fluid (from right to left), corresponding to regime B of Figure 7B (after the mobilization of most of flow paths) for M=o o. Gray denotes mobilized Bingham fluid, white denotes non-mobilized Bingham fluid and black is the displacing Newtonian fluid. The relationship to DLA patterns is apparent. 82 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. M log Ca DLA or Compact Invasion Percolation Yield Stress Dominant Figure 2.23: Phase diagram of the displacement of a Bingham fluid by a Newtonian fluid. The yield-stress dominant regime corresponds to the needle-like structure of the MTP pattern and may be affected also by the threshold distribution. The viscous dominant regime is of the DLA or the PL type depending on the values of the viscosity ratio (large or small, respectively). In the figure, the viscosity ratio is assumed large (hence the DLA-type regime shown). When the Bingham fluid displaces the Newtonian fluid, the yield-stress dominant regime is PL, as in [65]. 83 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 2.5 Summary and Conclusions Motivated by the need to understand better the flow and displacement of a fluid with yield-stress rheology, we described in this chapter 2-D pore-network simulations in which the effect of the microstructure was highlighted. First, we extended the IPM algorithm of Kharabaf and Yortsos [64] to account for viscous flow effects associated with the onset of flow. When viscous effects are important, the mobilization of the fluid, after the onset of flow, occurs in a layered-like fashion, with the result that fluid in many pore throats transverse to the flow is not mobilized. The fraction of the sites which belong to the open paths increases more slowly with the applied gradient, compared to the static case, to the case when viscous effects are not important, or when the threshold distribution is unimodal (homogeneous case). Using a qualitative analogy between the pressure gradient- open sites curve with the more conventional capillary pressure curve provides for a useful viewpoint, even though the two are not mechanistically similar. By the same token, a “relative” permeability concept can be also introduced to express the flow rate as a function of the open sites. We find that this relation is reasonable fitted with a cubic power. In all cases, the qualitative features of the quadratic behavior of the flow rate near the onset of flow can be captured by a model of parallel capillaries. When the conductances of the pore throats are widely varying, as expected in typical porous media, the viscous effects are much weaker. Then, the mobilization 84 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. of the fluid can be approximately described with the static IPM method. Under these conditions, the flow rate-pressure drop relation is effectively a straight line with an effective minimum threshold equal to the arithmetic mean of the threshold distribution. A sim ilar behavior applies for the problem in w h ich the fluid occupies the lattice only partly (the “two-phase” flow problem). Relative permeability and effective thresholds as a function of the fluid saturation can be derived, using the pore-network simulations in this case. Such curves can be useful for upscaling purposes. The algorithm was also used to classify the patterns that would arise in immiscible displacements involving a Bingham fluid. A phase diagram was constructed involving the capillary number and the yield stress number. Yield- stress dominant regimes were identified, where the displacement pattern is of the MTP (Minimum Threshold Path) type, with a 1-D structure that depends on the threshold distribution, if the Bingham fluid is being displaced, and of a compact pattern type, if the Bingham fluid is the displacing fluid. Capillary and viscous controlled patterns have the typical patterns as with Newtonian fluids. 85 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Chapter 3 FOAM GENERATION AND MOBILIZATION IN POROUS MEDIA 3.1 Introduction In the petroleum industry, foam is injected into porous media for two primary purposes: gas diversion to improve oil recovery [102, 113, 125] and acid diversion for matrix acid treatments [25, 46]. Foam is also used to direct the flow of remediation fluids in subsurface aquifer-remediation processes [53]. A number o f studies find that foam does not alter the relative permeability or viscosity of the liquid phase that makes up the foam [9, 31, 42, 111]. Foam does drastically reduce the mobility of gas, however. This reduction is inversely related to bubble size [36, 37, 54, 71], i.e. it is directly related to the number of liquid films, or lamellae, that separate gas bubbles, per unit volume of gas in the pore space. Lamellae resist the movement of gas by resisting movement themselves: this resistance arises from the drive by each lamella to minimize its surface area, according to its interfacial tension. As described below, lamellae form in pore 86 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. throats, Therefore, initiating movement of a lamella requires pushing it out of the pore throat, the position of minimum surface area, into the pore body, the position with larger surface area [98, 99, 132] The pressure drop required to mobilize a lamella, APmm, is inversely proportional to the minimum radius of curvature of the lamella as it is displaced, which is of the order of the radius o f the pore throat R: A Pm in ~ Ay IR (3.1) where y is the interfacial tension between liquid and gas. If lamellae block some throats, then the gas continues to flow as a Newtonian fluid, albeit with reduced relative permeability [26, 104] If enough lamellae are present to block gas flow completely, then the gas behaves as though it had a yield stress, because it cannot flow unless the pressure gradient is sufficient to mobilize lamellae along some pathway through the pore network. Much o f the gas may remain trapped even as some fraction of the gas flows [42, 6 6] The pressure gradient required to keep lamellae moving along this path may be less than that required to initiate flow [98, 99, 132]. Initiating flow requires mobilizing all the lamellae along the path from their initial positions in pore throats, where the resistance is greatest. Once lamellae are moving, their average resistance is less, because most of the time they are individually in positions o f lower resistance. In addition to the effective yield stress, the drag on moving lamellae imparts a shear- thinning viscosity to the foam [54,132], although this effect will be neglected here. Gas mobility then depends on the number density of lamellae in the gas phase. In turn, the latter depends on processes that create and destroy lamellae. In 87 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. particular, lamellae are destroyed by high capillary pressure, especially near the “limiting capillary pressure.” [10]. Lamella destruction plays a large part in controlling steady-state foam properties, but is not the focus of this chapter. 3.2 Lamella-Creation Mechanisms Lamellae are created by three different mechanisms [21, 31, 96, 101]: One mechanism is “leave behind,” i.e. the stabilization of liquid films or lenses in pore throats as gas invades the adjacent pore bodies through other throats, as illustrated in Figure 3.1. Depending on pore-throat geometry, leave-behind lenses may break, if surfactant is not present when the adjacent pore bodies are drained by gas. Although sometimes cited as a source of “weak” or ineffective foam [42], leave-behind can create a large number of lamellae. If leave-behind is the only lamella-creation mechanism, however, the gas has at least one continuous pathway for flow. A second mechanism is “lamella division”, where two or more lamellae are created from one. Each time a mobilized lamella passes a pore body with more than one pore throat unoccupied by liquid or another lamella, the lamella must either break or span both open throats, as illustrated in Figure 3.2. Lamella division could fill a large region of pore space with lamellae, starting with only a small number of moving lamellae, as illustrated in Figure 3.3. Models based on the concept of lamella division explain a variety of foam-generation data [44, 58, 60, 98, 123]. A 88 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. paradox afflicts this explanation on the pore-network scale, however: if lamellae must be moving to divide, and presumably they are moving downstream, then how are lamellae created upstream? If lamella division were the only mechanism operating, shouldn't there be a growing region near the inlet with few lamellae and high gas mobility? initial state: satur ated with liquid ‘ ‘ ‘ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ ♦ leave-behind lamellae gas invasion d d Figure 3.1: Schematic of lamella creation by the leave-behind mechanism. Gray diamonds represent sand grains; gaps between them represent pore bodies and throats (After Ref. [113]). Figure 3.2: Schematic of a single lamella-division event (After Ref. [103]). 89 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.3: Schematic of lamella division on the pore-network scale. In the first panel, one moving lamella produces 14 lamellae in place by division. In the second panel, a second moving lamella (taking a different path) produces 9 more in place, plus displaces an extra 2 moving ahead of it. In the third panel, a third moving lamella (taking a different path) produces 2 more in place, plus displaces an extra 5 moving ahead of it (After Ref. [113]). Figure 3.4: Schematic of snap-off in a pore throat. Black denotes pore-throat wall, gray the water and white the gas (After Ref. [113]). decreasing capillary pressure 90 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. A third mechanism is “snap-off’: lamellae are created in gas-filled pore throats (Figure 3.4) if the local capillary pressure falls to about half the capillary entry pressure o f the throat [101, 103]. (This depends on the geometry of the throat and the wettability of the medium, but the value of one-half is a reasonable representative value for 3D pore geometries.) The paradox of snap-off is as follows: If gas occupies the throat, then presumably capillary pressure is (or was) at the capillary entry pressure. How does capillary pressure at the throat fall to one- half this value to trigger snap-off? There are seven ways by which this can happen [103]. Three are relevant to steady flow in homogeneous media and to the discussion that follows: 1 Capillary pressure fluctuates as moving lamellae shift their positions and curvatures, altering the pressure drop between bubbles. 2 Since gas is nearly inviscid, pressure is nearly uniform along gas bubbles. Therefore, capillary pressure is lower at the rear of long bubbles than at the leading edge. We note that this same mechanism also applies to gas ganglia moving upwards against gravity. 3 When gas invades a wide liquid-filled pore body through a narrow pore throat, liquid draining from the pore body can sweep back into the throat and cause snap-off, as illustrated in Figure 3.5. This is called “Roof snap-off’ after Roofs study of oil-trapping in waterflooding [97]. 91 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.5: Schematic of Roof snap-off as gas invades a pore body. Black is pore wall, gray is water and white is gas. In Roof snap-off, alone among other snap-off mechanisms, the size of the pore body (relative to the throat) matters, because it is its radius that briefly controls local capillary pressure at the moment that gas fills the pore body. Some, but not all, pores then become “germination sites” for lamella creation by Roof snap-off [71, 96], but only when gas invades a pore body filled with liquid. In all other snap- off mechanisms, snap-off depends on the radius of the throat alone (and the local capillary pressure). 3.3 Models for Foam Generation A number o f foam-generation models are based on the assumption of repeated Roof snap-off, long after the pore body drains its liquid [26, 42, 48, 55, 65, 70, 71, 92 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 82, 111, 127, 133]. In particular, “break and reform” models [26, 55, 65, 111, 127, 133] assume that lamellae form by snap-off in pore throats, are mobilized, immediately break upon mobilization, and then, after a short delay, reform by snap- off. Rossen [103] reviews the evidence for this assumption and shows that it is problematic. These models do not explain why liquid would drain back to the throat at high capillary pressure, or how capillary pressure could be low enough to trigger liquid flow to the throats. Other models for foam generation are based on mobilization and division as the key to the creation of strong foam [58, 60, 96, 104, 123]. Rossen and Gauglitz [104] present a model to explain the minimum velocity (or pressure gradient) required to create strong foam in a porous medium already partially drained by gas. (Their experiments were designed to represent foam generation where gas is already present in the field and gas and liquid are co-injected to create foam.) They assumed that some population of throats is occupied by liquid before foam is created. Once surfactant is present and injection rates increased, these lenses drain to lamellae and can be mobilized, if pressure drop across the given throat is sufficient (cf. Eq. 3.1). The lamellae first mobilized would be those that block large stagnant gas clusters from flowing. The local pressure drop required to mobilize an individual lamella is related to the macroscopic pressure gradient |W y M m | and the length of the gas cluster Lc containing the lamella: APm in = VPm in Lc (3.2) 93 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Eq. 3.2 implicitly assumes that the pressure gradient is uniform across the pore network. Using percolation theory and a Bethe tree (also called a Caylee tree), one can derive a formula for the average cluster size as a function of the nearness to the percolation threshold for flow of the gas phase. The average cluster size diverges to infinity as the percolation threshold is approached, making foam generation easier. One can also relate the nearness to the percolation threshold to the injected liquid- gas volume ratio and the fluid viscosities, allowing one to predict the minimum pressure gradient \VPm,n\ for foam generation from the injected liquid fraction. Rossen and Gauglitz [104] show good agreement with data, after adjusting one parameter accounting for the number of lamellae that must be mobilized to trigger foam generation. The model predicts that \VPmm\ varies inversely with permeability, for porous media with scaleable pore geometry, like sandpacks and beadpacks, and the data of Gauglitz et al. [44] confirm this prediction. Rossen and Gauglitz [104] proposed a similar model for discontinuous-gas foam, i.e. one below the percolation threshold for gas, but still assuming (as in percolation theory) that the lamellae are placed randomly on the pore network. Like the version for continuous-gas foams, |VP| is assumed to be uniform across the pore network. For discontinuous-gas foams, the mobilization pressure gradient also approaches zero at the percolation threshold, because the average bubble size approaches infinity there. The results of this model are shown schematically in Figure 3.6. 94 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 150 100 Foam Mobilization Foam Plugging Continuous- Gas Foam 0.25 0.30 f Figure 3.6: Schematic of the results of the percolation model for foam generation by Rossen and Gauglitz. / is the fraction of pore throats through which gas can flow (not blocked by water-filled pores, lenses or lamellae). Here/=0.25 represents the percolation threshold, where the pressure gradient for foam mobilization |'W pm m \ approaches zero. The diagonal line on right side of figure represents foam generation during continuous injection of gas and liquid (starting from no-foam state of continuous gas flow). The diagonal line on the left side of the figure represents mobilization of discontinuous-gas foam (After Ref. [106]). The Rossen and Gauglitz [104] model makes several simplifying assumptions about flow in the pore network that Rossen et al. [106] later relaxed. An improved estimate of \VPm in\ for mobilizing discontinuous-gas foams in a Bethe tree network gives some curvature to the straight line to the left of the percolation threshold in Figure 3.6. Second, for continuous-gas foams, Rossen et al. [106] report a paradox. 95 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. While the average cluster of gas blocked by a lamella diverges in size to infinity at the percolation threshold, the number density of these clusters on the pore network approaches zero [100]. Thus, the number of clusters greater than any given size goes through a maximum and then to zero as the percolation threshold is approached. The implications of these two alterations in the theory are shown in Figure 3.7. The two parts of the model no longer agree at the percolation threshold. Rossen et al. [106] contend that the reason for the divergence is that the theory still assumes a uniform |VP| on the pore network near the percolation threshold, while it is well known that |VP| varies greatly in the network. Advancing the theory, and reconciling the two halves of the model, requires a study on a pore network, where the assumption of uniform IVi3! can be abandoned. Tanzil et al. [123] examined foam generation upon gas injection into a medium saturated with surfactant solution. They propose that foam generation depends on exceeding the critical pressure drop across the displacement front rather than on the pressure gradient. They propose a critical value of the capillary number, based on the pressure drop across the displacement front, as the threshold of foam generation. Experimental verification of this condition is difficult in fixed- injection-rate experiments, however, because the pressure drop increases upon foam generation. The measured pressure drop is as much a reflection as a cause of foam generation. 96 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 150 S S a 100 fl •P 4 T 3 « s - o a u s V ) (A a s . P h 50 Foam Mobiliz ation Foam Plugging <?.20 Foam Mobilization Continuous- Gas Foam 0.25 f 0.30 Figure 3.7: Schematic of results of the percolation model of Rossen et al. [106] modified from that of Rossen and Gauglitz [104], The two halves of the model diverge from each other at the percolation threshold (/= 0.25). A fully mechanistic foam simulator requires an accurate mathematical algorithm for lamella creation. Moreover, in some foam applications, foam generation is in doubt, while the issue of foam generation is crucial to process success. Therefore, it is important to test the various theories for foam generation at the pore-network level. This is the main motivation for the present study. We model on a pore network the three processes of lamella creation (leave-behind, lamella division, and snap-off), under various scenaria and flow conditions. The dynamics are complex. For instance, even upon steady flow, once strong foam forms (perhaps by 97 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. mobilization and division), gas mobility falls abruptly, gas drains additional pores, and additional lamellae are created by snap-off and leave-behind. Using a repetitive snap-off model we also probe the validity of the model o f Rossen and Gauglitz (cf. Figure 3.7). For comparison, we show results based on an assumption of repeated Roof snap-off, but by illustrating foam generation by other mechanisms we probe the assertion that repeated Roof snap-off may not be necessary for foam generation. Finally, we identify a new mechanism of snap-off near the pore entrance, which explains how lamellae are created near the entrance when foam is generated by lamella mobilization. 3.4 Pore-Network Model Development Modeling foam formation and propagation in a pore-network model must borrow from two different areas: the flow of two immiscible fluids, where capillary and viscous forces are important, and the flow of a fluid exhibiting a yield stress, for example of a Bingham plastic. Significant progress has been made in the first area in the last two decades using pore-networks (e.g. see recent review by Jackson and Blunt [16]). Modeling the onset and the subsequent flow o f a Bingham fluid, however, has been accomplished in Chapter 2. A key tool for solving the latter problem is the process of Invasion Percolation with Memory (IPM), introduced by Kharabaf and Yortsos [64] to model the onset of flow of a fluid with yield stress, and extended in Chapter 2 (also in Ref. [23]) to account for dynamic effects following the onset of flow. In the following, we will provide a brief summary of 98 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the methods used in the latter study, as they will be repeatedly used in the present analysis. Consider flow o f a yield-stress fluid and assume for simplicity the following flow rate-pressure drop relationship across a single pore f 1 A i - 1 i A p, k- / Ap( | when | Ap t |> — r. < 7, = 0 when | Ap t |< — (3.3) rt This model states that for flow to occur in pore i, the pressure difference across the pore must exceed a threshold. Following the onset o f flow, the flow rate is taken to vary in proportion to the departure of the pressure from its minimum, following a simple Poiseuille flow relationship. In the above, we have used dimensionless notation as follows 2jug _ R* D R ttQ , P = — P , r = — nyR 4 y R q = — 4 y Q , p =— p ’r = — (3-4) where capital letters denote dimensional and lower case letters denote dimensionless quantities, Q denotes volumetric flow rate, R* is a characteristic pore throat radius (taken later as the mean radius), and jug is gas viscosity (in the absence of foam). In this notation, which will also be used in the rest o f the chapter, the threshold in Eq. 3.3 is the dimensionless analog of Eq. 3.1, and it is inversely proportional to the pore radius. While, the flow rate is that which would result from the application of a pressure gradient equal to 4y/R*. 99 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Consider, now, mobilization and flow driven by a pressure difference applied across a network o f pores. Clearly, for flow of the yield stress fluid to occur, a minimum overall threshold must first be exceeded, A p , m n =m inV— (3 .5) < r i where the summation is over all throats over a connected path. Determining the onset of flow was analyzed in detail by Kharabaf and Yortsos [64] who showed that there exists a specific path, the Minimum Threshold Path (MTP), over which the corresponding sum of thresholds is minimum. Key to their analysis was the development o f a method, termed Invasion Percolation with Memory (IPM), which will also be used repeatedly below. Details on IPM are given in Ref. [64], Here we will simply note that for a uniform distribution of thresholds in the range [0 ,1], the minimum pressure gradient is about 0.3 for a 2-D square lattice [64], In the present case of a uniform distribution of radii in the range [0.5, 1.5], the corresponding value for a 2-D lattice is found to be about 0.96. The IPM algorithm defines the minimum threshold path, given a specified arrangement o f thresholds. The subsequent opening o f new paths, at higher pressure difference, cannot be strictly determined by a static application of IPM alone, however. Following the onset of flow, viscous effects become important, and the algorithm for opening new paths must be modified. Now, opening a new path is determined not only by the throat threshold l/r;, but also by the global velocity field, which affects the pressure distribution along the mobilized paths and 100 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. contributes additional resistance. Thus, even though the first path to be mobilized is indeed the MTP, identified using IPM, finding the subsequent fraction of the mobilized fluid requires a more elaborate approach. In chapter 2, we solved this problem by using updated dynamic thresholds, which depend on the flow rate. For example, the updated dimensionless threshold of a throat k in a flowing path now becomes where r t is the dimensionless threshold of a throat k and qk is the corresponding flow rate. The new paths to open at any incremental pressure above the minimum are determined using IPM, now with all thresholds in the flowing part updated following Eq. 3.6. The process works incrementally, by determining successively at small increments of the applied pressure drop, the corresponding minimum-energy (threshold) paths. In doing so, paths determined at the prior step always remain minimum paths at the next increment. Calculation details can be found in the previous chapter. To determine the pressure distribution due to viscous flow, and hence the flow rates needed to update the thresholds in Eq. 3.6, we assume incompressible fluids, and use Eq. 3.3 (or any other appropriate local expression) for the momentum balance. A mass balance for any site i along the flow path then gives z (3.7) j 101 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where is the volumetric flow rate from site i to its neighbor site j, and Z is the coordination number of the lattice (e.g. Z=4 for a square or Z=6 for a cubic lattice). The resulting equations for pressure were solved using a combination of Conjugate Gradient and Successive Over-Relaxation methods. Because of the inherent non- linearities of the problem, and the need for iterations, the algorithm is time- consuming, however. This limits the networks studied to relatively small sizes. For this reason, some finite-size effects may exist in the results to be obtained below. When the yield-stress fluid is the gas phase in the presence of foam lamellae, non-trivial complications arise due to the fact that the position of lamellae, which dictate the mobilization thresholds, is not fixed [98, 99, 132]. Rather, it constantly varies, both as a result of their movement, as well as a result of their generation and destruction, by the mechanisms described above. In turn, the rates o f generation depend on the pressure distribution, resulting in a complex, dynamic and non-linear problem. In the pore-network simulations to be shown below, lamella generation was incorporated into the model by snap-off, division and/or leave-behind mechanisms, as appropriate. When the gas phase is continuous, the calculation of the minimum pressure gradient is obtained from solving the fluid-flow problem in the continuous gas phase, along with the use of IPM, as described above. When the gas phase is discontinuous, it is obtained by applying IPM. Although quite important, the effect of the pore size distribution was probed only in two cases, one in which the size distribution was unimodal (pore size equal to 1, in the dimensionless notation used) and another in which it was uniform in the 102 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. range [0.5, 1.5]. The relevant dimensionless parameters are the capillary number, CasQw jU w /y, the viscosity ratio M ^uw //ug, the snap-off probability f s and the dimensionless capillary pressure, where appropriate, where subscript w refers to the aqueous liquid phase. The sensitivity of the results to these parameters is the objective of the numerical study. In the following sections we will consider various scenarios of foam generation under different conditions. In all cases, leave-behind and lamella division are standard mechanisms. The differences between scenarios concern how lamellae are generated, e.g. by repetitive Roof snap-off, by Roof snap-off during drainage only, or by capillary-pressure fluctuations. We will focus on foam states corresponding to incipient mobilization, namely always at conditions of minimum pressure gradient, as well as on states corresponding to a fixed applied pressure gradient, as discussed in detail below. 3.5 Lamella Generation during Displacement In this section we will focus on the formation of foam during the process of gas invasion and displacement of an initially liquid-occupied medium. We assume that the process is one during which the liquid flows at a constant rate, by occupying distinct pores in the form of a continuous phase (rather than residing in the comers of pores or as films along pore walls, which is a different case discussed later). For this reason, we expect that the foam generated in this case will be weak. This 103 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. scenario models the simultaneous steady-state flow of liquid and gas. The liquid flow rate is fixed and sets the capillary number. One o f the key objectives is to determine the mobility (or relative permeability) of the foam as a function of the liquid saturation, the capillary number and other parameters, as in a number of previous investigations. To simulate this process, we first use an Invasion Percolation (IP) algorithm (e.g. Feder [38]) to find a continuous path for the gas phase in the pore network. During this process, leave-behind lamellae, and snap-off lamellae at snap-off probability^, are created, as appropriate. The fraction f s represents the probability that a given throat-body combination has the right geometry for Roof snap-off [96]. For simplicity, these probabilities are assigned randomly to pore throats, i.e. independently of throat radius. Upon establishment of this state, namely right after breakthrough, a constant liquid velocity (hence a constant capillary number) is imposed in the liquid phase. Assuming spatially constant capillary pressure, namely that the pressure difference between gas and liquid is the same at both ends, the corresponding gas pressure profile can then be readily determined. If the overall pressure difference across lamellae is sufficiently high, lamellae in the gas phase are mobilized. These lamellae move, divide and create leave-behind lamellae, as appropriate, until a steady state is reached. When this state is reached, the gas flow rate corresponding to the pressure drop and the gas mobility are computed. In these and similar calculations, if several lamellae are moving, the first lamella to divide is the one requiring the minimum amount of time to move to a neighboring throat. 104 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. The lamella may be also divided into two or more, depending on the occupancy of throats and pore bodies by gas or liquid. The condition is also enforced that lamella generation by division occurs only in downstream throats. Then, the pressure drop in the gas phase is incrementally increased. The gas-liquid interface, corresponding to the maximum value of the dimensionless pressure difference f 1 \ max P * ~ P ' - T r is identified and advances into the liquid, displacing the liquid from the medium and reducing its saturation. In the process, new snap-off and/or leave-behind lamellae may be created. At the new pressure increment, the possibility of mobilization o f lamellae in the gas phase is sought, and then the previous process is repeated. We note that there is no destruction of lamellae in this study, which focuses on foam generation. Lamellae leave the region of interest by advection, but not by destruction. Figure 3.8 illustrates this scenario by showing patterns of the lamellae distribution for a 100x100, 2-D pore network with Ca=0.01, and at different values of the applied pressure (hence at different liquid saturations). In these illustrations, no lamellae were generated by leave-behind until gas breakthrough (top panel of Figure 3.8), although lamellae were generated by snap-off. However, this is not a necessary condition for this scenario. With the invasion of gas-liquid interfaces deep into the liquid phase, leave-behind and snap-off lamellae are generated. In the process, new snap-off lamellae may isolate the gas-liquid interface from the bulk 105 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 3.8: Sequence of patterns illustrating liquid displacement in a 2-D network and the creation of lamellae. The top panel shows the conditions at breakthrough, the middle panel conditions at a higher value of the applied pressure drop. The bottom panel is an enlargement of a region in the middle panel showing a number of lamellae generated by snap-off and leave-behind mechanisms. The liquid flows at a constant rate corresponding to Ca=10'2 . Gas is white, liquid and rock grains are black. Extensive lamella formation by leave-behind occurs as the liquid is being displaced. 106 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. gas. We assume that this type of interface will not be moved, unless the lamellae can become mobilized. Using the above approach, as the gas pressure difference increases, the liquid saturation decreases, at constant capillary number, and the relative mobility of the foam can be computed as a function of the liquid saturation. 0.6 - Liquid - Gas 0.5 0.4 0.2 0.0 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Sw Figure 3.9: Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability ^=0.1, M= 100 and Ca=10'3. Relative permeabilities (kr ) as a function of the applied pressure drop. The liquid saturation Sw falls asp rises. 107 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.5 Gas Leave-Behind Snap-Off Lamella-Division 0.4 0.3 l l 0.2 0.0 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 sw 0.025 0.020 0.015 > 0.010 0.005 0.000 ^ 0.00 0.03 0.06 0.09 0.12 0.15 IVpl Figure 3.10: Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability ^=0.1, M=100 and Ca= 10'3 . The fraction F of pore throats containing gas without lamellae, leave-behind lamellae, snap-off lamellae, and lamellae created by division, and the flow velocity-pressure drop relation, as a function of the applied pressure drop. The liquid saturation Sw falls as A p rises. 108 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.9 and Figure 3.10 show results obtained for a 3-D pore network (of size 20x20x20) and for the parameter values f s=0.l, M= 100, and Ca=0.001. Note that Ca as defined here depends on the liquid injection rate, not the pressure drop, and M reflects the v isco sity ratio betw een gas and liquid w ith no foam present. Thus a large value of M applies, even if foam is created and reduces the mobility of gas. Because of the small value of the snap-off probability in this case, the gas phase remains continuous throughout. The figures show that as the gas-phase pressure drop increases, the liquid saturation decreases, and the number of lamellae in the gas increases, mostly through the leave-behind mechanism. This significant number of leave-behind lamellae may be due to the larger coordination number for this 3-D network [101]. The creation of lamellae in the gas phase causes the relative permeability o f the gas to be considerably smaller than in the absence of lamellae, under the same conditions. In this case with a continuous gas path always available, the rheology o f the gas remains approximately Newtonian (in that the flow velocity remains roughly proportional to the pressure drop). An increase in the pressure drop leads to more lamellae being mobilized, but also to the displacement of the liquid and the potential opening of new flow paths. It must be noted that no shear- thinning mechanism was considered in the flow of the continuous gas phase, where unblocked by lamellae. Note that the velocity on the vertical axis is related to the dimensional superficial velocity V through the expression Tf v = —-^ -V (3.8) <k 109 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. where $ is porosity. Alternatively, the dimensionless velocity v is the relative velocity compared to the velocity the gas would have if it flowed in a porous medium consisting of identical pores of size R * under a pressure gradient equal to 4y/R*. 0.6 Liquid — - - G a s 0.5 0.4 0.2 0.0 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 SL Figure 3.11: Steady-state results corresponding to foam generation during the displacement o f liquid by gas, with snap-off probability f s=0.7, M= 100 and Ca= 10'3. Relative permeabilities (kT ), the fraction F o f pore throats containing gas without lamellae, leave- behind lamellae, snap-off lamellae, and lamellae created by division, and the flow velocity- pressure drop relation, as a function o f the applied pressure drop. Note the existence o f a critical gas saturation (and a critical pressure gradient), below which the gas relative permeability vanishes. 1 1 0 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.35 Gas Leave-Behind Snap-Off Lamella-Division 0.30 0.25 0.20 L l 0.15 0.10 0.05 0.00 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Sw 0.0020 0.0016 0.0012 > 0.0008 0.0004 0.0000 ■ — 0.00 0.04 0.08 0.12 0.16 IVpl Figure 3.12: Steady-state results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability f= 0.7, M= 100 and Ca=l0'i. The fraction F of pore throats containing gas without lamellae, leave-behind lamellae, snap-off lamellae, and lamellae created by division, and the flow velocity-pressure drop relation, as a function of the applied pressure drop. Note the existence of a critical gas saturation (and a critical pressure gradient), below which the gas relative permeability vanishes. I l l R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Under the same conditions, an increase in the snap-off probability fraction leads to a corresponding increase in the lamellae fraction and a decrease in the relative mobility of the gas (Figure 3.11 and Figure 3.12, with f s=0.7). The qualitative features of the foam generated are much more pronounced than previously. First, the relative permeability of the gas phase is much lower. Now, the gas phase is disconnected and cannot be mobilized until after the liquid saturation has been reduced significantly, and the liquid has retreated to narrow throats. In fact, there exists critical gas saturation, equal to about 0.65 in this case, below which gas flow does not occur. Properties of this critical value are currently under study. The trend of decreasing gas mobility with an increase in the snap-off probability is expected. Second, there exists a minimum pressure drop for the mobilization of the gas phase, as expected for a discontinuous gas, corresponding to conditions of a yield stress fluid. This value is a function of the capillary number, with higher capillary numbers leading to an earlier onset of mobilization. The resulting rheology is quite interesting. Due to the existence of many lamellae, a yield-stress behavior is exhibited, where a minimum pressure gradient is necessary for the onset of gas flow. Following the onset of mobilization, however, the system will ultimately reach a steady state, in which a continuous path becomes available. At that point, the flow rate increases rapidly, and the flow velocity-pressure drop relationship almost becomes discontinuous (as shown in Figure 3.12). Finally, the fraction of lamellae generated by snap-off is considerably larger than in the case of a low snap-off probability, as expected, although that by leave-behind is about the same. 112 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. We note, nonetheless, that without a mechanism of lamella creation near the inlet, lamellae must be mobilized until a continuous gas path forms there. The effect o f the capillary number is shown in Figure 3.13, under the same conditions as Figures 3.9, 3.10 and Figures 3.11, 3.12, but where the liquid flow rate has now increased by one order of magnitude. The most interesting feature of this figure is that the relative mobility of the gas is higher, at the same liquid saturation. This result is consistent with the interpretation above, namely the fact that at a higher pressure difference more lamellae are mobilized, with a resulting larger flow rate in the continuous path of the gas phase. This also reflects the paradox of foam generation by lamella mobilization: without a mechanism to recreate lamellae near the inlet, lamellae are advected downstream until a continuous-gas path forms near the inlet. The increase in the relative permeability also reflects the well-known beneficial relative-permeability effect of the capillary number for any immiscible displacement. We note that the two plots corresponding to the fraction of lamellae generated and the flow-rate pressure drop relationship for this case are similar to Figure 3.10 and Figure 3.12, and for this reason are not shown. 113 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.6 Liquid Gas 0.5 0.4 0.3 0.2 0.0 0.8 0.9 0.2 0.3 0.4 0.5 0.6 0.7 Sw 0.6 Liquid Gas 0.5 0.4 i_ 0.2 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Sw Figure 3.13: Steady-state relative permeability results corresponding to foam generation during the displacement of liquid by gas, with snap-off probability .£=0.1 (top panel) and £=0.7 (bottom panel), M=100 and Ca=10'2 . The increase in the liquid rate (capillary number) leads to increased relative permeabilities in both fluids. 114 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 3.6 Lamella Generation Due to Capillary-Pressure Fluctuations The second process we considered involves the generation of lamellae when the liquid saturation has been reduced to a small value, and liquid resides in pore comers and as films along pore walls. In this section, we assume that lamellae are generated by snap-off only initially. Subsequent lamella generation takes place when the capillary pressure in any pore throat drops below half of the entry value of that throat, as discussed above. A scenario of repeated Roof snap-off at germination sites, independently of pressure or flow rate, is discussed later on. The present mechanism requires calculation or estimation o f the liquid pressure field. While this was computed in the previous section using the pore network, here we will proceed with a simplification, namely by assuming a uniform liquid pressure gradient along the flow direction. For this, the liquid pressure at the exit is set to zero, without loss, that at the entrance to a fixed value, pio, and in-between it is linearly interpolated along the main flow direction. For a 2-D network of size a corresponding liquid capillary number can be roughly estimated, if one assumes that the liquid flows in the pore walls in the form o f macroscopic films of dimensionless thickness < 5 . Then, the capillary number is roughly equal to Ca=4adpio/L where a is a dimensionless parameter that relates permeability to the hydraulic diameter, and is assumed to be of the order of 10'2. In the following, we will assume <5=0.1, for illustration purposes. The gas-phase pressure downstream is set equal to a capillary pressure value, p c, and for the present purposes we will 115 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. assume that the gas viscosity is very small. Clearly, a sufficient condition for lamella generation is, in our dimensionless notation, LCa 1 1 P c - P i o = P c - - r ^ < - A = - = 0.1666 (3.9) ^ccS 4rm a x 6 while a sufficient condition for no lamella generation is LCa 1 1 c p ' - p '” = p ' - ^ > 4 ^ : = 2 = 0 -5 (310) If condition (3.10) applies, there is no lamella (or foam) generation by snap-off by the method of fluctuatingp c. When condition (3.10) is not satisfied, lamellae are generated, when the gas phase is continuous, mostly upstream. The upstream generation can be readily shown, as the capillary pressure at conditions of continuous gas phase is Pg-Pi= P c-P io + T P io (3-11) showing that it is easier to satisfy condition (3.9) upstream. Using the above interpretation, lamella generation will be favored at a low capillary pressure and at a high capillary number for the liquid. With this same interpretation, and for sufficiently large systems, condition (3.9) is always in effect. Following this scenario for lamella generation, we determined the state of the foam by proceeding in two different ways, either by operating always at the state of minimum pressure gradient (the state o f incipient mobilization) or by operating at a pressure gradient above the minimum. 116 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.4 0.3 0.0 i. 0 2000 4000 6000 8000 Step Figure 3.14: The minimum pressure gradient corresponding to the state of incipient mobilization, in which lamellae are generated by capillary pressure fluctuations. (f=0.3, p c =5A677, pio=7.6677, Ca=7.677x1c4 , L=40). The value fluctuates between two extreme states, one in which the gas has a continuous path (corresponding to a zero value of the minimum pressure gradient) and a state in which the gas is completely discontinuous. In- between states correspond to the displacement of lamellae from the blocked paths and the creation of a continuous path. The number of steps corresponds to the number of lamellae moved. 1 1 7 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.15: Lamella patterns corresponding to Figure 3.14. The system fluctuates between two extreme states, one in which the gas has a continuous path (shown above as a white line) and a state in which the gas is completely discontinuous. When the first state is reached, the capillary pressure drops and new lamellae are generated in the continuous path. These are advected downstream until a continuous gas path opens. For emphasis, in the top figure gas in the continuous path is white and the remaining gas gray; otherwise, gas is white, liquid and rock grains black. 118 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.40 0.35 0.30 0.25 0.15 0.10 0.05 0.00 10000 15000 20000 0 5000 Step Figure 3.16: The minimum pressure gradient corresponding to the state o f incipient mobilization, in which lamellae are generated by capillary pressure fluctuations. (£=0.7, p c =5.1677, pm=l.6677, Ca=7.677xl0"4, Z=40). The value fluctuates between two extreme states, one in which the gas has a continuous path (corresponding to a zero value o f the minimum pressure gradient) and a state in which the gas is completely discontinuous. In- between states correspond to the displacement o f lamellae from the blocked paths and the creation o f a continuous path. The number o f steps corresponds to the number o f lamellae moved. 1 1 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.17: Lamella patterns corresponding to Figure 3.16. The system fluctuates between two extreme states, one in which the gas has a continuous path (shown above as a white line) and a state in which the gas is completely discontinuous. When the first state is reached, the capillary pressure drops and new lamellae are generated in the continuous path. These are advected downstream until a continuous gas path opens. For emphasis, in the top figure gas in the continuous path is white and the remaining gas gray; otherwise, gas is white, liquid and rock grains black. 120 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.8 0.7 0.6 0.5 0.3 0.2 0.0 -10 ■ 8 ■ 6 ■ 4 • 2 0 -14 -12 Pc-PlO 0.5 0.4 - Q . > 0.2 - 0.1 - 0.0 1 1 1 1 -------------------- 0.0 0.2 0.4 0.6 0.8 1.0 fs Figure 3.18: The minimum pressure gradient at incipient mobilization corresponding to the maximum value of the fluctuating cycle between the states, as a function of the parameter P c-P io (for f s =03), at the top panel, and of the initial snap-off probability f s (for pc -pi0 =-2.5 corresponding to point A in the top panel) at the bottom panel. 121 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figures 3.14-3.18 show characteristic results for a 2-D pore network and for the scenario in which the system of lamellae is always at a state of incipient mobilization. At this state, we cannot strictly speaking define the rheology of the gas phase, as its flow rate is alw ays zero (state o f incipient m obilization). Figure 3.14 shows the value of the corresponding minimum pressure gradient (at incipient mobilization), which is shown to fluctuate and ultimately settle into a periodic behavior between two extreme states. Corresponding lamellae patterns at these two states are shown in Figure 3.15. These extremes correspond to an open (continuous) gas path, where lamellae are generated due to capillary pressure differences, and to a completely plugged state, where the local capillary pressure is high and no new lamellae are generated. In the latter state existing lamellae become mobilized upon consecutive increases in the pressure gradient (recall that the system is always at the incipient mobilization state), by means of which lamellae are advected downstream and eventually open a continuous gas path across the medium. This signals the rapid drop of the gas pressure and the subsequent generation of lamellae through the capillary-pressure mechanism. In reality, this change in capillary pressure, assumed to be instantaneous here, would take a short but finite period o f time, as liquid imbibes in response to the change in local gas pressure. Figures 3.16 and 3.17 show similar results for a different snap-off probability. The cycle is then repeated and involves a number o f intermediate states. As variable p c-pw decreases, the density of new lamellae generation increases, provided 122 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. that condition (3.9) is satisfied. A plot of the maximum value of the minimum pressure gradient (namely of the highest value of the quantity in Figure 3.14) as a function of p c-pio and of the initial snap-off probability fs, is shown in Figure 3.18. As expected, the maximum value of |V/?| increases as p c-pio decreases to negative values, ultimately reaching a limiting value (plateau) at sufficiently negative values of the difference p c-pio- This plateau is equal to about 0.96 and corresponds to the limit where all throats contain lamellae. The effect of the initial snap-off probability f s is not significant, however. It fluctuates around a mean value independent o ff s. Apparently, as long as there is significant potential for lamellae generation by capillary pressure fluctuations, the initial snap-off probability is not an important factor in the ultimate behavior of the system. Consider, next, the case in which the applied pressure gradient is above the minimum mobilization pressure, although a more appropriate scenario would be to operate at a constant gas flow rate. Contrary to the case of a Newtonian fluid, or even a Bingham plastic fluid, this does not lead to a steady state, however. Following the creation of a discontinuous gas phase, through the initial snap-off mechanism, no new lamellae can be generated, but lamellae are advected downstream (and leave the system), until a continuous path for the gas opens up, at which time the capillary pressure mechanism comes into effect and new lamellae form in the open path. This leads to the disconnection o f the gas phase, and the repetition of the cycle. We expect, therefore, the system to operate in a fluctuating state, consisting of a branch of low-pressure drop, where the gas phase is 123 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. continuous and the number of lamellae minimal, and a branch of high-pressure drop, where the gas phase is discontinuous and the number o f lamellae large, but slowly decreasing, as the lamellae are convected downstream. Ultimately, a continuous path opens up and leads to the generation of lamellae upstream, and the cycle continues. Working at constant gas flow rate is not computationally easy, however. For this reason, we simulated the sequence along the branch of high pressure drop by fixing the applied pressure drop to a constant value above the first minimum pressure gradient value and keeping it constant, thereafter. Figure 3.19 and Figure 3.21 show two different snapshots o f the lamella pattern at two different pressure gradients, which illustrates the mechanism discussed above. During this process, and until an open path develops, the flow rate increases with time (as shown in Figure 3.20 and Figure 3.22 for two different pressure gradients) in a fluctuating manner, until a continuous gas phase path opens. Then, new lamellae are generated and the cycle commences anew. We note that in Figures 3.20 and 3.22 time is made dimensionless by dividing by the time it takes for a passive tracer to be advected across the lattice, assuming a unimodal pore-size distribution, and based on the same applied pressure difference. Evidently, in this case, the time dependence of the flow rate on the applied pressure drop precludes the use of a constant mobility for the foam. 1 2 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.19: Lamellae patterns corresponding to two different states under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations (£=0.7, pc =S.\677, pi0 = 7.6677, |V/?|=0.17, Ca=6.8xl0'4 ). The system fluctuates between two extreme states, one in which the gas has a continuous path (shown in the bottom figure as a white line) and a state in which the gas is completely discontinuous (top figure, corresponding to the initial state, but which is also similar to one of the two extremes of the asymptotic state). Gas is gray or white, liquid and rock grains black. 125 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.0025 0.0020 - > 0.0000 0.0015 - 0.0010 - 0.0005 - 0 2 3 Time 4 5 6 Figure 3.20: The gas flow rate as a function of time for the process of Figure 3.19, namely under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations. One complete cycle is shown. At the end of the cycle, a continuous gas path opens, new lamellae are generated by capillary pressure and the process is repeated. The increasing flow rate as a function of time results from the mobilization and the advection downstream of lamellae in a pore network in which the gas phase is disconnected (/j=0.7,/>c =5.1677,p w = 7.6677, |Vp|=0.17, Ca=6.8xl0'4 ). 1 2 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.21: Lamellae patterns corresponding to two different states under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations (£=0.7, pc =5.\611, pi0 = 7.6677, |Vp|=0.275, Ca=l.lxl0'3 ). The system fluctuates between two extreme states, one in which the gas has a continuous path (shown in the bottom figure as a white line) and a state in which the gas is completely discontinuous (top figure, corresponding to the initial state, but which is also similar to one of the two extremes of the asymptotic state). Gas is gray or white, liquid and rock grains black. 127 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.008 0.006 > 0.004 0.002 0.000 1 1 1 1 ------------- 0 1 2 3 4 5 Time Figure 3.22: The gas flow rate as a function of time for the process of Figure 3.21, namely under the application of a constant pressure drop and for the scenario in which lamellae are generated by capillary pressure fluctuations. One complete cycle is shown. At the end of the cycle, a continuous gas path opens, new lamellae are generated by capillary pressure and the process is repeated. The increasing flow rate as a function of time results from the mobilization and the advection downstream of lamellae in a pore network in which the gas phase is disconnected ,p c =5.\611 ,p m = 7.6677, |Vp|=0.275, C a=l.lxl0'3 ). These figures illustrate two important results. First, snap-off o f lamellae near the inlet to the medium replenishes the supply of lamellae after mobilization and displacement o f lamellae downstream. Mobilization and division, combined with snap-off driven by periodic low pc near the inlet, can maintain foam at a sort of steady (albeit fluctuating) state. Repeated Roof snap-off, governed by pore throat- body aspect ratio, is therefore not necessary for steady foam regeneration. Second, 128 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the pressure drop fluctuates around this steady state. As is well known, fluctuations are an inherent part of the regeneration process for foam. The fluctuations represented here are related to, but not identical to, those cited by Rossen [99], Those fluctuations arise from the changes in shape and curvature of a fixed population of lamellae moving through the converging-diverging geometry of pores. The fluctuations represented here arise from the changing population of lamellae: creation of new lamellae by snap-off and subsequent displacement of these lamellae out of the region of interest. This process represents a new mechanism of snap-off and foam generation. 3.7 Repetitive Snap-off at Germination Sites The final case studied in this chapter involves repetitive snap-off o f lamellae in specific germination sites, randomly distributed in the pore network, according to snap-off probability f s. This scenario mimics foam lamella generation following the Roof snap-off criterion, and has been postulated as plausible in the previous literature [26, 42, 48, 55, 65, 69-71, 82, 111, 127, 133]. For example, it was used in a previous mechanistic pore-network study to model foam invasion in a porous medium [65]. Its key feature is that as lamellae are mobilized from germination sites, new lamellae are instantly generated to replace them. The model corresponds essentially to quasi-single-phase (gas and lamellae), as the liquid does not affect the process other than by supplying the small amount of liquid needed to form 129 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. lamellae. Given that lamellae are constantly generated, strong foam is expected. As discussed above, the state of the foam depends on the applied pressure gradient. Again, two different states were considered: (i) The foam state is determined by sequentially applying the minimum pressure gradient required to mobilize at least one lamella (state of incipient mobilization); (ii) The foam state is determined by the application of increments of higher pressure gradient. The second case will also result into a rheological characterization of this foam. As noted, the mechanism of repeated Roof snap-off after the pore body drains of liquid is problematic [103]. Assuming repeated Roof snap-off here allows comparison to the cases examined above, and, in one case, maps roughly on to the assumptions o f the foam-generation model of Rossen and Gauglitz [104]. 3.7.1 State o f Incipient Mobilization (Minimum Pressure Gradient) In this scenario, lamellae are initially distributed in the throats of the pore network according to probability f s. The gas phase may or may not be continuous, depending on the value of f s (percolation theory for a 2-D square lattice would require the condition/v <0.5 for a continuous gas). The algorithm to be implemented determines the minimum pressure gradient to mobilize at least one lamella, following which new lamellae are generated by lamella division mechanisms. In this case, and contrary to the previous, when a mobilized lamella leaves the 130 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. germination site, a new lamella is created to replace it at the same site. Then, the minimum pressure gradient to mobilize the new state is calculated and applied across the network, following which lamellae are mobilized, new lamellae are generated, as above, etc. Mobilized lamellae at the exit end leave the pore network; therefore we anticipate that an asymptotic state will be ultimately reached. One question is the nature of this state, namely whether it is steady, oscillatory or random. The pressure drops computed in this section do not represent pressure drop that would be measured in a laboratory experiment at fixed injection rate. Rather, they are the minimum pressure drop required to mobilize a lamella on the pore network. If the gas phase is continuous, and the actual pressure drop is less than this value, gas flows with no lamella mobilization. If the gas phase is discontinuous, and the actual pressure drop is less than this value, the gas phase would not flow. These results illustrate the dynamics of the minimum pressure drop for mobilization, but not the actual observable pressure drop. Typical results for this scenario are shown in Figures 3.23-3.28 corresponding to different snap-off probabilities f s for the case of a uniform pore-size distribution. Figure 3.23a shows the minimum pressure difference as a function of the state of the foam, with the corresponding numbers of lamellae generated by snap-off, lamella division and leave-behind shown in Figure 3.23b. They correspond to f s=0.2. The figures show that after about 200 steps, the system locks into a periodic state. While the number of snap-off lamellae rises only slightly above its initial 131 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. value, the number o f lamellae generated by division and leave-behind increases significantly, as a result of lamella mobilization. The minimum pressure gradient ultimately reaches a state fluctuating between about 0.17 and about 0.02, with a rather intricate pattern in between. In this case of low snap-off probability, the fluctuation was found to be periodic, with a period of 13 for the specific case of Figures 3.23-3.25. Figure 3.24 shows snapshots of the initial state and of one of the final states. While the generation of additional lamellae is evident, the gas phase is still continuous. The bottom panel in Figure 3.24 shows schematically the path being mobilized. It consists of a “lamella streamtube” (a schematic of which is shown in Figure 3.25) containing a train of lamellae, the sequential mobilization of which follows the pressure gradient fluctuation. Essentially, the system ultimately settles into a pattern in which all but the lamellae along the path depicted in the figure remain stationary. The fluctuating pressure gradient reflects the following: The mobilization path is “shielded” on all sides by lamellae that are stationary, or are generated by division as a mobilized lamella moves downstream (see also Figure 3.25). Thus, the flowing gas pressure is always evaluated at the same point. The fluctuating pressure drop is due to the fact that as lamellae move to new locations in the “lamella streamtube”, they encounter pore throats of different sizes. Figures 3.26-3.27 show corresponding results for a higher snap-off probability, f s=0.1, where in both the initial and the final states the gas phase is discontinuous. The patterns suggest the generation of strong foam, as the gas phase is discontinuous at all times. The density of snap-off lamellae and lamella by division 132 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. increases with the number of steps. Here it takes a much higher number of steps to reach the asymptotic state, compared to the case of low snap-off probability. As with the previous case, the asymptotic state is fluctuating, but here in an apparently random and chaotic manner. By comparison, the fluctuations in the previous case were periodic. It is of interest to note that many times the number o f the lamellae in the system is mobilized as this fluctuation regime is spanned. A chaotic analysis of the time series was conducted and it is found that the chaotic nature is weak. The fluctuation (Figure 3.26a) reveals a rather complex sequence of lamellae mobilization and division. Different values of the snap-off probability showed the same type of behavior, with lower values of f s leading to continuous gas phase and an oscillatory behavior, but with higher f s leading to discontinuous gas and a random fluctuation. The particular distribution of germination sites had an effect on the fluctuation patterns at the asymptotic state, although the average behavior remained similar. Figure 3.28 plots the minimum pressure gradient at the asymptotic state (i.e., sufficient to maintain mobilization indefinitely) as a function of 1 -F, where F is the final fraction of throats occupied by lamellae. For continuous-gas foam, and at low values of F (roughly below 0.5), if the actual pressure gradient is less than this value, then gas flows continuously, with no lamella mobilization and no generation of stronger foam. For discontinuous-gas foam, resulting from higher values of F, gas flow stops if pressure gradient is less than this minimum. 133 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.25 0.20 0.15 Q . > 0.10 - 0.05 - 0.00 1000 0.8 0.6 Gas Snap-Off Lamella-Division L l 0.4 0.2 0.0 800 1000 0 200 400 600 Step (b) Figure 3.23: Foam generation by repetitive snap-off at germination sites (/s=0.2, uniform pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization in panel (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is continuous in all states. 134 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6969657162^0087104092925 Figure 3.24: Foam generation by repetitive snap-off at germination sites (f=0.2, uniform pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states of the Figure 3.23a. The gas is continuous in all states. 135 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Y J V . j — \ J V . Figure 3.25: Foam generation by repetitive snap-off at germination sites (f=0.2, uniform pore-size distribution): Simplified schematic illustrating the concept of a lamella streamtube. The arrows in the latter indicate the mobilized lamellae and the direction in which they are mobilized. In the first panel, the moving lamella divides into two unoccupied downstream throats. In the second panel, the mobilized lamella moves to the throat on the left. In the third and consecutive panels (not shown here for simplicity), the lamella moves in the direction of the arrow and occupies consecutive pores without any new lamella generation. The pressure gradient fluctuates as a result of the fluctuating pore size. In the fourth and last panel, the lamella exits the system, thus recreating the top panel. This cycle is then repeated. 1 3 6 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0.8 0.6 CL > 0.4 0.2 0.0 2000 4000 6000 8000 10000 0 Step (a) Gas Snap-Off Lamella-Division 0.6 L L 0.4 0.2 0.0 2000 4000 6000 8000 10000 0 Step (b) Figure 3.26: Foam generation by repetitive snap-off at germination sites (£=0.7, uniform pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is discontinuous in all states. 137 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 3.27: Foam generation by repetitive snap-off at germination sites (/s=0.7, uniform pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states o f Figure 3.26a. Gas is white, liquid and rock grains are black. 138 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.8 0.6 Q. > 0.4 0.2 0.0 1.0 0.4 0.6 0.8 0.2 0.0 1-F Figure 3.28: The minimum pressure gradient as a function of the final fraction of the pore throats occupied by a lamella for the case of foam generation by repetitive snap-off at germination sites and uniform pore size distribution. Compare to the plot in Figures 3.6 and 3.7 (noting that F=\-f, and that here the results correspond to a network with different coordination number than that in Figures 3.6 and 3.7). Figure 3.28 thus maps approximately on to the theory o f Rossen and Gauglitz [104] without assuming uniform |VP| on the pore network. As Rossen et al. [106] conjectured, without this assumption the two halves of the theory smoothly merge near the percolation threshold. There are three differences between Figure 3.28 and the previous theory. First, Rossen and Gauglitz consider the onset of mobilization, 139 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. where here mobilization, division and repeated Roof snap-off lead to a sort of steady state where the minimum pressure gradient fluctuates in time around a fixed value. Second, the percolation threshold for the square network is 0.5, rather than 0.25 for the network in earlier studies (corresponding to a different coordination number). Third, the two branches of the theory are asymmetric; the minimum pressure gradient in the continuous-gas regime is much lower than that in the discontinuous-gas regime, compared to the earlier theory (Figures 3.6 and 3.7). Figures 3.29-3.32 show corresponding results for a slightly different case, where the lamellae moved out the system flow back into the system from the inlet. It is not surprising the flow patterns and also the pressure fluctuations are not much different. It may be caused by the assumption we have, i. e., any time there is only one lamella in a bond. If a lamella returns to a bond with a lamella in it, we’ll ignore a second one. Figures 3.33-3.36 reveal corresponding results for the case when the pore-size distribution is unimodal, namely of only a single pore size. Two key features are evident: While there is a randomly fluctuating initial rising part, as in the previous case of variable pore sizes, the asymptotic state in the case o f discontinuous gas is not fluctuating, but rather constant. A fluctuating periodic behavior was observed sometimes, for low values of the snap-off probability, although not in all cases. The final value was approximately the same to that for the variable pore size, with the deviation being more pronounced at low values of the snap-off probability f s. Indeed, it appears that at low values of the snap-off probability the minimum 140 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. pressure gradient was much higher than in the random case. Thus, the disorder in pore size manifests not only in altering the nature of the asymptotic state, from a fluctuating to a constant, but, in the case of low snap-off probability, also alters the value of the minimum pressure gradient. This trend towards higher values of the minimum pressure gradient is closer to the theory of Rossen and Gauglitz [104], which is based on an assumption of homogeneity. For completeness, we note that if we restrict the generation of snap-off lamellae only to the beginning of the process, the foam thus generated is weak; namely the gas phase becomes continuous, regardless of the value of the snap-off probability. Without allowing for snap-off near the inlet of the medium triggered by low pc, Roof snap-off appears to be necessary for continuous foam regeneration. The role of repeated Roof snap-off in the theory of this section is played by snap-off driven by low pc in the preceding section. The asymptotic states without lamella creation by either fluctuating pc or repeated Roof snap-off are oscillatory, with a well- defined period of oscillation, e.g. see Figure 3.37 for / s =0.2 and f s-0.7. The reason for this behavior is due to the fact that as lamellae leave the system they are not replaced by new snap-off lamellae and, although new lamellae are generated though lamella division, these are not sufficient to block the path o f the continuous gas phase and create strong foam. 141 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.25 0.20 0.15 C L > 0.10 0.05 0.00 1000 400 600 800 0 200 Step (a) 0.8 Gas Snap-Off Leave-Behind 0.6 Li. 0.4 0.2 0.0 200 400 600 800 1000 0 Step (b) Figure 3.29: Foam generation by repetitive snap-off at germination sites (£=0.2, uniform pore-size distribution, lamellae moved out flow back in from inlet): The minimum pressure gradient corresponding to the state of incipient mobilization (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is discontinuous in all states. 142 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.30: Foam generation by repetitive snap-off at germination sites (£=0.2, uniform pore-size distribution, lamellae moved out flow back from inlet): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.29a. Gas is white, liquid and rock grains are black. 143 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.6 CL > 0.4 0.2 0.0 10000 4000 6000 8000 0 2000 Step (a) 0.8 Gas Snap-Off Lamella-Division 0.6 LL 0.4 0.2 0.0 2000 4000 6000 8000 10000 0 Step (b) Figure 3.31: Foam generation by repetitive snap-off at germination sites (£=0.7, uniform pore-size distribution, lamellae moved out flow back in from inlet): The minimum pressure gradient corresponding to the state of incipient mobilization (a), and the fraction of pore throats occupied by continuous gas or lamellae of the various types (b). The gas is discontinuous in all states. 144 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.32: Foam generation by repetitive snap-off at germination sites (fs =Q.l, uniform pore-size distribution, lamellae moved out flow back from inlet): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.31a. Gas is white, liquid and rock grains are black. 145 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.4 0.3 > °'2 0 .1 - 0.0 1 1 1 1 ----------------- 0 100 200 300 400 500 Step (a) 0.8 — Gas - Snap-Off ■ Lamella-Division 0.6 LL 0.4 0.2 0.0 500 0 100 200 300 400 Step (b) Figure 3.33: Foam generation by repetitive snap-off at germination sites (f=0.2, unimodal pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization and the fraction of pore throats occupied by continuous gas or lamellae of the various types. The gas is continuous in all states. 146 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.34: Foam generation by repetitive snap-off at germination sites (£=0.2, unimodal pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.33a. The gas is continuous in all states. Gas is white, liquid and rock grains are black. 1 4 7 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 1.0 0.8 0.6 CL > 0.4 0.2 0.0 ' ' 1 '------------- 0 100 200 300 400 500 Step (a) 0.8 Gas Snap-Off Lamella-Division 0.6 LL 0.4 0.2 0.0 300 400 500 0 100 200 Step (b) Figure 3.35: Foam generation by repetitive snap-off at germination sites (£=0.7, unimodal pore-size distribution): The minimum pressure gradient corresponding to the state of incipient mobilization and the fraction of pore throats occupied by continuous gas or lamellae of the various types. The gas is discontinuous in all states. 148 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure 3.36: Foam generation by repetitive snap-off at germination sites (£=0.7, unimodal pore-size distribution): Two lamella patterns corresponding to the initial state and one of the final states of Figure 3.35a. Gas is white, liquid and rock grains are black. The gas is discontinuous in all states. 149 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.14 0.12 0.10 _ 0.08 Q . 0.06 0.04 0.02 0.00 300 400 500 600 0 100 200 Step Q _ > 2000 Figure 3.37: Foam generation by initial snap-off at specific germination sites, corresponding to f= 0.2 (top) and f= 0.7 (bottom). The minimum pressure gradient corresponding to the state of incipient mobilization is shown. The gas is continuous in all states, and very weak foam is observed. 1 5 0 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. ^ 3.7.2 Sequentially Increasing Pressure Gradient The final section considers a different approach to the formation of foam, in the case of repetitive snap-off, namely increasingly higher pressure gradients are applied. In this case, the flow rate-pressure drop relation, thus the foam rheology, can be determined. To understand the results to be obtained, consider, for example Figure 3.23a, and assume that initially a small pressure gradient is applied. Because of the initially continuous gas phase, a small pressure gradient leads to a corresponding gas flow rate, without mobilizing any lamellae. As a result, at a sufficiently small pressure drop, the flow rate-pressure drop relation is linear. As the threshold of mobilization is reached and slightly exceeded, lamellae become mobilized, new lamellae are generated, as before, and the mobility o f the gas phase now decreases. As long as the minimum pressure gradient is monotonically increasing, e.g. as in the early part of Figure 3.23a, the state of the foam would be the same as previously. However, Figure 3.23a also shows that at some point, the minimum pressure gradient decreases with time. At this point, in the scenario to be followed here, the pressure difference is set equal to that in the last time step. Clearly, under these conditions multiple lamellae, and not only one, will be mobilized. In this region, the gas phase flows as both a continuous and as a discontinuous phase. With this condition, the system will eventually reach an asymptotic state, following which the pressure gradient is increased again. On the other hand, if the snap-off probability is higher (e.g. as in Figure 3.26a), the gas phase is discontinuous from 151 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. the outset, hence a minimum pressure gradient must be applied to mobilize lamellae. Clearly, this value will be the same as the first value in Figure 3.26a. The subsequent sequence is then followed as before. The flow velocity-pressure drop relation resulting from this sequence is shown in Figures 3.38-3.39 for two different values of the snap-off probability (fs=0.3 and f s- 0.7). In the first case, Figure 3.38, the gas phase is continuous, and the fluid behaves initially as a Newtonian fluid, as all lamellae are stationary. As the pressure gradient increases, lamellae are being mobilized, and at higher values of pressure the gas flows both as a continuous and as a discontinuous phase. In that region, the rheology of the foam appears to be shear-thickening (again, we stress that we did not explicitly introduced a shear-thinning behavior for the flow of gas). The results for the case of higher snap-off probability are different. For a sufficiently high probability, an increasing number of lamellae must be mobilized before the system can flow (e.g. see Figure 3.39). As these are mobilized they create lamellae by division, which makes more difficult the further mobilization, unless the pressure is increased, when more lamellae become mobilized, and more lamellae created by division. The final result is that all throats are rapidly filled with lamellae, in which case the system becomes identical to that of the flow of a Bingham plastic, which is the problem discussed in Chapter 2. The corresponding flow velocity-pressure drop relationship is identical to that o f the flow of a fluid with yield stress. 152 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.25 0.00 1 1 1 ---------------- 0 500 1000 1500 2000 Step 0.035 0.030 0.025 0.020 > 0.015 0.010 0.005 0.000 0.00 0.03 0.06 0.09 0.12 0.15 0.18 IVpl Figure 3.38: Foam generation by repetitive snap-off at germination sites corresponding to f= 0.3 and a uniform pore-size distribution. The top panel shows the minimum pressure gradient at conditions of incipient mobilization for the specific f value. The bottom panel is a plot of the flow velocity-pressure drop relationship. Because the gas phase is continuous, the fluid behaves initially as Newtonian, with a shear-thickening behavior later, as lamellae are mobilized and additional lamellae are created by division. 153 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 2.0 > 0.8 - 0.4 - 0.0 1.6 - 1.2 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 IVpl Figure 3.39: Foam generation by repetitive snap-off at germination sites corresponding to f=0.1 and a uniform pore-size distribution. Because of the high snap-off probability, the system rapidly becomes fully occupied by lamellae (through snap-off and lamella division), hence it becomes identical to a Bingham plastic, as discussed in chapter 2. The corresponding flow velocity-pressure drop relationship is identical to that of the flow of a fluid with yield stress. In this chapter we used a pore-network model to study various scenaria for foam formation and propagation in porous media. All scenaria included mechanisms of lamella formation through snap-off, lamella division, following lamella mobilization, and lamella leave behind, following the displacement of gas-liquid interfaces. Three different processes were followed: One corresponding to a 3.8 Summary and Conclusions 1 5 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. drainage displacement, in which the liquid flows at steady-state and lamellae are generated by snap-off at advancing interfaces and by leave-behind; another in which the lamellae are generated by capillary pressure fluctuations, during the simultaneous flow of the liquid, with liquid flow in pore comers and in films along pore walls modeled schematically; and a third in which Roof snap-off occurs repeatedly at specific, randomly distributed germination sites. Because foam is always in dynamic equilibrium, a specific foam state depends on the path by which it was reached. Two different cases were examined in each scenario, if appropriate: a case of incipient mobilization, in which we considered the minimum pressure gradient to mobilize at least a single lamella, and a case in which a specific pressure gradient is applied. Clearly, the rheology o f the foam depends greatly on how these states were reached. Because of the intricate dynamic equilibrium involved in the flow of the gas and the liquid, and the attendant mobilization of a multitude of lamellae, foam properties were found to be fluctuating, in some cases with a well-defined period and in others in an apparently chaotic manner. When a steady state is reached, it does not involve mobilization or generation of lamellae (except for an ideal case of unimodal pore-size distribution), in which case the gas flows as either a Newtonian or a Bingham fluid, but with reduced mobility. In this relatively weak-foam case the relative permeability as a function of saturation can be computed. When lamellae form by capillary pressure fluctuations or by repetitive snap-off, the asymptotic state is dynamic. 155 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Our study shows that strong foam generation by lamella mobilization and division and capillary fluctuations is possible without the necessity of repetitive Roof snap-off. Minimum pressure gradient calculations using the latter led to results qualitative similar to the theory of Rossen and Gauglitz [104], without their assumption of a uniform \Vp\ in the pore network. The methodology proposed in this work can be used to explore additional properties of foam generation and propagation in porous media. 156 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 4 FLOW OF A BINGHAM PLASTIC AROUND A CYLINDER IN A POROUS MEDIUM USING A CONTINUUM MODEL 4.1 Introduction As discussed in Chapter 2, the flow of non-Newtonian fluids in porous media has been extensively investigated. Pore-network models have been widely used in modeling Newtonian fluid flow [16] and have also been applied in Chapter 2 and other literatures [92, 115, 121] to effectively capture the effects of micro geometrical structures of porous media on the flow o f non-Newtonian fluids. But macroscopic description of non-Newtonian flows in porous media from micro-scale is still a difficult task. A phenomenological description of non- Newtonian fluid flow in porous media requires some equivalent or apparent viscosities for the fluid to employ Darcy’s law. Experimental and theoretical 157 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. studies are needed to determine the rhelogical model and flow properties for a given non-Newtonian fluid in a specific porous material. Often used are also continuum models borrowed from the descriptions of N ew tonian fluid flo w in porous m edia to avoid introducing apparent viscosity. For example, a modified formulation [129, 130] of Darcy’s law can be used to study the flow o f Bingham plastic in porous media. This is an extension from the single capillary model to the flow in porous media. In this chapter, the continuum model is applied to a simple case, where the Bingham fluid flows around a circular cylinder. The corresponding Newtonian flow case can be found in Ref. [78]. For similar situations like the case studied below, several works in the literature can be found pertaining to flow in the bulk [1, 8, 15]. Adachi and Yoshioka [1] analyzed the creeping flow of a Bingham fluid past a circular cylinder by using variational principles. Beris et al. [8] solved the problem o f creeping motion of a sphere through a Bingham plastic numerically. Blackery and Mitsoulis [15] studied the creeping motion of a sphere in tubes filled with a Bingham plastic material. They computed numerical solutions for various Bingham numbers, and sphere and tube ratios. All of these studies show that there are regions around the cylinder or sphere where flow cannot happen. The emphasis of this chapter is to determine the yield surface which is between the yielded region and the un-yielded region adjacent and the flow field around the cylinder, but for the case of flow in a porous medium. 158 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. To study the creeping flow of Bingham plastic in porous media in this chapter, we employed the simplified extension to Darcy’s law to approximate this special flow. In the remaining, we will develop an asymptotic solution to the problem and determine approximately the yield surface around the circular cylinder and the corresponding velocity field of the fluid. 4.2 Problem Formulation For a Newtonian fluid, the macroscopic description of the flow in porous media is given by Darcy’s law. u = - — VP ■ (4.1) M To describe the Bingham fluid, we can construct an equation similar to (2.25) as a direct upscaling from single-capillary flow. The flow o f such a fluid can be described by the equations below u = 0, when G > |VP| u = - — VP P , when G < |VP| (4.2) where G is the minimum pressure gradient (for example, obtained from the flow of Bingham fluid in a single capillary with radius r [130]). Comparison between the average flow velocity by Buckingham [121] and (4.2) gives G = / d (4.3) 1 5 9 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where d is a length proportional to the radius of the capillary. So G is the pressure gradient corresponding to the yield stress of the fluid. Assuming that the Bingham fluid is flowing around a cylinder, far away from the cylinder, the v elo city w ill be uniform and denoted b y ux . The boundary conditions are then given by where R is the radius of the cylinder. To proceed, we will use the following dimensionless notation The problem expressed in dimensionless variables becomes rd = R , ur d = 0 ’ ur d = cos# and uM = - u x sin# (4.4) _ i 8 dr ^ |V^| / (4.6) r 56 ^ (4.7) V •« = 0 (4.8) with boundary conditions r = 1, ur = 0 r -» oo, ur = cos# r — » oo, ug = - s in # (4.9) 1 6 0 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Gk where e ------ M o o M 4.3 Results If £ is small enough, we can use a perturbation method to get an asymptotic solution to the problem. Let T C — 7C q + £7T[ + • • • Ur = Ur0 + SUr j + ■ • • u e = u eo + £Ugi H (4 .1 0 ) If we keep only the first two terms of the expansions of the above variables, we need to solve the two sets of the equations given below. ^-(™ro) + ^(w *o)=0 or 00 dr 1 dn. ueo ~ o r 50 r = \ , u r0=0 161 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. r — » o o, ur0 = cos 9 , u9 0 = -s in 9 (4.11) and d_ dr d n x v dr j + - d_ d9 r 1 d n x A vr _ a _ dr d n n + - d9 1 1 d n Q Vlv ^o| r d9 J U rl ~ ~ + d n x 1 d7r0 dr V ;rJ dr 1 d n x 1 1 d7V 0 um = ------------- r-------- r d9 |Vtt0 | r d9 r — > oo u r X w., = 0 uex= 0 (4.12) The first set of equations can be solved easily to get f 1N r + — v r j cos 9 (4.13) which is the same solution as for the Newtonian fluid flow case in Ref. [78]. The second set of (4.12) can be solved to get 1 6 2 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. nx = - [ d d J J - l o g i-n 4 4 ^- f ( V \ 2 > 1 + h . - 2 — cos(# - #0) (l + (r0r)'2 - 2(r0r)“ V r J [ deo f 7 - log .> - * 4 - 4^ - r f \ 2 r 1 + V ro j - 2 — cos(# - #0) (l + (r0r y 2 - 2{r0r)~l cos(# - #0 ))fdr0 (4 .1 4 ) where 2 cos#, f 1 2 O 0 /3 ^ ■ j ~ ro + 2 ro cos2#0 / = ' Vro 'o / 7 1 rn + — - - 2cos2#r i Hence, we can find the dimensionless pressure field n , 7 7 — 710 + S 7 U j while the flow velocities are given by U r = U rO + £ U r\ U 8 ~ U 8Q ^ 8 1 (4.15) where d7ro n 1 N , 1 ur0 = - — = (1— r )cos# or r 1 d7Tn U 80 ~ ur 1 = - r dQ d7Tx 1 d7T0 dnx dr (V /r01 dr dr '1 1 + "T V r > / sin# fi 1 2 2 f, 1 1 1 + — — — cos 2# V r 4 r 2 J I r J cos# 163 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. U g i -- — + 1 9 n x 1 1 9n0 _ 1 d7V x r 99 |Vtt0 | r 99 r 96 + 1 1 1 + —;--------- C 0 S 2 # 1 + — v r 4 r 2 r 2 sin 9 (4.16) where nj is given b y (4.14). The force on the cylinder can be obtained by Fx = j m xd9 = Jcos0(;ro + enx)d0 where and n x are given by (4.13) and (4.14), respectively and are functions of just 0 on the cylinder surface. Next, we would like to find the yield surface around the cylinder, which is defined by s = |V ;r|. This is obtained as shown below |V;r| = 9nn 9nx — - + e — 1 9r 9r \ 2 + 1 9 x n e 9n, \ 2 r 99 r 99 where n 0 and n x are substituted by (4.13) and (4.14), respectively. After some simplification, we have ( i+ ,> 2i-(l + e) r \ 1 1 2 1 + —---- — cos 20 r r \ ys ys y + £ ( xaS 1 9(f) \9 r j + 1 9 < j > \Jys~dQj - 2£(l + £)cos 0 V y* 9 < f > 9r + r ys \ r ^ — ' ys r V y*) ■ 2 s i n # - 1 - = £ 96 where < f > = n x- n ^ and r is a function of 9. 1 6 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Different a values give different configurations o f the yield surface. Figure 4.1 shows results with a varying in the range 0.1 to 0.7. In the figure, a increases in the direction of the arrow. The yield surface expands outwards and the non-yielding area becom in g larger as e increases. But one has to note that the perturbation method we used should only apply to the cases where a is small enough and higher order terms of a are not important. The figure here just shows the trend of the change. We can find that the result is consistent with the ones obtained by previous researchers even though different motion equations are used here. S urface s=0.1 e = 0.3 8=0.5 s=0.7 0.8 0.6 0.4 0.2 0.0 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 X Figure 4.1: The change o f the yield surface with different e values. The arrow in the figure shows the direction of increasing e. 165 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 4.4 Conclusions In this brief chapter, a simplified direct extension of Darcy’s law from the flow in a single-capillary to the flow of Bingham plastic in porous media is used to study a special case of flow around a circular cylinder. Assuming the dimensionless pressure gradient corresponding to yield stress in porous media is small, we used a perturbation method to solve the problem. The dimensionless pressure field, the flow velocities, the force on the cylinder and the yield surface were determined by a perturbation approach. The yield surface around the cylinder expands outward with the increase o f e , which is consistent with other researcher’s results. 1 6 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Chapter 5 CONCLUDING REMARKS 5.1 Summary and Conclusions In this dissertation, we studied the flow of fluids with yield stress in porous media, including Bingham plastics and foams. Our work is based on a new algorithm developed in this work and its applications in several problems. A general introduction and a brief literature review were given in Chapter 1. In Chapter 2, a new algorithm based on Invasion Percolation with Memory (IPM) was developed, which incorporates viscous effects caused by flow. The basic idea of this algorithm is to introduce a threshold for fluid in each bond and update it if necessary. This new threshold can be the original value of the fluid yield stress or the one with flow effects. To update the threshold for fluid in a bond belonging to flowing paths, the description of the fluid flow in a single capillary is required. The IPM algorithm is then employed to find the next flowing path. 1 6 7 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. The study in Chapter 2 showed that the macroscopic relation between the applied pressure gradient and flow rate is consistent with the observations made by Roux and Herrmann [108]. A quadratic regime is followed by a linear one, which demonstrates the gradually decreasing effect of yield stress on the flow. Comparisons with the results obtained by Kharabaf et al. [64] using static method reveal that a much higher pressure gradient is necessary to have the same amount of bonds belonging to the flowing paths and the flow patterns are essentially different. Different threshold distributions and conductances lead to similar results. When viscous effects are important, the flow after onset occurs in a layered-like fashion, with the result that fluid in many pore throats transverse to the flow is not mobilized. When the conductances are widely distributed, the viscous effects become weak. The flow rate and pressure gradient relation is effectively a straight line, and the minimum pressure gradient corresponding to the onset of the flow can be approximated by the arithmetic mean of the threshold distribution. A “relative” permeability introduced to express the flow rate as a function of the open sites. And the result can be fitted into a cubic relation. When the network is partially occupied by a fluid with yield stress, a relation between the effective threshold and the fluid saturation was obtained. The decrease of the effective threshold with increasing saturation is a result of the rising accessibility of small thresholds. This information is useful for further upscaling and a simple upscaling method was presented. 168 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. After applying the algorithm to the displacement of a Bingham fluid by a Newtonian fluid, patterns at different parameter values were classified. A phase diagram was constructed involving the capillary number and the yield stress number. If a Bingham plastic is displaced, the yield stress dominated region is of the MTP (Minimum Threshold Path) type, with a 1-D structure that depends on the threshold distribution. And if it is displacing a Newtonian fluid, the patterns are compact. Capillary and viscous controlled patterns have the typical patterns as with Newtonian fluids. In Chapter 3, the algorithm is applied to the study o f the formation and mobilization o f foams in porous media. In the study o f foam generation in drainage displacement of a liquid by a gas, foams are generated by snap off, leave-behind and lamella division. Steady states with no lamella generated or mobilized are identified by advancement of the gas-liquid interface and mobilization of lamellae. The relative permeabilities of the two phases are calculated as a function of the liquid phase saturation, which show that the presence of foams reduces the mobility of the gas phase. Large amount of leave-behind lamellae are created because of the porous medium multiple paths. For the case with relatively weak foam, the gas phase is always continuous and for stronger foam, the gas phase may become disconnected and lamellae have to be mobilized. The capillary number can affect the relative permeabilities of the two phases and the possibility of mobilizing lamellae. 16 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Foam generation by capillary pressure fluctuation is investigated schematically, where the liquid phase is mainly flowing in comers and thin films along pore walls. There are two cases studied, one with minimum incipient pressure drop and another w ith a constant pressure gradient. The m inim um incipient pressure gradient fluctuates between two extreme states, one with continuous gas path and foams generated by capillary pressure fluctuations; and another with blocked gas phase where foams may only be generated by lamella division. The amount of foam generated largely depends on the capillary pressure and the liquid pressure drop, and it is also affected to some extent by the initial state of the foam. The minimum pressure gradient at incipient mobilization corresponding to the maximum value of the fluctuating cycle between the states, changes from 0 where no lamellae can be generated by capillary pressure fluctuation to a maximum value, where lamellae are generated almost everywhere. Application of a pressure gradient larger than the first minimum pressure gradient leads to cycles between two extreme states similar to the first case. The flow rate o f the gas phase increases with time during the process that a continuous gas path is opening up. Then foam is generated due to capillary pressure fluctuations to block the gas path. We found that foams generated at the inlet of the medium driven by capillary pressure fluctuations and the foams generated by lamellae-division can replenish the supply of foams advected downstream or out of the system by mobilization. 170 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Thus, the system can maintain a periodic state. This observation means that repeated Roof snap-off is not necessary for steady foam generation. The last part in Chapter 3 involves repetitive snap-off of lamellae in specific germ ination sites, random ly distributed in the pore network, according to snap-off probability f s. The minimum incipient pressure gradient shows intricate patterns because o f lamella mobilization and generation. Depending on the snap-off probability, the system fluctuates either as periodic states or seemingly random ones. The results show that the minimum pressure gradient required to mobilize lamellae depends on the fraction of throats occupied by lamellae, which decreases when the fraction is roughly smaller than 0.5 and increases afterwards. This relation maps approximately on to the theory of Rossen and Gauglitz [104] without assuming uniform |VP| on the pore network. In the case with unimodal throat-size distribution, the nature of the asymptotic states and also the minimum pressure gradients are much different from the previous case. We also found that with one snap-off at the beginning and no further snap-off during mobilization process, weak foam is generated with continuous gas path in the network. Through the case with repetitive snap-off and increasing pressure gradient, we found that at small snap-off probability (0.3), the flow of the gas phase is Newtonian at the beginning and then demonstrates shear-thickening behavior later, because of lamellae generation by snap-off and lamellae-division. While very large 171 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. snap-off probability (especially close to 1) can result in a flow behavior similar to a Bingham fluid shown in Chapter 2. In Chapter 4, a continuum model was applied to a special case where Bignham fluid flows around a circular cylinder. An asymptotic solution to the problem is found and the yield surface around the cylinder is identified. 5.2 Recommendations for Future Work Various cases related to flow of fluids with yield stress are studied in this dissertation using the new algorithm presented in Chapter 2. Both the topic and the algorithm have potential extensions and applications for further work. • In the dissertation, the new algorithm was applied to the studies of relatively simple fluids. It can be employed to investigate more complex fluids with yield stress, with appropriate modification of the updating scheme for the thresholds. The incorporation of shear-thinning or shear-thickening properties may lead to interesting results. • There are different mechanisms proposed to describe the behavior of foams in porous media. As a mechanistic study approach, the method used in Chapter 3 of this dissertation can be extended to test other popular mechanisms, like break and reform. • For more realistic studies, the compressibility of the gas phase should be taken into account in flow and mobilization processes of foams. In various studies, it 172 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. is shown that the flow of lamellae in porous media is a complex process, which may have a non-linear relation between the pressure gradient and the flow rate of the foams (see, for example, Ref. [9]). Using such a relation after a lamella is mobilized may give a good fit to experimental results. • Foams in porous media are not stable, especially after they are mobilized. Destruction and coalescence may occur so that the amount of foam or lamellae can decrease and foam texture may change. A good description of foam destruction is necessary to relate experimental data to simulation results. • In studies o f multi-phase flow and transport in porous media, there are different descriptions of the pore space geometry proposed to make accurate predictions of macroscopic properties, as discussed by Blunt et al. in a recent review [16]. Employing a network model with varying coordination number, or more complex pore space shapes, or a random network where the connectivity is based on a representation of a real porous medium can be a choice of further study involving foam flow in porous media. 173 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. BIBLIOGRAPHY 1. Adachi, K. and Yoshioka, N.: “On Creeping Flow o f a Visco-Plastic Fluid Past a Circular Cylinder,” Chem. Eng. Sci., 28(1973), 215-226. 2. Apaydin, O.G. and Kovscek, A.R.: “Surfactant Concentration and End Effects on Foam Flow in Porous Media f TiPM, 43(2001), 511-536. 3. Aronson, A.S., Bergeron, V., Fagan, M.E. and Radke, C.J.: “Influence o f Disjoining Pressure on Foam Stability and Flow in Porous Media,” Colloids and Surfaces A: Physicochem. Eng. Aspects, 83(2), 1994, 109-120. 4. Barenblatt, G.I., Entov, V.M. and Ryzhik, V.M.: Theory of Fluid Flows through Natural Rocks, Kluwer Academic Publishers, Netherlands (1990) 5. Barnes, Howard: “The Yield Stress- a Review or ‘ izavm p e i’ - Everything Flows?” J. Non-Newtonian Fluid Mech., 81(1999), 133-178. 6. Bedrikovetsky, P: Principles of Gas and Oil Production, Kluwer Academic Publishers, Netherlands (1995). 7. Bergeron, V. and Radke, C.J.: “Equilibrium Measurements of Oscillatory Disjoining Pressures in Aqueous Foam Films,” Langmuir, 8(1992), 3020- 3026. 8. Beris, A.N., Tsampoulos, J.A., Armstrong, R.C. and Brown, R.A.: “Creeping Motion o f a Sphere through Bingham Plastic,” J. Fluid Mech., 158(1985), 219-244. 9. Bernard, G.G., Holm, L.W., and Jacobs, W.L.: “Effect o f Foam on Trapped Gas Saturation and on Permeability o f Porous Media to Water,” SPE J. (Dec. 1965), 195-300 1 7 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 10. Bertin, H.J., Quintard, M.Y., and Castanier, L.M.: “Development o f a Bubble- Population Correlation fo r Foam-Flow Modeling in Porous Media,” SPE J. (Dec. 1998), 356-362. 11. Bingham, E.C.: Fluid and Plasticity, McGraw-Hill, New York, 1922. 12. Bird, R.B., Armstrong, R.C. and Hassager, O.: Dynamics of Polymeric Liquids, vol. 1, Wiley, New York (1977). 13. Bird, R.B., Dai, G.C. and Yarusso, B.J.: “The Rheology and Flow o f Viscoplastic Materials,” Rev. Chem. Eng., 1(1982), 1-70. 14. Bird, R.B., Stewart, W.E., and Lightfoot, E.N.: Transport phenomena, Wiley, New York (1960). 15. Blackery, J. and Mitsoulis, E.: “ Creeping Motion o f a Sphere in Tubes Filled with a Bingham Plastic Material,” J. Non-Newtonian Fluid Mech., 86(1997), 229-252. 16. Blunt, M.J., Jackson, M.D., Piri, M. and Valvatne, P.H.: “Detailed Physics, Predictive Capabilities and Macroscopic Consequences fo r Pore-Network Models o f Multiphase Flow,” Adv. Water Resources, 25(2002), 1069-1089. 17. Brenblatt, G.I., Lions, J.L. and Papanicolaou, G.: Asymptotic Analysis for Periodic Structures, North-Holland, Amsterdam, 1978. 18. Bryan, S. and Blunt, M.J.: “Prediction o f Relative Permeability in Simple Porous Media,” Phys. Rev. A, 46 (1992), 2004-2011. 19. Carcoana, A., Applied Enhanced Oil Recovery, Prentice Hall, New Jersey, 1992. 20. Chambers, D. “Foams for Well Stimulation” in Foams: Fundamentals and Applications in the Petroleum Industry edited by Schramm, L.L., ACS Advances in Chemistry Series No. 242, Am. Chem. Soc., Washington, DC (1994). 21. Chambers, K.T., and Radke, C.J.: “Capillary Phenomena in Foam Flow Through Porous Media,” in Interfacial Phenomena in Oil Recovery, N. R. Morrow, (ed.), Marcel Dekker, New York, 1990. 175 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 22. Chauveteau, G.: “ Rodlike Polymer-Solution Flow through Fine Pores - Influence o f Pore-Size on Rheological Behavior,” J. Rheol., 26(2), 1982, 111- 142. 23. Chen, M., Rossen, W.R., and Yortsos, Y.C.: “Flow and Displacement o f Fluids with Yield Stress in Porous Media,” accepted by Chem. Eng. Sci.. 24. Chen, M., Yortsos, Y.C. and Rossen, W.R.:‘V 1 Pore-Network Study o f the Mechanisms o f Foam Generation,” paper SPE 90939, presented at the 2004 SPE ATCE, Houston, TX, 26-29 Sept. 2004. 25. Cheng, L. et al.: “Simulating Foam Processes at High- and Low- Foam qualities,” paper SPE 59287 presented at the 2000 SPE/DOE Improved Oil Recovery Symposium, Tulsa (3-5 April 2000). 26. Cheng, L., Kam, S.I., Delshad, M., and Rossen, W.R.: “Simulation o f Dynamic Foam-Acid Diversion Processes,” paper SPE 68916 presented at the 2001 SPE European Formation Damage Symposium, The Hague, 21-22 May. 27. Chou, S.I.: “Percolation Theory o f Foam in Porous Media,” paper SPE 20239 presented at the 1990 SPE/DOE Symposium on EOR, Tulsa, OK, April 22-25. 28. Churchill, S.W.: Viscous Flow, Butterworths, USA (1988). 29. Cohen, D., Patzek, T.W., and Radke, C.J.: “Two-Dimensional Network Simulation o f Diffusion-Driven Coarsening o f Foam Inside a Porous Medium,” J. Colloid and Interface Sci., 179(1996), 357-373. 30. Cohen, D., Patzek, T.W., and Radke, C.J.: “Onset o f Mobilization and the Fraction o f Trapped Foam in Porous Media,” TiPM, 28(1997), 253-284. 31. Cox, S.J., Neethling, S., Rossen, W.R., Schleifenbaum, W., Schmidt- Wellenburg, P., and Cilliers, J.J., “ A Theory o f the Effective Yield Stress o f Foam in Porous Media: The Motion o f a Soap Film Traversing a Three- Dimensional Pore,” Colloids and Surfaces A: Physicochem. Eng. Aspects (2003), submitted. 32. de Vries, A.S. and Wit, K.: “ Rheology o f Gas/Water Foam in the Quality Range Relevant to Steam Foam,” SPERE, May 1990, 185-192. 33. Du, C., and Yortsos, Y.C.: “ A Numerical Study o f the Critical Gas Saturation in a Porous Medium,” TiPM (1999), 35, 205-225. 1 7 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 34. Dullien, F.A.L.: Porous Media. Fluid Transport and Pore Structure, Academic Press, New York (1979). 35. Ettinger, R.A. and Radke, C.J.: “Influence o f Texture on Steady Foam Flow in Berea Sandstone,” SPERE, February 1992, 83-90. 36. Falls, A.H., Hirasaki, G.J., Patzek, T.W., Gauglitz, D.A., Miller, D.D., and Ratulowski, T.: “Development o f a Mechnistic Foam Simulator: The Population Balance and Generation by Snap-Off ” SPERE, August 1988, 884- 892. 37. Falls, A.H., Musters, J.J., and Ratulowski, J.: “The Apparent Viscosity o f Foams in Homogeneous Bead Packs,” SPERE, May 1989, 155-164. 38. Feder, J.: Fractals, Plenum, New York (1988). 39. Fergui, O. and Quintard, M.: “Transient Aqueous Foam Flow in Porous Media: Experiments and Modeling,” J. Petr. Sci. Eng., 20(1998), 9-29. 40. Flumerfelt, R.W. and Prieditis, H.: “Mobility o f Foam in Porous Media,” in Surfactant-Based Mobility Control: Progress in Miscible-Flood Enhanced Oil Recovery Smith, D. H. Ed., ACS Symposium Ser. 373. AM. Chem. Soc., Washington, D. C. 1988. 41. Flumerfelt, R.W., Chen, H-L., Ruengphrathuengsuka, W., Hsu, W-F., and Prieditis, J.: “Network Analysis o f Capillary and Trapping Behavior o f Foams in Porous Media,” SPEFE, March 1992, 25-33. 42. Friedmann, F., Chen, W.H., and Gauglitz, P.A.: “Experimental and Simulation Study o f High-Temperature Foam Displacement in Porous Media,” SPERE, Feb. 1991, 37-45. 43. Gauglitz, P.A., Friedmann, F., Kam, S.I., and Rossen, W.R.: “Foam Generation in Homogeneous Porous Media,” Chem. Eng. Sci., 57(2002), 4037-4052. 44. Gauglitz, P.A., Friedmann, F., Kam, S.I., and Rossen, W.R.: “Foam Generation in Porous Media,” paper SPE 75177 presented at the 2002 SPE/DOE Symposium on Improved Oil Recovery, Tulsa, OK, 13-17 April. 45. Gauglitz, P. A., Laurent, C.M.St, and Radke, C.J.: “Experimental Determination o f Gas-Bubble Breakup in a Constricted Cylindrical Capillary,” Ind. Eng. Chem. Res. 27(1988), 1282-1291. 177 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 46. Gdanski, R.D.: “Experience and Research Show Best Designs fo r Foam- Diverted Acidizing," Oil and Gas Journal (Sept. 6,1993), 85. 47. Gillis, J.V. PhD Dissertation 1990 University of California at Berkeley. 48. Gillis, J.V., Radke, C.J.: “ A Dual Gas Tracer Technique fo r Determining Trapped Gas Saturation During Steady Foam Flow in Porous Media," paper SPE 20519 presented at the 1990 SPE ATCE, New Orleans, Sept. 23-26. 49. Hanssen, J.E.: “Foam as a Gas-Blocking Agent in Petroleum Reserviors II: Mechanisms o f Gas Blockage by Foam," J. Petroleum Science and Technology, 10(1993), 135-156. 50. Hazlett, R.D. and Furr, M.J.: “ Percolation Model fo r Permeability Reduction in Porous Media by Continuous-Gas Foams," Ind. Eng. Chem. Res., 39 (2000), 2709-2716. 51. Herschel, W.H. and Bulkley, R., Proc. Am. Assoc. Test Mater., 26(1926), 621. 52. Hirasaki, G.J. and Clifford, P.J.: “ Analysis o f Factors Influencing the Mobility and Adsorption in the Flow o f Polymer Solution Through Porous Media," SPE J., 1974,337-346. 53. Hirasaki, G.J., Jackson, R.E., Jin, M., Lawson, J.B., Londergan, J., Meinardus, H., Miller, C.A., Pope, G.A., Szaffanski, R., and Tanzil, D.: “Field Demonstration o f the Surfactant/Foam Process fo r Remediation o f a Heterogeneous Aquifer Contaminated with DNAPL," in NAPL Removal: Surfactants, Foams, and Microemulsions, S. Fiorenza, C. A. Miller, C. L. Oubre, and C. H. Ward, (eds.), Lewis Publishers, Boca Raton (2000). 54. Hirasaki, G.J. and Lawson, J.B.: ‘ '“ 'Mechanisms o f Foam Flow in Porous Media: Apparent Viscosity in Smooth Capillaries," SPE J., April 1985, 176- 190. 55. Holm, L.W.: “The Mechanisms of Gas and Liquid Flow Through Porous Media in the Presence of Foam,” SPE J., Dec. 1968, 359-369. 56. Huh, D.G., and Handy, L.L.: “Comparison o f Steady- and Unsteady-state Flow o f Gas and Foaming Solution in Porous Media," SPERE, 4(1989), 77- 84. 178 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 57. Jimenez, A.I. and Radke, C.J.: “Dynamics Stability o f Foam Lamellae Flowing Through a Periodically Constricted Pore,” in Oil-Field Chemistry: Enhanced Recovery and Production Stimulation, Borchardt, J. K. and Yen, T. F. Eds. ACS Symposium Series 396; American Chemistry Society: Washington DC 1989. 58. Kam, S.I., Li, Q., Nguyen, Q.P. and Rossen, W.R.: “Dynamic Simulations With an Improved Model fo r Foam Generation,” paper SPE 90938, presented at the 2004 SPE ATCE, Houston, TX, 26-29 Sept. 2004. 59. Kam, S.I. and Rossen, W.R.: “ A Model fo r Foam Generation in Homogeneous Media,” paper SPE 77698, ATCE 2002. 60. Kam, S.I., and Rossen, W.R.: “ A Model fo r Foam Generation in Homogeneous Porous Media,” SPE J., 8 (Dec. 2003), 417-425. 61. Katz, A.J., and Thompson, A.H.: “Quantitative Prediction o f Permeability in Porous Rock,” Phys. Rev. B, 34(1986), 8179-8181. 62. Kharabaf, H. Ph. D. dissertation, 1996, University of Southern California. 63. Kharabaf, H., Shah, C., and Yortsos, Y.C.: “Flow and Displacement o f Bingham Plastics in Porous Media,” 7th UNITAR International Conference on Heavy Crude and Tar Sands, Beijing, China (Oct. 27-30, 1998). 64. Kharabaf, H. and Yortsos, Y.C.: ''‘ 'Invasion Percolation with Memory,” Phys. Rev. E, 55(1997), 7177-7191. 65. Kharabaf, H. and Yortsos, Y.C.: “ Pore-Network Model fo r Foam Formation and Propagation in Porous Media,” SPE J., 1998, 42-53. 66. Khatib, Z.I., Hirasaki, G.J., and Falls, A.H.: “Effects o f Capillary Pressure on Coalescence and Phase Mobiliies in Foams Flowing through Porous Media,” SPERE, August 1988, 919-926. 67. Kornev, K.G., Neimark, A.V., and Rozhkov, A.N.: “Foam in Porous Media: Thermodynamic and Hydrodynamic Peculiarities,” Advances in Colloid and Interface Science, 82(1999), 127-187. 68. Kovscek, A.R. and Bertin, H.J.: “Foam Mobility in Heterogeneous Porous Media (I: Scaling Concepts, II: Experimental Observations),” TiPM, 52(2003), 17-49. 1 7 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 69. Kovscek, A.R., Patzek, T.W., and Radke, C.J.: “ A Mechanistic Population Balance Model fo r Transient and Steady-State Foam Flow in Boise Sandstone ." Chem. Eng. Sci., 50(23), 1995, 3783-3799. 70. Kovscek, A.R., Patzek, T.W., and Radke, C.J.: “Mechanistic Foam Flow Simulation in Heterogeneous and Multidimensional Porous Media," SPE J., Dec. 1997, 511-526. 71. Kovscek, A.R. and Radke, C.J.: “Fundamentals o f Foam Transport in Porous Media," in Foams: Fundamentals and Applications in the Petroleum Industry edited by Schramm, L. L., ACS Advances in Chemistry Series No. 242, Am. Chem. Soc., Washington, DC (1994). 72. Kovscek, A.R. and Radke, C.J.: “Gas Bubble Snap-Off under Pressure-Driven Flow in Constricted Noncircular Capillaries,” Colloids and Surfaces A: Physicochem. Eng. Aspects, 117(1996), 55-76. 73. Kovscek, A.R. and Radke, C .J.: “Pressure-Driven Capillary Snap-Off o f Gas Bubbles at Low Wetting-Liquid Content," Colloids and Surfaces A: Physicochem. Eng. Aspects, 212 (2003), 99-108. 74. Kovscek, A.R., Tretheway, D.C., Persoff, P. and Radke, C.J.: “Foam flow through a transparent rough-walled rock fracture," J. Petr. Sci. Eng., 13 (1995), 75-86. 75. Kraynik, A.M and Hansen, M.G.: “Foam and Emulsion Rheology: A Quasistatic Model fo r Large Deformations o f Spatially- Periodic Cells," J. Rheol., 30(1986), 409-439. 76. Laidlaw, W.G. and Wilson, W.G.: “ A lattice Model o f Foam Flow in Porous Media: A Percolation Approach." TiPM, 11(1993), 139-159. 77. Lake, L.W., Enhanced Oil Recovery, Prentice Hall, New Jersey, 1989. 78. Lamb, H., Sir Hydrodynamics. Dover, New York, 1945. 79. Lenormand, R.: “Liquids in Porous Media," J. Phys.: Condens. Matter, 2 (1990), 79-88. 80. Lopez, X., Valvatne, P.H. and Blunt, M.J.: “Predictive Network Modeling o f Single-Phase Non-Newtonian Flow in Porous Media," J. Colloid and Interface Science, 264(2003), 256-265. 81. Mills, C.C., Rheology of Disperse Systems, Pergamon Press, Oxford, 1959. 180 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 82. Myers, T.J. and Radke, C.J.: “Transient Foam Displacement in the Presence o f Residual Oil: Experiment and Simulation Using a Population-Balance Model," Ind. Eng. Chem. Res., 39(2000), 2725-2741. 83. Nguyen, Q.N.: PhD dissertation, Technical University of Delft (2004). 84. Nguyen, Q.P., Currie, P.K., and Zitha, P.L.J., “Determination o f Foam Induced Fluid Partitioning in Porous Media Using X-ray Computed Tomography,” SPE 80245 presented at the 2003 SPE International Symposium on Oilfield Chemistry held in Houston, Texas, U.S.A. (5-8 February 2003). 85. Nguyen, Q.P., Currie, P.K., and Zitha, P.L.J.: “Motion o f Foam Films in Diverging-Converging Channels,” Journal of Colloid and Interface Science, 271(2004), 473— 484. 86. Oldroyd, J.D ,:“ A Rational Formulation o f the Equations o f Plastic Flow fo r a Bingham Solid”, Proc. Camb. Phil. Soc., 43(1947), 100-105. 87. 0ren, P.-E. and Bakke, S.: “ Process Based Reconstruction o f Sandstones and Prediction o f Transport Properties,” TiPM, 46 (2002), 311-343. 88. 0ren, P.-E., Bakke, S. and Amtzen, O.J.: ‘ ‘ ‘ 'Extending Predictive Capabilities to Network Models,” SPE J., 3 (1998), 324-336. 89. Owete, O.S., and Brigham, W.E.: “Flow Behavior o f Foam: A porous Micromodel Study,” SPERE, 1987, 315-323. 90. Papanastasiou, T.C.: “Flows o f materials with Yield,” Journal of Rheology, 31(5), 1987, 385-404. 91. Patzek, T.W.: “ Verification o f a Complete Pore Network Simulator o f Drainage and Imbibition,” SPE J., 6 (2), 2001, 144-156. 92. Pearson, J.R.A and Tardy, P.M.J.: “Models fo r Flow o f Non-Newtonian and Complex Fluids through Porous Media,” J. Non-Newtonian Fluid Mech., 2087(2002), 1-27. 93. Princen, H.M.: “Rheology o f Foams and Highly Concentrated Emulsions I: Elastic Properties and Yields Stress o f a Cylindrical Model System,” J. Colloid Interface Sci., 91(1983), 160-175. 181 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 94. Princen, H.M., Aroson, M.P. and Moser, J.C.: “Highly Concentrated Emulsions II Real Systems: The effect o f Film Thickness and Contact Angle on the Volume Fraction in Creamed Emulsions," J. Colloid Interface Sci., 75(1980), 247. 95. Radke, C.J., and Gillis, J.V.: “A Dual Gas Tracer Technique fo r Determining Trapped Gas Saturation during Steady Foam Flow in Porous Media f paper SPE 20519 presented at the 1990 SPE Annual Technical Conference and Exhibition, New Orleans (23-26 September 1990). 96. Ransohoff, T.C., and Radke, C.J.: “Mechanisms o f Foam Generation in Glass-BeadPacks," SPERE, May 1988, 573-585. 97. Roof, J.G.: “Snap-Off o f Oil Droplets in Water-Wet Pores," SPE J., 10(1), 1970,85-90. 98. Rossen, W.R.: “Minimum Pressure Gradient fo r Foam Flow in Porous Media: Effect o f Interactions with Stationary Lamellae," J. Colloid Interface Sci., 139(1990), 457-468. 99. Rossen, W.R.: “Theory o f Mobilization Pressure Gradient o f Flowing o f Foams in Porous Media. I. Incompressible Foam. II. Effect o f Compressibility. III. Asymmetric Lamella Shapes," J. Colloid and Interface Science, 136(1990), 1-53. 100. Rossen, W.R.: “Size Distribution o f Blocked Clusters on a Bethe Tree," J. Phys. A: Math Gen., 24(1991), 5155-5161. 101. Rossen, W.R.: “Enhanced Oil Recovery” in Foams: theory, measurements, and applications edited by Prud'homme, R. K. and Khan, S. A. New York, Marcel Dekker, Inc., 1996 102. Rossen, W.R.: “Snap-off in Constricted Tubes and Porous Media," Colloids and Surfaces A: Physicochemical and Engineering Aspects, 166(2000), 101— 107. 103. Rossen, W.R.: “ A Critical Review o f R oof Snap-Off as a Mechanism o f Steady-State Foam Generation in Homogeneous Porous Media,” Colloids and Surfaces A: Physicochem. Eng. Aspects, 225(2003), 1-24. 104. Rossen, W.R. and Gauglitz, P.A.: “Percolation Theory o f Creation and Mobilization o f Foams in Porous Media," AIChE J., 36(8), 1990, 1176-1188. 1 8 2 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 105. Rossen, W.R. and Mamun, C.K.: “Minimal Path fo r Transport in Networks,” Phys. Rev. B, 47(18), 1993, 11815-11825. 106. Rossen, W.R., Shi, J. and Zellinger, S.C.: “Percolation Modeling o f Foam Generation in Porous Media,” AIChE J., 40(6), 1994, 1082-1084. 107. Rossen, W.R. and Wang, M.-W., “Modeling Foams fo r Acid Diversion,” SPE J (1999), 92-100. 108. Roux, S., and Herrmann, H.J.: “ Disorder-Induced Nonlinear Conductivity,” Europhys. Lett., 4(11), 1987, 1227-1231. 109. Sahimi, M.: “Nonlinear Trnsport Processes in Disordered Media,” AIChE J., 39(3), 1993, 369-385. 110. Sanchez, J.M. and Hazlett, R.D. “ Foam Flow through An Oil-Wet Porous Medium,” SPERE, 7(1992), 91-97. 111. Sanchez, J.M., and Schechter, R.S.: “Surfactant Effects on the Two-Phase Flow o f Steam-Water and Nitrogen-Water through Permeable Media,” J. Petr. Sci. Eng. 3, 185-199, (1989). 112. Sanchez, J.M., Shcechter, R.S. and Monsalye, A.: “The Effect o f Trace Quantities o f Surfactant on Nitrogen/Water Relative Permeabilities,” paper SPE 15446 Proc, Ann. Mtg of the Soc. Pet. Eng., New Orleans, LA, 1986 (October). 113. Schramm, L.L. and Wassmuth, F.: “Foams: Basic Principles,” in Foams: Fundamentals and Applications in the Petroleum Industry edited by Schramm, L.L., ACS Advances in Chemistry Series No. 242, Am. Chem. Soc., Washington, DC (1994). 114. Schramm, L.L. (ed.): Foams: Fundamentals and Applications in the Petroleum Industry, ACS Advances in Chemistry Series No. 242, Am. Chem. Soc., Washington, DC (1994). 115. Shah, C., PhD Dissertation, University of Southern California (1994). 116. Shah, C., Kharabaf, H., and Yortsos, Y.C.: “ Immiscible Displacements Involving Power-Law Fluids in Porous Media,” Proc. Seventh UNITAR International Conference on Heavy Crude and Tar Sands, Beijing, China (Oct. 27-30, 1998). 183 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 117. Shah, C. and Yortsos, Y.C.: “ Aspects o f Flow o f Power-Law Fluids in Porous Media,” AIChE J., 41(5), 1995, 1099-1122. 118. Sheffield, R.E. and Metzner, A.B.: “Flows o f Nonlinear Fluids Through Porous Media,” AIChE J., 22, 736(1976). 119. Siddiqui, S., Talabani, S., Saleh, S.T., and Islam, M.R.: “Foam Flow in Low- Permeability Berea Sandstone Cores: a Laboratory Investigation,” J. Petr. Sci. Eng., 36(2002), 133- 148. 120. Skelland, A.H.P. Non-Newtonian Flow and Heat Transfer, John Wiley & Sons Inc., New York City (1967). 121. Sorbie, K.S., Clifford, P.J. and Jones, E.R.W.: “The Rheology o f Pseudoplastic Fluids in Porous-Media Using Network Modeling,” J. Colloid Interface Sci., 130(2), 1989, 508-534. 122. Tanner, R.I. Engineering Rheology Clarendon press Oxford, 1988. 123. Tanzil, D., Hirasaki, G.J. and Miller, C.A: “Conditions fo r Foam Generation in Homogeneous Porous Media,” paper SPE 75176, SPE/DOE Improved Oil Recovery Symposium, Tulsa, Oklahoma, 13-17 April 2002. 124. Taylor, K.C. and Nasr-El-Din, H.A.: “ Water-Soluble Hydrophobically Associating Polymers fo r Improved oil Recovery: A Literature Review, ” J. Pet. Sci. Eng., 19(1998), 265-280. 125. Terdre, N.: “ Foam-Assisted Injection Trials Could Spread to Other North Sea Fields,” Offshore 63 (Aug. 2003), 70. 126. Tian, J.-P., and Yao, K. -L.: “ Immiscible Displacements o f Two-Phase Non- Newtonian Fluids in Porous Media,” Physics Letters A, 261(1999), 174-178. 127. Veeningen, D., Zitha, P. L. J., and van Kruijsdijk, C.P.J.W.: “Understanding Foam Flow Physics: The Role o f Permeability,” paper SPE 38197 presented at the 1997 European Formation Damage Symposium, The Hague, June 2-3. 128. Vorwerk, J. and Brunn, P.O.: “Shearing Effects fo r the Flow o f Surfactant and Polymer Solutions through a Packed Bed o f Spheres,” J. Non-Newtonian Fluid Mech., 51(1994), 79-95. 129. Wu, Y.S., Pruess, K. and Witherspoon, P.A.: “ Displacement o f a Newtonian Fluid by a Non-Newtonian Fluid in a Porous Medium,” TiPM, 6(1991), 115- 142. 184 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 130. Wu, Y.S., Pruess, K. and Witherspoon, P.A.: “Flow and Displacement o f Bingham Non-Newtonian Fluids in Porous Media,” SPERE, 1992, 369-376. 131. Wu, Y.-S. and Pruess, K.: “ A Numerical Method fo r Simulating Non- Newtonian Fluid Flow and Displacement in Porous Media,” Adv. Water Resources, 21(1998), 351-362. 132. Xu, Q. and Rossen, W.R.: ‘ ‘ ‘ ‘ Effective Viscosity o f Foam in Periodically Constricted Tubes,” Colloids and Surfaces A: Physicochem. Eng. Aspects, 216(2003), 175-194. 133. Yang, S.H. and Reed. R.L.: ‘ ‘ M obility Control Using CO2 Foams,” paper SPE 19689, Proc. Ann. Mtg. O f the Soc. Pet. Eng., San Antonio, 1989. 134. Yortsos, Y.C., Xu, B., and Salin, D.: “Phase Diagram o f Fully-Developed Drainage in Porous Media,” Phys. Rev. Lett., 79(23), 1997, 4581-4584. 185 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Appendix A FLOW OF FLUID WITH YIELD STRESS IN PARALLEL CAPILLARIES Consider capillaries of dimensionless radius r distributed according to a pdf a ( r ) , containing a yield-stress fluid of the rheological behavior of equation (2.6). If the capillaries are in parallel, it is not difficult to show following overall flow behavior i i < q >=< Ap > r 3a(r)dr - r 3r(r)a(r)dr where Ap > A p ^ (A.l) A p Ap where Apm in = — and brackets denote average. Near the onset of flow, we can rm max take Ap - Apm in (1 + x) where x is small, and expand (A l) to further obtain a (rt=r '' 'm a il <q >» x 2 (A.2) Equation (A2) indicates that a quadratic dependence between flow rate and departure of the pressure drop from the critical is consistent with flow in a bundle 1 8 6 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. of parallel capillaries. In the same configuration, the fraction of open capillaries is given by the equation Xpc = Jp * a (r )d r = £ f( r ) d r (A.3) A p which near the onset of flow becomes (A.4) The linear relation near the onset of flow is consistent with the argument of Roux and Herrmann. In the opposite end, if the capillaries are in series, then it is not difficult to derive the linear relationship, < Ap >-< q > ^j+ < t > (A.5) which can model the linear part of the flow rate-pressure drop relationship. 1 8 7 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Appendix B SIMULTANEOUS TWO-PHASE FLOW IN A PORE NETWORK Two-phase flow in porous media is an important phenomenon, which involves a large number of factors, such as capillary, viscous, and gravity forces, the viscosities of the two fluids and the interfacial tension separating them. Pore-scale network models, based on capillary-dominated displacement mechanisms, provide explicit calculations of relative permeabilities, as well as capillary pressures and saturations. Numerical network models with various pore-space representations have been employed to investigate many phenomena related to multi-phase flow in porous media. In chapter 3, we discussed the generation and flow of foams in displacement process, which involves two-phase flow of a wetting phase, liquid and a non-wetting phase, gas. A simple network model was employed to represent the interactions between the two phases. Here we use a similar approach to 188 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. investigate co-current and countercurrent two-phase flow of Newtonian fluids in porous media. A 3-D pore and bond network was used. Periodic boundary conditions were applied except the inlet and the outlet. The follow in g assum ptions are m ade in the algorithm: 1. At the outlet, the pressure of the wetting phase is zero. (Pn w P'S)inlet outlet ' 3. There is no mass interchange at the interface o f the two fluids. 4. All of bonds have the same length. 5. The two fluids are all incompressible. 1. Co-Current Flow For a fixed total flow rate of the two phases, i.e. constant capillary number, the flow rate fraction of the wetting phase is f. The flow in a single capillary follows Poiseuille’s Law 4 T T X ^ ^ 0 = 7—7 A P (B -!) 8 /j I Viscosity ratio of the two phases is given by M, which is M - . At any Aw time, only one o f the two phases can be in a bond and it occupies the whole bond. For each of the two phases, at each pore i , the equation 1 8 9 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 0 < B-2) j is obeyed, where qtj is the volumetric flow rate from pore i to its neighbor pore j , and Z is the coordination number, here Z = 6 . According to Darcy’s Law, for a porous media full of one single phase, Ak Qf = VPF (B.3) A And for two-phase flow, Q = - ^ ^ V P . (B.4) Hence, the relative permeability k r is given by Q A P f kr = - ----- (B.5) Qf a p In a drainage process, the interface between non-wetting phase and wetting phase advances toward the wetting phase if the following is satisfied. 2 Y P n w ~ P w > — (B.6) r The network is initially fully filled with wetting phase and Invasion Percolation is employed until breakthrough to find a continuous path for the non-wetting phase. Interfaces are moved so that the flow rate of the non-wetting phase increases until it reaches the value determined by f The interface in the bond with max f \ _ 2 y P inw P iw V r9 j advances if necessary, where i and j are the two adjacent 1 9 0 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. pores. During the process, the wetting phase may be trapped which will affect the connectivity, thus permeability of the non-wetting phase. So in the invasion of the non-wetting phase, we assume that the wetting phase left alone in a bond without direct connection to bulk wetting-phase will finally escape through microscopic wetting film and the bond will be occupied by non-wetting phase as shown in Figure B .l. Similar algorithm can be used for imbibition process subject to the following 2 V changes, a) The criteria for moving the interface is that p n w - p w < — ; b) Choose r the interface with min f \ P inw P jw r~ V V J The effect of capillary number and viscosity ratio on the flow is studied. For the case Cci — > 0 where capillary forces dominates, the advancement of the interface is totally determined by bond radius and interfacial tension, and the process becomes equivalent to Invasion Percolation that large bonds are invaded first. The relative permeability of the two phases of an IP process is shown in Figure B.2 (solid line), which will be compared to the cases with larger Ca later. Because the viscous forces may be neglected, the viscous ratio doesn’t affect the relative permeabilities. From the figure, it is obvious that wetting phase has a residual saturation about 0.3. The invasion process will generate some blobs with small permeability surrounded by invading phase. We assume these trapped blobs will not be displaced and thus result in the high residual saturation. 191 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure B.2 also reveals the variation of relative permeability o f the two phases with increasing capillary number. For relatively small capillary numbers (Ca=10'4, dotted line in Figure B.2), the permeability curve is very close to the IP case, which indicates that the capillary forces are dominant. When capillary number increases further, the curves move, especially after capillary number is larger than 10'2 as shown in the figure. The residual saturation of the wetting phase decreases from about 0.3 for the IP case, to about 0.2 for capillary number 0.5(dash-dotted line in Figure B.2). This is the result of the increasing viscous effects, which make it possible for small bonds to be invaded and thus less wetting phase is trapped. Another observation is that with increasing capillary number, the permeability curves tend to be linear, especially the permeability of the non-wetting phase (see, Ca=0.5). Figure B.3 demonstrates that the pressure gradient of the system changes with wetting phase saturation. Because of the competition between saturation, thus relative permeability and fraction of flow rate, the necessary pressure gradient increases at the beginning (region A), reaches some largest values and then decreases (region B). The pressure gradient is obtained from the wetting phase flow rate. For wetting phase, the relative permeability and flow rate decrease at the same time, but from the pressure gradient curve, it is obvious that in region A, permeability decrease more rapidly than flow rate fraction of wetting phase, while in region B, the case is converse. And the position of the peak of the curve changes with capillary number. 192 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. At small capillary numbers, the wetting phase fraction appear to be s-shaped, while for large capillary numbers, the curves are approaching straight lines in Figure B.4. The extrapolation is that with infinite capillary number, negligible capillary forces and large dominant viscous forces, the curve will approach the diagonal with the residual saturation of wetting phase vanishing. Figures B.5-B.6 show the result of different viscosity ratios at a constant capillary number, Ca=10'4. The fractions of wetting phase flow rate as functions of the wetting phase saturation for different viscosity ratios are plot in Figure B. 5. With increasing M values, the curves change from an almost concave line to strong convex lines. At same wetting phase saturation, / increases with M, which is consistent with fractional flow result found in literatures. The viscosity ratio doesn’t have strong effect on permeability curves at this specific capillary value. For different M values studied, the preference of large bonds in invasion doesn’t change. This also explains the fact that the residual saturations of the wetting phase at different M values are almost the same. But M value affects the pressure gradient trends as shown in Figure B.5. The peak value of the pressure gradient decreases with increasing M and the position of the peak value moves towards small wetting phase saturation. For same wetting phase saturation, the permeability values keeps constant regardless the different viscosities, so the pressure gradient is proportional to f/M and f/M decreases with M (Figure B.5). Another trend shown in these figures is that the monotonic part in 193 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. small saturation region becomes smaller with rising M and it finally vanishes for M -0.5. For imbibition processes, the corresponding results are shown in Figures B.7- B .ll. The effect of capillary number on relative permeabilities is similar to the drainage case and large capillary number leads to small non-wetting phase residual saturation because of the increasing viscous forces (Figure B.7). The pressure gradient values increase with capillary number, while the positions of the peak don’t change obviously. Comparison of Figure B .8 and Figure B. 3 shows that the pressure gradients are higher in the imbibition case than in drainage under same conditions. This is because in imbibition, small bonds are chosen in initial IP process and also in further invasion when capillary number is relatively small. The decrease of residual saturation of the non-wetting phase in Figures B.7 and B.9 is a direct result of the fact that at larger capillary numbers, more large bonds are invaded and less are trapped. Figure B. 10 shows very similar trends of the fractions o f the wetting phase as a function of wetting saturation at different M values. Because the weak effect of the viscosity ratio on relative permeability, the result is not shown. We can find in Figure B .l 1 that the peaks of the pressure gradient moves towards small saturation as in drainage case, but here large M values lead to high pressure gradients, which is contrary to the drainage case. Different from the drainage case, in imbibition, the pressure gradient is determined by non-wetting phase flow. With small changes in relative permeabilities at different M values, the pressure gradient is nearly 194 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. proportional to fM which increases with M. The weak viscous effects also explain the small difference between the residual saturations of the non-wetting phase at different Ms shown in Figure 10, which is similar to the drainage case, where the wetting phase residual saturation keeps almost constant. 2. Countercurrent Flow Heat pipe technology has been developing since the mid-1960s. They are widely utilized in spacecraft application, in large industrial application, or as heat sinks to cool electronic components and packages. Heat pipe operates on a closed two- phase cycle and only pure liquid and vapor are present within the container. A heat pipe has three distinct regions: an evaporator, a condenser and an adiabatic region. The wick in a heat pipe can be viewed as a typical porous medium. At the bottom, the liquid is vaporized below by a constant heat flux source, and at the upper end, the liquid flows down. What takes place in heat pipe is a typical countercurrent two-phase flow. At steady state, we can assume that there is a mass conservation in heat pipe, i.e., the mass of the vapor flows up is equal to the mass o f the liquid flows down. We can use similar algorithm in co-current flow to simulate this process, under the same assumptions. In this system, the flow rates of vapor and liquid are all constant. If the total flow rate o f the two phases is denoted by Qt; the flow rates of the two phases are given by 195 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Gv Q t \ + Pl/ and Q t \ + Pv, Pi With Darcy’s law for the two phases &=■ Akk„. Pi dP; dz L + P tg where i can be / or v and we can get where, /3 PlPv PlPv dPv = dPl dz dz Qh =Qvq{r*)2 p vLv, CO - P v Q h ^rv^rl Lvp vkApgA krl + p kr Flows of vapor and liquid in a single capillary d are transformed to / A 7irA a dP L + { r ^ P v _ L dz kco Ap j + p v P i and < h i = - 4 nr-- V Pv 00 Pi dPij + (r *)2 p , ___ 1 dz kco Ap | + p v Pi (B.7) (B.8) (B.9) (B.10) (B .ll) (B-12) (B.13) (B-14) 1 9 6 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Because of the flow direction, the gravity effect on the move of interface should be considered for flows in vertical throats. Two different cases are demonstrated in Figure B. 12. For case a, P, - P 4 + (p v + p l) g h > ^ - r and 2 y p4 - pl -(P v + P i)sh > - r for case b, where we assume the interface is at the mid-point o f the throat and h is kco F f changing variable, like the flow rate fraction / was in co-current flow. With the half of the throat length. To simulate this process, - 7 —— is the independent kco increase of - — — , the capillary number should vary correspondingly, too. (r*)2 The effect of density ratio and viscosity ratio on the drainage and imbibition processes is investigated, see Figure B.13. For different viscosity ratios and density ratios, the highest value of ^ varies. The shapes of ^ as a function of (r *) (r *) vapor saturation in drainage processes are quite similar for different parameters used. This shows that for a value of co, depending on the path to approach the state, there may be two different steady states which corresponding to different vapor saturations. This is as predicted theoretically. 1 9 7 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure B.l: Snapshots of invasion of the non-wetting phase a) Before, and b) after a step. Black part is non-wetting phase, white part is wetting phase. 198 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.8 0.6 Ca=10 Ca=10" Ca=10': Ca=10': Ca=0.5 Ca=0.5 0.4 0.2 0.0 0.0 0.2 0.4 0.6 0.8 1.0 Sw Figure B.2: Relative permeabilities of the two phases as functions of the wetting phase saturation, for drainage case, M= 1 and different capillary numbers. Solid line, IP; dotted line, Ca= 10'4 ; short dash, Ca~ 10'2 ; and dash-dot line, Ca=0.5. 1 9 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.0012 0.0010 0.0008 Q . > 0.0006 0.0004 0.0002 1.0 0.0 0.2 0.4 0.6 0.8 Sw 0.07 Q . 0.05 0.04 0.03 0.0 0 .2 0.4 0 .6 0 .8 1 .0 Sw 2 .8 2.4 2 .0 0.0 0 .2 0.4 0 .6 0 .8 1 .0 Sw Figure B.3: The necessary pressure gradient of the system as a function of the wetting phase saturation, for drainage case, M= 1 and different capillary numbers. Top panel, Ca=10"4 ; middle panel, Ca-10'2 ; and bottom panel, Ca~0.5. 200 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Ca=10 Ca=10': Ca=0.5 0.8 0.6 0.4 0.2 0.0 0.0 0.2 0.4 0.6 0.8 1.0 Figure B.4: Change of fraction of wetting phase flow rate with wetting phase saturation, for drainage case, M= 1 and different capillary numbers. Solid line, Ca=10"4 , dash line, Ca-10'2 , and dotted line, Ca-0.5. 201 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 1 .0 M=0.5 M=1 M=2 M=5 0.8 0.6 0.4 / / 0.2 0.0 1.0 0.0 0.2 0.4 0.6 0.8 Figure B.5: Change of fraction of wetting phase flow rate with wetting phase saturation, for drainage case, Ca-\0~A and different viscosity ratio M values. Solid line, M= 0.5, medium dash, M -1, short dash, M=2 and dotted line, M=5. 202 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0 .0 0 1 6 0.0012 - [> 0.0008 - 0.0004 - 0.0000 0.0 0.2 0.4 0.6 0.8 1.0 Figure B.6 : The necessary pressure gradient of the system as a function of wetting phase saturation, for drainage case, Ca=10'4 and different viscosity ratio M values. Solid line, M= 0.5; medium dash, M= 1; short dash, M=2; and dotted line M= 5. M=0.5 M=1 M=2 M=5 2 0 3 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 1.0 Ca=10‘ Ca=10“ Ca=10'“ Ca=10': Ca=0.5 Ca=0.5 0.8 0.6 0.4 0.2 0.0 0.0 0.2 0.4 0.6 0.8 1.0 Figure B.7: Relative permeabilities of the two phases as functions of the wetting phase saturation for imbibition case, M- 1 and different capillary numbers. Solid line, Ca=10'4 ; short dash, C a-10'2 ; and dotted line, Ca=0.5. 2 0 4 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0.004 0.003 [ > 0.002 0.001 0 .0 0 0 0.6 0.8 1.0 0.0 0 .2 0.4 Sw 0.24 0.20 0.16 0.08 0.04 0.00 0.4 0 .6 0 .8 1 .0 0 .0 0.2 Sw 4.0 3.2 Q . 2 .8 2.4 2 .0 0 .6 0 .8 1 .0 0 .0 0 .2 0.4 S w Figure B.8: The necessary pressure gradient of the system as a function of wetting phase saturation, for imbibition case, M= 1 and different capillary numbers. Top panel, C a-10'4 ; middle panel, Ca= 10-2 ; and bottom panel, Ca=0.5. 205 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Ca=10 Ca=10': Ca=0.5 0.8 0.6 0.4 0.2 0.0 1.0 0.4 0 . 6 0.8 0.2 0.0 Figure B.9: Fraction of wetting phase flow rate as a function of wetting phase saturation, for imbibition case, M= 1 and different capillary numbers. Solid line, Ca= 10"4 ; short dash, Ca= 102 ; and dotted line, Ca=0.5. 2 0 6 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 1.0 M = 2 0 M =5 M=1 M = 0 .2 M = 0 .0 5 0.8 0.6 0 .4 0.2 0.0 0.2 0 .4 0.6 0.8 1.0 0.0 Figure B .1 0 : Fraction of wetting phase flow rate vs. wetting phase saturation, for imbibition case, Q z = 1 0 '4 and different viscosity ratio M values. Solid line, M= 2 0 ; long dash, M= 5; short dash, M= 1; dash-dot, M=0.2; and dotted line, M= 0 .0 5 . 2 0 7 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0 .0 1 6 M=20 M=5 M =1 M=0.2 M=0.05 0.012 0.004 0.000 0.8 1.0 0.2 0.4 0.6 0.0 Sw Figure B .ll: The necessary pressure gradient vs. wetting phase saturation, for imbibition case, Ca-10'4 and different viscosity ratio M values. Solid line, M= 20; long dash, M= 5; short dash, M= 1; dash-dot, M= 0.2; and dotted line, M= 0.05. 2 0 8 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. Figure B.12: Different flow directions of the moving interface, downwards in (a) and upwards in (b). 2 0 9 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission. 0 .0 0 0 o.c 0.4 0.6 S v M=1 pv/pi=0.5 0.020 0.015 0.010 0.005 0 .0 0 0 • Drainage o Imbibition • Drainage o Imbibition 0.030 0.025 t Q. 0.040 0.035 M=1 p y p p o .i o,4 o.e S v 0.040 0.035 • Drainage 0.030 0.025 | 0.020 0.015 0.010 0.005 0 .0 0 0 1.0 0.6 0 . 8 0.0 0.2 0.4 0.040 r 0.035 • 0.030 ■ 0.025 0.020 0.015 0.010 0.005 0,000 M=0.5 p / P p O .0 5 Drainage Imbibition 8o Sv 0.4 0.6 Sv Figure B.13: kco/r*2 as a function of vapor saturation at different viscosity and density ratios. 210 R eproduced with perm ission of the copyright owner. Further reproduction prohibited without perm ission.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Dynamics of combustion fronts in porous media

PDF

Effect of non-monotonic mobility at various scales in porous media

PDF

A study on relative permeabilities in three-phase flow in porous media

$Aspects of fracture heterogeneity and fault transmissivity in pressure transient modeling of naturally fractured reservoirs$

PDF

Aspects of fracture heterogeneity and fault transmissivity in pressure transient modeling of naturally fractured reservoirs

PDF

Characterization of permeability fields between horizontal wells using a hybrid of cross -hole imagery and repeat interference test

PDF

Characterizations of permeability fields and their correlations structure from pressure transients.

PDF

Analysis of retroviral production and transduction systems

PDF

Evaluation and analysis of microscale flow and transport during groundwater remediation

PDF

Pore Network Simulations Of Immiscible Displacements In Heterogeneous And Anisotropic Porous Media

PDF

Application of transverse flow equilibrium in miscible displacements

PDF

Characterization of plastics, carbon black filled rubbers and biodegradable polymers by ultrasonic techniques

PDF

Computers simulation of particle transport processes in porous media

$Fractal theory and groundwater flow in fractured media$

PDF

Fractal theory and groundwater flow in fractured media

PDF

Direct numerical simulation and stability analyses of three-dimensional swirling jets and wakes exhibiting vortex breakdown

PDF

Characterization of spatially correlated permeability fields using performance data

PDF

A percolation biofilm-growth model for biomass clogging in biofilters

PDF

A study of porous media displacements controlled by distributed thresholds

PDF

A method for automating delineation of reservoir compartments and lateral connectivity from subsurface geophysical logs

PDF

Creep deformation in superplastic yttria -stabilized zirconia ceramics and yttria-stabilized zirconia /alumina composites

PDF

Mobility Of Polymer Solutions In Porous Media

Asset Metadata

Creator Chen, Min (author)

Core Title Flow of fluids with yield stress in porous media

School Graduate School

Degree Doctor of Philosophy

Degree Program Chemical Engineering

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag engineering, chemical,OAI-PMH Harvest

Language English

Contributor Digitized by ProQuest (provenance)

Advisor Yortsos, Yanis C. (committee chair), Ershaghi, Iraj (committee member), Sammis, Charles G. (committee member)

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c16-566737

Unique identifier UC11340886

Identifier 3155392.pdf (filename),usctheses-c16-566737 (legacy record id)

Legacy Identifier 3155392.pdf

Dmrecord 566737

Document Type Dissertation

Rights Chen, Min

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the au...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus, Los Angeles, California 90089, USA