THE UNIVERSITY OF MICHIGAN INDUSTRY PROGRAM OF THE COLLEGE OF ENGINEERING COUNTERCURRENT GRAVITY SEGREGATION IN POROUS MEDIA tft tC James E. Briggs A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in the University of Michigan Department of Chemical and Metallurgical Engineering 1963 January, 1963 IP-604

,_e ti t r ~~ o.t3e - J (It f\A2e>

ACKNOWLEDGMENTS The author wishes to express his sincere appreciation to the many individuals and organization who have contributed to this research. Thanks are due particularly to: Professor Donald L. Katz, chairman of the doctoral committee for his interest, advice, and assistance in the preparation of this thesis. Doctor K. H. Coats, committee member, for initially guiding this work and for his continuting interest and advice. Professor L. E. Brownell, B. A. Galler, M. R. Tek and Mr. M. Sheffield, committee members, for their interest and helpful suggestions. Messrs. R. L. Nielsen, H. B. Kristinsson, M. C. Miller, S. C. Jones, and L. Yemin for their interest and valuable discussions throughout the course of the investigation. Mr. H. G. Olson of the Phoenix Memorial Project and the office and shop personnel of the Chemical and Metallurgical Engineering Department for their help in securing laboratory space, materials, and constructing equipment. The National Science Foundation for three years of support in the form of fellowships, and the Michigan Gas Association for their support during the last stages of the work and for financial aid in obtaining equipment. My wife, Anne, and my parents, for their encouragement. Messrs. R. E. Carroll, D. L. Danford and the personnel of the Industry Program of the College of Engineering for their efficient and accurate rendering of the final form of this dissertation.

TABLE OF CONTENTS Page ACKNOWIEDGEMENTS.........................,... ii LIST OF TABLES..................................................... v LIST OF FGURES............................................... vi LIST OF APPENDICES....................................... ix NOMENCLATURE................................... x ABSTRACT.. x,,,,,.,........,. xiii I. INTRODUCTION............... oesdsoo*,oooo o 1 II. SCOPE OF THISRESEARCH 3 III. LITERATURE SURVEY.....o*~o........... 5 A. Saturation Measurements.......... 5 1. Use of Neutrons in Measuring Saturations.. 6 2. X-Rays.......... 6 3. External Gamma Ray Sources. 7 4. Gas Phase Tracers............ 8 5. Liquid Phase Tracers 8 B. Mathematics of Multiphase Flow in Porous Media............ 11 IV. MATHEMATICAL DEVELOPMENTS.................................. * 16 A. Flow Through Porous Media............................. 16 1, Darcy's Law...................... 0 a 16 2. Capillary Pressure..................... 18 B. Mathematical Description of the Problem..20 C. Difference Equations...........................0.....,.0 23 V. EXPERIMENTAL PROGRAM............................*,.....*.. 29 A. Systems Studied.................... 29 B. Choice of Isotope............................... 30 iii

TABLE OF CONTENTS (CONT'D) Page C. Experimental Apparatus and Materials....................... 31 1. Preparation of Beads and Packing....................... 31 2. Detection and Counting System.......................... 32 3. Radioactive Solutions.................................. 34 D. Measurement of Liquid Saturations.......................... 34 E. Preliminary Experiments and Calibrations...5............... 39 1. Permeability................................... 5 39 2. Porosity... 40 3. Density of Glass Beads... 40 4. Detector Plateau............... 41 5. View Factors....... 41 6. Self Absorption of Aqueous Solutions................... 43 7. Halflife of 131...................................... 45 8. Absorption of Gamma Rays by Lucite..................... 45 9. Coincidence Losses...................................,. 49 F. Experimental Procedure..................................... 49 VI. EXPERIMENTAL RESULTS AND DISCUSSION.. 54 A. Capillary Pressure and Relative Permeability Curves........ 54 B. Transient Saturation Distributions......................... 62 C. Comparison of Data and Theory and Discussion of Results.... 73 D. The Need for More Complete Capillary Pressure Data......... 83 E. The Case of a Column Open at the Bottom........ 91 VII. RESULTS OFFIELD CALCULATIONS.................................. 96 VIII. SUMMARY AND CONCLUSIONS........................................ 103 BIBLIOGRAPHY..........................................105 APPENDICES................................. e..... *109 iv

LIST OF TABLES Table Page I Dimensions and Properties of Experimental Columns........ 31 II Accuracy of the Method for Measuring Saturations......... 39 III Values of Parameters for Experimental Runs............. 70 IV Sample Data Sheet....................................... 115 V Time Interpolated Saturation Values from Run No. 3...... 123 VI Time Interpolated Saturation Values from Run No. 4...... 124 VII Time Interpolated Saturation Values from Run No. 5...... 125 VIII Time Interpolated Saturation Values from Run No. 6...... 126 IX Time Interpolated Saturation Values from Run No. 7...... 127 X Time Interpolated Saturation Values from Run No. 8...... 128 XI Time Interpolated Saturation Values from Run No.10...... 129 XII Preliminary Readings for Run No. 5...................... 130 XIII Preliminary Readings for Run No. 6...................... 131

LIST OF FIGURES Figure Page la. Typical Relative Permeability Curves..................... 19 lb. Typical Capillary Pressure Curves......................... 19 2a, Difference Representation................................. 25 2b, Finite Difference Numbering Scheme for Columns............ 25 3. Sketch of Column and Pistons Used in Experiments.......... 33 4. Shielding Arrangement..................................... 35 5. Numbering Scheme for Determination of View Factors, Ai... 37 6. Detector Plateau.......................................... 42 7. View Factors (Ai)................................... o... 44 8. Self Absorption of Aqueous NaI131 Solution..............O 46 9o Halflife of I1310.......................................... 47 10. Absorption of I131 Gamma Rays by Lucite.................O 48 11o Coincidence Loss Measurement....O........................ 50 12. Photograph of Equipment........................ O O........ 52 13. Capillary Pressure Data, Water-Air System................. 55 14. Capillary Pressure Data, 43% Glycerin-Water Solution-Air System.................................................... oooooo o o ooooo c 56 15. Capillary Pressure Data, n-Heptane-Air System..o..o....... 57 16, Drainage Capillary Pressure Curve - Water-Air System, St. Peter Sandstone....................................... 59 17. Relative Permeabilities, 5.3% Connate Water............... 60 18. Relative Permeabilities, 10.8% Connate Water.............. 61 19. Experimental Run No. 3, Water-Air......................... 63 20. Experimental Run No. 4, Water-Air..., O.......o...O.OOOO. 64 vi

LIST OF FIGURES (CONT'D) Figure Page 21. Experimental Run No. 5, Water-Air........................ 65 22. Experimental Run No. 6, Water-Air........................ 66 23. Experimental Run No. 7, 43% Glycerin Solution-Air........ 67 24. Experimental Run No. 8, 43% Glycerin Solution-Air........ 68 25. Experimental Run No.10, n-Heptane-Air.... e............... 69 26. Effect of Average Saturation on Length of Run............ 72 27. Effect of Physical Properties on Length of Run........... 74 28. Composite Column......................................... 76 29. Simulation of Run 3, Using Imbibition Capillary Pressure Data..................................................... 79 30. Simulation of Run 4, Using Imbibition Capillary Pressure Data..................................................... 80 31. Simulation of Run 5, Using Imbibition Capillary Pressure Data..................................................... 81 32. Simulation of Run 6, Using Imbibition Capillary Pressure Data..................................................... 82 33a. Typical Capillary Pressure Curves........................ 84 33b. Hypothetical Path Through the Transition Region.......... 84 34. Hypothetical Network of Capillary Pressure Curves........ 87 35. Superposition of Capillary Pressure Curves on Final Distribution of Run 7................................... 90 36. Countercurrent Gravity Segregation in a Column Open at One End.................................................. 93 37. Hypothetical Saturation Distributions in a Column Open at One End.................................................. 94 38. Disappearance of Bulge and Transition Zone............... 95 39. Field Cases Treated...........................0........ X*97 vii

LIST OF FIGURES (CONT'D) Figure Page 40. Average Saturation as a Function of Time for the Various Field Cases................................... 99 41. Saturation Distributions - Field Case I................. 100 42. Saturation Distributions - Field Case II................ 101 43. Saturation Distributions - Field Case IV................ 102 44. Interfacial Tension, n-Heptane-Air.................... 119 45. Density of n-Heptane...................................... 120 46. Viscosity, n-Heptane, 43% Glycerin Solution............. 121

LIST OF APPENDICES Appendix Page A Stability Analysis.................. 109 B Calculational Procedure1..1..,,,,........... l13 C Data Processing Methods........................,.,. 114 D Physical Properties of Materials,................ 117 E Tables of Data........... 122 F History of Initial Distributions................. 130 ix

NOMENCLATURE Although a single letter tl, in some cases, several different meanings, the context in which each appears will serve to avoid confusion. a constant defined by Equations (23a), (23b), (23c) A cross-sectional area (cm2) or A geometrical view factor (dimensionless) C constant defined by Equation (25) or C constant defined by Equation (24) or C constant defined by Equation (57) D constant defined by Equationl (25) E error in numerical solution f flowing fraction (dimensionless), flow rate of one phase divided by sum of all phase flow rates or f the j-component error E defined by Equation (A-4) g = 9.66 x 10-4 (atm cm2/gm), gravitational acceleration or g the n-component of error E, defined by Equation (A-4) position or counting index J Leverett J-function (dimensionless); defined by Equation (6) k = 0.00565 (hours- 1), radioactive decay constant or k relative permeability (dimensionless) or k constant defined by Equation (24d) K single phase permeability (darcys) M coefficient defined by Equation (21)

n counting index or n total number of zones as in section IVD P pressure (atmo.) PC capillary pressure (atmo) Q flow rate (cc/seco) R variable defined by Equation (18) Ra true activity (counts/minute) Rm measured activity (counts,"minute) R activity of an infinitely thick solution (counts/minute) S liquid phase saturation (dimensionless) Sc residual wetting phase saturation (dimensionless) Sgr residual non-wetting phase saturation (dimensionless), measured in terms of wetting phase saturation Savg average liquid phase saturation avg S normalized saturation defined by Equation (54) Si aS/6Pc, reciprocal slope of capillary pressure curve (atm-1) t time (seconds) v velocity (cm/sec) VB bulk volume (cc) VP pore volume (cc) x thickness of solution (inches) X constant defined by Equation (58) y vertical coordinate Zry y-direction dimension or increment xi

Greek Letters angle of tilt from horizontal 7 interfacial tension (dynes/cm) IL denotes difference Q contact angle x constant defined by Equation (A-5) viscosity (centipoise) p density (gm/cc) 0 porosity (dimensionless) = P + pgy, flow potential (atm) Subscripts 1,2 denote particular phases avg. denotes average property g denotes gas phase i dummy index or position index j position index denotes liquid phase n time index or counting index nw denotes non-wetting phase T denotes total w denotes wetting phase xii

ABSTRACT The purpose of this study was to investigate the factors influencing the rate of vertical, countercurrent gravity segregation of two phases in porous media. This process is important in underground gas storage operations in which the injected gas rises to displace the water above it in the formation. Experiments were performed using a closed vertical column packed with glass beads. An initial distribution was established in which the liquid phase was aobe the gas phase. As the two phases exchanged places, their movement was followed using a quantitative method based on measuring the gamma radiation from I131 tracer in the liquid. The liquids used were water, a 43% by weight solution of glycerin in water, and n-heptane. The gas phase was air. Use of the various liquid phases allowed the determination of the effects of viscosity, density difference, and interfacial tension on the rate of the segregation process. The length of time required for completion of the process was found to increase with decreasing average liquid saturation. For a constant average saturation the length of time was found to increase for incrzeasing values of P/y(Ap). An attempt was made to simulate the experimental results by numerically solving the equation describing the segregation process. A fair degree of success was obtained using imbibition capillary pressure data for the calculations but the principal conclusion was that xiii

static capillary pressure data as ordinarily obtained can not be used to describe a process in which both drainage and imbibition occur at the same point at various times. In view of this fact, a method is presented for obtaining and using the necessary data. Numerical calculations were also performed using typical reservoir properties in order to study the effect of permeability variations on the rate of drainage of water from a column of porous media closed at the top but communicating with other porous media at the bottom. This study showed that thin, low permeability layers can greatly slow down the rate of drainage. xiv

I. INTRODUCTION Underground storage of natural gas has grown steadily in importance in recent years. It is an excellent answer to the problem of maintaining reserves at a location near the market and it allows more economical year round use of pipeline facilities. It is natural therefore, that more and more research is being done on problems concerning storage operations. A great deal of the research on natural gas storage has been concerned with finding ways of calculating injection-production schedules. This refers to the case in which for a specified future production schedule it is necessary to determine an economic injection schedule which will provide the amounts and pressures required. Other work such as determination of properties of the underground formation by core analyses and surface tests has been done in order to gather the information necessary for the calculations of field performance. Most of the work described above has been based on the theory of single phase flow in porous media, i.e., at any given point there is either gas or water flowing but not both. At the same time, research in the oil industry has been concerned with applying the concepts of two-phase or multi-phase flow in porous media to the many problems associated with the various stages of oil recovery. In many cases the problems of natural gas storage and of oil recovery are analogous. As an example, the injection of water into an oil bearing formation in order to force the oil out is very similar to the injection of gas into a water filled reservoir.

-2Although there are differences caused by the properties of the fluids involved, these operations should both conform to the theories of twophase flow in porous media. With the above in mind, the work described in this dissertation was performed in order to test the theory as applied to the case of countercurrent gravity drainage in a gas storage reservoir. This process in which injected gas rises to displace the water initially occupying the formation, is one which includes both capillary and gravity forces and thus provides a severe test of the concepts of two-phase flow in porous media.

II. SCOPE OF THIS RESEARCH The purpose of this research is to obtain some knowledge of the variables which influence the rate of gravity segregation in a vertical column of porous media. There is a practical aspect to the problem in that the rate of gravity segregation can influence the performance of gas storage reservoirs as well as oil reservoirs. When gas storage is started in an aquifer (an underground formation containing water but no gas initially) it has been found(42) that sometimes the injected gas will find a high permeability layer and spread out in a pancake or disk having a thickness of perhaps one or two feet and a radius of 600 to 700 feet. After the pancake has formed, gas is present below a water saturated zone and it is to be expected that the two phases will exchange places because of the difference in densities. It has been observed(25) that this exchange does take place but very slowly (perhaps as long as 3 to 5 years). It is desirable that the gas should rise to the top and displace the water downward as quickly as possible so that a "bubble" will form at the top of the structure and the storage operation can proceed by enlarging the bubble by means of gas injection. The work described here is divided into two parts. The first consists of performing experiments which are similar to the gravity drainage problem described abov-e and which provide a severe test of the mathematics of two phase flow in porous media. The experimental work, as described in a later section, consists of measurements of -3

saturation distributions as a function of time and position in a vertical column of glass beads closed at both ends. At the start of a run the gas phase is at the bottom of the column and the liquid phase, containing a radioactive tracer, is in the upper part. As the two phases exchange positions, measurement of radioactivity at different points of the column allows liquid saturations to be determined. This arrangement corresponds to the situation in a shut-in reservoir where an amount of gas has been injected and then the wells are closed. The second part of this research is concerned with mathematically predicting the results obtained experimentally and also with carrying out calculations using physical parameters typical of actual reservoirs. In contrast to most of the work previously reported in the literature, the analysis to be used here will include both capillary and gravity forces. The calculations using reservoir properties will deal not only with the case of a shut-in reservoir but also with the case of a reservoir where injection is continuous. The purpose of these calculations is primarily to investigate the effect of permeability variations on the rate of gravity drainage.

III. LITERATURE SURVEY The literature concerning flow through porous media is vast. For the purposes of this report the two topics of principle concern are the various methods of measuring transient saturation distributions and the mathematics of two phase or multiphase flow in porous media. A. Saturation Measurements An essential part of the research described here is the measurement of saturation distributions during unsteady state flow in porous media using a radioactive tracer in the liquid phase. The purpose of this section is to summarize some of the work on saturation measurement which has been reported in the literature. Methods for measuring saturations which have received, consideration include light transmission(36), magnetic susceptibility(52) electrical resistivity, X-ray absorption, neutron diffraction, and radioactive tracers.(24) The resistivity method has use in the determination of saturation profiles but needs checking by independent means and is also not applicable to oil-gas systems. The X-ray method has been used successfully, but requires elaborate equipment and calibration. At present, gamma ray absorption is practical only for very large cores.(l8) The neutron method does not appear adaptable to measurement of saturation variation in short cores since it is a scattering method. -5

-61. Use of Neutrons in Measuring Saturations Neutrons are very penetrating and their scattering by substances containing hydrogen, such as oil, is qualitatively different from that caused by other materials. A convenient source of neutrons is an intimate mixture of radium and beryllium. The alpha particles from the radium and its disintegration products cause the ejection of neutrons from the nuclei of the beryllium atoms, Such a neutron source also produces gamma radiation which may disturb the apparatus used to detect the neutrons. The number of slow neutrons emerging from an object irradiated with fast neutrons gives a rough index of the quantity of hydrogen it contains. This is the means used to determine saturations. One way to measure neutrons is with an ionization chamber lined with boron or lithium or filled with boron triflouride gas. This is the most sensitive method although a high gamma background can cause trouble. Another is to induce beta activity in a suitable material such as silver, indium, or rhodium and use a Geiger-Muller counter to measure the beta radiation.(5) It should be noted that the neutron method is only good for near steady-state measurements since quite a bit of time is nenessary for each reading. Also, the method can not easily distinguish water from oil. 2. X-Rays In order to use x-rays to measure saturations, powerful absorbers used to accentuate the effect of one component must be

-7available and soluble in that component alone. As an example, cerium sulphate or cerium chloride make the detection of water possible, In using x-rays, two sets of data are usually taken(4'12'28'33). one for the sample being tested and the other for a standard core run at the same time. This is because of the possibility of fluctuations in beam intensity and also because the experimenters used white x-rays rather than monochromatic x-rays and sometimes obtained non-linear calibration curves. The reason for using white x-rays is that filtering to obtain monochromatic rays cuts down the intensity of the beam. The problem of fluctuations in voltage is covered by Morgan, McDowell, and Doty(33), who explain that the output of an x-ray tube varies in an approximately linear manner with the tube current, but to the second or higher power of the anode voltage. To overcome this a second beam can be taken off at a small angle to the main beam and used as a control. X-rays have been questioned as to the uniqueness of saturation distributions obtained during cycling of flow.(15530) 3. External Gamma Ray Sources There has not been much published about this technique because it generally requires that a thick column be used in order to get an appreciable absorption of gamma rays. Isb:n et al.(22) and Hwa and Beckmann(2) describe experiments on liquid holdup in packed beds (large particles as compared to porous media) in which 13 millicuries of Hg205 were used to measure holdup (or saturation) in a 6 inch

-8diameter column and Thulium-170 was used for a 12 inch tube. But, once again, these diameters are fairly large, thus having more void space and therefore more liquid in a given volume than there would be in a porous bed. Recently, Nielsen( 36) has used an external source of visible light with photo cells to detect the light transmitted through a rectangular 1/2 inch thick model. An advantage to this method is the absence of potentially dangerous radioactive materialsc 4. Gas Phase Tracers In the work under consideration in this dissertation, it would have been advantageous to have had a suitable gas phase tracer. This is because the gas is the non-wetting phase, i,e,, the gas does not in general touch the solid material because the liquid phase is preferentially adsorbed onto the surface and forms a covering film. Thus the problem of adsorption of the tracer onto the porous media would be greatly reducedo A gas phase tracer must have a low solubility in the liquid phase and must be a gamma emitter with a half-life long enough to allow the experiments to be carried out. Kniebes, Burket, and Staats(27) have used argon-41 to measure natural gas flow in a pipe but the 109 minute half-life is a serious drawbacko Clayton(9) also measured gas flows in a pipe but used krypton-85. This isotope has a long halflife but a very low gamma activity. 5. Liquid Phase Tracers As mentioned above, a tracer in the liquid phase is prone to adsorption onto the solid material. But, since there is a large number

-9of isotopes and compounds to choose from, much more work has been done with the liquid phase than with gas tracers. Watkins and Mardock(49) give a very useful list of the characteristics which a water phase trancer should have: a) low adsorption on solid material, b) high solubility in water, c) low or negligible solubility in oil or other phases, d) a wide range of solubility of compounds that may be formed by chemical reactions with ions present in the solid material or water, e) high detectibility in low concentrations by portable apparatus, f) general availability at low cost, g) non-hazardous under ordinary working conditions. Russell, Morgan, and Muskat (41) did experiments on the mobility of interstitial water, i.e., the movement of that water which forms the irreducible saturation and is not removed by flushing with another phase. They used radioactive vanadium produced by bombarding titanium metal in a cyclotron with deuterons. Two radioactive tracers were tested as a means of determining core saturations in multiphase flow studies by Josendahl, Sandiford, and Wilson. (24) Cesium-134 chloride was tried as a water phase tracer, but complications in its use in low permeability cores resulted from adsorption of cesium by the core. Iodo-131 benzene proved very satisfactory as an oil phase tracer.

-10Hull(20) used cesium-134 chloride as a water phase tracer but was not working with porous media but rather with flow through pipes and therefore did not have the adsorption problems mentioned above. Iodine-131 has been used quite frequently for studies of the type under consideration here. Stanley (46) used I13l to study the effectiveness of sand filters in the removal of suspended particles from water. For the particles, he used a suspension of ferric oxide. The isotope was inserted by mutual coagulation between negatively charged colloidal AgI and positively charged hydrous ferric oxide. Counting at various depths was accomplished by the use of a geiger tube having a high efficiency in gamma counting, mounted in a specially designed shield. Unfortunately, no details are given concerning the amounts of tracer used or the counting rates obtained. Watkins and Mardock(49) found that radioactive iodine is a satisfactory tracer material because it has most of the characteristics listed previously for an ideal water phase tracer, namely: a) loss of iodine, as iodide, through adsorption on solids in low compared to other potential tracers, b) iodide remains preferentially soluble in water and insoluble in oil, c) the possibility of an insoluble iodine compound being formed and precipitated is low, d) cost is low and I131 is readily available,

-11lAn extensive study of the application of radioactive tracer techniques to flow through sands was carried out by Coomber and Tiratsoo. They followed the movement of oil and water through a sand pack using first ethylene dibromide and then elemental iodine. The dibromide was made radioactive by irradiation but apparently there was some molecular breakdown and recombination leaving free bromine-82 in the liquid. Since this caused the radioactivity to be distributed in the oil and water phases and also since the bromine was easily vaporized and inhaled, their attention turned to I131. With this element there was no need for complicated syntheses of organic compounds, because elemental iodine is considerably more soluble in hydrocarbon oil than in water. For counting, a G.M. end window tube was employed. Josendal et al.( ) object to the use of elemental iodine because the presence of reducing substances such as hydrogen sulfide in either the oil or water, or the presence in the water of complexing agents would greatly affect the distribution of activity in the phases. B. Mathematics of Multiphase Flow in Porous Media Since one of the main purposes of the work reported here is to find a way of mathematically describing gravity counterflow segregation, the principal references cited relate to the mathematics of two phase flow in porous media. In general, references to other aspects of reservoir engineering are cited in the appropriate sections of this dissertation.

A mathematical description of immiscible two phase flow in porous media consists of combining Darcy's Law, continuity equations, capillary pressure relations, and if necessary, an equation of state for the phases. One example of an equation for one-dimensional, incompressible flow is the "fractional flow equation" of Buckley and Leverett(6'29), given here as Equations (1) and (2). Kk1 A [PC Ss fg g(1) kgg l+ kg9 Sg = 7 (2) Another example is Equation (19) of section IVB which is the one used in this work. The above equations and various forms of them have been solved by many investigators for various cases under study. Before the advent of high speed computers, linearizing assumptions were necessary in order to solve the equations. One of the most common assumptions is that capillary and gravity forces are unimportant. This allows Equation (1) to be written as: f 1 (3) g k2g 1+ kgetl Equations (2) and (3) can be solved graphically but the resulting answers are triple-valued, i.e., at some points three different saturation values are predicted to exist at the same time,(6) and it is

-13necessary to be careful in choosing the correct values.(7,50) After the use of computers became wide spread, it was possible to include either capillary or gravity forces or both in the solutions of the equations. Some of the investigators who solved the equations for various flow situations and the assumptions made by them are discussed below. Jones-Parra and Calhoun(23) solved the fractional flow equation for a linear displacement without the effects of capillary or gravity forces. Cardwell and Parsons(8) neglected capillary effects and solved the case of gravity drainage from a column open at both top and bottom. Terwilliger et al.(48) also solved the drainage problem and achieved excellent agreement with a concurrent experimental investigation. They state that their work shows that static capillary pressure and relative permeability data can be used to calculate transient behavior. West, Garvin, and Sheldon(51) treated horizontal linear and radial flow in solution gas drive reservoirs neglecting capillary forces. Blair, Douglas, and Wagner(3) solved the one-dimensional equations numerically for linear two phase flow of incompressible and immiscible fluids including capillary forces but neglecting gravity. Blair(2) also included capillary forces to predict saturation and pressure distributions in linear and radial systems open to water at one end (i.e., countercurrent water imbibition). McEwen(32) used the fractional flow equation including capillary effects but no gravity forces to calculate the performance of a water flood in a linear system. He found that at high rates the solution approached that obtained by neglecting capillary forces. Hovanessian and Fayers solved the one dimensional equation

-14including gravity and capillary forces with injection and production. They concluded that gravity has a pronounced effect on both saturation and pressure distributions. The effect of linearizing assumptions (neglecting capillary and/or gravity forces) was studied by Fayers and Sheldon. (14) They conclude that the inclusion of the capillary term eliminates the triple valued saturation behavior. They used a Lagrangian form (S is an independent variable) of the equations which would allow triple values to be computed but found them only at high flow rates and even then could eliminate them by reducing the time step size, They also point out some interesting compatibility relations which must be taken into account. Several authors have treated the problem under consideration here, i.e., one-dimensional gravity segregation of incompressible, immiscible fluids, Gilmour(16) derives essentially the same equation as that presented in section IVB but makes no attempt to solve it, Martin() neglects capillary forces and derives equations for closed systems where the fluids are initially uniformly distributed in the column or where the phases are completely separated with the heavy fluid on top. He imposes the no-flow condition at the ends of the column by requiring the saturations at these points to be such that the relative permeability of one of the phases is zero. The author feels that rather than setting saturations at the ends, it is better to reflect the flow potentials across the boundary, thus leaving the saturations free. Templeton(47) did not attempt a complete solution of the equations but did perform experiments very similar to those described

in this report. Perhaps the main difference between his work and this is that his initial distributions were nearly uniform with an average liquid saturation of 50 to 60% whereas the experiments described here had semi-segregated initial distributions and a wider range of average saturations~ As will be pointed out later, this author feels that certain aspects of the problem show up more clearly in this work than in Templeton's~ The discussion above has been limited to one dimensional case and also immiscible fluids. It should be noted that Douglas, Peaceman, and Rachford(l3) have published an excellent paper on two dimensional flow including capillary and gravity effects. Their paper also serves as an example of the application of the alternating direction implicit method of numerically solving two or three dimensional equations~ Nielsen(36) provides another example of the solution of a two-dimensional problem along the same lines as Douglas et al. Examples of solutions of the problem of miscible flow are very few in number. An excellent one is that of Peaceman and Rachford (38)

IV MATHEMATICAL DEVELOPMENTS The purpose of this section is to present the derivation of the equation for countercurrent gravity drainage. Two cases are considered; that of a column of porous material closed off at both top and bottom and that of a column closed at the top but communicating with other porous media at the bottom. As an introduction the basic concepts of two phase flow in porous media are presented. Then follows the derivation, and finally the transformation of the differential equation to a difference equation. A. Flow Through Porous Media The purpose of this section is to summarize briefly some of the concepts used in describing two phase flow through porous media. lo Darcy's Law The equation relating potential difference to flow rate for a single phase in porous media is Darcy's Law.* Q = K grad. (4) A where Q = flow rate (cm3/sec) A = cross sectional area (cm2) p = viscosity (centipoise) = flow potential (atm) = P + pgy * Written for incompressible fluids

-17It should be noted that in general: = P + pgy + pv2 + electrical forces, etco but, since in this work P and pgy are the only important forces, all other terms have been dropped and thus T is always taken to be P + pgy. The constant of proportionality K is called the permeability and is a function of the medium and not of the flowing fluid. The unit of permeability is the darcy (d) or millidarcy (1 md = 00ld) and is defined as follows: A cube of a porous medium 1 cmo on edge will have a permeability of 1 darcy if water flows between the front and back faces at a rate Q of 1 cm3/seco under a pressure drop of 1 atm. at a temperature of 68 F, where the viscosity is 1 centipoiseo(26) When there are two or more fluids flowing through a porous material, the above equation is modified to take into account the fact that only a portion of the original cross-sectional area is available for flow of a given phase, The modification consists of introducing relative permeabilities into Darcy's Law~ As an example, the following equations can be written for a system in which there is two phase flow. Q1 gradl (5a)..... grad 1 A -- = - 2 grad ~2 (5b) A ~2

The relative permeabilities must be determined experimentally for the media and fluids in use although several methods have been proposed for correlatirng them or deriving them from other properties of the system(35'5 ) A typical set of relative permeability curves is shown in Figure la. Several references(43''4453) support the contention that relative permeabilities are a fumction of the saturation onlyo Although it is known(37) that the relative permeabilities are dependent to some extent on which phase is driving the other, the curves used in this study are assumed to be without any hysteresis effect. 2. Capillary Pressure Because of the minute size of the pores and channels in which flow occurs, capillary forces can be a very important factor in the flow of more than one phase through porous media. In the majority of cases, one of the phases wets the solid more than the other, i.e., the wetting phase forms a film over the solid and the second phase is kept out of contact with the solid. Because of the curvature of the interface between the two phases, there exists a pressure difference between them. This pressure difference is defined as the capillary pressure Pco Pc = P non-wetting - P wetting The curves in Figure lb are typical. Leverett(29) proposed a useful correlation for capillary pressure curves in various unconsolidated porous materials. His dimensionless "J" function is defined as:

r - kw kw O 1 Figure la. Typical Relative Permeability Curves. w C,) I. 0 WETTING PHASE SATURATION Figure lb. Typical Capillary Pressbilite Curves. wL xI

J(S) = c (K)12 (6) 7 cosG X where Pc = capillary pressure K = permeability 0 = porosity 7 = surface tension 0 = contact angle The contact angle is defined as the angle measured through the wetting phase which the interface between the two fluids will assume at the point of contact with a solid. It has been found that even for a single porous material and a given pair of fluids, there exists a difference between capillary pressure curves obtained when the wetting phase saturation is increasing (imbibition) and when the saturation is decreasing (drainage)o This hysteresis effect can have a great effect on the flow behavior of a system as can be seen later in this report. B. Mathematical Description of the Problem The problem to be solved is that of two phase counter-current flow in porous media. The flow takes place in a vertical column and is considered to be one-dimensional. The fluids are assumed to be imcompressible and immiscible and to have constant physical properties such as viscosity and density. The porous media may be non-uniform with respect to porosity and permeability.

-21Darcy's law may be written as QW -KA ay (7a) kw W knw ~nw Qnw 6k -w (7b) lnw aY where the subscripts w and nw refer respectively to the wetting and nonwetting phases. A continuity or material balance equation may be written for each phase as: aQw = _,A S (8a) ay at OQnw = _A Snw (8b) ay at Using the following definitions, Pc = Pnw - Pw (9) Tw = + PwgY (10) ~nw Pnw + Pnwgy (11) Pc can be written as Pc = nw - ow + pogy (12) Combining (7b) and (12) results in nw KA kn (P + - gy) (13) Pnw 6Y Since we are dealing with counter-current flow, it follows that: Qw + Qnw = Q = 0 (14)

-22and therefore (7a) and (13) may be added to produce: kw a$w knw'"Pc cw g \ - KA tw ay -KA -- 0 +a )= (15) or knw APWY)) * ~n(16) -w (PC - Apgy)) (16) kw+ knw Now, using Equations (7a), (8a), (10), and the definition: S = 1 - Snw = Sw Equation (17) results. a [ Kkwknw _ (Pc )gy) = (17) ay -kwilnw + knwPw ay at For convenience, the variable R is defined as R = 1/2 (P - Apgy) (18) and Equation (17) becomes a F Kkwknw aR] =- ()9 aR ay Lkwinw + knww y- (19) Equation (19) is second order and thus requires two boundary conditions and an initial condition. The initial condition is supplied as the known saturation distribution at time zero. For the case of a column closed at both ends, a reflection of the potentials at the ends allows no flow across these boundaries. Thus:

-236y J, L Ka L J or (o w)L = = (P - gy), 2 (a TY nw 0) L 3y 0L0o, L or (-) = 0,L The difference equation used to represent the above equation and boundary conditions is presented in the following section while a description of the numerical procedure for obtaining solutions and an analysis of the stability of the difference equation are given in Appendices A and B. For the field calculations discussed in section VII the boundary condition at the bottom of the column must be such as to maintain the saturation at a constant value. Therefore we may write: R(y = 0) = X C. Difference Equations In order to lessen difficulties in stability and convergence, an implicit difference form was chosen to represent the differential equation under consideration. The equation is: a F Kknwkw j=R S' -R (19) ay kw.lnw + knww aY at where S' = aS/aPc.

-24Following the procedure of Douglas, Peaceman, and Rachford(l3) the term Kkwknw/*winw + knwpw) is treated as an "interblock" coefficient and evaluated between blocks. The difference form can be written as Mj+2,n+2 (Rj+l,n+l j,n+) (Rj,n+l - Rjl,n+) (Rn+l - Rj 1) tJ-2 (Rj n+l - Rj,n) (20) At where: ~Mj+ ( n+ay)2 =Kw w [ K+wKk +~ kwknw ]n (21) 2 12 2 2(zy)2 _ Lkwnw + knwwjl knwPlw + kwxlnw-j n4 Figure 2a shows where the variables are evaluated. Equation (20) may be rearranged to yield: Mj+1 n+ Rj+l,n+l + (-Mj+,n+ - Mj n++ + S.n+) Rj n+l 2) 2 2 2 2 2 At + Mj_. L R _ S,n+ Rj,+ (22) -,n+ j-1, n+l j,n At Thus a tri-diagonal matrix is generated. A typical equation may be written as: aj,j-_l Rji + aj,j Rj + aj,jl Rj+l Rj kj (23) where aj,j_l = Mj,ni (24a) jj = -Mj+2 n+2. _ Mj_, n+ n+2 (24b) aj,n+l = Mj+.,n+ n. (24c) k = AS,,n+ RO (24d)

-25WI X DENOTES COEFFICIENT I —n+ 0 C) ILI j-1 j+l DISTANCE COORDINATE (y) Figure 2a. Difference Representation. j =jm I I j+1 VERTICAL COLUMN — j+1/2 —. j-l L.- j=2,y=O Figure 2b. Finite Difference Numbering Scheme for Columns.

-25A simple technique for inverting a tri-diagonal matrix has been given by Richtmeyer (40) The following presentation of the method follows notes prepared by K. Ho Coats(lO) for a course at the University of Michigan. It is assumed that there exists a relation of the form: Rj_1 = CjRj + Dj (25) Insertion of Rj_1 from (25) into (23) and rearrangement leads to: a. j+l k. - DR. j, _ R+ + J -iDj (26) j,j + aj,j-iCj j+ a, iComparison of (26) to (25) shows that C. aii+1 (27) aj + aj-l D kJ = Ikj - a.j_-lD (28) 3+1 aj,j + aj,j-_Cj Thus if any Cj, Dj can be found, all others may be calculated using (27) and (28). In the closed column, at the bottom (see numbering scheme in Figure 2b), R1 = R3 since ( =R) = O R3 -R1 (29) o 2zy Thus Equation (23) may be writtn with j = 2 as a2,1R1 + a2,2 R2 + a2,3 R3 = k2 (30)

or (a2, 1 + a2, ) R3 + a2,2R2 = k2 (31) Solving for R2 R2 k2 a&2,1 +a2 R C3R3 + D3 (32) a2,2 a2,2 Thus C3 = a2,1 + 2, (33) a2,2 D3 = k2/a2,2 () and C4... Cjm, D4..o Djm may be calculated using Equations (27) and (28). At the top of the column, Equation (23) can be written with j = jm am,jm-l Rjm-l + ajm,jm Rjm + ajm,jm+l Rjm+l = kjm (35) and since in the case of the closed column Rjm+l = Rjm_l, (a + ajm.j+ aj,jm R. =k. (ajm,jm-l + ajm,jm+l) Rjm-1 + ajm, jm =jm (36) also R.jml = C. R. + D. (37) jm-l jm Jm jm Therefore (ajm,jm-l+ ajm,jm+l)(CjmRjm + Djm) + ajmjmRjm = kjm (38) or R - kjm - Djm[ajm, jm-l + ajm, jm+l] jm a + [a + a. ]C m jm,jm [jm, jm-l + jm,jm+l ]C jm

-28Thus from the known values of Cj, Dj and the coefficients aj.j and kj, Rjm can be found and then using Equation (25) all the Rj can be calculated. For the field calculations, the boundary condition at the bottom of the column requires that R2 = X. Thus Equation (23) is written with j = 3 as: a312R2 + a3,3R3 + a3,4R4 = k3 (40) 3 a3,3 a3,2 Therefore C4 = - a3,4/a3,3 (42) k3 - a3-1 2X (4 i 4 ra3me(43) All other relations remain the same.

V. EXPERIMENTAL PROGRAM The six sub-sections presented here are concerned with all the factors involved in measuring transient saturation distributions by means of a radioactive tracer in the liquid phase. The apparatus is described and the equations for calculating saturations from radioactivity measurements are derived in detail. All preliminary experiments and results are presented except for capillary pressure curves which are discussed in section VIA. The procedure followed in performing the experiments on countercurrent gravity drainage in a closed column is given. Ao Systems Studied In order to obtain some insight as to how various parameters affect gravity segregation, several different fluid systems were investigated. They were: 1, Water - air 2, 43% glycerin solution - air 30 Heptane - air. The glycerin solution provided an easy means of changing the viscosity of the liquid phase while leaving other properties such as capillary forces and density essentially unchanged from those of water. Heptane was employed as the liquid phase in one run in an attempt to show the effects of a marked change in capillary forces. -29

-30B. Choice of Isotope In all of the systems studied, the radioactive isotope used was I13l as NaI in aqueous solutions or 12 in the run using heptane. The choice of what isotope to use was simplified somewhat by results previously reported in the literature. Having in mind the problems mentioned by previous experimenters it seemed that the criterion to be met was a half-life of less than one month but longer than a day since each run took up to a week. Also, the isotope had to give off gamma rays and be soluble in some form in the liquid phase. The half-life stipulation reduced the number of possible isotopes to about twenty, ranging from Cadmium - 115 (53h.) to Cerium-141 (32.5 d.). Setting an acceptable price at less than five dollars per mc restricted the choices to the following: Cd.115 Re186 Au198 i131 Hg197 Bal4o Mo99 Rb86 Sb122 Cr51 Au199 Ce141 Any of these isotopes would probably have been worth trying but since good results have been reported with I131 and since it is readily available and cheap, the experiment was designed for 1131. Another reason for choosing 1131 was that Archibald(1) reports 131 tests made to determine contamination of glassware by I13 The results showed that radioactive iodine does not contaminate glassware but does

-31adsorb on metals - especially brass. Thus it was felt that adsorption 131 of I on the glass beads used in these experiments would be negligible, C. Experimental Apparatus and Materials The apparatus used in the experimental work consisted of the column in which flow took place, the detection and counting system, and auxiliary equipment for introducing fluids into the column. Two separate columns were used, the first being discarded after two experimental runs for the reasons given below. The dimensions and properties of the columns are shown in Table I. TABLE I DIMENSIONS AND PROPERTIES OF EXPERIMENTAL COLUMNS Column Length Diameter K(Darcys) 1 341/16" 2" 11.84 0.389 2 297/8" 2" 11.79 0.394 The glass beads used were obtained from the Minnesota Mining and Manufacturing Company and are designated by them as No. 996. They are spherical with an average diameter of 113 microns with 90% of the beads within 15 microns of the average. 1. Preparation of Beads and Packing The glass beads were treated by soaking tnem in a ~uyl niu~ aqueous solution for several days and then washing with distilled water and drying at about 2000C. for four or five days. This helped to insure that the beads would be strongly water-wet. A lucite tube (2" I.D. 1/4" wall thickness) was fitted with end plates and tie rods and a

-32piston with a nickel screen to retain the beads was inserted at one end, The column was filled halfway with water and the beads were slowly poured in the top while the column was vigorously tapped at the pack-water interface, After the pack had risen to the desired height, a piston was inserted into the top end and the bolts on the tie rods were tightened until there was no further compaction of the pack. The pistons were then sealed to the lucite tube with Duco cement. The details of the column and pistons are shown in Figure 3. One of the reasons for building two columns was that the first one had pistons with 1/8" diameter bores and no grooves for distributing the fluids and this caused great difficulty in filling and emptying the column because of plugging in the screens and piston bores, The second reason was that the column was too long and thus there were over two inches at each end which were not covered by the traverse of the detector. 2. Detection and Counting System The detector employed in this study was a model DS5-2 scintillation detector purchased from the Nuclear - Chicago Corporation. The detector utilized a 1" x 1" thallium activated sodium iodide crystal and a 10 stage photomultiplier circuit to provide pulses to drive a binary scaler. The scaler, manufactured by Nuclear - Chicago Corporation, recorded the pulses and displayed the total by means of a mechanical counter and interpolation lights. The positioning of column and detector is shown in Figure 12. d ac the detector was oalways adjusted to be 1/2" from the face of the column and facing the centerline, A total of ten three inch high zones could

ENDPLATES (1/4"thick brass, 3 1/2" diameter) TIE RODS (1/4" diameter steel rods) LUCITE COLUMN (2" I. D., 2 1/2"O.D.) 1/4" BOLTS STAINLESS STEEL TOGGLE VALVE 1/4"DIAMETER HOLE THREADED FOR VALVE \FACE OF PISTON SHOWING DISTRIBUTOR NETWORK Figure 3. Sketch of Column and Pistons Used in Experiments.

be covered by a traverse, Figure 4 shows the detailed geometry of the shielding and crystal with respect to a 3" zone, This shielding arrangement was decided on after experimenting with various schemes such as a lead plug with 1/8" diameter holes to collimate the radiation and/or various amounts of shielding around the detector. The shielding used kept the "noise" level down to about 20% of the reading. 3. Radioactive Solutions As mentioned previously, the tracer used in these experiments was I1310 The tracer was supplied as a virtually carrier free aqueous solution of NaI3 by the isotope unit of the University of Michigan Hospitalo In the water-air and glycerin solution-air runs, the liquid was made radioactive by merely adding an appropriate amount of NaI151 solution. In the heptane-air system it was necessary to oxidize the iodide to iodine using concentrated HN03 and then extract the heptane. Because of serious adsorption problems it was also necessary to add a large amount of carrier (inactive NaI) to the solution being oxidized so that essentially all of the adsorbed iodine was inactive. Do Measurement of Liquid Saturations The purpose of this section is to derive the method used for relating radioactivity measurements to saturations. The problem that arises in the use of a radioactive tracer in experiments of this type is that one is attempting to measure the activity of a small part of a large radioactive mass. Thus there is interference or noise caused by radiation from parts of the column other than the one being measured.

-35-'- - OH) wc, 00 I~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~I w~~~~~~~~~wc 1 — >wii[ cn O Cl)o~~ w-J_w w w cn 1.. Cj) -J V)~~~~~~~~~~~~~~~~~~~~~~~~~~~/ U) on < c ) _I VI Z -J I 0z -i —I

This problem can be partially overcome by the use of heavy lead shielding but a point is soon reached where a significant reduction in noise requires such a large amount of shielding that the equipment becomes too heavy and bulky for efficient operation. It was found experimentally (section VE-6) that for a single zone of the column the relation: Rai = CSi (44) is valid. When a number of zones are considered, the interference with the reading at zone j is accounted for by assuming that the measured radiation at point j (see Figure 5) can be represented as: n n Rm = AiRai = C L AiSi (45) i=l i=l where Rm. = measured activity at point j (counts/minute) a Rai = actual activity at point i Ai = coefficients or view factors depending on distance from detector to point i (Aj = 1, Aj_1 =Aj+l, etco) C = constant factor in Equation (44) S. = saturation at point i equals fraction of void space in zone i filled by liquid. It can be seen that since there are n zones, the above equation can be written n times resulting in a system of linear equations in n unknowns, i.e,,

-37i=n COLUMN - i=j DETECTOR| I I i=3 i=2 i=1 Figure 5. Numbering Scheme for Determination of View Factors, Ai.

-38CA1S1 + CA2S2 +............ + CASn = Rm CA2s1 + CA1s2 +................ + CA 2 1 1 2 + Co + CA IS = Rm n-in 2 (46) CAnS1 + CAn 1S2+ + CASn = Rmn Therefore, if Rmj, Ai, and C are known, Si can be determined by solving the above system of equations thus resulting in saturation measurements as a function of position and the time of the readings. The determination of the Ai is discussed in section E5. The constant C was determined for each run by carefully taking readings on an unknown saturation distribution and solving for the n quantities Rai or CSi in Equation (46). Use is then made of the relation n n / Rai = C Si = nCSavg (47) i=l i=l or Ra. C = (48) n S avg. It should be mentioned that since I 31 has a half-life of approximately 8 days, the decay during a run can be appreciable. Therefore, the actual form of Equation (45) is n -kt Rmj = C AiSi e (49) i=l where k has been experimentally determined to be 0.00365 hours-1 Table II indicates the accuracy of the method. It lists calculated saturations versus height for a column filled almost completely with liquid. The low reading at the top is due to the imcomplete filling.

-39TABLE II ACCURACY OF THE METHOD FOR MEASURING SATURATIONS Height (cm.) Saturation Per Cent Deviation From 1.0 3.81 O.9963 - 0.37% 1.43 1.0089 0.89 19.05 0.9963 - 0.37 26.67 1.0017 0.17 34.29 1.0036 0.36 41.91 1.0168 1.68 49.53 1.0382 3.82 57.15 0.9973 - 0.27 64.77 1.0147 1.47 72.39 0.9261 - 7.39 E. Preliminary Experirents and Calibrations The series of experiments described below were carried out to obtain necessary data on physical properties of the columns, the performance of the counting and detection system, and the information necessary to convert radioactivity measurements to saturations. 1. Permeability The permeabilities of the glass bead packs were measured by evacuating the column and then filling with water. A Ruska constant rate pump was then use= to flow water through the packing at several different rates. The pressure drop between taps near the ends of the columns was measured u$ing manometers.

-4oThe temperature of the effluent was measured so that the proper water viscosity could be used and the permeability was then calculated by straight forward application of Darcy's lawo In order to test the homogeneity of the glass bead packs a special column was prepared which had five pressure taps along its length. Typical results on a pack gave the following resilts (listed from one end of the pack to the other in order): 11.805, 12.23, 12.35, 11.53 with an overall measurement of 11.78 darcys. 2. Porosity The measurement of porosity consisted of evacuating the column and filling it volumetrically using the Ruska pump to measure the amount of water injected. Knowing the pore volume, the porosity is calculated as: Vp VP= - - (50) VB where Vp = pore volume V = bulk volume of pack 30 Density of Glass Beads A pyncnometer was weighed empty to establish the tare weight, then weighed when full of distilled water to determine its volume. It was then dried and filled with beads and weighed. Water was then added with stirring to force any trapped air out until the water filled the bottle to the same mark as before, The bottle was then weighed again, From these weighings and the density of water at the temperatures

of the weighings, the weight and volume and thus the density of the beads were found. The value of the density was determined as 2.48 gm/cm3. 4, Detector Plateau The power which drives the scintillation detector comes from an adjustible high voltage supply (0 to 2500 V. ) contained in the scaler. The response of the detector to a given source of radioactivity is not related linearly to the voltage supplied. Thus the relationship of voltage to counting rate must be determined experimentally. (The relationship is unique for a given isotope.) The results of this calibration are shown in Figure 6. As can be seen, the most stable counting region is the interval from 1330 volts to 1630 volts. The detector is always used in this range or plateau, in this case at 1430 volts. 5. View Factors As noted in section VD, there was a need for the values of the view factors which account for the fact that radiation arrives at the detector from several or all of the zones of the column and not just from the one directly in front of the detector. The method used to obtain these data is as follows: a) Construct ten lucite cylinders 2" diameter, 3" high and 1/4" wall thickness, fill them with glass beads and various amounts of radioactivity and seal the ends with thin plastic caps. (When all ten blocks

32000 30000 28000 26000 OPERATING VOLTAGE PLATEAU 24000 W I22000 20000 Z 18000 0 16000 14000 12000 10000 J I I I I I I 900 1000 1100 1200 1300 1400 1500 1600 1700 1800 PHOTOMULTIPLIER VOLTAGE (DC) Figure 6. Detector Plateau.

-43are mounted on top of one another, they represent a complete column with a known saturation distribution.) b) Take a reading on each block to determine the relative strength (or saturation) of each one. c) Take readings on various combinations of blocks. The view factors may then be calculated as: Rn Rn-l (51) An (51) Rn where Rn = reading for n blocks piled up n-l R = reading for n-1 blocks Rn = reading for nth block when placed directly in front of detector. Figure 7 shows the results of the experiment. At least six different determinations were made for each view factor. The data are plotted to show the range of values obtained and the average value for each Aio The large range at position 10 is due to the fact that as more blocks are piled up, the value of Ai is found as the small difference of large numbers and the values are difficult to determine precisely. The values used in subsequent calculations are: A1 = 1l0, A2 = 0.093, A3 = 0.00839, A4.. Alo = 0.0 6. Self Absorption of Aqueous Solutions In order to determine whether a linear relationship exists between activity and saturation the self absorption of an aqueous

.9.8 7\ 1VIEW FACTORS. ~~~I ~VS. POSITION.6 p t ( DISTANCE BETWEEN POSITIONS = 3") t I a:.4 - AVERAGE VALUE RANGE 3 >.3.2 O -I.I I I 2 3 4 5 6 7 8 9 10 POSITION (NUMBERED AWAY FROM DETECTOR) Figure 7. View Factors (Ai).

-45solution of NaI131 was measured. This was done by placing solutions of various thicknesses in front of the detector and taking readings. As can be seen from Figure 8, self absorption is not a problem for solution depths of less than one inch, and does not become appreciable until the thickness reaches 2 inches or more. Since the effective thickness of the aqueous phase can not exceed about one inch in the column, it may be assumed that the relation Rai = CSi is valid. It may be of interest to note that the data can be correlated by the equation: Rm = _ e-1502x (52) R 00 where Rm/RN = ratio of reading to that of an infinitely thick solution. x = depth of solution (inches). 7. Halflife of Il31 An experimental halflife for I131 was determined by placing a source in front of the detector and following its decay for almost one halflife. The data are shown in Figure 9. The value determined for the halflife was 190 hours as compared to a published value of 193 hours. (17) 8. Absorption of Gamma Rays by Lucite Figure 10 shows data determined for the absorption of gamma rays by lucite. The method used was to place various thicknesses of lucite between a source and the detector.

0 0 m In~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*6 0 00 C3 z I-, X oC t

0 oJ 0 Ch N # 0o xxOD 0'C O Q O O 31nNIW ti3d SlNlOOZ

/ 14 CL' / 0 (03dZ 03 a, co mX.. - 0(3~3.nI 8 SNO0 3IfnNIIN t13d S.lNNl03

-499. Coincidence Losses Coincidence losses are caused by the fact that there is a minimum time interval between events which can be registered by a counting system. The time, following an event in a detector, during which the detector gives no pulse or can not give a response large enough to be registered is called the resolving time. The method used to ascertain the magnitude of coincidence losses for the counting system used in these experiments was to place a source which just barely jammed the mechanical counter in front of the detector and then follow the decay of the source. After a few days the source (I131) decayed to the point where the counting rates were low enough that any coincidence loss was negligible. The data were then extrapolated back to the high initial rates. Any coincidence loss would then appear as an increasing deviation of the data from the extrapolation line at high rates. Figure 11 shows the data which all fall on the same line thus indicating no appreciable coincidence loss at any rates up to the point of mechanical jamming of the scaler. F. Experimental Procedure The procedure followed during an experimental run varied somewhat but consisted essentially of three major steps: 1. Establish an initial distribution. 2. Invert column. 3. Take readings during movement of the phases. The main differences between runs on any given liquid phase were the method of reaching the initial distribution, the initial

."500 0 rE) 0 N cm 0 - l)a cU C)~~~~~a 0 u) 0 0 0iZ PC x c) 0 Ill.- 0 lw~~~~~0 )a w S ( C) O~~~b -H C)~~~~~a _o _ _o _o _ o _o _ J(~~~' 0) aD I'~ ~D R') 31rIN II, b13d $1N1'0 0

-51distribution, and the average liquid saturation. After recording the background readings for all the zones, a measured amount of the radioactive phase was allowed to flow into the column which had previously been evacuated by means of an aspirator. In most cases the column was then closed off and inverted and then opened at the top. This insured that all portions of the pack contained at least connate water or the residual liquid saturation. After allowing the liquid phase to settle to the bottom and wet the column, the system was closed, readings were taken to establish the initial saturation distribution and the column was inverted to start the run. In a few cases (runs 7 and 10) the column was filled with almost one pore volume and then allowed to drain into a constant head device. After the drainage had ceased, the amount that came out was measured and readings were taken to determine the drainage capillary pressure curve. This distribution was then used as the initial distribution for these runs and the column was inverted. The inversion of a column consisted of releasing it from its stand, lifting it out, turning it, and replacing it on the stand. The operation consumed about 1/2 minute and the time of turning was noted and used as the starting point of the run. The procedure for traversing the column and taking readings was to start with the detector in the top notch (Figure 12) and, as soon as possible after the inversion of the column, take a reading. The detector was then moved down one notch and another reading was taken. This procedure was repeated until the reading at the bottom notch had been taken. Then, depending on how fast the fluids were

-52Figure 12. Photograph of Equipment.

-53moving, readings were taken again starting at the top notch. Each reading consisted of turning the scaler on and allowing counts to accumulate for a definite period of time (usually 0.6 minutes or 1 minute). The total amount was recorded along with the time and duration of the count. Thus a reading resulted in a certain number of counts per minute (c.p.m.). The time required to move the detector from one notch to the next and align it with the column was about 1/2 minute. Thus a series of ten 0.6 minute counts required about eleven minutes while a series of 1 minute counts required about fifteen minutes to complete. After a run was completed, the column was flushed with distilled water and dried by passing air through it and evacuating with an aspirator.

VI. EXPERIMENTAL RESULTS AND DISCUSSION A great deal of data is tabulated in Appendix E. The results presented here are shown graphically and for the sake of clarity, only enough data are plotted to give a clear picture of the results and support the discussion of them. The capillary pressure curves are included here because of their particular relevance to the discussion. A. Capillary Pressure and Relative Permeability Curves The imbibition and drainage capillary pressure curves for the water-air system in glass beads were measured several times. The method used was to pack a column as described in section VC, perform a drainage or imbibition and then remove the bottom piston and endplates and dig out the packing. Each inch of packing was deposited in a weighed dish and covered. The packing was weighed, dried, and weighed again. Then, knowing the weight of packing, of water, and the porosity of the pack and the density of the beads, the saturation of each section was computedo The capillary pressure was then calculated for each sample as: Pc = Apgy (53) where y is the vertical distance from the constant head inlet or outlet to the midpoint of the segment. The data for water-air are shown in Figure 13. The capillary pressure curves for the other systems studied are shown in Figures 14 and 15 and were obtained by performing an imbibition -54

-55-.17.16.15 X K=12.52 D, =0.394 O K 11.78 D, ~:0.398.14 A Kl 11.20 D, 4:0.397 V REF. (36).13.12.11.10.09 Ir.08 w 0.07 - o DRAINAGE c~ 0.06.05.~~~04 ~~0.03.02.01 USED FOR CALCULATIONS ACTUA 0 0.1 0.2 a3 Q4 0.5 0.6 0.7 0.8 0.9 1.0 WETTING PHASE SATURATION, Sw Figure 13. Capillary Pressure Data, Water-Air System.

-56-.10.09 x.08 K=11.79 D x -=0.394 T= 840 F DRAINAGE.06 - X.0 6 aL.04 c-.03.02-.x 0.1.2.3.4.5.6.7.8.9 1.0 WETTING PHASE SATURATION, Sw Figure 14, Capillary Pressure Data, 4131 Glycerin-Water SolutionAir System.

-57-.08.0 7.0 6 K 11.79 D. -- 0.394.0 5 T-800F _ X.04 x re' w.0 3 x - x,..0C 2_x DRAINAGE xa.0 I 0.0-.01 -.021 1 I I I I I I I 0.1.2 3.4.5.6.7.8.9 1.0 WETTING PHASE SATURATION, SHEPTANE Figure 15. Capillary Pressure Data, n-Eeptane-Air System.

-58or drainage using a radioactive solution of the liquid phase and the counting and calculation procedures described previously. An imbibition capillary pressure curve for the heptane-air system is not shown because the data obtained for it were considered to be unreliable. The reason for this seems to be instability of the radioactive heptane - iodine solution. The imbibition required approximately one week to stablilize as compared to about two days for the drainage. It appears that over a period of a week the solution is unstable in contact with the glass beads. Visual observation showed a definite lightening of the red color of the solution at the bottom of the column where the saturation would be expected to be greatest. This is thought to be due to the transformation of the iodine to some form such as HI and subsequent movement of this compound in the liquid. The capillary pressure curve used in the field calculations is shown in Figure 16. It represents data from a core sample taken from a gas storage reservoir. Since calculations showed that the solution of Equation (19) is relatively insensitive to errors in relative permeability as compared to capillary pressure data, it was decided to use the correlation proposed by Corey as given by Wyllie and Gardner.(54) Thus the relative permeability curves were calculated from the relationships = (s - SC)/(Sgr- (54) knw = (1 - l) (1 - 2 ) (55) kw=S (56)

-591.9 1.8 1.7 1.6 1.5 1.4 1.3 K =322 md 1.2 *=0.158 I i. I _ CONNATE WATER S=10.8% c3 w.9 >-.8 J.7 o.6.5.4.3.2.1 0 -. i I I I 0.1.2.3.4.5.6.7.8.9 1.0 WETTING PHASE SATURATION Sw Figure 16. Drainage Capillary Pressure Curve - Water-Air System, St. Peter Sandstone.

1.0 RELATIVE PERMEABILITY CURVES BASED ON COREY APPROXIMATION CONNATE WATER S 5.3 % 0., kg kw IQ w 0.4 r 0.3 0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 WETTING PHASE SATURATION, Sw Figure 17. Relative Permeabilities, 5.3% Connate Water.

1.0 RELATIVE PERMEABILITY CURVES BASED ON COREY APPROXIMATION CONNATE WATER S= 10.8 % as 0.7>? 2 ~ ~ k k9k 0. 0L w. 4Q3 0. 0.1 0 0.1 0.2 0.3 0.4 Q5 0.6 0.7 0.8 0.9 1.0 WETTING PHASE SATURATION, SW Figure 18. Relative Permeabilities, 10.8% Connate Water.

-62The curves are shown in Figures 17 and 18 for different Sc and S = 1. gr Bo Transient Saturation Distributions The experimental results are presented as values of liquid phase saturation as a function of time and position. Although the data actually consist of radiation measurements at various positions and times, it was felt that the saturation values would be most useful to the reader and that their tabulation would avoid the necessity of performing the complicated data processing necessary to relate readings to saturations. A total of seven experimental runs were made. Runs 3, 4, 5, and 6 are for the water-air system while runs 7 and 8 are for an aqueous solution 43 per cent by weight glycerin. Run 10 is for the system n-heptane-air. The non-consecutive numbering is due to the fact that runs 1, 2, and 9 were abandoned soon after their start because of either errors in filling the column (runs 1 and 2) or obvious adsorption problems (run 9 - no inactive carrier was added and thus the I131 adsorbed onto the packing). Figures 19, 20, 21, 22, 23, 24, and 25 show graphically several of the saturation distributions obtained for each run. Many more of the distributions are presented in tabular form in Appendix E. Table III lists the values of the various system parameters. The shapes of the distributions, which in some cases are somewhat unusual, will be discussed in the following section. The discussion contained here will be confined to the effect of the physical

90 K=11.84 D., 4=0.389 /w 0.855,.L9 0.0186 38.6 HRS. SAV.O. 498 80 4 IHR..J9 MIN. X 0 MIN. 70 I -DRAINAG x 60 0 x>:~~ o W X ~~~~~TRANSITION X 40 30 x fz r IMBIBITION 20 0 -30 I I I.1.2.3.4..6.7.8.9 1.0 WETTING PHASE SATURATION,(Sw). Figure 19. Experimental Run No. 3, Water-Air.

90 K=11.84 D., #-0.389.LW' 0.768, Ag= 0.0186 SAVG. 0.453 80 107.6 HRS. 90 MIN. 30 MIN. 9(~~ XSO~~ E0IMIN. 70-' 6 0 oV D R AI N A G E 60- / 150.I soIx I: TRANSITION ~40- X' 30 ~/20~~~ t ~ ~~IMBIBITION 20-I -? 10. I I I I I, I 0.1.2.3.4..7.86.9 1.0 WETTING PHASE SATURATION,(Sw)-. Figure 20, Experimental Run No. 4, Water-Air.

-6590 K:11.79 D., = 0. 394 FLw =0.818,.Lg0O.0186 so SSAv.x 0. 699 50 MIN. 10 MIN. O MIN. 70 X \ DRAINAGE),40 1^1I VTRANSITION t(3 ~~~~~~ X b~IBI XBI TION 30 20 / i 0.1.2.3.4.8.7.8.9 1.0 WETTING PHASE SATURATION,(Sw)-,. Figure 2I. Experimental Run No. 5, Water-Air.

90 K = 11.79 -_ - = 855 U/Lg =. 0186 80C _ -.394 80 - S =. vg 7 24.7 hrs. (; 50 min. 10 m. 0 min. 70 60 _ DRAINAGE. / 22.~ R N 6 tXX I 40/ I / -9 I O I / /0 X Er 0fr~/ 0.1.2.3.4..6.7.8.9 1.0 Figure 22. Eperimental Run No. 6, Water-Air.

-6790 K = 11.79: O. 394 Lz = 3.3 80 / Lg =.0186 Svg=.645 48 hrs 300min 80 min O min_ 7'0 DRAINAGE ~4 0 I 40 TRANSITION "It" 3 0 IMBIBITION 0 O.1.2.3.4.6..8.9 a._ WETTING PHASE SATURATION,(Sw)-.. Figure 23. Experimental Run No. 7, 43* Glycerin Solution-Air.

-6890 K = 11.79: 0.394 /4c = 3.3 8S0 Lg:.0186 Sovg=.882 240 min m 30mini x min 70 DRAINAGE 60 50 TRANSITION:40 30i x I IMBIBITION I -I I 0.1.2.3.4.5.6.7.8.9 1.0 WETTING PHASE SATURATION,(Sw) — Figure 24. Experimental Run No. 8, 43% Glycerin Solution-Air.

-6990 K 11.79, 7 =0. 394,w. 392 Ag -.0186 80 Savg =. 7 34.8 HR 140 MIN. 70MIN. O MIN. 70 60 DRAINAGE 4o l 4 50 4 0 TRANSITION 30 IMBIBITION 20 I I I I I 1 i. 0.1.2.3.4.5.7.8.9 1t.0 WETTING PHASE SATURATION,(Sw)-.., Figure 25. Experimental Run No. 10, n-Heptane-Air.

-70-.H C-)' H =L bO s.. O O O O H Q> 1 o Lfl~ Lf\ OIN U) r 0 0' O O O O O O O D Cf)'H'H O 0 0 0 0 H H ~ 03. o O CO \O \O \ \ --.RT * co co co co co cC Co H'H u0 o - 0 0 0 0 0 0 H *rl O Cd H > CO L- r cO CO Lr\ C0U H E O 0' L"0 0, 0 0 H F~ OO. tD a) O U) CO hCO 0 Ln C C0 U] z W1 X zt o0 - CO 0 olEl Cl) co \ C oH o o o e cii 00 H0 0o 0- Lrx O -- O A 4 CO ON CO CO CO CO - U) 1-4 C43: 1;3: ~3: ~~3: 0 co 0 co a Hi > OC OCl >!u 4u rC J C 0 t- C

-71properties of the liquid phases on the speed with which the segregation process is completed. The establishment of a criterion for measuring the rate of countercurrent gravity segregation in a closed system is rather difficult since a final equilibrium distribution can not be predicted a priori with any degree of certitude. However, one can inspect the data and find for each run a time at which the movement of the phases is essentially finished. This time value is somewhat uncertain but does indicate the relative speeds of the various runs. It is listed in Table III for each run. Figure 26 shows the variation of time for a run with the average saturation. It should be noted that the initial distributions for the runs were all at least similar, and if they had not been then this correlation would not be as good as indicated. The effect of fluid properties can be seen by comparing the times for runs 5, 6, and 10 which all had an average saturation of 70%. The glycerin solution can be included in the comparison by reading the time for a 70% Savg. run from Figure 26. This indicates a time of 13 hours. The three physical properties which could be varied were the liquid viscosity, the density difference between liquid and gas phase, and the surface tension. The water-air runs (numbers 5 and 6) and the glycerin solution-air system had essentially the same surface tension and density difference but differed markedly in liquid phase viscosity. As expected, the high viscosity of the glycerin solution results in the

-720 0 - ~0~~ 2 0 I- -..J o 0 (n o _z w t~~~~~~~C U~~~ -I Ctl) 0 z z uJ 0~~ z~~ 4 ~~~~~~~~D0 oo iI. 10 0 4-) w O a) 00~~D~W 10 q~~' jf) N BAY~~~~~~~~~~~~~~~k s/I~~~~~~~~~~~~s co r — (D N ~ ~ - ~ *GAIIS/ I~~~~~~~~c

-73time for a run being much greater than the time for the water-air runs. Run No. 10 differs from the other runs with respect to all physical properties. Since it is reasonable to expect that the time for a run increases with increasing values of viscosity and decreases with decreasing values of the density difference, one would expect that run No. 10 would be completed faster than runs 5 and 6 because the viscosity of heptane is lower than that of water and while the density of heptane is lower than that of water the difference is not so great as to slow the movement to the extent shown by the data. Thus it seems that in run 10 the effect of changing the surface tension is the dominant feature along with the effect of the transition zone as described in the next section. The effect of surface tension is not obvious since in the experimental system under consideration here, capillary forces should aid gravity in imbibing the liquid phase into the bottom section of the column while acting against gravity in the top section by retaining the liquid. As mentioned, run 10 was slower than the water-air runs and thus it seems that the net effect of lowering surface tension is to increase the time for a run. In other words, the imbibition process in the bottom of the column must be a controlling factor on the rate of segregation. Some degree of correlation of physical properties is indicated in Figure 27 which is a plot of kw/y(Ap) versus time for a run. C. Comparison of Data and Theory and Discussion of Results The use of existing theory to predict the experimentally observed results has proven to be a difficult task. This section

-?7k-74 — - 0 0 tr)~~~ u, I W - w i cCI P 0, CD U. 0 o > c'z zz 0~~~~~~~~~~~~C p Co I.t CHF~ ~ 0) IL, P4 _oco - o r) ci c CC)~~~~C *r4 o r~ aU. r 0 0 CH Wr 0) Ic 8~~~~~~~~~~~~~~I 0 4w N~~~~~~~~~~~N0 90 Q0 0 0 6ooooo d d 0 d (dV)X

-75is devoted to a discussion of the discrepancies between observed and calculated behavior. The experimental data possess at least two notable features. The first is that there is a zone comprised of approximately the middle third of the column in which the liquid saturations both increase and decrease at various stages of the process. The second feature, which is related to the first, is that the distributions have an unexpected "bulge" in them. Examples of this can be seen to various extents in the distributions of all the runs. The method used to simulate mathematically the experimental runs was to solve numerically the difference equation presented in section IVC. It is obvious that the final distribution calculated in this manner will conform to the capillary pressure curve used as data in the calculations. This is because at steady state we have: R = constant or Pc = Apgy + C (57) which is precisely the equation of the capillary pressure curve. With this in mind, it seems that because of the bulge which remains in some of the final distributions, the major point of difficulty is the problem of what to use as capillary pressure data. Careful consideration of the process and the data shows that a column can be described as a composite of three sections. In the top zone, there is a continual decrease in liquid saturation, in the middle the saturations increase and decrease at various times, and

TOP DRAINAGE TRANSITION IMBIBITION BOTTOM Figure 28. Composite Column.

-77in the bottom section there is a continual increase in saturation. Thus a column may be represented by three zones as in Figure 28. On Figures 19, 20, 21, 22, 23, 24, and 25 which show some of the saturation profiles obtained in the experimental runs, the three zones have been marked off and labelled. The extent of the transition zone was found by looking through the data and finding the points at which the saturation increased at some times and decreased at other times. As an approximation, it appears that for the runs with S of 70% or less, the transition zone extended from avg. o.38 (Y/L) to 0.65 (Y/L) while for the run having an average saturation of 88% (run no. 8) the zone included most of the column and extended from 0.15 (Y/L) to 0.85 (Y/L). Since, as noted in the literature survey (section III), good results have been obtained using static capillary pressure curves to calculate transient behavior, it is reasonable that one would attempt to do the same in this case. Three obvious alternatives are to use the experimentally determined final distribution as the Pc curve or, if either drainage or imbibition is the controlling process, to use the drainage or imbibition capillary pressure curve respectively, or, finally, to attempt to use a drainage curve for the upper half of the column and an imbibition curve for the bottom. The choice that Templeton(47) made in attempting to do calculations for a similar experiment was to use the final distribution. However, he found, as we did, that this is inadequate for predicting transient behavior although it naturally causes the correct final

-78distribution to be computed. Also, if the final distribution is considered as a necessary part of the data for the calculations, then they will have limited usefulness for predicting anything other than the run which produced that distribution. It should be. mentioned that while Templeton ascribed the inadequacy of this method to either the inability of static capillary pressure data to predict transient behavior or to errors in relative permeabilities, we found, as did Blair(2) and Nielsen(36) that the calculations are somewhat insensitive to changes in relative permeabilities while they are extremely dependent on the capillary pressure curves used. This is due, in great part, to the fact that Equation (19) contains the reciprocal slope of the Pc curve as a coefficient (6S/6Pc). The second choice is to use either the imbibition or the drainage curve. As we pointed out in the preceding section, imbibition seems to be the dominant process in these runs. Thus, by using the imbibition Pc curve, it may be possible to at least approximate the experimental data. Figures 29, 30, 31, and 32 show the results of this method for runs 3, 4, 5, and 6 which were all water-air systems. As can be seen, the results of these calculations do come close to the data but are still a rough approximation. An undesirable feature in the use of the imbibition curve is that it must be arbitrarily extended to S = 1 in order to be used in runs where saturations greater than residual gas occur. The curve for water-air in Figure 13 indicates the extension used in the calculations. Also, it proved to be impossible to even approximate the data for the glycerin solution runs or the heptane run. (The numerical solution predicted faster saturation changes than the data indicated).

-79-'0 i HR 1OMIN. 4HR 6 10 OMIN. 80 X X o X I HR 0 4HR 70 X X I0X \ 0 60 II 0 CX. soI A 0.6.7.o.9 1 I0 0.1.2.3.4.6.7.8.9 1.0 WETTING PHASE SATURATION,(SIw). Figure 29. Simulation of Run 3, Using Imbibition Capillary Pressure Data

90 25H I HR 10 MIN 80 00 X 70 o o El 0 60 0!~05ox 40 2.5 HR ~X0 E0l I0 WETTING PHASE SATUR AT I ON,(Sw)- Figure 30. SimujLation of Run 4, Using Tmbibition Capillary Pressure Data.

-8190 X-6 MIN Q_16 MIN 80 _ \x 70 50t 6 x ox 30 20\ x 1o0 V do_ / Ao xtA /iX / O.1.2.3.4.5.6.7.8.9.0 WETTING PHASE SATURATION,(Sw)*. Figure 31. Simulation of Run 5, Using Imbibition Capillary Pressure Data.

-8290 X 6 MIN O 16MIN El I HR 80 IHR 16MIN 6MIN _\ S \ \ XX 70 N) 60 2~4~0 I~\ 100 30 - fII 20 0E O.1.2.3.4.5.6.7.8.9 1.0 WETTING PHASE SATURATION,(SW).. Figure 32. Simulation of Run 6, Using Imbibition Capillary Pressure Data.

-83The use of two capillary pressure curves - drainage for the top and imbibition for the bottom-introduces an artificial discontinuity at the point in the column where the two sections meet and the results calculated in this manner are erratic and not representative of the data. Thus, while it is felt that the data are accurate and that the measured distributions do represent those existing in the system, we must conclude that capillary pressure data as conventionally measured are not adequate for the task of predicting the distributions which occur in a system undergoing successive imbibition and drainage. In connection with the case of the closed column, it is interesting to consider the case in which a very small amount of gas is present at the bottom while above it the pores are completely saturated with liquid. As indicated by the data presented here, as gas rises in the column some of it is trapped and the liquid saturations below the rising gas reach a maximum which is less than one. Thus if the amount of gas is small enough, it will all be trapped before reaching the top of the column and a final distribution as indicated in Figure 37a will be reached. This unusual situation is resolved by the work of Leverett(29) who found that after a long period of time, diffusion of the trapped gas through the liquid will eventually allow it to rise to the top and the liquid sati-rations in the lower part of the column increase to one. D. The Need for More Complete Capillary Pressure Data From the discussion above it is obvious that a more complete understanding of capillary forces and their effects is necessary. A typical set of Pc curves is shown in Figure 33a. The upper curve is

Lr) C) Lr DRAINAGE 0I IMBIBITION 0 i WETTING PHASE SATURATION Figure 33a. Typical Capillary Pressure Curves. L,J C) u) Li a: DRAINAGE CI c< | IMBIBITION 0 1 WETTI NG PHASE SATURATION Figure 33b. EIypothetical Path Through the Transition Region.

-85produced by desaturating as porous medium originally at 100o saturation and measuring the pressure required to produced the displacement. Similarly, the lower curve is produced by saturating the medium starting at the residual non-wetting phase saturation. As the arrows indicate, each curve is valid only for saturation changes in one direction. The question raised by the results presented here is what happens when saturation cycling occurs. If a saturation value at some point in the column has been reached by a drainage process, it can be represented by a point in the drainage curve. If the saturation is then increased, an imbibition is occurring and the point can not move along the drainage curve, Likewise, it can not move along the usual imbibition curve since the only way to get from one curve to another without a jump is to follow the arrow up the drainage curve and then back down the the imbibition curve. Therefore, it will move along some intermediate path. Figure 33b shows a hypothetical case of a point starting on the drainage curve and successively imbibing and draining. The principle to be inferred from the above is that the behavior of an element of porous material which is undergoing successive drainage and imbibition is dependent on the past history of the element. (In view of this dependence, a more complete history of the initial saturation distribution is given in Appendix F.) It is also implied that in order to predict transient behavior of a system which contains a transition zone such as the one under consideration in this dissertation, a great deal of experimental work must be done to determine the paths lying in the region between the drainage and imbibition curves.

-86A series of capillary pressure cells could be made and all completely saturated with liquid. They then would all be desaturated to various extents so that each represented a different point on the drainage curve. Then an imbibition would be performed on each cell, measuring the pressure after each saturation change. This would produce a family of imbibition curves. A similar procedure could then be followed starting at points on the various imbibition curves and proceding by drainage to determine a family of drainage curves. We suggest that the result would be a network of curves such as that shown in Figure 34. In order to use these data in calctlations it would be necessary to formulate some parametrical representation of the network. One possibility would be to define each imbibition curve by the point at which it leaves the limiting (conventional) drainage curve. Similarly each drainage curve could be characterized by the point at which it crosses some chosen imbibition curve, Thus each point in the network would be defined as the intersection of a drainage and an imbibition curve. For the actual calculations there are several possibilities as to what could be done. Using the equation presented in this work (section IVB), the saturation changes could be noted after each step and then used to decide whether an imbibition or drainage is occurring at each point. Then for the succeeding time step, the appropriate curve would be used in the calculations. Another possibility would be to check the saturation changes after each iteration within the time step. However this procedure would probably be the less stable

DRAI NAGE IMBIBITION C,) w 0 \\ \ 0 WETTING PHASE SATURATION Figure 34. Hypothetical Network of Capillary Pressure Curves.

-88of the two. A third choice would be to use the equation presented here in conjunction with a similar equation for the flow potential distribution [for examples see Douglas et al(13) or Nielsen(36)]. The advantage here is that the potential distribution at the beginning of a time step could be used as an accurate method of predicting the direction of saturation change during the time step. The disadvantage is that in this method it is necessary to solve two simultaneous partial differential equations. The question of how much of the past history of a system under consideration it is necessary to know is a difficult one. In many cases it is no problem. As an example, in a reservoir containing oil, gas, and water, it is assumed that the initial distribution conforms to the imbibition capillary pressure curve(29, 34), In other cases where the previous history has been one of continual increase or decrease in saturation the assumption can be made that the principal or limiting Pc curves can be used to obtain capillary pressure values at the time of interest. Probably the only case presenting serious difficulties is the one in which a number of saturation cycles have occurred in the past. Here, there would be guess work involved and the more known about the cycles the more accurate the calculations. On the basis of the mechanism postulated here it is possible to give a reasonable explanation for the existence of a distribution such as that observed as the final distribution of run no. 7. If the capillary pressure curves for the glycerin solution-air system, as shown in Figure 14 are converted to a plot of height versus saturation by use of the relation

-89Pc = Apgy + C they can then be superimposed on a plot of the final distribution of run no. 7. This is done by matching one of the data points in the top of the column to the Pc value corresponding to the same saturation on the drainage curve. The result is shown in Figure 35. The solid line is made up of the drainage curve in the upper part and the imbibition curve in the lower, This represents a situation which is at least semi stable in that between any two points, Yl, Y2, the capillary pressure Pc is balanced by the gravitational force Apg(y2 - Y1). The fact that the Pc curves fall so close to the data indicates that the process was one of drainage in the top part and imbibition at the bottom, It should also be noted that there is still one point in the transition region between the two Pc curves. The question as to whether this type of distribution will remain stable over an extended period of time is open to conjecture. Runs no, 10 (heptane-air) and no. 6 (water-air) which had average saturations (70%) near that of run no. 7 (64.5%) exhibited the same type of bulge for a while and then it disappeared. This may mean that over a period of time, all the points will finally pass through the transition zone and come to rest on the imbibition curve. Thus the effect of the transition zone may be principally that of greatly increasing the time required for the segregation process to be completed. This is borne out by the fact that the numerical solutions using the imbibition curve for runs 7 and 10 predicted a much shorter duration for each run than was actually observed. The calculations indicate that the speed with which saturations change depends on the slope of the Pc curve. In other words,

90 80 X - DATA (RUN NO.7) 70 60 x DRAINAGE E 0 L300 IMBIBITION i O I 1. 2.. -....6.I.8.9 1.0 WETTING PHASE SATURATION,(Sw) — Figure 35, Superposition of Capillary Pressure Curves on Final Distribution of Run 7.

-91the amount by which a given saturation is changed during a time interval is a function of the slope of the capillary pressure curve at that point. If the slope is increased (this corresponds to a decrease in absolute value of S' = as ) then the rate of saturation change is decreased. Thus, if, as indicated in Figure 33b, the saturation of a point in the transition zone follows a path leading from the drainage curve to the imbibition curve or vice-versa, the slope of the path will be greater than that of either the drainage of imbibition curves at the same saturation, and the change in saturation with time will be smaller. This argument would explain the presence of the bulges in the distributions since the saturations in the transition zone change much more slowly than those in the drainage or imbibition regions. It should be mentioned that Templeton's data(47) do not exhibit, at least to any great extent, the bulges reported here. This can be explained by the fact that he started his runs with an almost uniform saturation distribution, e.g., 50% saturation everywhere in the column, and thus does not seem to have produced transition zones anywhere near as large as those found in this work. A further observation based on the data is that the chance of obtaining a "final"distribution exhibiting a large bulge, is greater for higher average saturations and also is greater the more completely segregated the phases are initially. E. The Case of a Column Open at the Bottom At this point it is isnstructive to consider the case of a column closed at the top but open at the bottom. In this case, gas is

-92supplied to the column while water drains out and is carried away. Also, the average saturation continually decreases. This situation is indicated in Figure 36. This is more like the situation in a gas storage reservoir than is the completely closed system. The rate at which gas is supplied to the column would determine whether or not a transition zone could occur. If gas were supplied at a high enough rate (or pressure), the water would not accumulate within the column and thus the process would be one of drainage alone and the saturation distributions during the displacement would be expected to follow the sequence shown in Figure 37. In the event that a distribution with a bulge in it as in Figure 22 happened to occur during the process, it would not be stable since as the water leaves the bottom of the column the saturations below and in the bulge would decrease and it would eventually disappear. This sequence is indicated in Figure 38. In this case, the final distribution would show no bulge since at low saturations the imbibition and drainage curves merge into one.

TOP GAS /. \ WATE R Figure 36. Countercurrent Gravity Segregation in a Column Open at One End.

-94. IU Xa b c d LLI 0 0 0 1 0 S Sw Sw Sw Sw TIME - I! e f g 0 1 0 1 0 S SW Sw Sw TIME -- Figure 37. Hypothetical Saturation Distributions in a Column Open at One End.

-95y y y o S 0 S 1 0 s TIME -- 0 S 1 0 S 1 TIME Figure 38. Disappearance of BulgE: and Transition Zone.

VII. RESULTS OF FIELD CALCULATIONS The calculations using typical reservoir parameters were designed to simulate a countercurrent gravity drainage in a vertical column open at the bottom. The situation is sketched in Figure 36. Although it is realized that the actual movement of gas injected into a reservoir is a three dimensional problem and is much more complex than the case treated here, we feel that some useful qualitative information can be gained by treating a simplified version of the problem. Thus the flow of water and gas is assumed to be one-dimensional - occurring only in the vertical direction. The assumption is also made that at the bottom of the column the water saturation has been decreased to some value (in this case 1.5 Sc) and that it remains at this value during the remainder of the process. Thus, as noted in section IVC, the boundary condition at the bottom is: S2 = 1.5 Sc or Pc = constant (58) or R2 = X and with this taken into account, the calculations are carried out as described in Appendix B. The principal variation in the calculations was the permeability distribution in the column. Four cases were run, three of them the same except for permeabilities, and the fourth differing from the others in

-97=0.158 36 K =200 md 21 K = 200 md K| 50 md K =200md y,(ft) K=5Omd K = 10md K= 50 md K=500md K=200md K=500md CASE I CASE II CASE m CASE'I Figure 39. Field Cases Treated.

-98that it was a different height. The four cases are shown in Figure 39. The initial saturation distribution in all cases is: 0 < y < 4 ft, S = 0.162 y > 4 ft, S = 1.0 The capillary pressure data for St. Peter Sandstone(25) in Figure 15 were scaled to the various permeabilities by use of the J-function (Equation (6)), which can be rewritten as: Pc2 Pcl ( K )(59) Pc2 = K when y, G, and 0 are not variables. The results of the calculations are shown in Figure 40 as a plot of the average saturation in the column versus time. For the 30 ft columns the initial average saturation was 0.919 while the 21 ft column has an initial average saturation of 0.86. The most prominent feature shown by this figure is the great effect of the permeability variations on the length of time before any appreciable decrease in average saturation occurs. This indicates that even a thin layer of low (but not zero) permeability media can greatly affect the displacement process. Thus the occurrence of a layer of this type in a reservoir may account for the length of time necessary for a "bubble" to be formed at the top of the structure. Some of the calculated distributions are shown in Figures 41, 42, and 43 for cases I, II, and IV. The distributions for case III are similar to case I but occur at later times. The calculations were not carried out to completion of the drainage because of the amount of computer time involved.

-99-. O N E E E U 2 0 00 0 0 Co w E E E Q w 8 o o 00 Wj Iii 800 - 0~~~~~~ U t +r~~ I I 4 I g E E E IL I L 0 "r bw 0 0 -- o~~~~~~~~~~~~~~~~ I4~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~I z 10, la vv 30 00 0 1 r~~. v'V o~~~~~~~~_ o ~ Co (O 0 0, 1-H 4-~ ).. >->- >-'U 10'I)~~~0%3 O 1(o~ ~~~~~r 0 II I I N _ I I 9 0 0 t. 0 0 t IE) N -: oO~~~~~~~~~~~~~ NOILYIAUNLVS 30Vt:13AY

-10038 36 34 32 30 28 - FIELD CASE I =0. 158 26 24 22 20 18 IL 1412- K 200 md --- LAST CALCULATED TIME = 27.34 Days I6_ INITIAL K 50 md / K =500md 0.1.2.3.4.5.6.7.8.9 1. WETTING PHASE SATURATION, Sw Figure 41. Saturation Distributions - Field Case I.

-10138 10 Weeks 36 4 weeks 34 32 FIELD CASE It UNIFORM K=200 md 30 = 0.158 28 26 24 22 - 20 s L I I I Week I6 14 128 6 2 I I I I I WETTING PHASE SATURATION, S. ----—:-e Figure 42. Saturation Distributions - Field Case II,

-10222 20 18 FIELD CASE =0.158 16 14 -L 12 _LAST CALCULATED 12;~ TIME = 23.1 Days -'- K 200md -10 z 8 INITIAL K 50 md 4 =500 md 0 0.1.2.3.4.5.6.7.8.9 1.0 WETTING PHASE SATURATION, SW - _ Figure 43. Saturation Distributions - Field Case IV.

VIII. SUMMARY AND CONCLUSIONS Summary 1. A quantitative me:thod has been developed for the use of radioactive tracers in the measurements of saturations. 2. Saturation data have been experimentally obtained for the countercurrent flow of two phases in a vertical, closed column of porous material. 3. A mechanism has been postulated to account for the experimental results and some success has resulted from the use of a single capillary pressure curve as data for a numerical solution of the equation governing the process. 4. Experiments have been proposed which would lead to a more complete understanding of the role of capillary forces in processes undergoing saturation cycling. 5. Calculations have been performed using reservoir properties for the case of a column of porous media closed at the top but continually supplied with gas at the bottom. The effect of permeability variations has been studied to find their effect on the time required for the gas to displace the water initially in place. Conclusions 1o The use of a radioactive liquid phase tracer is an accurate method for measuring saturations in a one dimensional system. The method presented here for taking into consideration the radioactivity in parts of the system other than the one under consideration results in saturation measurements at least as quantitative as any other method. -103

2. The use of the ordinary drainage and imbibition capillary pressure curves to calculate transient behavior in a system undergoing both drainage and imbibition leads to erroneous results. Capillary effects in such a system need further experimental work before a complete solution of the problem is possible. Also, because of the great influence of the capillary pressure data upon the solutions obtained to the equations for flow in porous media, it is necessary to have very good data for any calculations involving capillary forces, 3. A thin, low permeability layer can have a great effect on the speed of countercurrent drainage. This may indicate one of the principle reasons for the slow bubble growth observed in some gas storage fields.

BIBLIOGRAPHY 1. Archibald, R. S., "Radioactive Tracers in Flow Tests," Journal Boston Soc. Civil Eng., 37_, (1950), 49. 2. Blair, P. M., "Calculation of Oil Displacement by Countercurrent Water Imbibition," Presented at the Fourth Biennial North Texas Section AIME, Secondary Recovery Symposium, Witchita Falls, Texas, May 2/3, 1960. 3. Blair, P. M., Douglas, J. M., Jr., and Wagner, R. J., "Calculation of Linear Water Flood Behavior Including the Effects of Capillary Pressure," Trans. AIME, 213, (1958), 96. 4. Boyer, R. L., Morgan, F., and Muskat, M., "A New Method for Measurement of Oil Saturation in Cores, " Trans. AIME, 170, (1947), 15. 5. Brunner, E., and Mardock, E. S., "A Neutron Method for Measuring Saturations in Laboratory Flow Experiments,?" Trans. AIME, 165, (1946), 133. 6. Buckley, S. E., and Leverett, M. C., "Mechanism of Fluid Displacement in Sands", Trans. AIME, 146, (1942), 107. 7. Cardwell, W. T., Jr., "The Meaning of the Triple Value in Non Capillary Buckley-Leverett Theory", Trans. AIME, 216, (1959), 271. 8. Cardwell, W. T., Jr., and Parsons, R. L., "Gravity Drainage Theory", Trans. AIME, 179, (1949), 199. 9. Clayton, C. G., "Precise Tracer Measurements of Liquid and Gas Flows," Nucleonics, 18, (July, 1960), 96. 10. Coats, K. H., "Selected Topics in Numerical Analysis," notes prepared for a course in the Chem. Eng. Dept., University of Michigan. 11. Coomber, S. E., and Tiratsoo, E. N., "The Application of Radioactive Tracer Technique to the Study of the Movement of Oil in Sands,?" Inst. of Petr., 36, No. 321, (1950), 543. 12. Dombrowski, H. S., "Residual Equilibrium Saturation of Porous Media," Ph.D. Thesis, University of Michigan (1952). 13. Douglas, Jim, Jr., Peaceman, D. W., Rachford, H.H., Jr., "A Method for Calculating Multi-Dimensional Immiscible Displacement," Trans. AIME, 216, (1959), 297. 14. Fayers, F. J., and Sheldon J. W., "The Effect of Capillary Pressure and Gravity on Two-Phase Fluid Flow in a Porous Medium, Trans. AIME, 216, 147 (1959). -105

15. Geffen, T. M., and Gladfelter, R. E., "A Note on the X-Ray Absorption Method of Determining Fluid Saturation in Cores," Trans. AIME, 195, (1952), 322. 16. Gilmour, H., "Gas and Oil Segregation in Oilfield Reservoirs," J. Inst. Pet., England, 46, No. 444 (Dec. 1960). 17. Hodgman, C. D., ed., Handbook of Chemistry and Physics, 40th Edition, Chemical Rubber Publishing Co., Cleveland, Ohio, 1958-59. 18. Holmgren, CO R., "Some Results of Gas and Water Drives on a Long Core," Trans. AIME, 179, (1949), 103. 19, Hovanessian, S. A., and Fayers, F. JO, "Linear Water Flood with Gravity and Capillary Effects," SPE Journal, 1, (1961), 32. 20. Hull, D. Eo, "Total-Count Technique in the Refinery," Ind. Eng. Chem., 50, (Febo 1958), 149. 21. Hwa, C. S., and Beckmann, R. B., "Radiological Study of Liquid Holdup and Flow Distribution in Packed Gas-Absorption Columns," J. AIChE, 6, (1960), 359. 22. Isbin, H.S., Sher, N. C., and Eddy, K. C,, "Void Fractions in TwoPhase Steam-Water Flow," Jouro AIChE, 3, (1957), 136o 23. Jones-Parra, J., and Calhoun, J. C,, Jr., "Computation of a Linear Flood by the Stabilized Zone Method," Trans. AIME, 198, (1953), 335. 24. Josendahl, V. A., Sandiford, B. Bo, and Wilson, J. W., "Improved Multiphase Flow Studies Employing Radioactive Tracers, " Trans. AIME, 195, (1952), 65. 25. Katz, D. L., (Private Communication). 26. Katz, D. L., et alo, Handbook of Natural Gas Engineering, McGraw-Hill Co. Inc., New York, 1959. 27. Kniebes, D. V., Burket, P. V., and Staats, W. R., "Argon-41 Measures Natural-Gas Flow," Nucleonics, 18, (July, 1960), 142. 28. Laird, A. D. K., and Putnam, J. A., "Fluid Saturation in Porous Media by X-Ray Technique," Trans. AIME, 192, (1951), 275. 29. Leverett, M. C., "Capillary Behavior in Porous Solids," Trans. AIME, 142, (1941), 125, 30, Lipson, L. B., "Theoretical Note on Linear Absorption Methods or Determination of Fluid Saturation in Porous Media," Trans. AIME, 192, (1951), 375.

-10731. Martin, J. CO, "Some Mathematical Aspects of Two-Phase Flow with Applications to Flooding and Gravity Segregation Problems," Producers Monthly, 22, No. 6, (April, 1958), 22. 52. McEwen, C. R,, "A Numerical Solution of the Linear Displacement Equation with Capillary Pressure," Trans. AIME, 216, (1959), 412. 33. Morgan, F., McDowell, J. M., and Doty, E. C., "Improvements in the X-Ray Saturation Technique of Studying Fluid Flow," Trans. AIME, 189, (1950). 183. 54. Muskat, M., "Calculation of Initial Fluid Distributions in Oil Reservoirs," Trans. AIME, 179, (1949), 119. 35. Naar, J., and Henderson, J, H., "An Imbibition Model; Its Application to Flow Behavior and the Prediction of Oil Recovery," SPE Journal, 1, No, 2, (June, 1961), 61. 36. Nielsen, R. L., "On the Flow of Two Immiscible Incompressible Fluids in Porous Media," Ph.D. Thesis, University of Michigan, (1962). 37. Osoba, J. So, Richardson, J. G., Kerver, J. K., Hafford, J. A., and Blair, P. M., "Laboratory Measurements of Relative Permeability," Trans. AIME, 192, (1951), 47~ 38. Peaceman, D. Wo, and Rachford, H. H., Jr., "Numerical Calculation of Multidimensional Miscible Displacement," Presented at 35th Annual Fall Meeting of the SPE of the AIME, Denver, Colorado, Oct. 2/5, 1960. 39, Radioisotopes Catalog, Oak Ridge National Laboratory, Oak Ridge, Tennessee. 40, Richtmeyer, R. D., Difference Methods for Initial-Value Problems, Interscience Publishers, Inc., New York, 1957. 41. Russell, R. G., Morgan, F., and Muskat, M., "Some Experiments on the Mobility of Interstitial Waters," Trans. AIME, 170, (1947), 51. 42. Rzepczynski, W. M., Katz, D. L., Tek, M. R., Coats, K. H., "The Mt. Simon Gas Storage Reservoir in the Herscher Field," Paper Presented at Research Conference on Underground Storage of Natural Gas, July 16/17, 1959, University of Michigan. 43. Sandberg, C. R., Gournay, L. S., and Sipper, R. F., "The Effect of Fluid-Flow Rate and Viscosity on Laboratory Determinations of Oil," Trans. AIME, 213, (1958), 36. 44. Scheidegger, Ao E., Physics of Flow Through Porous Media, MacMillan Company, New York, 1957.

-10o845~ "Selected Values of Physical and Thermodynamic Properties of Hydrocarbons and Related Compound.," API Res. Project 44, Carnegie Press, Pittsburgh, 1953. 46. Stanley, D. R., "Sand Filtration Studied with Radiotracers," Am. Soc. C. E. Proc., 81, No. 592, (Jan. 1955). 47. Templeton, E. E., Nielsen, R. F., and Stahl, C. D., "A Study of Gravity Counterflow Segregation," SPE Journal, 2, No. 2, (June, 1962), 185. 48. Terwilliger, P. L., Wilsey, L. E., Hall, H. N., Bridges, P. M., and Morse, R. A., "An Experimental and Theoretical Investigation of Gravity Drainage Performance," Trans. AIME, 192, (1951), 285. 49. Watkins, J. W., and Mardock, E. S., "Use of Radioactive Iodine as a Tracer in Water-Flooding Operations," Trans. AIME, 201, (1954), 209. 50. Welge, H. J., "Simplified Method for Computing Oil Recovery by Gas or Water Drive," Trans. AIME, 195, (1952), 91. 51. West, W. J., Garvin, W. W., and Sheldon, J. W., "Solution of the Equations of Unsteady State Two-Phase Flow in Oil Reservoirs," Trans. AIME, 201, (1954), 217. 52. Whalen, J. W., "A Magnetic Susceptibility Method for the Determination of Liquid Saturation in Porous Materials," Trans. AIME, 201, (1954), 203. 53. Wykoff, R. D., and Botset, H. G., "The Flow of Gas-Liquid Mixtures Through Unconsolidated Sands," Physics, 7, (1936), 325. 54. Wyllie, M. R. J., and Gardner, G. H. F., "The Generalized KozenyCarman Equatiorn," World Oil, (March and April, 1958).

APPENDIX A STABILITY ANALYSIS* OF EQUATION (22) For the purpose of stability analysis, the equation must be considered to be linear since no adequate theory exists for the stability of non-linear equations. The procedure commonly used is to assume all coefficients to be constant and then when a criterion for stability is established the coefficients are assigned their most adverse values. This will usually insure that the stability criterion is valid. Thus with M and S' constant, the difference equation can be written as: M (R 2R n+l + R n+l) (Rjn+l R (x)2 j-l, - 3, j-l,,,n) (A-l) or upon rearrangement 1 1 R+ l (2 +-R + R= R? (A-n(R Rj+l,n+ (2 + a) Rj,n+l + Rj_l,n+l a j,n (A-2) where MAt a = - Mt OS'(Ax)2 Let Rn represent the exact solution of (A-2) and + E. j,n 0Jn denote the actual or numerical solution. The error Ej here is due to rjn round-off error and stability denotes the case where: This analysis is based mainly on notes by K. H. Coats.(lO) -10o9

-110Ej,+l <1 Ej,n that is, the error should decrease from one time step to the next. Since both solutions Rj and R* + Ej must satisfy Equation (A-2), the J,n j,n J,n equation can be written twice (once for each solution) and subtraction then yields Ej+l,n+l - (2 + 1) Ej,n+l + Ejn+l =(A-3) a ~ j-l+ a j,n The substitution E. = f jg gives gn+l fj+l 2 - 1+ -l = - (A-4) gn fj a fj a or: 1 g fj+ 2 1 jl (A-5) a gn+l f; a f and 1 n gn = ( aA ) = aX The stability criterion is thus: laX < 1 Kl(A-6) The solution to the equation: fj+l - (2 + X)fj + fl = (A-7)

-111which results from (A-5), is obtained by letting 2+ - X = 2 cos y (A-8) a so that solutions involving sine and cosine will result. Therefore: X = 2 + - 2cos y (A-9) and (A-6) becomes 1 + 2a (1 - cos <1 (A-10) since y will take on only real values, cos y will be in the range - 1 _ cos y < 1. Equation (A-10) is satisfied for any positive a and thus there is no limit on the size of a as long as it is positive. Since S' is always negative it follows that -MAt a = Os' (z~)~ is always positive. The special case of y = 0 or cos y = 1 should be considered. One of the eigenvalues ym determined by the boundary conditions may be zero. One of the eigenfunctionsin the Fourier expansion of Ej,n E.n = n (a m)n f j,n m j,m will then have a time dependence of (1)n and an error introduced at any point neither grows nor decays with increasing time. This case is not

-112serious if, as is usually the case, the round-off errors are random in sign. Thus, from the above analysis it is to be expected that the difference form employed in this work will be stable and convergent especially for judicious choices of Ayx and At. However, in practice it was found that the rate convergence does depend on the choice of time step size, i.e., At must be kept small in order to converge in a reasonable number of iterations. The maximum At's which can be used can be found only by trial. It is assumed that the non-linearity of the equations and computer round off cause this.

APPENDIX B DESCRIPTION OF NUMERICAL SOLUTION A complete program was written in order to perform the calculations necessary to solve the algorithms derived in section IVC. The program was used on I.B.M. 704, 709 and 7090 computers at the University of Michigan Computing Center. As can be seen from the equation and difference form, the necessary input to the program consists of tables of capillary pressure versus saturation, relative permeabilities versus saturation, porosity and permeability as functions of position, the initial distribution, and physical properties such as viscosities and densities. Perhaps the most important part of the numerical solution is the iterative procedure used to obtain convergence during a time step. After experimentation with various ways of time centering the coefficients, the scheme described here was evolved. The two coefficients involved in the iteration are the term S' and the interblock coefficient M. For S', the procedure of Douglas Peaceman, and Rachford(l1) was used in which the calculated value of R is averaged with the value at the beginning of the time step and this average value is used to evaluate the coefficient. The coefficient M was calculated using the saturation values at the beginning of the time step and then again using the saturation resulting from one iteration. These two values were then averaged and this average value was then used for all successive iterations. The iterative routine was considered to have converged when the maximum per cent difference between values of R at any point on successive iterations was less than some tolerance-between 0.1% and 0.001%.

APPENDIX C DATA PROCESSING METHODS Table IV is a sample data sheet showing the various facts recorded during two scans of the column. The processing of data from a run follows the format outlined below. (All calculations were performed by IBM 709 and 7090 computers,) 1. All readings are punched on cards in sets of 10 with an average time assigned to each set. 2. The set comprising the readings on the initial distribution is used to calculate the constant C of Equations (44) and (45) and also the initial distributions by means of the method presented in section VD, 3. Successive sets of readings are converted to net counts per minute and used to calculate the saturation distributions at each average time value. 4. Using the output from the first data processing program described above, the uncorrected readings expressed in counts per minute are plotted versus time of reading for each of the ten points covered by a scan. After connecting the data points with as little smoothing as possible, the readings are interpolated at constant time values. 5. The values obtained by interpolation are used along with the constant C as input to the second data processing program which converts the interpolated values to saturation distributions. The -114

-115TABLE IV SAMPLE DATA SHEET Experiment: Run No. 7 Date: 8/3/62 Scaler No. FPL - 2 Voltage: 1500 Base Time = 13:06 Temperature = 850F Position Time Time Total cpm Background Net cpm of Count Count (cpm) (min) 10 19:20:5 O.6 6168 10280 5144 5136 9 21.6 o.6 6326 10543 2190 8353 8 22.9 o.6 9523 15872 1495 14377 7 24.1 o.6 19168 31947 1185 30762 6 25.1 o.6 25612 42687 805 41882 5 26.4 o.6 21328 35547 703 34844 4 27.6 o.6 19552 32587 1014 31573 3 28.75 0.6 22336 37227 773 36454 29.8 o.6 22543 37572 918 36654 1 30.9 0o6 19314 32190 2425 29765

interpolation procedure is only used to cover periods of time when rapid changes in saturations are taking place. For the periods when the phases are moving very slowly, the distributions computed by the first data processing program are considered to be valid since not much change occurs between the first and last readings of a set.

APPENDIX D PHYSICAL PROPERTIES OF MATERIALS Data for the physical properties of the materials used in this research were obtained either experimentally or from reliable sources. Each material is listed below along with the source of the data, and the figures on the following pages present some of the data graphically. 1. Water Distilled water from the engineering building supply was used for runs involving water and also for flushing the columns. Data for water and also for air were obtained from the Handbook of Chemistry and Physics.(17) 2. Glycerin Commercially pure glycerin was used to prepare the glycerinwater solution. The specific gravity of the solution was measured with a hydrometer as 1.105 gm/ml and the viscosity was determined with an Ostwald viscometer calibrated with an oil supplied by the National Bureau of Standards. The viscosity of the solution is plotted in Figure 46. The surface tension of the solution (43% glycerin by weight) was measured with a tensiometer and found to be approximately 90% that of water. 3. Iodinel31 The properties of I131 are taken from (39) and are as follows: Half Life: 8.05 days Radiations (MEV): Beta- 0.250 (2.8%), 0.335 (9.3%), o.6o08 (87.2%), 0.815 (0.7%) Gamma- o.080 (2.2%), 0.163 (0.7%), 0.284(5.3~), o.364 (8o ), 0.637 ( 9% ), 0.722(3%) -117

Chemical form: Carrier free NaI in basic sodium sulfite solution pH: 7-11 Concentration: greater than.1 mc/ml Total solids: less than 2 mg/mc Heavy metals: (as Pb) less than 2 micrograms/mc Purity: greater than 99o9% (exclusive of I133) 1133: less than 2% 4. Heptane Physical data for n-heptane are shown in Figures 44, 45 and 46. The data were taken from API Research Project 41o(45)

-11923 E 22 DATA FROM (45) E CT z X u 21 IZ -J 19 I s ILl X 170 10 20 30 40 50 TEMPERATURE, T(C) Figure 44. Interfacial Tension, n-heptane -air.

-120-.700 DATA FROM (45) 0 E 2t.685 QI..680 _ 15 20 25 TEMPERATURE, T (~C) Figure 45. Density of n-heptane.

-121-.54 3.8.52 \ X N-HEPTANE 0 GLYCERINE SOLUTION 3.6.48 N-HEPTANE DATA FROM (45) 3.4.46 3&2.44 X -3.0 "t.42 x 2.8.40 2.6.38 x 2.4.36 2.2.34 X 10 20 30 40 TEMPERATURE, T (~C) Figure 46. Viscosity, n-heptane, 43% Glycerin Solution.

APPENDIX E TABLES OF DATA -122

-123I ( C UIN n I,; It- II D M NI IJ4If I'UliIt -r1LA0 444. N N Ap' In i LN I *1 I N 4 t-~ 4 I 11' -4 4 C I4 n.I. U n 10 3lC jo 10'0 IS 10 t I. I I,.,1' I.-I o 1 r- I ";I 1,I, It N I I M I I W I II~~~~~~~~~~~~~~~~~~~ II I 0 LA U. I L I 0 4 i 14 A IP A 0' I u. -' Ln 1'Irn' N t 4 co CD A4 In' 0 I'W C:-' I- a' lo 11 it 1'A P.- IA lA C4 I1;.4~~~~~~~I I(f I?Ic Q U,\, I.I I I n I c IN gn * g t * v I'''. I U)I I( n co U, Igo I I4 CI' I IJ aD I~0,91 J 10 I I~o I o, Im g I V I I I" I~I Ic Im I I Is S In I I I, en I Ie r — O'O 10 1, C4 I c' N r' A 10 It M 1 D 10 C'1 t L U t U L'1 4 4 0 C3 I I I I I I~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~C 6 m 9 uc I~ )~ I I' I~~~I l0 Iu'~ ~ ~ ~ ~ -)"0 L I ia' f O ~a c~~~o ~ L 10 10a P- P1 It lI 0: 1 I i,..- 1N110 10'17' N'0h 110 ~o a'T CU L O It IN* a' NC C'It 0 1 00000 Is Ile Icn a) Ivl (s. ~~~10 LA I 1 00 LA N -'0 I 0 0 0 C) 0 119 -4 co it co I cI O' I, I L I $I 4 N o 111 L0 a I f I m U% cgr o I',U-% 9 lr CN I I- I~ I ~ I ~ I I 4 I I I I~ I I I,- le-dI ~ II Um P1 ) I't IT 10 P 110 It 1CLACD1 Itc a' N w U 0 ~ r IU yIHN:'-''~ 10 41' 1a 0 I 4 i-AI U-0 N I"4 I: 01' I P-' CO a' 11 L K*4'-4 p. 10.:;0'.' 41 Ia 4.'44 00000 ~ ~ I It\'''' 10:I-:0- t'9` LA U' It 1~ w n -'0 I, * * Ia' I II'. In.f I.0 I U- I 4n' ~ ~~~ ~ ~ ~:1 Il * It l (P' I *r Lo co'I i (j 00% 1,4 io%; it;( I v 1 i I;a' -I N I t LA LA 10 1 1- 1 ZO 401,~~~~~~~~'... ~I t AL 110I Iz i: 3''4 LA LAI lI l0 19 CO L N ('a'1.- - N 1 41:o Ia aco~L LA to, CO U-m I t 1w (4'N I Ln I U, Q r- cv w~~~~~~~~~~~' I I It L LA 10j Ic- I I cu I ~ n Ir 9 rC 9 N I IC IVI 9 9 I co IN I I enI 0 ~ ~ O LA I. L'.o 1 I't I Ln U'I% I II I I 4 v t j II I N 14 10If'- - 1:' -4 inJ L co I'0 v - IN LA IN g CO',~~~~ 110 ~~~~~~~~~ ~~~~ 0:w'4'~~~~~~~~~~~~~~~ 10 w~~~~~% (P a' 1- 110L'4 rIN 1-4 Ir'41 0 I a'I IN N % 1.44'Ia L10 LA N CO IT w lI lt CO' o a''-'4 1 ~~~~~~~~ I''i LC I e 10O' LA It LA LA- 1 10 COIt EH co 1, rn co I' 1' I I I 1 4 ir c L 1I -- t I I ~~~~~~IL Ic ~ L' — - I It 1 IN'I:' Iaco I It LI I~<30 N 4 rs'L 1 0 a'PIa' uW-!''- I IN 1e I 4 ~~~~~~~~~~~~~~~~~~C 0 NC 111 11.4.a It- p'- 1 It LA'' LA 10o1w lI I. q IS I.. I u.,I I I I~~ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~' r'i I.I C' y I 4 u-4K~~~~~.rt —O\ 1 1 1 II~~~~~~~~~~P- I~ I - I I It~~~~~~~~~~~~~~~~~~~~~~~~~I3 0.'0,Nc~I~'" IiI'' ir co~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 110 N P1 LA r ~ ~~' l 1': L'- 10 p. "'lt=I 1 I.-.KII' I~0 AI- N L N I'I' 4 - 1 N U " 4 I t-4 10% IN 1 0 0'N 10'It CZ ~ 11 A I 10w I ~ t 0 a gI~I~~I-1~1 IN. I-4. I I'' *''I I I I Iff'I I4 U'\ ID I t* I~ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ I 0 I~~~ 1001010010 J 1 01001010 ololo..ooco,o'' I~ ~ ~ ~~~~ I~ I I 1 iI 0, 1.0 10 0 0 10 0 1 0 0 10 0' I' I, I I I~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~, I-','

.124I~I'O co u.:i.. "co I- N in. p -o in' " " oI~..11 -' 1 iN It t tin. NO i or. * -'.4 1' - I' r 1co 1I.on 1 0! ~ ~ ~~~I I -t %O N I:ji co N *' 0 % 0Y lLf in 100 M N LA i Sf- I 4N N D I e I$ 001D - I4' O c4 c' I -vf CO0 01- CD.- n In Ol in 0 P.M Oli N P M 0 No LN Ny 0 01 cm N0 toOAf- i co110',C0 in i" CY - 00 "4 0o I to ot't i'. It'0 " co N P P. i 4 4 M ".co t *i;t lI 66c~~~~~~~~'' I' in iN in Q'M' Uf.0 in M' MI N P. co CD'I" Lf'% %O Io 0.'i M oO i I 01 4'lIC " C o -,Y %tto), N C'. N O M; P. io -0 to tin 011 P. N'.4 P.- I I I~~~~~~~~~~Y1 tCYc 0 in "4 lIt'It:00 100 In IM 0. 0 in It *. m ".4 rn C cq 0\u-f\\0 i. P. P. co "g 0 P. tot N o U 0 4 i 4 I t M N'N 0 c o t02 g0 I N P- % 0010'"4 It lotiM lIt "4101 N 01 i-4 "4 ~~ ~~~~~~~~~~ ~~~~~~~0 - It CO. "4 P. " P 0 g t~~~~~~~~~~k o. VOt. tot I 00 It.-. 0' in N CY It P. iM tn 0 1't1 1 I''1 M 0 t I 0 M N' inOM 00 N co0010 co loN N 01 " A I I I I I ~~I, I g ~~.4 g I I I I~~~~~~-. g CY l NI IM NM o I I'Y CY SOt,00 C" 1-t -t "O 5." It I. *- ~~~~~~~~~~~~~~ It ~~~~~~'in lo% t If"* CY E-1 le~~~~~~~~~~~I1H 1 I I IN co I, I. I4 eP A' St co O O NI'M I I c o r4 ~ ~ I I I I I I tn CD I LA Co~~~~I lot l'~~~~~~~~~~ N t 00IU 100 t UN ICOti ItN tin - j l01'P50.00 10 4.M-~No M M 00 1'.'~~~~~ iM 0 i I 1 1 I %O P. P. N in SM'I t ja. N 100 H 0I IN'00 01 P.l1010 in 101 110 00''4 N, 1-41-4 ( "4 M it" 4. It P. P 00~ t P.jo It'tn im n ct ti 10 l N N 1* %* I' 1* co N1 0 (~o-.0 I. 0 10 1. % m D 0 0 0 0 I..10 1 0 1 A I'O CIY la CD I I 1 1 Ic I% co e.4 CDI I IdI' P. % 014 H 4 tot lIt- N 100 P. "4 fin9 1014 01 IN n1 0o IN Cy IL - in "4 r- in t00 4 P. It 10' 0,0 m Ic c' 000'IN 0T1 P.W 1-41-4 C m cp 6 IN 4.I00 -" It It'It in a', tin IM N I'. JY14~P:., 1 0 0 1 Iin im Il CY N in ti.M o M'41 0I 0 MI 001IN 4 "4 I 1 1 50 0 F- I I In I- I I 1 1 1~~~~~~~~~~~~* n InCy t ul -4I IIY I In P

2LE VII TIME lmROIATED SATURATION VALUES FRCM RUN NO. 5 Heliht (a) Initial Distribution Satton 3.81 41.91 0.2028 O. 8336 11. 43 49.53 O. 2680 O. 9056 19.05 57.15 0. 4920 0. 9787 26.67 64.77 0.6728 o.985s4 34.29 72.39 0.7733 o.8778'rime (mitn.) Height~. e. _. 2 FOs......... Height (m.) 2 4 6 8 10 13 ~ —----------- ~ — - - -- -------- - 3.81000.24309.27980.31292.34273.37936.45774 11.43000.32640.38346.43603.47968.54107.57029 19.05000.54987.59874.64112.67361.69571.70120 26o67000.70528.72518.73621.74063.73576.73560 34.29000.78565.78202.76933.75888s.75177.74370 41.91000.83723.81823.79091.78251.77378.78645 49.53000.86491.80423.80783.80613.80324.80078 57.15000.94243.91829.88454.85591.82933. 80411 64.77000.88097.86488.83675.81238.78569.75061 72.39000.85538.81645.77515.73780.69423.64156 K~~~~~~~~~~~~~~~~~ ~~ ~ ~~~ ~ ~ ~ ~~~~~~~~~~~x.t ~.)..... i o2 o~.......... ^...5.0 i.8b0ght 1........16. 20 25 30 40 50 -.3.816i00.54826.66610.79551.814C2.82867.86770 11.43000- -.600-5 ---.63811l-.68 1.40.8409 8622 19.05000.70816.71483.71990.73217.75945.82339 26.67000.73749.73606.73918.75161.77302.79133 -~34.29000 -.74078 473443.74230.75844.77582.17817 41.91000.75602.74825.74032.74246.75625.16708 49.53000.80420.79098.76971.75620.78222.79939 - ---------- 51150.78412. —16049. —13108 -- -- -- -- - -743 -~ -658-.05 64!I'sIn000^.71611l f -w,0 Ar I-1.6652 6933I505 r.44824.3715

- co cc in c 0 I(% J%', I'-'C - F4 10 LA N 1o N H lo IA'!,,o - I N,.,,. I., iU iA I r_ lo 4' r.- - co co t LI No!~~~~~~~~~~~~~~~~~~~~i r.~~ N 4' r~ -.-0;'4'' L m 00 CM f" " co p- W 00 LA, P.- 0-'D I"'D 0 O. I U I',,: I I I.10 eOILA 0-4 - p.. r0N leo 0 IN,..N. cry 1~~~ 0 N ff1 ~~- CO CO 0 je'''IN I 0' N r Oh' ff 0 0 )r. Oh N U 1' H Ie It' / N N LA' N'. l~~~~ 1. - 0 0%.01.1~~~~~~~~~~~~~~~~~~~~t LA 00010o c l* je * I. 0 I*1 1 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~L P.1,n IN.1eI i1 0 0 0 0 0 40 1.4 -t UN co -r %O Ln 0% NO~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~0 IN 0 r''' 0,, Ul% 0 N0Na,,, N C i in 1n jN 0% m r- *0 Ul% en ~~~O 0, w U-N~.'' 1** I, o'C / I "o ff -' Oh I CC4 L 0%'O I c o m co'c. 1 I~~~~~~~~~~~~' I l a' LA LA 0\' N N i L'C LAi~~~~ in I~ Qr mInn.- * I leo MO'jff 1a'o 1 1ND P- 4' Oh 0 Y% I 4 L4 co'co co N % co r- LA 4' 4 I~~~~~~~~~~~~~ i ~ ~~~~~~~~ 0 i "C oI N 0 112'p.. IN co ino'Cfr I 0 %a m CY 0 11 w 0 0 0 I eN C uh N Co N CM OhP' LA Nv'C NO LA n N m 4' 1I'CY~I f" o O1 0 %Dt o w c co 1I * *- n e I I I I,:, Ii I- t I i i I ri f"6 - 40 0-, C.O 0 CO II-.4 0O 0 ff1 Ica'4 co LA 0" 0 I- kn tn 1,0 u, I a-. I'4'.4 0% 10.-. LAO 0 en.~~~~~~~~~~~~~~~~~~~ N 0~% 04 lIr N4.0 C0 N L ~\ON I''' LA' L''4' f':',~~~~~~~~! 10 &A 1.4 1 co IN. U., 101 W I~ 1~ N I 4 0"4.:.I,o ~0'I co /w,l:.o Ic o. I 0 N'C 0 O I' - I~ 14'l loL'O, l ~'C 0i, 0 L ~~~~1 I U l.a''~ hl. *.. L: [, C), I I I I I~~~I U'*,i i:'11 l~l~~~~~~~~ I I1 H Cq K%10 10 1,o 0 a.10 c 10 P V It 1, len enIA0 O ICat % - ln LA -f C". - -t I N~ cm 1,4c 04 *t IN Oh L.-o p. 1.1 1UN NO lco UN UN4 I.. o l.. 0 * 0 a *l LA II 10% lo - 0 " 11 0 1010 01 0! 0 toh ~ Ia' fr INlo l 4 O I 1 14' LA'C If'-'4~~~~~~~~~~~i L'I N 4 4'lA n ~' 6 I I n P I' I [ C 11 1 fill — 4 —I LJI 10: CY en:- 1I -1N C

Ij n r o m! 4 N.r m jt- joj UN - W% G IN 0%m,.c ~. o ra UN W4 %V-. % - 0- 0o rA V. CM'C00A "410 r C co0 LI-"4 I- I1, t 0 0 0%.A iN 0 V. VA "4 "4.,I', I 14 o (In -4 1 cJi N %o uI i u co W'I in1' o P- ^ ~ Uo.4'o *Z fl rN co b 4 o1'0 ^0 t ~ ro 0 lo r N t ts Ln CD co ax4 co t c 1. rl i', I i I I I i~~: o.:!: i~, I r' 1 1- > 1 ~, o1 ^ l 1~ 8 0\4 0 " | "!N IN - 0 0 j NI| L"rA V co 1 4 co 104 10. 0Qc VA P. & o~ 50.0o 14 4 LVo 0.0 0 cN I 0 in C, a, 10 (7 mo? W I- o in O - co P4 xI^' 0 0It10' I i IPrI, I- z P I; ls CI), ~O CO c' - 4 l C 1 U-%' I jo..%! 10.0 4. 0% 0 "4 1r jo M A j - 0c co P-~ VAi.0 ILA 0%%- tn inI Q~ ~~~~0 "4 I~ 4 iID ^ 1 nl S ~1 II I i I I I0 "4 "4 L I 1 f ( 0 V "4 I % *! im LA 0 o!.! o! 0-4 1 w g - OND.~4'f'.1. * * *~~~~~~% 1*. 00 N I 0 j I~~~~~~~~. 1 tLn I4 I I 1%0 e 4 i'VA 0 1 *i*1*0 0 10 co %iVAIN 0LA1 co N'.00 VA. 1LA UN 0 m t-~ 0% W\ P4 i4 CD0 N "4 0% 14 LA 10 N V A V- I,% 0~ t 0 V.- 0~~~~~~~~~~~~~~~~~L I'Or P I IC I iN VA N 14 VA 0 10.~~~~~~~~~~~ 90 10 a"4 co 0 N rN 0 CO N 1V.10 co Vn Nl 0 IAitAo N t' ( %' N coI I,' P-. - o.. 9 01 0 co 0 0%en CY, 410 CY rA m u UN It D c t 0 Is 10 1 1 01% 0 010 4.0 0~~~ n. ~~~ 10 0 4 0% N 0 0I0 10 t in 0% 4 C11 0 0 co " -4.- IQ VA 4 C 0 4 14 co 0I 0 VA NO 410 O,1 0Q jo I VA "4n 0% C to Qr v, &I 0 w0 * 1*' 04 1r. N N 0}'~~~~ "4jN.V.. 00%jO 010%~~~~~~~~1 n * n101 cli ~ ~ ~ ~ ~ ~ ~ II f- % L 0 "41.0,0 c 3 Ni 0% IN 10- 0 1.0 VA.0 w c 4O % 1'.0 1.0 co — CI N 0 P-.0 1o0 in V.. 41 UN tn -t.4.r~~0 "4'4 O 10 ~ j 0 i n o a'0 0%'0 0 0 10 0% 0 1e S o 0 0 0je

-128TABLE X TIE INTERPOLATED SATURATION VALUES FROM RUN NO. 8 Height (cam) Initial Distribution Sattion 3.81 1.91 0.8650.91 11-3 ~4953 0.9409 09721 19.05 5715 o.8833.9433 26.67 6.77 0.9126 0972 34.29 7239 o. 9230 o.6o6 Time (min. ) Height (c) 10 20 30 40 50 6o 3.81000.58030.66449.72930.74185.74540 74235 11.43000.90685.90036.90384 92844.9503.9700 19.05000.93042.93443.93726.94482.95668.96929 26.67000.92427.92659.92766.93341.94196.94497 34.29000.93786.93589 93128.93067.92968.93033 41.91000.93196.92100.91151.90997.92279.93393 49.53000.92866.90202.88770.88135.88075.8c3203 57.15000.88688.87491.87153.87083.86936.86709 64.77000.92634.92028.91817.92422.93612.93869 72.39000.87399.85276.81696.76667.69347.64359 HeighLt (am.) 80 100 120 P4o 180 2140 3.81000.74101.74362.74587.74766.750 51.75401 11.43000.98432.98204.97654.96953.96459.96746 19.05000.96146.97059.97454.97844.98132.98345 26.67000.95060.95812.95811.94786.94338.94-540 34.29000.93706.93532.93150.93048.93391.93740. 41.91000.9431-6.93426.93677.94405.94220.92827 49.53000.89485.90504.91277.91813.91865.91520 57.15000.86480.86501.86184.86217.86145.86137

-L?-9-: ~;, i i I': I I i: i t I I'' ~: ] o i.o,o 1,0 =,,0 ~~,...-,., ~ ~.,~ o,,.-,..,,-. ~.t I',*, t%l I~' ~0 i,~! o I,,, I,.;. ~. o o o ~. ~.' o,. j. i ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~., i':,..=1,,: I I'! I*~ ~: *..-~ I'-, C0', I',- 0",~ CO **..I O~,=4 ~ I',= ~0' ~ ~',. i tn ~n,... Ln.o,'m,,o m o ~ m.=,,o,o ao,o m o ~ ~ i0' eo ~.1' i~o.o,4' --..~ tn o~.4' aD m eo,o:i~,!~ ~ ~'om m,o ~ o'om,.-,,o,o,0 m ~ ~ m o~,-, ~ ~ O~ 0.0~!0', ~ O~ O~ 0', ~ eO O, ~ 0,=~ t.-4.~%1'0 ~0 ~! t"*' O~.~0 O~ ~0,-.~ 4"',0 ~,,1' O ~,,t' ~0,,1'' ~ ~ ~ ~ ~''i' I'.....'' o i:'"", i' ~1o~ o, l....' j, ~: l,o,o o,,o,,',,,, ~ 0, o o ~ ~ —, ~ ~.-, ~ ~ i ~ r,,. 0 0~!0~ O' ~' r,,.. ~0 0~ oc,-, ~ ~ ~ i ~ ~ ~ ~ ~ ~ ~. ~ } i'' ~' -~' i i ~ i ~.,' ~ u..-, o.,' ~,.- e ~,',, ~ ~., c;c;ddd, ~1 ~ CO 0~ 0:0 ~ I,~ ~1' I',,.. CO ~, O~ aO,41",,.=1'.. ~ ~!",., ~ t"~ P,- ~ ~ tr~,4*,,O ~ ~ 0~ I~, r' —.a 0 0 ~ a~ 0,, o O, o,. ~,4' ~ ~ co 0', o', a~,0 r~ $ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~' v,~_ O~ ~ ~-~ O~ O~ ~ 0 0 ~ P' ~!~. ~ b'~,0 I~' ~r~ O~ ~.~ 0~ _-:~ ~ j I",.. ~:C~,43 CO O.O ~,:~ ~ ~ 0',,0 r~ O~ ~-"' —~' 0 0,-4 I~. ~0 O~ 0, O~ O~ r~ I' —.co O~ O~ aO.~0 ~ i ~ ~ ~ ~ ~ ~ ~ -e ~ ~ ~ ~ i ~ ~ ~ ~,'.' i'' f-4,.4 i ~ ~,..-4 I~. 0 ~ a0 ~~ ~0,,1' O~ ~ ~ ~ ~ ~ ~ ~' Ii~ ~, ~0,0 O' cO aO ~.~ 0.-4 ~!~. i~..~ m'l ~,0! Oe O ~ ~ 00 0', O O 0'~ ~, rr~ ~ ~ I"',* aD C") 00 00 ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~. ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ i cu':'i -. ~,~ ~ ~ g,0,.- o,o.,-.,-,o o O ~,,0!0~ O O,O i t,-4 ~"1:t"~'r-, co ~', O a0 o~........ ~{ I~'' t ~ * ~ ~ i* ~' ~ * ~ I' i' I i i i i. i~: ~ oooo'ooooo: o, o o o,ooooooo o o o o o:!: o o o o o o o o;~ v O. c.~,,=4 ~,-~,'~ ~,'~ ~7',.-, m'~ ~ ~ 0" ~ I.~ P- O~..,4 ff'~ t,~ P- 0~ ~......j ~,:: ~,,......;:.~. ~,=4 N ~ -,1',,1" ~n,,o ~',...,,=., ~ ~ ~.,,t un,,o e-.. I~

APPENDIX F HISTORY OF INITIAL DISTRIBUTIONS The fact that the saturation history of every point in the column is important was not realized until after runs 3 and 4 were completed. For this reason, all that can be said about them is that the initial distribution was reached by imbibing the desired amount of liquid into the evacuated column and then inverting it. The initial distribution for run 5 was reached in the same way as those of runs 3 and 4, but the following data was taken before inverting the column. TABLE XII PRELIMINARY READINGS FOR RUN NO. 5 Position Height Reading j ~~~~~Y (cm) Rmj(cpm) 1 5~~~ ~~.81 125,666 2 11.45 155,875 5 19.05 158,7355 4 26.67 154,575 5 54.29 154,~507 6 41.91 1355,155 7 49.55 127,~276 8 57.15 159,~941 9 64.77 47,~962 10r -725 2,624 -1 --

-131Using the procedure of section VD these data could be converted to saturations, thus a complete history is provided for the initial distribut ion. Similarly, for run 6 the data listed below was taken. It should be noted that in this case the column was not inverted since all portions were wet during the filling operation. TABLE XIII PRELIMINARY READINGS FOR RUN 6 Position Height Reading Y (cm) Rmj (cpm) 1 5~~~~.81 4ooo6 2 11.45 4632,5 5 19-05 46077 4 26.67 45840 5 34.29 435565 6 41.91 41009 7 49.53 358194 8 57.15 59150 9 64.77 24815 10 72.139 11852 Note: C 40,970; Savg = 07 The initial distributions for runs 7 and 10 were reached solely

-132of run 8 although the initial distribution for this run was reached by imbibing into a vacuum and inverting. UIVERSITY OF MICHIGAN 3 9015 02526 1192