THE UNIVERSITY OF MICHIGAN COLLEGE OF ENGINEERING Department of Chemical and Metallurgical Engineering Progress Report EFFECT OF ULTRASONIC WAVES ON BIOLOGICAL MASS TRANSPORT H. Scott Fogler ORA-,Pro ject. 0319'l supported by: U. S. PUBLIC HEALTH SERVICE DEPARTMENT OF HEALTH, EDUCATION, AND WELFARE NATIONAL INSTITUTES OF HEALTH GRANT NO. HE-10549-02 BETHESDA, MARYLAND 20014 administered through: OFFICE OF RESEARCH ADMINISTRATION ANN ARBOR January 1970

I-' J

TABLE OF CONTENTS Section Summary... I Ultrasonic Gas AbsorptionM.L. Cadwell and H.S. Fogler-Paper accepted for publication in the Chemical Engineering Progress Symposium Series. This work is also to be presented at the 67th National A.I.Ch.E. Meeting in May, 1970............................... II Effect of Ultrasonic Waves on Membrane Transport............................................ III Acoustically Induced Facilitated Diffusional Transport in Membrane Ducts. This paper has been submitted for presentation at the 68th National A.I.Ch.E. Meeting in August, 1970........................... IV Acoustic Cavitation in Viscoelastic Fluids-H.S. Fogler and J.D. GoddardPaper presented at the 62nd Annual A.I.Ch.E. Meeting in Washington,D.C.Accepted for publication in the Physics of Fluids....................................... V Personnel......................................... VI iii

SECTION I Summary by H. S. Fogler I

Summary Page Grant Number: HE 10549-02 Effect of Ultrasonic Waves on Biological Mass Transport. Principal Investigator: Hugh Scott Fogler Sponsoring Institution: Department of Chemical and Metallurgical Engineering University of Michigan, Ann Arbor, Michigan 48104 Period Covered: Sept. 1, 1967 - Dec. 31, 1969 Date of Report: Jan. 31, 1970 Summary The following reports describe the preliminary results of the effect of ultrasonic waves on mass transfer. The mass transfer studies conducted under NIH Grant Number HE 10549 were in two areas: 1) gas absorption and 2) membrane transport. The studies on gas absorption were divided into two areas: large scale bulk streaming and our recent findings on high speed streaming in thin films. In the bulk streaming experiments, increases in the gas absorption rate of up to 800% were observed in the ultrasonic runs over the runs in which no ultrasonics were applied. However, one of the most interesting and potentially useful things we found in our study was that of the small scale streaming. This streaming was extremely rapid and was induced in thin films which were fractions of a centimeter in height. The fluid velocities in this film were greater than any previously reported streaming velocities which were either experimentally observed or theoretically predicted. These intense streaming currents will bring about an even greater increase in the mass transfer rate than was observed in the bulk streaming experiment mentioned above. In the membrane transport area, both theoretical and experimental investigations were also undertaken to determine the effect of ultrasonic waves on mass transfer. 3

In the theoretical analysis by Fogler and Lund, first results show that significant increases in the rate of mass transport through membrane ducts should be found with the application of ultrasonics. Experimental studies do indeed reveal that ultrasonic waves do increase the rate of transport through membranes. Steady flow experiments utilizing Amicon PM-30 and XM-100 membranes showed increases in mass transfer rates of 300% and 400% respectively, when ultrasonic waves were applied to the system. With the application of ultrasonics to biological systems, one must take care that the wave conditions are properly adjusted so as not to induce any harmful side effects to the tissue and body fluids by acoustic cavitation. Since it is known that various biological fluids behave as viscoelastic liquids, a brief and preliminary investigation was undertaken on cavitation in viscoelastic fluids to determine whether the degrading effects of cavitation could be accelerated or retarded in this fluid type. Preliminary results show that the elasticity can significantly retard the collapse process and certain situations produce damped oscillatory motion of the cavity rather than the usual damaging catastrophic collapse observed in purely viscous liquids. This suggests that under proper control, the damaging effects of cavitation could be minimized in biological fluids, and possibly bring about an additional enhancement in the transport rates as a result of the stirring produced by the oscillatory bubble motion. There was no cavitation present in either the membrane or gas absorption experiments reported above, as the observed increases were for systems operated below the cavitation threshold. 4

SECTION II Ultrasonic Gas Absorption The following paper represents a portion of the work by Mr. M. L. Cadwell and Prof. H. S. Fogler which has been accepted for publication in the Chemical Engineering Progress Symposium Series. This work is also to be presented at the 67th National A.I.Ch.E. meeting in May. 5

ULTRASONIC GAS ABSORPTION M. L. Cadwell and H. S. Fogler Division of Chemical Engineering University of Michigan Ann Arbor, Michigan 48104 7

I. Introduction The intent of this research was to study the influence of ultrasonics on the absorption of a gas contacted directly with a liquid. It was originally anticipated that ultrasonically induced "acoustic streaming" or "microstreaming" currents would be the primary mechanisms enhancing the mass transport process. A brief derivation and review of the equations governing acoustic streaming is presented in the first part of this phase of this report, and solutions to the streaming equations are presented for a few simple geometries. Initial results on the enhancement of gas-liquid mass transport are presented for streaming occurring in a cylindrical tube in which a high frequency (800 kHZ) ultrasonic transducer was aimed perpendicular to the liquid surface, thereby inducing gross circular streaming patterns within the liquid. The final phase of this investigation involved visual observations of acoustically induced convection in various geometrical configurations at frequencies of 800, 175, and 20 kHZ. Emphasis in this study was on how the convective patterns established might enhance absorption rates at the gas liquid interface. It is believed that the results of these observations fall into two distinct categories. At high frequency (800 kHZ) bulk acoustic streaming occurred which was similar to that observed 8

by Ghabrial and Richardson (1954) and described mathematically by Piercy and Lamb (1954), Nyborg (1964), and Eckart (1949). At low frequency (20 kHZ), when the liquid depth between the transducer and gas-liquid interface was much smaller than the sonic wavelength, high intensity vortex patterns were observed. These patterns apparently can not be explained in terms of simplified acoustic streaming theory, since the velocities observed were an order of magnitude greater than that predicted from purely elementary considerations. The Equations Governing Acoustic Streaming The following is an abbreviated development of the equations describing ultrasonic streaming in a fluid. Starting with the 3 basic equations of fluid mechanics Equation of Continuity aP + V(pU) = 0 (1) Equation of Motion p 4 P[It + (U.V)U] = -p + 4+K V(VU) - IVxvxU (2) Equation of State Vp = C2Vp + RV-P (3) Markham (1952) has shown that the inclusion of the term 9

RP- in the above equation can adequately account for heat transfer effects in many situations. It is assumed that the fluid is being perturbed ultrasonically with an angular frequency w = 27f, and that there is no movement of the fluid by any other external forces. As usual, we assume that the velocity, pressure and'density can be represented by U = Uo + U1 + U2 + (4) P = Po + P1 + P2 + (5) P = Po + P1 + P2 + (6) where Ui =a(E) Ui,, P. = a(E) P and pi = a() p - ~ I "I i-li i-1 Substituting the above expressions into equations (1), (2), and (3) and collecting terms of like orders of magnitude, noting that U0 = 0 since there are no applied external forces other than ultrasonic perturbations, one arrives at the zero, first, and second order perturbation equations: Zero Order Uo = 0 (7) VP0 = 0 (8) Vpo = 0 (9) 10

First Order aP + po(V U1) = (10) at Po = -VP1 + [V+K V(VUi ) - iVxVxUi (11) VP1 = C2Vpl + R-P (12) Second Order aPa + V*piUl + poV'U2 = 0 (13) t rU2 Ui, Po T + (Ui.V)Ui + Pt = -VP2 + I 4+KlV(VU2) - PVXVXU2 (14) VP2 = C2Vp2 + RVaP2 (15) Both pressure and density may be eliminated from equation (11) to yield a2U1 au1 p0 = poC2V(VUi) + +K+RpoV V - VxVxU (16) at2 t ~ Multiplying equation (10) by U1, and adding it to equation (14): aU2 a (17) Poat + o0{(U1'V)U1 + U1(V'U1)} + (p ) (17) 11

=- VP2 + L +KV (VU2) - 1VXVxU2 Equations (16) and (17) in more simplified form provide the basis for ultrasonic acoustic streaming. Since U1 represents the velocity due to sonic perturbations, it is sinusoidal with respect to time; it can therefore be represented by: U1 = Ui(r)ejwt (18) where U (r) is a velocity vector dependent only on spacial position. Because equations 10, 11, 12, which relate pi and P1 to U1 are linear, it is seen that pi and P1 also vary sinusoidally with time. Substitution of equation (18) into (16), and dividing by X2po/2 2C + 2j f43 +K+poRJ] V(V*U1) + 2U1 -= 2 VxVxUl (19) 2 opo 3 3 Po Po~ The following parameters frequently appear in the literature on ultrasonics: k = w/C, defined as the wave number; C = 2 4 K Rp+ +, defined as the attenuation coefficient; poC 3 + - and 3 = /wpo/2p, recognized as the reciprocal of the acoustic boundary layer thickness. Substitution of these terms into equation (20) 2 + 4 K + Rpo oV(V-U1) + 2U1 = VxVxU1 Because, under most conditions a << k << B: 12

2 j F!+ 22 K + Rpo 2 k2 3 ~ - (k-jea) 2 Therefore, equation (21) takes the form (k-a)2 V(V*U1) + 2U1 = VxVxU (21) The above solution can be divided into two terms, (Eckart, 1948) U1 = Ulc + U1R such that VxUir = 0 and V2Uir = - (k-ja)2Uir (21a) V-Uic = 0 and V2Ulc = 2jB2Ulc (21b) Procedures suggested for solving equation (21) or equations similar to them have been reviewed extensively (Nyborg, 1964, Rayleigh, 1945, Schlicting, 1960, and Eckart, 1948). Once U1 is obtained, pi and P1 can be obtained from equations (7) and (8); one can then proceed to solve Equation (19) for U2. In practice, this is a formidable task, for the underlined terms in equation (17) are nonlinear. A simplified approach to the solution of equation (19) is to time average each individual term over a few sonic cycles. The nonlinear nature of the term (U1'V)U1 + Ui(V*Ul) suggest that a solution to U2 is of the form 13

U2 = U2(r) eLjwt + U2(r) (22) where U2(r) and U2 (r) are dependent only on spacial position. Substituting the above expression into Equation (19) and performing the necessary time averaging (indicated by an overbar) po{(U1'V)U1 + UI(V.UI)} + U (VpUu) =VP2 + [4 +K)V(V7U2) - PVxVxU2 (22) From the solutions to the first order equations, it can be shown that at (Piu) = and it has been shown that V(V'U2) = (V2U2) o[2L) (Rayleigh, 1945); thus, for most systems: V(V*U2) <<< V2U2 Therefore equation (22) becomes Po{(Ui'V)Ui + Ui (V.Ui)} = -VP2 + UV2U2 (23) This equation is known as the basic streaming equation. Several examples from the literature will be reviewed here to 14

illustrate the types of solutions which can result upon application of equations (23) and (21). In the first example, an ultrasonic beam of radius Rb is propagated along the axis of a very long circular cylinder whose length is L and radius R. The solution to equation (21) is assumed to be irrotational. Thus, V2U1 = -(k-ja)2U1 If the length of the cylinder is sufficiently long, aL >> 1, and the first order solution is: U1 = Ae cos(wt-kx) r < Rb U1 = O R < r < Ro The discontinuity in first order velocity profile arises from the assumption that U1 is irrotational. The streaming equation becomes, for radial symmetry: ~2 U2 1 i2] 3P + p<aA< r2 r r ] = + A2 r Rb mI r-" J a+ 2 Rb < r < Ro The boundary conditions are that U2 = O at r = Ro 15

B2 = 0 at r 0 Br and if the ends are capped so that there is no net bulk flow across any cylindrical cross section, Ro U2 rdr = 0 0 The solution is: cpoA2Rb2 J R F 2 R 2 U2 = n R- + j 1 -r l r 1t for r < Rb and p09A2 Rb2 j r ( r I,21(Rb 2 U2= - R - 1- R for Rb < r < Ro This solution was derived by Piercy and Lamb (1954), and represents a slight variation from that originally derived by Eckart (1948). Nyborg (1964) presents an analogous solution for the case of propagation of a beam of sound with Zb, between two parallel plates separated by a distance h. Defining Z = 0 at the edge of the lower plate, his solution is: 16

U2 po- A2 e-2ax (h-Zb)Z + Z(Z-h) 2-31 + j 2 for Zb/2 < Z < h/2 These solutions are characterized by a bulk, non-circular (except at the ends, x = 0 and x = k) "fluid motion. A somewhat different pattern results, however, if the first order solution along the plane mid-way between the plates is that of a standing wave, U1 = A coswt cos kx and if the sound beam propagates along the entire distance between the plates. Under these conditions, both a horizontal and vertical velocity exists: Uo2sin2kx -y U2 = -.8C eY (4 sin~y + 2 cossy + e ) 3 2( 22 () (s + 3cosy + V =............ + 2 B(h-y) - - 3hil - ^V The fluid patterns are circular, and repeat themselves at integral ~-wavelengths. This solution was originally derived by Rayleigh. 17

List of Symbols Used in Equations Governing Acoustic Streaming A Amplitude of velocity wave C Speed of sound K Ditational viscosity L Length of cylinder P Pressure Po Zero order pressure P1 First order pressure P2 Second order pressure R Radius of cylinder Rb Radius of transducer U Velocity vector Uo Zero order velocity vector U1 First order velocity U1(r) First order velocity, a function of position only Uic Incompressible part of first order velocity Ulk Irrotational part of first order velocity U1 First order velocity in axial direction or parallel to plates U2 Second order velocity vector U2 (r) Second order velocity vector which is a function of position only U2 Streaming velocity vector, time-averaged value of U U2 Streaming velocity along axis of cylinder V2 Streaming velocity component perpendicular to flat plates Z Coordinate perpendicular to plate Zb Thickness of sound beam propagating between flat plates 18

f Frequency h Distance between flat plates j Complex number, /Vk Wave number r Position vector r Radial coordinate x Coordinate along axis of cylinder t Time a Attenuation coefficient 3 Reciprocal of acoustic boundary layer thickness e Perturbation parameter p Density Po Zero order density Pi First order frequency P2 Second order frequency w Radial frequency 19

Experimental Gas Absorption Studies Experimental Apparatus The apparatus used to measure the amount of gas absorbed into the liquid in the absorption cell is shown in the accompanying diagram (Figure 1). Three vertical manometers, M4, M3, and M1, measure the pressure of tanks 1, 2, and the absorption cell itself, respectively. The pressure regulator and vacuum pump allow operation of the system over a pressure range from a few millimeters of mercury absolute to about 2 atmospheres. The slant manometer uses water as the working fluid and measures the difference in pressure between tanks 1 and 2. The gas used in these initial studies was C02, although in later studies it is anticipated that other gases will be used. Tanks 1 and 2 can be isolated from each other by closing valves VS1 and V10; in addition, if valves V16 and V13 are open, tank 2 is connected directly to the absorption cell, and thus serves as a source of gas for the cell. The number of moles of gas absorbed during a run can then be calculated from the change in height of the slant manometer. Assuming, therefore, that (1) the pressure in tank 1, tank 2, and the absorption cell is/the same, (2) the pressure in tank 1 is always uniform, (3) the pressure in tank 2 is also uniform, and equal to that in the absorption cell at any time, (4) the system operates isothermally, and (5) the density of the manometer fluid is much greater than that of the gas, the following equation gives a relationship between amount of gas absorbed, nL, and change in height of the slant manometer in terms of system properties: h: RT (1) sincAPo A 1 (2 1) PmV1 ( Amh\ Amh) 2n

VALVE~ ^^^;l ~ ^(~TANK 2 ) r TANK I QIU TANK I~ ~....ABSORPT ION ( ^ ANT SWITCH GAS MANOMETER I CYLINDER 1 (~^^ ^ I I ^j MANOMETER MANOMETERS VACUUM MANOMETERANT GAS MANOMETERS Figure 1. Gas absorption measuring system.

where: nL = amount of gas absorbed in cell h = change in slant manometer height Po = initial system pressure a = angle slant manometer makes with horizontal Am = slant manometer cross sectional area V1 = volume of tank 1 V2 = volume of tank 2 P = density of manometer fluid g = acceleration of gravity R = gas constant T = absolute temperature For the system shown in the accompanying diagram (Figure 1), V1 = 600 cm3, V2 = 8,200 cm3, Am = 0.785x10 cm2 p = 1.0 g/cm3, and sina = 0.1. Because the manometer is designed so that h is always less than 100 cm, V1/V2, Amh/2V1, and Amh/2V2 are all less than unity. Substituting the above numerical values into Equation (1) and taking Po = 2 atm. and T = 298"K, Equation (1) becomes approximately h = 3.79x105 nL It is seen that the system is quite sensitive to very small amounts of gas absorbed. The absorption cell used in this study was cylindrical, 12" long and 1-7/8" inside diameter, and was constructed of 3/8" thick transparent plexiglass. The cell was designed so that a Macrosonics 800 kHz transducer could be bolted to the cell bottom. Thus, ultrasonic energy propagated upward along the axis of the absorption cell; however, the diameter of the radiating diaphragm of the transducer is 1.45", and therefore the ultrasonic beam did not fill the tube entirely. An O-ring seal was placed 22

between the transducer and bottom of the cell to prevent leaks. Three openings were placed near the top of the cell. One of these openings was connected to manometer M1 to measure the pressure in the cell, another to the gas absorption measuring apparatus, and the third to a line leading to the vacuum pump. Ultrasonic power to drive the transducer-was supplied by a Macrosonics 500-2 Ultrasonic Broad Band generator whose maximum output is 500 watts. Procedure At the start of a run, the absorption cell was filled with liquid to approximately the desired height and inserted into the system. (If the fluid was distilled water, it had been boiled for approximately 1 hour before insertion into the system to drive off any dissolved C02 vapor.) The cell was then evacuated for 1 hour to draw off any gases that might have dissolved in the fluid as a result of exposing the fluid to air during its transfer into the absorption cell. The height of the liquid above the transducer was then measured with a cathetometer. Once the cell was inserted into the system, gas was fed from the cylinder through the humidifier into tanks 1 and 2. Humidification of the gas was essential when using water as the absorption fluid to minimize evaporation effects. At this time valves V3, V4, V5, V6, V7, V8, V9, VS1, V16, and V10 were open, while valves V12, V13, V14, and V15 were closed. When the system reached desired operating pressure, valves V4 and V3 were closed, and the system was allowed to stand until pressure equilibrium was attained. While this was occurring, the ultrasonic generator was adjusted to give the desired output to the transducer. Measurements of the room temperature 23

3.0HW |^~~~ / I~* /! With Ultrasonics ao / V Without Ultrasonics w O m ) 2.0 / o / 00 (Lf) LL Z / I ^~%~ ~ From Equation (2) 0.0 1.0 2.0 3.0 4.0 /( TIME) IN (MIN.)/2 Figure 2. Representative comparison of CO2 absorption in water with and without ultrasonics.

H Qa * Wit h Ultrasonics ^ x 2.0 * v Without Ultrasonics c d. (In NZ 0 r ^^ / LL 0^ 1.0 (f) (.) (I)Q UJ /n 0 C\1 0c'J 1.0 2.0 3.0 4.0 /~~~~~~~~~~~~/ ^(IME) IN (MIN.)^ Figure 3. Comparison of C02 absorption in glycerol with and without ultrasonics.

and pressure differences in all manometers were taken when equilibrium was attained. The actual absorption run was initiated by opening valve V13 to allow gas to enter the absorption cell. When the pressure in the absorption cell was equal to that in tank 1, valves V10 and then VS1 were were closed, thereby isolating tank 2 from tank 1, which was connected to the absorption cell. The timerwas started when valve VS1 was closed. The change in slant manometer height was recorded as a function of time. Measured values of h were corrected for system leaks by operating the system without the absorption cell and following the change in manometer height under conditions nearly identical to that in which the actual absorption run was performed. Preliminary Results Some results, which might be considered representative of the best data taken thus far, are shown in Figures (2) and (3). Most of the data taken with ultrasonics lie between the two curves depicted in these accompanying figures. Two systems were studied: CO2 absorption into water and CO2 absorption into glycerol. The number of moles of gas absorbed was calculated according to the procedure described above and then corrected to 25~C and 1 atmosphere pressure. All data obtained was compared to the amount of gas dissolved if only pure diffusion of gas into the liquid were occurring. The amount of gas absorbed for this latter case (i.e., no ultrasonics) is given by the equation AP aPw 4DABt (2) a; =?I —i (2) where: Wa = number of moles of gas absorbed at time t t = time 26

A = absorption cell cross-sectional area P = partial pressure of absorbing gas above the liquid H = Henry's law constant P = density of the liquid Mw = molecular weight of the liquid DAB= diffusion coefficient of dissolved gas in the liquid w = mass density of the fluid Thus, results are plotted as the number of moles absorbed versus the square root of time. It must be pointed out that there is as yet no firm theoretical basis for plotting the data labeled "with ultrasonics" in this way. This method, however, does serve as a convenient basis for comparative purposes. With reference to the data shown in Figure 2 for CO2 absorption into water, the run labeled "without ultrasonics" was carried out at a pressure of 72.3 cm Hg (CO2 partial pressure of 69.6 cm Hg) and temperature of 81~F with the cell filled with fluid to a depth of 22.4 cm above the transducer. The run labeled "with ultrasonics" was carried out at a pressure of 75.35 cm Hg (CO2 partial pressure of 73.6 cm Hg) and temperature of 67~F with the cell filled with water to a height of 9.43 cm. The plate voltage and current of the ultrasonic generator during this run were 0.52 volts and 260 milliamps, respectively. In the C02-glycerol runs shown in Figure 2. the run labeled "without ultrasonics" was carried out at 74.9 cm Hg pressure and 73~F temperature with a liquid height of 9.7 cm, while the run labeled "with ultrasonics" was carried out at 115 cm Hg and 75~F with a liquid height of 9.72 cm. The generator plate voltage and current in this run were.80 volts and 300 mi 1 liamps, respectively. At this time, these results can only be regarded as preliminary. There are inherent experimental difficulties in the operation of any 27

unsteady state system; this is especially true for gas-liquid systems. As seen in both figures, extrapolation of data "with ultrasonics" to zero moles absorbed indicates that there is a response time lag which is either operational or inherent to the system. Efforts have been made to reduce this lag through slight modifications of the system and through improvements in operating procedure, but as yet no trend can be observed of how this time lag varies with system parameters (in particular, height of the liquid above the transducer and system operating pressure). Nevertheless, this time lag can be estimated by extrapolation of the curve of CO2 absorption into water "without ultrasonics", and comparison of this curve with Equation (2) using values of diffusion coefficient and Henry's law constant from the literature. On this basis the time scale for this run was adjusted to the left. The dotted line in this figure is that predicted by Equation (2). At the end of the curve for CO2 absorption into water without ultrasonics, an abrupt increase in the number of moles of gas absorbed was seen to occur. This increase might be caused by natural convection which has been observed by other investigators in similar gas-liquid systems in which the density of the solution (solute + solvent) is greater than that of the liquid, causing an unstable density gradient to be established in the liquid. The occurrence of natural convection is governed by the Rayleigh number, which in inherently unsteady state systems is time dependent. The runs made in the latter portion of this study used glycerol as the absorbing fluid; because of its high viscosity, natural convection effects should therefore be less important than in the C02-water system. Another difficulty encountered in operation of this system was due to the slant manometer, which is very responsive to only slight upsets 28

within the system. This manometer had to be refilled frequently or replaced. Although care was taken to fill this manometer with clean, distilled water, after continued use dirt did accumulate at certain points within the capillary tube of the manometer, resulting in erratic manometer response, or manometer fluid separation. Data taken under these conditions were very difficult to interpret. However, a number of runs were carried out without the above difficulties. An estimate of the increase in the mass transfer rate can be made by comparing the slopes of the curves shown in Figure (1) and (2). On the basis of this comparison, it can be seen that the increase in mass transfer rate with the application of ultrasonics is on the order of 800% for both liquid systems. While no firm conclusions have yet been drawn as to what the cause of this enhancement is, these data were obtained without the effects of cavitation and atomization, and therefore the mechanism causing these increases is most likely convective streaming. It appears, therefore, that surface renewal concepts might be particularly appropriate for modeling this phenomena mathematically. Preliminary Visual Studies of Streaming Qualitative visual observations of acoustically induced convection were made at three frequencies, 20, 175, and 800 kHz in various geometrical configurations. Power to the transducers was furnished by a 500-1 Macrosonics Broadband Generator, capable of delivering a 500 w. input over the frequency range of 10 to 1000 kHz. Glycerol was used as the fluid; because of its low vapor pressure (.0025 mm Hg) and its high viscosity (1440 cp) at 250C, cavitation effects should be small. The fluid movement was followed by the use of tracer particles of either 29

alumina or cigarette ashes. In some cases, the observed convection patterns were photographed with a motion picture camera. High Frequency Studies. Two geometrical systems were used with 800 kHz insonation. In the first, a geometry identical to that used in the gas absorption studies was employed. A cylinder 9" long and 1-7/8" in diameter was filled to a depth of about 4-1/2" with glycerol. Insonation was accomplished by fastening an 800 kHz transducer to the bottom of the cylinder such that the vibrating diaphragm of the transducer was in direct contact with the fluid. The diameter of this diaphragm was 1.45", and therefore the sound beam did not fill the tube. When the liquid was insonated, intense streaming currents were observed throughout the liquid. These currents appeared to be circular, symmetrically located about the axis of the cylinder, and ordered - that is, tracer particles within the fluid appeared to move along well-defined paths, which they repeated after a certain period of time. The most intense streaming appeared along the axis of the cylinder; velocities here were estimated to be on the order of 10 cm/sec. The magnitude of these currents increased with an increase in the power supplied to the transducer. There was no evidence that cavitation was present during insonation in this experiment. A sketch of the observed streaming patterns appears below. Tra nsdeeer Figure 4. Approximate Fluid Patterns in 800 kHz Cell 50

In the second study at 800 kHz, a cylinder whose diameter was about 5-1/2" was filled to a depth of about 1" with glycerol. An 800 kHz transducer identical to that used in the study described above was fastened to the bottom of the cell about 2" from the axis of symmetry of the cylinder. Enough power was supplied to the fluid to produce streaming, but not enough to make the interface unstable. The Induced fluid motion brought about by insonation was recorded on standard motion picture film at a speed of 18 fps. Streaming speeds near the transducer were estimated to be 2 cm/sec, about 5 times smaller than that observed in the first 800 kHz study. Fluid motion was most prominent near the transducer, and there were pockets of stagnation in areas away from the transducer. Again, no cavitation phenomena were observed. A sketch of the patterns observed appears below.:to 4 or-F e e Figure 5. Streaming Patterns in a 5-1/2" Diameter Cell at 800 kHz 531

Low Frequency Study. Observations at 20 kHz were performed in three different geometries. In the first about 50 ml of glycerol were added to a 150 ml beaker, and the tip of a sonic gun transducer was immersed about 1/2 cm from the fluid air interface. Insonation was applied and the streaming pattern from the side appeared as if two vortex rings had formed between the bottom of the beaker and the tip of the transducer. Photographs were taken of the streaming patterns from the side, but because of extensive dispersion of the particles, an estimate of the streaming velocity could not be obtained.:" ~'~ Tr-a~s~ve er en fe r Aame Figure 7. Streaming Patterns at 20 kHz in a Beaker In the second low frequency experiment, the box used for the 175 kHz study was filled to a depth of 3/4" and insonated indirectly by placing the tip of the transducer gun perpendicular to the bottom of the stainless steel bottom of the box. Moving pictures of the ultrasonically induced fluid motion were taken from both the side and the top view of the box. A sketch of the streaming motion is shown in the following figure. Surface velocities, measured by following tracer particles at a distance of about 1" from the axis defined by the center of the tip of 32

Intermediate Frequency Study. A 175 kHz transducer with a radiating diaphragm 59.5 mm in diameter was placed in a rectangular box 5" wide, 8" long, and 5" deep. The box was filled with glycerol until the top of the transducer housing was slightly over an inch below the face of the transducer. The transducer was placed approximately in the center of the box. The streaming patterns formed upon-insonation are sketched below; particle speeds were greater than those in the last 800 kHz study, and were on the order of 5 cm/sec at a point about 1 cm below the gas liquid interface. This motion, however, was confined to that area directly above and slightly to the sides of the transducer. Streams of cavitation bubbles were emitted at points of rapid fluid motion; there was also evidence of the existence of resonating bubbles existing on the transducer surface which also might have contributed to the fluid motion.:5 +Y~I Trnsdveer Figure 6. Streaming at 175 kHz in a Rectangular Box Moving pictures of this motion did not come out as well as expected because of mechanical difficulties encountered in operating the camera, and because of very rapid particle dispersion upon insonation. 33

/ OIP.i VI I onlfIof Trahsdveer SiDE IVJ-W TOPVIE0 V Figure 8. Streaming Patterns in a Box at 20 kHz the transducer, were on the order of 3 to 6 cm/sec and increased rapidly as these particles approached the point of insonation. Also, small rapid vortices, on the order of 3/16" in diameter, existed around this point. There was, however, activity around the point of insonation indicating that cavitation was present. In the last set of observations carried out at 20 kHz, a cylindrical tube,l 11/16"in diameter, fitted with a stainless steel bottom, was filled to a depth of 1-3/8" with glycerol and again insonated indirectly by placing the tip of the transducer on the stainless steel bottom. Streaming patterns were symmetrical about the point of insonation, as shown in the sketch below and in the photographs taken. The velocity along the axis of symmetry was estimated to be about 1.0 cm/sec, although this figure may be low due to the extensive amount of dispersion of tracer particles. The presence of cavitation could not be discerned from the photographs. 5^~fv~ecing P;'tte s~Batton T1v) 5fqnlss Sm4ec Tw Hsd \cer I Figure 9. Streaming Patterns in a Cylinder at 20 kHz - Moderate Liquid Depth 34

The liquid height in the cell was then reduced to 1/2". Exceptionally rapid mixing was observed when the fluid was insonated and streaming was especially rapid a few millimeters from the point of insonation. A sketch of the streaming patterns appears below. rrahsivd er Figure 10. Rapid Streaming Patterns at 20 kHz in a Cylinder - Small Liquid Depth An attempt was made to take moving pictures of the ensuing motion, but because of the difficulty of simultaneously operating the camera and turning on the sonic power, and as a result of very rapid particle dispersion, it was not possible to obtain a good film strip suitable for giving and estimate of the magnitude of the streaming velocity. The pictures that were taken, however, indicate the presence of cavitation. When the transducer was moved off center and forced upon the stainless steel bottom of the cell, very rapid streaming velocities were encountered, in the pattern similar to Figure Pou,+ OPl h.o,, af ~o L ria ttes Figure 11. Top View of Rapid Streaming at 20 kHz in a Cylinder - Small Liquid Depth This pattern may be due, however, to the "fountain" effect encountered in ultrasonics at high intensities, i.e., interfacial instabilities. 55

Discussion of Observed Patterns These preliminary studies indicate that with the equipment available in the Sonochemical Laboratory, low frequency insonation induces larger convective velocities than high frequency insonation. This observation is especially true for the case of streaming occurring in liquids whose distance between the transducer and the gas-liquid interface is much less than the wavelength of sound. Because of the wavelength/associated with 800 kHz insonation is only.19 cm. this observation is difficult if not impossible to verify for high frequencies. While the studies described here can at best be termed exploratory, some comparison and comment on the modes of induced convection normally attributed to ultrasonic insonation can be made. Order of magnitude estimates of the streaming velocity can be made from elementary streaming theory reviewed previously. For a traveling wave, an order of magnitude estimate of streaming velocity is (Piercy and Lamb, 1954) 2w2 u = 4-+ R2 "2 3 u u/ R And for a standing wave, the maximum streaming speed should be (Rayleigh, 1945) 3u1 u2 16c The first order velocity u1 is equal to the product of particle displacement amplitude/and the frequency. If the value for s is taken as 10 microns at 20 kHz, and 0.12 microns at 400 kHz, and if it is assumed that K/u << 1 and Rpo/u << 1, for the system used in this study Piercy's

and Lamb's (1954) analysis gives a streaming speed of 0.6 cm/sec at 20 kHz, and 1085 cm/sec at 800 kHz; using Rayleigh's expression, u2 is.002 cm/sec for both high and low frequency. These theoretical estimates are at least one to three orders of magnitude lower than that observed at low frequency, and, for the case of a traveling wave, two or three orders of magnitude higher. The experimental observations do not agree at all with the theoretical prediction for the traveling wave case that an increase in frequency results in an increase in streaming velocity. Several investigators (Eckart, 1949, and Westervert, 1960) have advanced the explanation that at sufficiently high ultrasonic energy level inputs, a transition occurs from laminar to turbulent streaming. In turbulent streaming, particle velocities and vortex sizes far exceed that predicted by the classical equation governing acoustic streaming. One study (Westervelt, 1960) in infinitely long cylinders subjected to oscillation perpendicular to the axis of the cylinder indicates that the transition from laminar to turbulent streaming occurs at a "streaming Reynolds number," s2w/v, equal to unity, as long as the wave displacement amplitude (s), to cylinder radius, R, is sufficiently small. For the same system, others (Andres and Ingard, 1953) have shown theoretically that the classical solution offered by Schlichting (1932) applies for higher streaming Reynolds numbers. An estimate of the streaming Reynolds number as defined by Westervelt can be made since it is independent of system geometry. Glycerol, the fluid used in this study, has a viscosity of 1440 cp and density of 1.2 g/cm3 at 250C. With a particle displacement amplitude of 20 microns at 20 kHz, the streaming Reynolds number is about.001; at 800 kHz, with 37

a particle displacement amplitude of 0.12 microns, the streaming Reynolds number is about 10. According to Westervelt, then, both low and high frequency regimes studies here should be described by the "classical" equations of acoustic streaming. While this appears to be the case for high frequency studies, it does not appear so for the low frequency cases. This analysis assumes that the particle displacement amplitude, s, is small compared with some (as yet undefined) length of the system. However, the cylinder radii, liquid heights, and wave lengths (7.5 cm at 20 kHz,.1875 cm at 800 kHz) are all greater than s, and hence his criteria appear to be satisfied. These observations were performed at sonic power inputs of about 300 - 500 watts. Assuming that the energy efficiency is about 25%, the average sonic energy intensity with respect to the cross-sectional area is 16-26 watts/cm2 Since the threshold cavitation intensity for castor oil, whose viscosity is about 630 cp at 250C, or about 1/2 that of glycerol, is 5.3 watts/cm2 at 20 kHz, it is conceivable that cavitation was present in the observations made at low frequency. This might be especially true near the point of insonation where the sonic intensity would be a maximum, possibly two or three times as great as the above estimated values. On the other hand, because the viscosity of glycerol is several times that of water, whose cavitation threshold is 250 watts/cm2 at 800 kHz, no cavitation should be present in high frequency insonation. This is in substantial agreement with what was observed. What role cavitation plays in the observed convection process has not been ascertained at this time. Although agreement with simplified theory is poor, it appears that 38

at high frequency a form of acoustic streaming similar to that observed by others predominates as the convective mechanism; however, the very high speeds and rapid vortices observed at low frequency, do not appear to be explainable right now in terms of elementary streaming theory. Recent Discoveries on Rapid Streaming in Small Liquid Depths As mentioned previously, intense streaming vortices were observed when the liquid level, h, was reduced such that h/X << 1. These rapid streaming vortices cannot be explained in terms of previously derived streaming solutions which, in most cases, neglect the depth of the liquid above the transducer. In this regime, it appears that steep velocity gradients are induced because of the low liquid level. Because of the short distance between the transducer and the gas-liquid interface, these gradients induce exceptionally rapid streaming vortices which extend from the point of insonation up to the interface; thus, this type of streaming appears to be particularly applicable to gas-liquid mass transfer, especially when analyzed from a surface renewal viewpoint. To our knowledge, this type of streaming has not been uncovered before, either experimentally or theoretically. Since the efficiency of most mass transfer equipment, and especially apparatus like the artificial lung and kidney machines, is increased considerably as the surface area to volume ratio is decreased, the form of streaming described here will be studied in greater depth. Experimental work will be pursued to define and specify the conditions under which this unique streaming occurs. Photographic studies will be conducted to establish the optimum streaming pattern with respect to mass transfer. Our preliminary photographs indicate that the important parameters in this high-speed streaming are the ratio of liquid height/sonic 39

wavelength, type of fluid, ultrasonic intensity, point of insonation, and cavitation intensity. It is also believed that these patterns can be adequately explained mathematically if the appropriate approximations are made to the streaming equations. Literature Ci tati ons Andres, J. M., and Ingard, U., J. Acoust. Soc. Am., 25, 928, (1953). Eckart, C., Phys. Rev., 73, 68, (1948). Markham, J. J., Phys. Rev., 86, 497, (1952). Nyborg, W. L., "Acoustic Streaming," Ch. 11 in Physical Acoustics, Vol. II-B, Muson, W. (ed), Academic Press, New York, N.Y., (1964). Piercy, J. E., and Lamb, J., Proc. Roy. Soc., A226, 43, (1954). Rayleigh, Lord., Theory of Sound, Vol. II, Dover Press, New York, N.Y., (1945). Schlichting, H., Physic Z., 33, 327, (1932). Lo

SECTION III Effect of Ultrasonic Waves on Membrane Transport (To be submitted for publication) 4'

EFFECT OF ULTRASONIC WAVES ON MEMBRANE TRANSPORT H. S. Fogler, W. Franklin, V. K. Verma and M. L. Cadwell Division of Chemical Engineering University of Michigan Ann Arbor, Michigan 48104 43

I. Introduction The resistance to the transfer of a solute species across a membrane is due to the concentration boundary layers which exist on either side of the membrane, and to the membrane itself. Membrane transport is the principal mode of mass transfer in nearly all vital organs of the body, but because of the extremely high ratioJsurface area for mass transfer to volume of most of these organs, artificial duplication of these organs has been very expensive, both in terms of operatic and investment. Artificial machines designed so far have been of limited utility due to their large size. Efforts have been made to reduce the size and inherent complexities of operation. The size of these machines can be reduced by reducing the resistance to membrane transfer. In steady flow systems, the boundary layer resistances on either side of the membranes can be significantly reduced by increasing the through-put of material, and in batch operated systems, these resistances can be reduced by agitation. However, the induced fluid-mechanical shear rates brought on by these increases almost always result in permanent damage to cells transported in vital fluids. In addition, if the principal resistance to mass transfer is in the membrane, increasing throughput in steady flow systems and increasing agitation in unsteady systems are only marginally effective in increasing the efficiency of the artificial machine. Under proper conditions, the application of ultrasonics can significantly reduce all three resistances to mass transfer. Ultrasonic microstreaming has been observed at solid-fluid interfaces; this microstreaming

should significantly reduce the resistances to mass transfer outside the membrane. Also, since the ultrasonic wave can be propagated into the micropores of the membrane, this mechanism should also enhance mass transport within the membrane itself. At sufficiently low intensity levels, both of these reductions can be achieved with a minimum of cellular or membrane damage. In summary, the application of ultrasonic waves to membranes has shown significant increases in the mass transport rate. In addition to enhancing transport rates in artificial membrane devices, ultrasonics could also enhance internal cellular mass transport in the body fluids and cells. II. Experiments Utilizing a Non-flow System A. Apparatus Two pyrex glass half cells, each of 300 ml capacity, were separated by a membrane and clamped together with slip-on flanges. An O-ring seal was inserted to prevent leakage. The solutions in both cells were stirred by the propeller type variable speed stirrers introduced into position as shown in Figure (A). Samples were drawn at suitable intervals (1/2 or 1 minute) through sample holes, and the samples were returned after conc. measurement to maintain a constant volume of solution in the cells for.0~~~~~~~~~Tr. d ce o m~o^1/_ _~^-^ ^^ JTjsdu^r Mewbvene/ T; \HClf C Figure A 45

each run. Since sodium chloride solutions were used, the concentration of the solution was monitored by a Beckman conductivity bridge in conjunction with standard conductivity cells. Ultrasonics was applied directly to the solution, in a direction perpendicular to the membrane, by a transducer operated at its natural frequency of 800 kHz, and the intensity of the ultrasonic waves was controlled indirectly by monitoring the power output from a Macrosonic 500-1 generator. B. Mathematical Analysis Pseudo-steady state mass balance for the solute in the two half cells gives dC1 dC2 VI d V2 dt KA(C1 - 2) (1) where: V1, V2 = solution volumes, constant C1, C2 = solute concentrations A = area for transfer K = overall mass transfer coefficient Integrating Equation (1) with initial conditions of: at t = 0, C = C1 and C2 = C2, we have: 0 0 V1C + V2C2 = V1Cl + V2C2 = a (a constant) (2) 0 0 Eliminating C1 from Equations (1) and (2), we have: dC2 - V2C2 V 2( KA(V (-02) (3) 2 dt Equation (3) can be integrated between concentration C* at some reference time t*, and concentration C at time t to give: In (3 - C) = In (3 - C*) - vK(At) (4) 46

where: V1C + V C B = + V1 V + V2, constant V1 + V v = ~VV, constant At =(t - t*) Equation (4) shows that the plot of In (3 - C) vs. time will be a straight line with a slope of vK, giving the value of K, since v is known from initial conditions. C. Results Obtained in the Non-flow System The purpose of these experiments was to determine the difference in mass transfer coefficient with ultrasonics and without ultrasonics at the same constant stirrer speed. This was necessary in order to determine the magnitude of the effect of ultrasonics on the boundary layer and membrane resistances. At low stirrer speeds, the boundary layer resistance is expected to be predominant, whereas at high stirrer speeds the membrane resistance to mass transfer would be controlling. We may note that the membrane resistance will essentially be unaffected by the stirrer speed which can only be effective in reducing the boundary layer thickness. The results (in spite of difficulties as mentioned in Section D) had the following trends: (1) there was an increase of 15-200% in the overall transfer coefficient; (2) this increase was high at low stirrer speeds and low at high stirrer speeds; so much so that at a stirrer speed of over 1300 rpm the increase was less than 15%; (3) these results indicate that ultrasonics is more effective in reducing the boundary resistance than the membrane resistance. 47

From conclusion (3) above, it is obvious that in systems using low velocities near the membrane, resulting in a larger boundary layer resistance, the application, of ultrasonics can be very effective in increasing the mass transfer rate. D. It must be noted that the use of the stirrers, and the occurrence of entrained bubbles in the system, interfered somewhat with the ultrasonic field, thereby reducing its effectiveness. Also, although efforts were made to control the temperature during each run, isothermal operation was difficult to achieve and maintain. Thus, the data at very high stirring speeds or high intensity ultrasonics was difficult to reproduce. Nevertheless, the trends outlined above are clearly discernable. It was felt that many of these difficulties could be circumvented in switching to a steady flow membrane transport system. III. Experiments Utilizing a Continuous Flow System A. Apparatus The apparatus depicted in Figure 1 was used to study steady-state mass transport across a membrane. The diffusion cell designed to handle either co-current or counter-current flow is shown in Figure 2. The cell was constructed in a sandwich manner of transparent plexiglass plates, 6" high, 7/8" thick, and 24" long, which are separated in the middle by two 1/16" thick metal sheets. A rectangular channel, 1" high, 1/16" wide, and 18" long, was milled in the interior portion of each plexiglass plate. At both ends the channels were beveled to eliminate entrance effects, and holes were drilled to accommodate both feed and discharge of material. Thus, the 48

~-__.!- - 11 CONSTANT HEAD - ~TANKS SALT SOLUTION WATER ROTAMETERS:I: I I I MANOMETER -VALVES TO GENERATOR~ ~ TRANSDUCER ___ / MEMBRANE -1 —--—. — ";= —o - a-.... --.._....._..... Side View \ COLLECTING VESSELS _ Figure 1. Apparatus for study of steady-state mass transport across a membrane.

TO GENERATOR TRANSDUCER BRASS PLATES PLEXIGLASS -MEMBRANE \ EXPLODED VIEW Figure 2. Membrane diffusion cell. 50

metal sheets which were placed between the two plexiglass plates separated the high and low concentration streams. In the middle of each metal plate, a 1" diameter hole was bored to allow diffusion. Amicon Corp. Diaflow ultrafiltration membranes (type XM-100,PM-30) were fitted and held between each plate by bolting together the entire cell. The bottom plexiglass plate was designed to accommodate a Macrosonics 800 kHz stainless steel transducer. The radiating diaphragm was placed in a direct line with the membrane. Therefore, ultrasonic waves were propagated perpendicular to the cross-sectional area of the membrane, and perpendicular to the direction of channel flow. The total experimental apparatus and its arrangement are shown in Figure 1. Two five gallon tanks located 2.5 feet above the cell were used to provide an approximately constant rate of flow for both the high and low concentration streams. Each stream was fed from the storage tanks through a rotameter and into the cell. A tap was inserted in each stream between the rotameters and the cell and connected to a mercury manometer to measure the pressure difference between the high and low concentration streams. Both streams exited into collection tanks and could be tapped for sample analysis. A Beckman type electrical conductivity cell was used to determine the concentrations of each stream. The transducer was driven by a power supply from a Macrosonics 500-1 ultrasonic generator capable of delivering 0-500 watts of power. The power delivered to the transducer was measured by a Tectronix oscilloscope. Voltage 51

and current probes were attached to the line between the generator and the transducer and fed into the oscilloscope. Thus the oscilloscope displayed the voltage and current amplitudes and the phase angle between the voltage and current. B. Procedure Before a series of runs, a 23% weight salt solution was prepared from Merck reagent grade salt and distilled water,,and inserted into the constant head tank. The low concentration constant head storage tank was filled to the same liquid level with distilled water. Flow was initiated by opening the valves on the rotameters. Steady state was assumed to have been reached when constant flow was observed. Samples were taken from the discharge of the diffusion cell and their concentrations determined with the conductivity bridge. The volumetric flow rates were noted and the temperature was taken from the solutions in the collection tanks. Runs with the ultrasonics on both the high concentration side and the low concentration side were performed. During these runs, the voltage and current amplitudes, the wave length, and the phase angle were read off the oscilloscope display. C. Mathematical Analysis For the well developed steady state flow of two streams exchanging solute across a permeable membrane of constant width, the diffusional flux can be defined as dJ= k. dA.AC, (1) where: dJ = diffusional flux across an element of area dA of the membrane AC = difference in bulk concentrations,(C1-C2) k = overall mass transfer coefficient as defined by Eqn.l. 52

The mass balance for the solute transferred gives V2 2 1 c U 2 Ci C222 dJ = V dCl- V2dC2 (2) where: V1,V2 = flow rates on sides 1 and 2 dC1, dC2 = differential change in bulk concentration during the flow over the area dA, and the sign is plus for counter current and negative for co-current. Integration of Eqn. 2 shows that for constant V1 and V2, the plots of J vs. C and J vs. AC are straight lines; thus, d(C) = AC2-AC1 (3) dJ J a~ to att Uto Eliminating dJ from Eqns. 1 and 3 and integrating over the total area (hence length) for constant k, we get Jo = mass transferred/time = k ATACln (4) ~Cln = log mean concentration difference 53

100 NO ULTRASONICS AMICON XM-IO MEMBRANE - NO ULTRASONICS, AMICON PM-100 MEMBRANE 80 KHZ ULTRASONICS, AMICON XM- IO0 MEMBRANE o -V 800 KHZ ULTRASONICS, AMICON PM-10 MEMBRANE 0 800 KHZ ULTRASONICS, AMICON PM-30 MEMBRANE x 50z Ur UC) 0J50 F I l S 25 1 2 5 10 20 50 100 WATER FLOW RATE, cc/min Figure 3. Mass transfer coefficient vs. water flow rates with and without the application of ultrasonics.

= 2-1 1 n (AC2/tCl) Thus, knowing JO andCln1, which can be computed by noting the solute balance on either of the two streams, and knowing AT, the value of the mass transfer coefficient can be determined. D. Results Obtained in the Continuous Flow Systems A typical set of data taken with the constant flow membrane cell is shown in Figure 3. In these runs, the flow rate of the stream whose concentration was high in salt content("salt side" of the cell) was held constant and the flow rate of the distilled water stream ("water side" of the cell) was varied from 4 - 30 cc/min. Under these conditions, flow on both sides of the membrane was laminar. However, the salt side concentration and flow rate was arranged such that the principal resistances to mass transfer occurred in the membrane and on the water side. If under these conditions, the membrane resistance is also much smaller than the water side resistance, a log-log plot of the mass transfer coefficient vs. flow rate should be linear with a slope of 1/3. The fact that the mass transfer coefficient for the runs taken with the PM-30 membrane is always lower than that fot the XM-100 membrane at the same conditions indicate that there was membrane resistance to mass transfer, and therefore the lines drawn through the data points in Figure 3 are first approximations.For the data shown, the water side of the cell was subjected to ultrasonics. During the runs carried out with ultrasonic irradiation, the power input to the transducer was held constant, with the exception of the data taken using the PM-30 membrane, for which the power level varied due to some fluctuations in the power output of the generator. 55

It is seen from the data using the XM-100 membrane that the application of ultrasonics increases the mass transfer coefficient between 100-300% for all conditions in this study; for the PM-30 membrane, this increase is of the order of 350%. The fact that the increase differed with the type of membrane used suggests that ultrasonics not only favorably affects transport of solute from a solid surface to the bulk, b-ut also perhaps transport through the membrane itself. This latter result is of special significance since there is no mechanical means of enhancing mass transport through a membrane. It should also be noted that,at the power levels used thus far in this study, no membrane destruction has been observed when ultrasonics is employed. Further work is anticipated with the constant flow membrane cell. Efforts will be made to delineate further the effects ultrasonics has on the three principal resistances to mass transfer. The recently initiated study of the effects of ultrasonic power level applied under constant flow conditions will be continued with some emphasis given to examination of the properties of the membrane before and after insonation. Ultrasonic irradiation will be applied to the side of the cell in which the least resistance to mass transfer is encountered to ascertain whether the ultrasonic waves propagate through the membrane. Experiments with and without ultrasonics will be carried out at flow rates in which the primary resistance to mass transport is in the membrane itself. 56

SECTION IV Acoustically Induced Facilitated Diffusional Transport in Membrane Ducts This paper has been submitted for presentation at the 68th National A.I.Ch.E. meeting in August, 1970. 57

ACOUSTICALLY INDUCED FACILITATED DIFFUSIONAL TRANSPORT IN MEMBRANE DUCTS Scott Fogler and Kasper Lund Division of Chemical Engineering University of Michigan Ann Arbor, Michigan 48104 59

ABSTRACT The present work considers enhancement of mass transport in membrane ducts by superimposing a convective transport induced by ultrasonics on a diffusional transport. As a result of the non-linearities in the Navier-Stokes equation, time independent secondary flow will be produced when an acoustic wave is passed through a medium. This secondary flow which is in the form of vortices is commonly called acoustic streaming. Between adjacent vortices or cells molecular diffusion is the only means of transport; however, within each cell, mass transport is primarily by convection. A frequency range of 5-800 kcps in liquids and gases was investigated for different values of acoustic and physical variables. Initial results show increases in the rate of mass transfer of the order of 150% above the normal diffusional flux. Ultrasonics may be applied to increase mass transfer through membrane systems (i.e. dialysis) and to increase the efficiency of very active catalytic systems. 60

I. INTRODUCTION One is most always looking for new methods of increasing the rates of heat and mass transfer in various transport operations. In particular, we are interested in this study in the enhancement in mass transfer through membrane type ducts. The use of ultrasonic waves as a means to achieve this goal has been receiving ever increasing attention in recent years. Sonic and ultrasonic waves have been found to increase heat and mass transfer rates by several orders of magnitude in a number of situations (Fogler (1966), Richardson(1967a), Boucher (1959), Nyborg (1965)). The advantageous effect of acoustic waves on various transport phenomena and chemical reactions has been discussed in some detail by Fogler (1967). In the present work attention will be focused on the enhancement of mass transport in gas and liquid diffusion resulting from acoustic streaming. Various studies have been conducted on the effects of acoustic- and micro-streaming on heat transfer. The majority of these studies have emphasized streaming on the exterior of vibrating objects such as streaming around a cylinder or sphere. In certain limited situations heat transfer studies have been conducted on streaming inside ducts and tubes. Complete reviews of acoustic streaming have been given by Nyborg (1965) and by Richardson (1967b). As previously mentioned, in this study we shall consider cases in which the convective streaming transport is superimposed on the diffusive transport in diffusion in ducts and tubes in 61

membrane pores. The streaming cells induced by acoustic waves should greatly enhance the rate mass transport through the system. II. THEORY Acoustic streaming is a secondary flow which produces time independent vortices when an acoustic wave is passed through the medium. The formation of these vortices or cells inside ducts, tubes and pores can increase the rate of mass transfer through the enclosures. The figure below shows a sketch of these cells in a wide membrane duct. Between adjacent cells molecular diffusion is the only means of mass transport; however, within each cell transport is primarily by convection. Before engaging in discussion of the coupled transport processes, it would be beneficial to give a brief development of the streaming equations. Starting from the continuity and NavierStokes equation + V- (pu) = 0 (1) ^(pu) + p( -V)u + uV- (pu) = -Vp + [u' + u] VV-u - iVxVxu (2) at' J3 Since there is no external force other than the sound field, the solutions to these equations will take the form 62

U = U1 + ~ U 2 P-P = P1 + P2 (3) 2 2 Po l= E + p P2 u1 and wl are time dependent velocities which directly arise from the sinusoidal movement of the transducer surface. The streaming velocities u2 and w2 consist of a time independent and a time dependent term. In this investigation we are only interested in the time independent part of the streaming velocities. The first order velocities give rise to increasing momentum with time in the system; dynamic balance is restored by viscous stresses from the second order velocities and second order pressure gradient. Upon substitution of equation (3) into equation (1) the first order approximation becomes a1 2 4 Po at = c - RVpl + (p + 3)VV-ul - PVxVxu1 (4) and the time average of the second order velocity u2 is p <(u1 V)u1 + u1V-u1> =-VP2 - pVxVxu2 (5) It is precisely the second order velocity u2 which gives rise to the circular streaming cells. A. Solutions to Streaming Equations for Transport in Membrane Ducts The solution to Equation (4) in terms of first order velocity approximation was given by Rayleigh (1895, 1945) to be u = A cos kx[cos wt - eZ cos(wt - Sz)] (6) V1 = -Ak sin(kx) (cos(wt - i) - e zcos(wt - z -) (7) 63

Upon taking the curl of both sides of Equation (5), one obtains pVx (V u2) =-VxF (8) where -F = pO<(ulV)u1 +1(V-ul)> (9) With Vu2 approximately zero, one can substitute for u2 in terms of the stream function. V2 (V2i) = VF (10) Upon substituting Equations 6, 7 and 9 into Equation 10, one obtains V42 = - 1- [Bp kA2sin 2kx] [2c + S-e 2n (11) V2 = dz The solution to this equation in terms of the x and z components of the streaming velocity U2 = 3A sin 2kx[e2Bz + 2S-1 + 6nl(1-nl)] (12) (2) 3kA2 -2az W2 = - 8c- cos 2kx[e2z + 2(S+C)-3 + 2Shnl(l-n1) (1-2nl)] (13) B. Coupling of the Streaming and Diffusional Transport Equations If the fluid in the duct consists of a mixture of species and their concentrations are different at the two ends of the duct, mass will be transferred through the duct partly by diffusion and partly by convection. For dilute solutions, the mass transfer equation can be coupled with the acoustical induced velocity field. The unsteady mass balance is + V - VC = DV2C 6k

with the boundary conditions C = C0 at X = O C = C at X = L aC = 0 at Z = 0 and Z = h/2 aZ The velocity, V, can be taken as the mass average velocity and D as the binary diffusion coefficient for a dilute solution. Further, the length of the duct L is taken as an integer multiple of r/2k. We now substitute for V to obtain DC ac DCc a2C ~~2C ~C + (U1 + U2) - + (W 1 + W)- = at- 1 2 (Z1 2 D( 2 ax az Time averaging the equation over several cycles and neglecting terms like ac <U DC1 > and <W Z > the equation simplifies to <C> = C aC aC ~ZC a 2 U aX + W D( + Method of solution The partial derivatives can be approximated by the following finite difference formulas: 2 C. - 2C + C. 2 C - 2C + C aC C ij- 2Cij Ci,j+l C _ i-l,j, j + Ci+lj ax (AX) a2 (AZ) C C a_ Ci,j+l i,j ac Ci+l,j Ci-,j ax 2AX Z 2AZ Let AX = AZ Boundary conditions: C.i = C at X = 0 and C. = C at X = L 1,1 o zi,m+l 1 65

C2, = Clj at Z 0 and Cn+lj = Cnj at Z = h/2 2,j lj n+1,3 n,3 i = 1,2,3,.o..., m+l and j = 1,2,3.....,n+l Substituting into the differential equation and rearranging yields Ci = C. + Q[(l+U. * AX/2) Ci + (1-Ui AX/2)C, i3j 4 ij i,j +l ijij-1 + (l+Wi. AX/2)Ci+ j + (l-W i- AX/2C j -1 4Cij] ~irJ-il~3j i-l3j where Q is an accelerating factor 1<Q<2. The system of linear equations is solved by a Gauss-Seidel iteration scheme (Carnahan, et al. 1969). III. RESULTS The dimensionless concentration field is plotted in Figure 1 for the case of a frequency, f, of 20 kcps, a maximum displacement -3 of the transducer surface, s, of 7 x 10 cm and a duct height, h, of 2 cm. The increase in mass transfer for this case is 96% and the plot clearly shows that the acoustic streaming pattern strongly modifies the concentration field which would be present when only diffusional mass transfer was occurring. Near the two ends of the duct, large concentration gradients are present and diffusional mass transport will be the most important means of transport at these points. In the middle of the duct concentration gradients are smaller and convection will be the primary means of transport. In Table 1 a few representative results are given. In Figure 2 the increase in mass transfer is plotted versus the Peclet number, 22 transport by convection 3kh 2A2 Pe. Pe transport by diffusion Pe 8D c where c is the speed of sound, k the wave number and A the first order velocity amplitude. 66

TABLE 1 Fluid System Frequency Number of % Increase kcps''vortex' cells''. N2-O2 20 1 8.1 N2-02 20 2 3.8 N2 -2 20 3 5.2 N2-O2 20 1 96. 02 in water 20 1 59. 02 in water 100 1 145. 02 in water 800 1 59. It is seen that at each frequency the Peclet number has to reach a certain magnitude before any significant increase in mass transfer can be expected. The plot also indicates that much higher increases in mass transfer can be expected but this could not be investigated because of stability problems with the numerical scheme. The present numerical technique is being modified so as to extend our results into other regions. The results show that by the application of ultrasonics to mass transport processes, substantial increases in the rate of mass transfer can be expected. Application of ultrasonics will be very beneficial for mass transfer through membranes when a substantial resistance to mass transfer lies in fluid boundary layers and for solid-fluid reactions (catalytic or non-catalytic) where the rate of mass transfer is limiting. Discussion Our initial investigation using the accepted ultrasonic streaming equations shows increases of up to 145% in mass transfer can be 67

.50 o~~.43~0 t t ~ X kx Figure 1.. Concentration profiles. Co - C Figure. Concentration profiles.68 68

1000 w co0m N2 -02 100 z / I Cell i-. /; ^ 5 kcps o L ap 10 kcps A 10W/~ ~20 kcps o 10 40 kcps 0 IL I I I I I 20 40 60 80 100 % INCREASE Figure 2. Increase in mass transfer as a function of the Peclet number. 69

expected when ultrasonic convective transport is superimposed on a diffusional transport. In the derivation of the ultrasonic streaming equations Nyborg made certain approximations which limit the region of validity of the equations; the results from the initial investigation will, because of these approximations, only be semiquantitative, but we would expect that the trends and the magnitude of the mass transfer increase will be unaffected. To get accurate results predicting the increase in mass transfer as a function of physical properties of the fluid as well as the acoustic variables, it will be necessary to derive more general solutions of ultrasonic streaming equations without making some of the very limiting approximations made by Nyborg. This theoretical work is currently being done. The new streaming equations will then be coupled with the mass transport equation and solved numerically. IV. SUMMARY Upon the passage of an acoustic wave through a duct, one can observe the formation of time independent vortices. This phenomena, which is commonly called acoustic streaming, can enhance mass transport by superimposing a convective transport on a diffusive transport. The differential mass transport equation was coupled with the second order time independent streaming equations in a rectangular membrane duct and solved by finite difference techniques. A frequency range of 5-800 kcps in liquids and gases was investigated for different acoustic and physical variables. The acoustic streaming strongly modifies the concentration field which would be present 70

when only diffusional mass transfer takes place. Regions with large concentration gradients are formed near the two ends of the duct and only small concentration gradients are present in the middle of the duct. At each frequency the Peclet number, and thereby the first order velocity amplitude, has to reach a certain magnitude before any significant increase in mass transfer can be expected. At high Peclet numbers the increase in mass transfer is proportional to the logarithm of the Peclet number. Preliminary results show that with the application of ultrasonics increases are up to 150% above the normal diffusive transport. V. BIBLIOGRAPHY Boucher, R.M. G., "Ultrasonics in Chemical Processing," Chem. Engr. 66 (Sept. 1959). Carnahan, B., Luther, H.A., and Wilkes, J.O., "Applied Numerical Methods," Wiley, New York (1969). Fogler, H.S. and Timmerhaus, K.D., "Effect of Ultrasonic Waves on Mass Transfer Rates of Selected Fluids," AIChE 12, 90 (1966). Fogler, H.S., "Applications and Research in Sonochemical Engineering," Sound and Vibration 1, No. 8, 18 (1967). Nyborg, W.L. in Physical Acoustic Vol. II B, pp 265-331. W.P. Mason, ed., Academic Press, New York (1965). Richardson, P.D., "Heat Transfer from a Circular Cylinder by Acoustic Streaming," J. Fluid Mechanics 30, 337 (1967a) Richardson, P.D., "Effects of Sound and Vibration on Heat Transfer," Appl. Mech. Rev. 20, No. 3, 201 (1967b). 71

SECTION V Acoustic Cavitation in Viscoelastic Fluids This paper by Prof. H. S. Fogler and Prof. J. D. Goddard was presented at the 62nd National A.I.Ch.E. meeting in Washington, D. C., and has been accepted for publication in the Physics of Fluids. T3

ACOUSTIC CAVITATION IN VISCOELASTIC FLUIDS As previously mentioned the overall objective of this research is to utilize ultrasonic waves to accelerate kinetic and mass transport phenomena in biological type systems. Preliminary results, both theoretical and experimental, have shown that substantial increases in mass transfer can be brought about by ultrasonics. With the application of acoustic waves to biological systems, one must be sure that the wave conditions are properly adjusted so as not to induce any harmful side effects to the tissue and body fluids. In some cases, acoustic cavitation could degrade these materials if it is induced beyond a certain level. It is known that various biological fluids behave as viscoelastic liquids [Trans. of the Soc. Rheology 9, Part 1, p. 448 (1965)]. Consequently, a brief and preliminary investigation was undertaken on cavitation in viscoelastic fluids to determine whether the degrading effects of cavitation could be accelerated or retarded in this fluid type. As a result of the complexity of the problem, only a few limiting cases were studied. These cases were chosen such that if they showed the cavitation process was unaffected by the elastic effects in the liquid, then the other situations would in all probability show the same results. The paper which follows was presented at the 62nd Annual A.I.Ch.E. meeting and has also been accepted for publication in the Physics of Fluids and concerns the collapse of a spherical void in a viscoelastic fluid. These preliminary results show that the elasticity can significantly retard the collapse process and certain situations produce damped oscillatory motion of the cavity rather than the usual catastrophic collapse observed in purely viscous liquids. Since first results show viscoelasticity could 75

-2quite significantly retard material and fluid degradation, we have outlined a program for further study in this area in our recent proposal to NIH for continued support on this project. 76

ON THE STABILITY AND COLLAPSE OF SPHERICAL CAVITIES IN VISCOELASTIC FLUIDS H.S. Fogler and J.D. Goddard University of Michigan, Ann Arbor, Michigan Preprint 63A Presented at the Session on Fundamental Research in Fluid Mechanics - Part II SIXTY-SECOND ANNUAL MEETING Washington, D. C. November 16-20, 1969 AMERICAN INSTITUTE OF CHEMICAL ENGINEERS 345 East 47 Street, New York, New York 10017 75 cents 77

ABSTRACT An analysis is given of the collapse of a spherical cavity in a large body of an incompressible viscoelastic liquid. Proceeding from a linear rheological model for the liquid, one obtains a non-linear integro-differential equation for the motion of the cavity. Analytical solutions are derived for certain limiting values of the parameters governing collapse, and some numerical solutions are presented for various other values. As one of the more interesting results of this work, it is found that elasticity in the liquid can significantly retard the collapse of a void and produce prolonged, oscillatory motion whenever the relaxation time of the fluid is moderately large in comparison to the Rayleigh collapse time. This is in sharp contrast to the catastrophic collapse which will always occur for voids in purely viscous liquids. Both numerical and approximate analytical solutions are presented to demonstrate the damping effect of liquid viscosity on the cavity motion. 78

1. INTRODUCTION The term cavitation usually refers to the phenomenon of growth and collapse of flow-induced voids or vapor bubbles in liquids. The effects resulting from cavitation are known to produce metal erosion, luminescence, and increases in various chemical reaction rates. In the previous works on this subject, attention has been restricted mainly to classical liquids. The earliest theoretical treatment is apparently that of Lord Rayleigh ( 1 ), who considered the collapse of a spherical void in an inviscid liquid. In later theoretical works, attempts have been made to account for viscous effects in both the bubble phase and in the surrounding liquid, and most of the analyses have dealt with Newtonian (Flynn, 2; Plesset, 3; Fogler, 4 ) or purely viscous fluids (Yang, 5 ). An interesting question arises as to the effects that elasticity might have on cavitation in viscoelastic liquids. In other contexts, it has been observed that the presence of elasticity, such as that produced by addition of small amounts of high polymers, can drastically change the flow behavior of liquids. Hence, one might well inquire as to the possible and, perhaps, beneficial effects of viscoelasticity on bubble collapse, such as suppression or reduction in the intensity of cavitation. An analysis of bubble growth in viscoelastic fluids has already been given by Street ( 6 ), but, because of the applications contemplated in his analysis, inertial effects were neglected. It is precisely these effects, 79

however, that tend to predominate in the collapse phenomena usually associated with cavitation. This provides part of the motivation for the present work, in which we shall focus our attention primarily on the collapse of spherical voids, (i.e., regions containing no gas) in an idealized viscoelastic fluid. We recall that previous studies have shown that collapsing cavities which contain permanent gases will generally always rebound short of actual collapse, such that the cavity radius never actually decreases to zero. On the other hand, a void will generally always collapse to zero radius at least in purely viscous fluids. It is therefore interesting to reconsider this question of rebound versus complete collapse for the case of a void in a viscoelastic fluid. 2. EQUATIONS OF MOTION We wish to treat here the motion of a spherical bubble contained in a large body of an incompressible liquid. Initially at time t = 0 the system is at rest, with a bubble radius R and a uniform pressure Po. It has been previously (Flynn, 2; Plesset, 3 ) shown that the equation for the spherically symmetric motion for a bubble,in which there is no condensation or evaporation of fluid, can be reduced to oo 8 P 0 Po 1 r (1 80

divergence of the deviatoric or "extra" stress for the liquid phase, and r = R(t) is the radial position of the bubbleliquid interface, with P and Po denoting the pressure in the liquid at r = R(t) and r = o, respectively. The dots denote derivatives with respect to t, and p is the liquid density. Irrespective of the fluid rheology, the radial velocity at any radial position r in the liquid is required by continuity, incompressibility, and the assumed symmetry to be *2 RR u = -2 (2) r By the usual force balance at the bubble-liquid interface in the cavity, the term P in (1) can be expressed in terms of surface tension and the radial stresses as 2c (3) T +P - P +T + ) rr,g g Trr,, R where g refers to any gas which may be present in the cavity, and Z refers to the liquid phase. Since neither surface elasticity nor viscosity are considered in this analysis, the surface tension force is given by the static surface tension a. * As in Fogler ( 4 ) we adopt here the sign convention of Bird, Stewart, and Lightfoot for the stress tensor: The symbol T (or T.i ) denotes the deviatoric stress tensor reckoned as a comp essive stress (7). 81

For bubbles containing an ideal gas in a uniform state the gas-phase stress at the bubble surface is equal to the pressure alone; hence Trr = 0 (Fogler, 4 ), and rr,g the liquid-phase interfacial pressure is given by P=P- -T " (4) P Pg R Trr, 2 g R r = R Furthermore, the term (V * T) which occurs in (1) can be written in terms of three normal stresses as (V T)3 2 (T + (T) (. T) rr rr r' ) and, since the sum of these deviatoric stresses is by definition zero, one can express the ~ and 9 stresses in terms of the radial stress as + T = - Trr (6 which with (5) yields aT T (V T) = rr rr (7) V r ar r Then, upon substituting (4) and (7) into (1), one obtains the equation 00 l +3 R2 =Pg Po 3 T rrdr R ~ R + R: R fO? the bubble radius R(t). In order to complete the description of motion, we must now relate the liquid-phase radial stress T to the bubble motion. r r8 82

In the case of a general viscoelastic fluid exhibiting long-range memory effects, the stresses will depend on the past history of strain or rate of strain. For the simple, radially symmetric flow field considered in.Eqn. (2), the strain consists merely of an unsteady simple extension. Hence, we expect that for an isotropic material the instantaneous radial stress T (t) can be expressed as a functional on the rr past history of the radial strain rate e (t'), O t' _ t. rr Here, as in the following analysis, t = 0 corresponds to the beginning of the collapse process, where we assume the liquid to be in a completely "relaxed" state of purely hydrostatic stress. As with other analyses involving viscoelastic fluids, we must now postulate a relation between the strain and the kinematic history of the motion to be considered, and, for this purpose, we adopt the usual material coordinates. Thus, we let r' denote the position at past time t', 0 _ t' t, of a particle which is at position r at the present time t, so that,with the velocity field given by Eqn. (2), we have (r')3 = r3 + R3(t') - R3(t). (9) Now, at any position (r,t) the radial deformation rate is given by Du 2RR err (r,t) R (10) rr - - 310 and, therefore by (9) and (10), the history of the deformation rate is determined by 83

e(t',r') = 2R(t') R2(t') e (11).)rr \3 3 3 r + R (t') - R (t) For the present work, we shall employ a rather simple, linear viscoelastic fluid model, in which the normal radial stress is related to the corresponding strain rate by t T (t)= -2 N(t-t') er(t')dt' (12) 0 where N(t) is a "memory" function or relaxation modulus. On combining (11) and (12) we have t T 4r N(t-t') R(t') R2(t') dt', ( rr r + R(t) - R3(t) 0 and the integral in equation (8) becomes oo t oo r 3 3 3 5 T dr =4 L N(t-t') R(t') R2(tt)dt'dr r(r + R (t') - R3(t)) R R t N(t-t') R(t') R (t')ln(R(t')/R(t))dt' R(t') - R(t) o Under these restrictions the complete equation governing the collapse of a cavity is the non-linear integro-differential equation: R * R + 3 R2 = -_ _~ - 2aR R + 2 p Rp 12 N(t-t') R(t') R (t')ln(R(t')/R(t))dt' (15) R3(t') - R3(t) ~o 84

For the purpose of the analysis to follow we shall adopt an elementary form of the relaxation modulus N, consisting of linearly viscous or "Newtonian" contribution and "Maxwellian" contribution, as follows: N(t) = 6 (t) + Get/, (16) where 6 denotes the delta function, p a constant viscosity, X a relaxation time, and G an elastic modulus. In terms of dimensionless variables, (15) becomes 7 32 g o 2 P-0Po NWe NRe t* 12N EL ( (t*-tl ) 11 ln(- /) -N N' exp( N n~ /3 3 dtl (17) Re NDe 3 3 0 with p(0) = 1 and p(0) = 0, where NDe = t (A Deborah number ), G t NEL = c (An elastic number), EL P-W PR2 NRe = t- (A Reynolds number), Re Pt PO pR3 NWe = _ (A Weber number), t do and R = R/Ro, t* = t/tc, C1 = (tl Also, t = R 4 is a characteristic ("Rayleigh") collapse c o-v r0 time, with P being the initial pressure. In this manner one can readily identify the relevant physical parameters characterizing the collapse process. In view of the number of parameters, even in this relatively simple model, one is practically forced to consider some special limiting cases where certain effects may be assumed * Reiner ( 8 ), Metzner et al. ( 9 ) 85

to predominate. Thus, we focus our attention first and foremost on fluids with long relaxation times, corresponding to NDe+ o. Here, as in the remainder of the analysis, we shall consider only voids,such that P = 0 in (17). 3. COLLAPSE CRITERIA AT LARGE DEBORAH NUMBERS 3.1 Large Reynolds Number To begin with, we treat the case where both the Deborah and Reynolds numbers are large. In this limit, NDe+, NRe+, the fluid behaves essentially as a purely elastic material, and one obtains in effect a "conservative" dynamical process characterized by an energy integral. Considering first the case where surface tension is negligible, N We, and reverting to dimensional variables, one has for the equation of motion, P R R + 3 R2 = G (18) 2 p P where R 1 l in (R1/R) dR ~~~~~~~~~G =12G. 1= G s 1 GGds 1Ro () On multiplying equation (18) by 2R p dR and integrating, we obtain R p 23 =2 R3 - R3) - 2 GR2 dR. (19) R - Po 2o R One will immediately recognize that this equation is an energy integral, with the left-hand side representing 86

the total kinetic energy of the liquid, which is expressed as the difference between the stored elastic energy and the work done by the ambient pressure. Rebound short of collapse is therefore possible and will occur at a "rebound" radius R, which is the root of the equation R Po(Ro - R) 2 GR2 dR = 0, (20) R 0 corresponding to zero kinetic energy in (19). With the substitution into (18) y = l/s, z = (R/Ro)3 the integral in (20) can be written as R 3 z x H= 2~ GR dR 8R 9 G s(y dy dx, (21) R 1 1 0 which, after changing the order of integration, can be expressed as the infinite series Thus (20) becomes 0 4 I y (lZ)n (In (l/z) 2 GI 3 | n2 Z 2 (-z) (23) 0 n=l n which provides the criterion for rebound, giving the rebound radius R = R E Rozl/3 as a function of P /GO 87

In the "marginal" case, that is, rebound at R = 0, we have co 3 2 = 2.1932... (24) n=l from (23), and therefore the condition for collapse without rebound is P 2 G > 29 (25) G 9 0 whereas, for rebound short of collapse, R* > 0, we must have P 2 0 < 2 (26) G 9 0 In the latter case, equation (23) provides us with a plot of rebound radius R*/R versus the ratio of initial pressure to elastic modulus, Po/Go, which is displayed as the lower curve in Figure 1. If we consider the case of a finite Weber number, where surface tension is included in the equation of motion, the criterion for collapse is no longer independent of the initial bubble radius. By an analysis similar to that above, one can show that the condition now becomes P + o R 2 o 2 7 0G > (27) 0 instead of (25). From this relation one sees that surface tension effects will tend to be important only in small bubbles. 88

We should like next to determine the importance of viscous retardation on the collapse process, corresponding to a finite Reynolds number in (17). Whenever NRe and NDe are finite, the system is no longer conservative and, hence, does not in general admit an energy integral like (19). We are thus forced to treat (17) with numerical techniques, as will be discussed below. First, however, it is worthwhile to note that for infinite Deborah numbers a cavity is characterized by a certain "equilibrium" radius Req, as determined by the static balance between pressure, surface tension, and elastic forces. One can easily derive an expression for this radius, and, considering the case of negligible surface tension NWe = o, one finds from (20) that the We condition of static equilibrium is 00 n=l to a static situation. For the purposes of obtaining the numerical solutions, a finite-difference technique was employed to treat eqn. (17). In particular, a modified Milne "four-point predictor" formula was used, and the numerical solutions thus obtained were compared for accuracy with existing numerical solutions for bubble collapse in ordinary liquids(Flynn, 2 Fogler 4 ). In all cases, the solutions were the same. 89

3.2 Finite Reynolds Number - Viscous Damping Figure 2 gives a plot of cavity radius versus time for infinite NDe and NWe. The cavity is seen to oscillate about an equilibrium radius subject to viscous damping which increases with NRe. Under these circumstances, one would expect to observe "critical" damping below some threshold value, NRe, say, which we shall refer to here as the "critic cal" number. To obtain an estimate of this number, we shall make use of the techniques of linear stability theory. Thus, letting 4e represent the dimensionless equilibrium radius and 4 a small perturbation about this radius, we have = e + 4' 4 << We (29) Then substituting Equation (29) into the equation of motion (17) and neglecting terms of the second order in 4, we obtain the corresponding linearized equation for a collapsing void, which in the case of infinite Weber numbers becomes e - i + 4.. + 12 [G(*e) + at- ] = 0 (30) e NRe~e e Since 12 G(*e ) = 1 by Equation (18), the preceding equation becomes 3l(e +r 44G n~... A N 4*e) 4o -e)1, (31) NRe 2 - Po (e) 2(1 - e or, simply, 4' + b4' + c = 0, (32) 90

where b and c are constants. In the usual way, it can be seen that the oscillation of the cavity will be critically 2 damped whenever b = 4c. With the appropriate values of these constants from Equation (31), this criterion becomes 2 No e N N (33) NRe = Re G- 2 3 /3 c o e ln(l/~ ) e e which on rearrangement and making use of the definition in (17) becomes 3 1 - NRe NEL =2 ~ 3 (34) c ie ln(l/e) Since the equilibrium radius corresponding to e is determined by Po/Go, we may express the critical Reynolds number as given by Equation (34) in terms of Po/G or, alternatively, in terms of e'. In the latter case one obtains a plot of the critical Reynolds number as a function of the equilibrium radius as shown in Figure 3. A physical interpretation can be given to the shape of the curve in the following way. Near *e = 1 where the elastic force is, relatively speaking, not very large, a greater viscous force is required to damp oscillations as the cavity approaches its equilibrium radius. However, when ie is only slightly less than unity (e.g., We =.7 as in e e Figure 2), the elastic force, which increases rapidly in a non-linear way, exerts a greater degree of retardation on the motion, and consequently a smaller viscous force is necessary for critical damping. 91

Owing to the method of derivation, the present expression for the critical Reynolds number, at which cavities move from their initial radius to their equilibrium radius on a critically damped path, can be regarded as strictly valid only for cavities in which We is close to unity. For cavities with equilibrium radii close to zero, the departure from equilibrium p at the initial state = 1 is effectively much greater than the equilibrium ratio We' and hence the above linearization technique cannot provide an adequate description of the cavity motion from i = 1 to 4 = We One notes, however, that for an equilibrium radius ratio of.74, the critical Reynolds number obtained from Figure 3 is 1.25, and from Figure 2 it is observed that for this value of the Reynolds number the cavity does indeed approach equilibrium in a critically damped way. Thus, the linearization is evidently valid in this range. In general, one might also use Figure 3 to determine the Reynolds number at which critical damping occurs when the equilibrium radius is shifted from some value We to a second value, e, say, by a change in the total pressure. A In this case, one would use ie in Figure 3 to find NRe c This prediction of the damped shift would be valid irrespecA tive of the magnitude of e', provided only thatle - Wel << ~e For large but finite Deborah numbers, Figure 4 shows the numerically computed motion of the cavity. One observes complete collapse for Po/G = 100, with a collapse time very nearly equal to the Rayleigh collapse time for an inviscid, non-elastic liquid. Furthermore, it is evident that for the 92

case NDe = 1000 shown there, the motion on the first few cycles is effectively the same as for NDe = a. Also, it can be observed that the Reynolds number has a significant effect on the initial motion only when it is numerically on the order of magnitude of ten or less. Because of the greater energy dissipation at the lower Reynolds numbers, it appears that the rebound radius decreases with the increasing Reynolds numbers. 4. COLLAPSE AT SMALL DEBORAH NUMBERS While it is evident that for any finite Deborah number a void must eventually collapse to zero radius, it is nonetheless of interest to investigate how collapse is delayed by the elasticity of the fluid. In particular, we may consider the first cycle of motion, as in Figure 5. For a given Po/Go, the rebound radius on the first cycle decreases with decreasing Deborah number as shown there. If the fluid is'lnviscid" (Ne = o) the "critical" Deborah number at which Re the cavity collapses completely on the first cycle is.51 for a Po/Go ratio of 1.43, whereas for a finite Reynolds number the cavity no longer collapses on the first cycle at NDe =.51, but instead rebounds as shown in the figure. For various cases, the numerical solutions were carried out for several cycles of the motion, and some of the results are shown in Figure 6. One observes in this figure that for a Deborah number of 2/3, the cavity collapses in approximately three major cycles. One also notes that the maximum radius reached after each rebound decreases in an almost linear fashion for the first few oscillations when NDe, 1. The "modulation" within the later cycles 93

and the exact radius values in final stages of collapse are uncertain at this time, since numerical integration difficulties were encountered at long times. (The longest time shown represents some thirty to forty minutes of IBM 360 computation time for a single run). 5. CONCLUSIONS The results of the preceeding analysis indicate that elastic effects may well have a strong influence on cavitation in viscoelastic liquids. We should certainly expect such effects to occur at high Deborah numbers X/tc, where the relaxation time X of the fluid is long compared to the classical Rayleigh collapse time tc. In particular, for the Maxwellian liquid considered here, the present analysis shows that in the limit of large Deborah numbers, X/t+ c, a spherical void may either collapse or undergo oscillations about an equilibrium radius, depending on whether the ratio of ambient pressure to the elastic modulus of the fluid exceeds a definite, critical value. The presence of viscosity in the fluid tends to damp the oscillations, and a critical-damping phenomenon occurs for Reynolds numbers below a certain value. Even for finite and moderate Deborah numbers, X/tc = 0 (1), the ultimate collapse of a void is delayed for several cycles of expansion and contraction. Although we have not considered in detail the possible effects of gases or vapours in the collapsing cavity, we should not expect such effects to greatly alter the role 94

of liquid elasticity in the collapse process. In fact, one might reasonably anticipate that the combined effects of volume elasticity in the gas and shape elasticity in the liquid would reinforce one another in such a way as to retard or completely suppress the collapse of bubbles. From the results of previous studies of gas-filled bubbles in Newtonian fluids, we might also expect that, in many instances, the effects of liquid elasticity would be important at a much earlier stage in the collapse process. In such cases, the build-up of the liquid-phase momentum, which gives rise to catastrophic collapse, would be greatly suppressed. In addition to any experimental work which may be suggested by the present study, it would also be of some interest to investigate theoretically the hydrodynamic stability of the spherically-symmetric motion of cavities collapsing in viscoelastic liquids. While one might be tempted to employ a somewhat more refined theological model for the liquid, this would probably lead to rather difficult analytical and computational problems, without necessarily providing much additional insight on the physics of the collapse phenomenon. 95

ACKNOWLEDGEMENT A portion of this work was supported by NIH Grant HE-10549. 96

err radial strain rate, Sec-1 GO an elastic modulus, dyne/cm G defined by Eqn. 18 H defined by Eqn. 21 NDe a Deborah Number De NEL an Elastic Number NRe a Reynolds Number NWe a Weber Number N(t) memory function, gm/cm/sec n integer 2 P pressure, dyne/cm R radius of the cavity wall, cm. R1 cavity radius at some previous time, cm. RO initial cavity radius, cm. r radial coordinate, cm. u velocity, cm/sec. t time, sec. t modified collapse time, sec. x dummy variable z (R/R0)3 Subscripts c critical g gas L liquid e equilibrium radius 1 refers to a previous time Greek Symbols'p -dimemsionless cavity radius R/RO 2 T.. normal stress, dyne/cm surface tension, dyne/cm a 97

p liquid density, gm/cc X relaxation time, sec. p viscosity, gm/cm/sec. 6 (T) delta function 98

REFERENCES 1. Ra.yleigh, Lord, 1917, Phil. Mag. 34, 94. 2. Flynn, H. G., 1964, Physical Acoustics W. P. Mason Ed., Academic Press, New York. 3. Solomon, L. P. and M. S. Plesset, 1967, Report No. 85-38, Division of Engineering and Applied Science, California Institute of Technology. 4. Fogler, IH. S., 1969, ChemL Engrg. Sci. 24, 1043 (1969) 5. Yang, W. J. and H. C. Yeh, 1966, AIChE J. 12, 927. 6. Street, J. R., 1968, Trans. Soc. Rheol. 12, 103. 7. Bird, R. B., W. E. Stewart and E. N. Lightfoot, 1960, "Transport Phenomena", John Wiley, New York. 8. Reiner, M., 1964, Physics Today 17, 62. 9. Metzner, A. B., J. L. White and M. M. Denn, 1966, AIChE J. 12, 863. 99

Po NDe = D Go L.3 N'o 1 C Nwe = I I I I i 0 0.2 0.4 0.6 0.8 1.0 R4/Re Figure 1. The initial rebound radius and equilibrium radius as a function of P /G. The middle curve was computed from Eq. (17). 100

:/Go = 10/7 NDe = Ny/e = 00 N =1.25 ^Re^ 40 H-i 1.0 H 0.8 0.6 0.4 Re =10 0.2 0 1.0 2.0 3.0 4.0 5.0 6.0 Figure 2. The influence of the Reynolds number on the damping of sably oscillaing caviies. Figure 2. The inf luence of the Reynolds number on the damping of stably oscillating cavitis

12 NDe aO) 60 10\ NWe= 4 e 0 0.2 0.4 0.6 0.8 1.0 Je = Req / Ro Figure 3. The "critical" Reynolds number as a function of the equilibrium radius. 102

NDe = 1000 N = 140 we 1.0 ~G 100 2 | e \ O } AXNOR/Po/ G o 2.8 4 NR e 1 000 1~ PoN/Go 100/ Npe =-100 o}oo.2 0.2.4.6.8 1.0 1.2 1.4 1.6 1.8 t* Figure 4. The effect of the Reynolds number and the P /G ratio on the collapse of a cavity.

P/Go = 10/7 NWe NRe = 1.0 0.6 NDe CD.4 N 51 NDe=.5 De.2, I NDe =.51 NReN =100 De 8NReOO NDe =.55 0.2.4.6.8 1.0 1.2 1.4 t Figure 5. The effect of the Deborah number on the initial motion of a cavity.

NWe =CO NRe =OD P/Go=10/7 0.8 NiDe I. 0 0.6 0.4- -De. 667 0.2 0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 t* Figure 6. The effect of the Deborah number on the cavity motion.

SECTION VI PERSONNEL 107

UNIVERSITY OF MICHIGAN ii11Jii 111111ii lllill1 l 3 9015 02826 3419

PERSONNEL H. S. Fogler is an Assistant Professor in the Division of Chemical Engineering. He has worked on the application of ultrasonic waves to chemical reaction and physical transport systems for the past eight years. M. L. Cadwell is in his last year of study for the Ph.D. degree. His thesis topic is ultrasonic gas absorption. V. K. Verma is a third year graduate student and K. Lund a second year graduate student in the Division of Chemical Engineering. W. Franklin is a senior anticipating graduation in April and graduate school in September. 109