Appl. Sci. Res. 33
January 1978
THE BEHAVIOR OF VISCOELASTIC SQUEEZE FILMS SUBJECT TO NORMAL OSCILLATIONS, I N C L U D I N G THE EFFECT OF FLUID INERTIA JOHN
A. TICHY
Dept. of Mech. Eng., Aeron. Eng. and Mech. Rensselaer Polytechnic Inst. Troy, New York, U.S.A. Abstract
The problem of the squeeze film flow of a viscoelastic fluid between parallel, circular disks is analyzed. The upper disk is subject to small, axial oscillations. Lodge's "rubber-like liquid" is used as the viscoelastic fluid model, and fluid inertia forces are included. An exact solution to the equations of motion is obtained involving in-phase and out-of-phase components of velocity field and load, with respect to the plate velocity. Peculiar resonance phenomena in the load amplitude are exhibited at high Deborah number. At certain combinations of Reynolds number and Deborah number, the in-phase and/or out-of-phase velocity field components may attain an unusual circulating type of motion in which the flow reverses direction across the film. In the low Deborah number limit, and in the low Reynolds number limit, the results of this study reduce to those obtained by other workers. Nomenclature a, b
F f, 9 h K m, n
p, q r R Re
S IJ s t v
parameters, defined by equations (29), (30) load arbitrary functions, equation (2) film thickness memory function arbitrary functions, equation (9) functions, equation (19); velocity field components, equation (23) radial coordinate disk radius Reynolds number deviatoric stress tensor past time = t -- z present time fluid velocity -
501
-
502 V c~,fl )/J, 7,j e 0, 0~ p, #~ ~ p 4~ @ ~u co
JOHN A. TICHY upper disk velocity auxiliary functions, equations (27), (28) metric or strain tensor dimensionless oscillation amplitude relaxation time Deborah number viscosity convected coordinate density time angular coordinate phase angle stream function oscillation frequency
Subscripts mean value radial, angular, axial component sine (in-phase) component; cosine (out-of-phase) component convected coordinate
m r, ¢, z s, c o
Superscripts * amplitude ' derivative = d/d~0 dimensionless parameter, equations (14) and (22) § 1. Introduction
The behavior of a thin film of viscoelastic liquid that is subject to a n o r m a l squeezing m o t i o n between parallel, circular opposing surfaces is considered in this report. This geometry, the simplest form of the so-called squeeze film, has been widely studied by m a n y workers since Stefan in 1874 [ 1]. Stefan considered the behaviour of a N e w t o n i a n fluid and disregarded fluid inertia forces in the equations of motion. In m o r e recent years, researchers have extended the work of Stefan by including inertia effects, and by using constitutive relations other than the N e w t o n i a n fluid model. In this study, an analytic solution to the equations of motion, including inertia terms, is presented for an interesting and useful class of squeeze film flows using a simple integral constitutive equation to describe the behavior of the viscoelastic fluid. The squeeze film is of practical interest to lubrication and to rheometry. It is representative of all unsteady h y d r o d y n a m i c lubrication. Knowledge of the changes in pressure, as the separation between bearing surfaces
VISCOELASTIC SQUEEZE FILMS
503
varies, is desirable for satisfactory bearing design. The parallel surface case, while not often used as an actual bearing configuration, is instructive for study of the mathematics and mechanics of unsteady lubrication. The parallel surface squeeze film has also been used as an instrument to measure viscosity and other fluid material properties; see for example, the discussion of Leider and Bird [2]. The class of squeeze film flows to be analyzed is the case of parallel circular coaxial disks, one of which is subject to small sinusoidal oscillations in a direction normal to its surface. These conditions may exist in lubrication when bearings are subject to periodic load variations. The inclusion of inertia effects is desirable due to recent tendencies to increased machine speeds and low viscosity lubricants. Modern multigrade lubricants are known to be viscoelastic, thus warranting the use of nonNewtonian fluid models in bearing analyses. The normally oscillating squeeze film, including the effect of fluid inertia forces, was studied for the Newtonian case by Kuhn and Yates [3], and Terrill [4]. The latter achieved an exact solution to the Navier-Stokes equations, which can be obtained as a limiting case (low Deborah number) of the results presented herein. Terrill's solution shows the existence of inphase and out-of-phase components of fluid velocity and pressure (or load), while the solution of lubrication theory requires that these variables be exactly in phase with the plate velocity. Of particular interest to this study is a paper by Kramer [5], who was concerned with the case of flow between parallel circular disks subject to any (not necessarily sinusoidal) normal motion. He showed that if the problem is formulated in a system of convected coordinates, using for the constitutive equation Lodge's "rubber-like liquid" [6], a very simple form for the equations of motion results. Kramer's technique was successfully used by Grimm [-7] to include inertia effects in a study of Newtonian fluid films. Kramer also considered the small sinusoidal motion case and shows a further reduced form for the governing equation. He finds for his solution that the fluid and plate velocities are exactly in phase, regardless of the values of the parameters in the constitutive equation, although the load experiences a phase shift. He does not include inertia effects in this solution. The result of Kramer is also obtained as a limiting case of this analysis (low Reynolds number).
504
JOHN
A. TICHY
§ 2. Analysis F o r m u l a t i o n . Consider the configuration described above of two parallel, circular, axially aligned disks, the upper of which is subjected to small sinusoidal oscillations in a direction normal to its surface. The plate motion can be described by
h = hm(1 - e cos coz), V = ~
e ~ 1,
= hmcoe sin cot,
(1)
where h is the film thickness at any time ~, h mis the mean film thickness, co is the oscillation frequency, e is the dimensionless oscillation amplitude, and V is the plate velocity. Kramer showed for this geometry, using the axisymmetric assumption and a linear variation of radial velocity in the radial direction, that the particle pathlines are of the form r = rog(Zo, 7:)
~ = (a o
z = f(Zo, z)
(2)
The cylindrical spatial coordinates r, q5and z are functions of the convected cylindrical coordinates r o, qS0, z o and time, see Fig. 1. When T = 0, the convected and spatial coordinates of a particle are identical, i.e., 9(Zo, O) = 1 and f ( z o, O) = z o. The integral constitutive equation used is that of Lodge's "rubber-like liquid," in which the contravariant deviatoric stress tensor S ~ at present
T=O I/,," / / / /
I
/////////i
Y'//
/ / /
I' j'//
;-ol ,) / /
/ /
/ /
/ ,/ / / /
~
,tl / i/
/ / / / / A
"~.
z=t ( z o r)
//
//
1"~O
/ ///4Materlal
I
/ / / / / / / /
h
I
/
/
/
/ / /
////J
iz I///////////I/~
R
plens Zo=constant
in determed s t a t e
+/~i F positive as shown ~/J/
deformed e'~ate (r~,O)
/
/ / / /
~
/f
/d
Fig. 1. Squeeze film geometry.
VISCOELASTICSQUEEZEFILMS
505
time t is determined from the contravariant metric tensor yij and a memory function K according to g ( t - r)[TiJ(z) - 7iJ(t)] dr.
SiJ(t) =
(3)
The metric t e n s o r yij can be found from the relation between the spatial and convected coordinates, equations (2), and the condition ds2(r) -- Yij('c)d~id~ j = dr 2 + (rdqb) 2 + dr2; ?ijY jk = (~ik,
(4)
employing the usual summation convention on repeated indices, with ~i representing the convected coordinates in summation notation, and ~i k the Kronecker delta. Lodge's model is known to portray, at least qualitatively, the behavior of many viscoelastic fluids in a variety of experimental conditions. Kramer goes on to show that the constitutive equation (3), the equations of motion, and continuity can be combined to give the following expression:
P 020 (Zo, t) -~ O(Z o, t) Or 2 -
~2Oezg(z°' t) r [ K(t- r)g4(Zo,V)dvg(Zo, t) .j_oo ff
3
-oo K ( t - z ) g
02g 4F(t) (Zo, r)~Z2o (Zo, Z ) d z - ~zR4
(5)
where F is the load on the plates defined as positive in tension (upward), R is the plate diameter, and p is the fluid density. If the deformation is assumed to be of small amplitude about the undeformed state, that is (6)
g(Zo, T) = 1 + eg~(zo, z)
then -pea-r2
(Zo, t ) + e
~g~
K(t-Qdr-
~z 2 (zo, t ) oo
-
retaining only terms of order e.
e
f
t
O2ge oo K ( t - r) ~
4F(t) dr = rcR4
(7)
506
JOHN A. TICHY
The memory function is now represented as the sum of I exponentials l
K(s) = Y~ ~2. e -~/°~ i=l
(8)
i
where the 0~ are relaxation times, the #i are viscosities and the past time s = t - z. Periodic forms for g~ and F are now assumed, g~ = - m(Zo)[COS0)r - 1] + n(Zo) sin 0):
(9)
and
F = e[F~ sin cot + Fo cos 0)t]
(i0)
where m and n are arbitrary functions of order 1. Note that the condition g(z o, 0) = 1 is satisfied by (9). Equations (8), (9) and (10) are now substituted in the governing equation (7). The coefficients of the sine and cosine terms are equated and rearranged after the integration has been performed to give d2n 0)2 i #iOi i= 1 1 + ((DO/)2 -Jr-
e p0)2n + ~
d2m -~-~0)
{
dam0)2 i
e -P0)2m - dz~oo
1 fli } 4F~e i= =~1 1 -~ ((DO/)2- = nR 4
(11)
~iOi
i=t 1 + (coO/)2 + d2n
1
+ d- o
kq
1+ (
}
4Fe
oy
A single relaxation time parameter is defined by fl iOi
0 = i= 1 1 + (0)0i)2 l
(12)
i~=1 1 + (0)0i)2 Similarly, a viscosity parameter can be obtained: l = [1 + (600)2] E ~/i i= I I -~ (0)01) 2
(13)
507
VISCOELASTIC SQUEEZE FILMS
Note that if the memory function K has only one exponential term, 0 = 01, and/2 =/21. The following dimensionless parameters are also introduced:
zo O=coO, 5 O - h ( t ) ,
Re=
pcoh 2 /2 , F -
4fhZm rcR4/2co.
(14)
The symbol 0can be interpreted as a Deborah number (at least in the case when K(s) is a single exponential); that is, a ratio of a characteristic time of the material 0, to a characteristic time of the process 1/co. Similarly the Reynolds number is indicated by Re. Substituting (12), (13)and (14)into (11), rearranging, and denoting d/dff0 by the prime symbol ('), yields
e{Ren(~o)+l@[rn"(~o)+On"(~o)]}=eFs
(15) e - Re m(#o) + ~
[n"(g0) - 0m"(5o)] = eFc
Only terms to order e are retained, which requires that the higher order terms from the coordinate transformation between z 0 and ~o be omitted. The symbols m(£o) and n(io) now represent the functions obtained from m(Zo) and n(zo) by making the substitution z o = ~oh. The equation of continuity for an incompressible material is reported by Kramer as
C~Zo
= (g(zo, ~))-4
(16)
1 -
(17)
which for equation (6) becomes
af
-
2eg~
8z o retaining terms of order e. Assume the following form for f:
f ( z o, r) = z o - ehm{P(Zo)[COS cot - 1] - q(Zo) sin coz}
(18)
noting that f(Zo, 0) = z 0. From equations (9), (17) and (18) it follows that
m(zo ) _
h m dp 2 dz 0
h m dq n(z°) = - 2 - dz-~-
(19)
508
J O H N A. T I C H Y
Substituting equations (19) into (15), retaining terms of order e, and raising the order of differentiation by one to eliminate the unknown F's gives as the differential equations to be solved: 1
Re q"(2o) + T ~ -Rep
t!
(Zo) + ~ -
[Pw(Z°) + OqW(2°)] = 0
(20)
1
[q (Zo) - Op'V(2o)] = 0 IV
-
As in the development of equation (15), only terms to order e are retained; and the functions p(2o) and q(5o) are obtained from p(Zo) and q(Zo) by the substitution z 0 = 2oh. Finally, since e ~ 0, it can be removed from the lefthand side of the equations (19) to yield the result as shown. Boundary conditions on the p and q are obtained from the velocity components v,, % and v~ which are determined from the particle pathlines according to
v~(r, z, t)
~'c]~_t_
v4' =
~=t' Vz(Z, t) = c3z
(21)
The material coordinates ro, q~0, and z o are regarded as constants in performing the differentiation of (21). Substituting among (2), (6), (9), (18) and (19) into (21), introducing the dimensionless variables Vr
Uz
r0
v= = - - ,
~r - - (.oh m ,
ro = ~-,
and t-= cot,
(22)
coh m
transforming the coordinate z o to go, and retaining terms of order e, gives foR e[p'(2o) sin t-+ q'(2o) cos t-] 2h m
f, -
(23)
~ = e[p(5o) sin t-+ q(2o) cos }-]. The boundary conditions on v, and vz are z = O ( 2 o =0):
vr=O
t)z=O
z = h ( 2 o = l):
v,=O
vz = (dh/dz)~= t = h,~coesin cot
which from equations (22) and (23) gives
(24)
509
VISCOELASTIC SQUEEZE FILMS
p(0) = q(0) = p'(0) = q'(0) = p'(1) = q'(1) = q(1) = 0
(25) p(1) = 1 as the b o u n d a r y c o n d i t i o n s on p a n d q. Solution. The s o l u t i o n to e q u a t i o n s (20) subject to the b o u n d a r y c o n d i t i o n s (25) m a y be f o u n d b y c o n v e n t i o n a l m e a n s to be:
P(eo) =
~(1)~(2 o) + fl(1)fl(2o) ~2(1)+~20)
'
(26)
~(1)~(2o) - ~(1)~(~o) q(/o) = ~2(1)+fl2(1) where • (2o) = cosh a(1 - 20) cos b ( l - Zo) - cosh a 2 o cos b2 o + a2 o sinh a cos b
-b2ocoshasinb-coshacosb+
l,
(27)
fl(2o) = sinh a(1 - 20) sin b(1 - 20) - sinh a 2 o sin b2 o + a 2 o cosh a sin b + + b2 o sinh a cos b - sinh a sin b a =
--[(1
b =
--[(1
e
(28)
+ 02), _
(29)
+ 02) ~ + O]
(30)
The velocity d i s t r i b u t i o n can n o w be readily e v a l u a t e d from e q u a t i o n s (22) a n d (23). It is also necessary to e v a l u a t e the l o a d s on the plates, F, in response to the i m p o s e d plate m o t i o n . Substituting e q u a t i o n s (19) into (15), a n d retaining terms of o r d e r e as usual, p r o d u c e s
e Req,(~o) + ~
Ep,,,(2o) + Oq,,,(2o)]
2
2
{_Rep,(2o)
1
= efts )
(31)
+ 1-T0r [ q " G ) - 0p'(eo)]>) = ~Fo
as expressions for the dimensionless c o m p o n e n t s of the load. At 2 o = 0 (or 2 o = 1), p' = 0; thus the l o a d s are m o r e easily e v a l u a t e d from
510
JOHN A. TICHY
Fs =
-1 [p"(o) + Oq"(o)] 2(1 + 02 )
(32) -1 F,: --.= [-q"(0) - 0p"(0)] 2(1 + 0 2) The stream function will prove to be of interest to provide a pictorial representation of the flow behavior. In cylindrical coordinates, for the axisymmetric case, the stream function ~, is defined by v~ -
1 aq/ r c~z
1 aV/ vz = - -- - r ~r
(33)
Defining a dimensionless stream function ~, and transforming to the dimensionless variables of (14)and (22) gives 1 R
Or)
17r = r0 hm (~Zo' gz =
1 aV) fo gf0 ' q = coR;hm ,
(34)
by retaining order e terms. The dimensionless stream function then becomes v) = - ~ -e f0e[P(io) sin }-+
q(io) cos }-],
(35)
requiring that 0 = 0 at io = 0, which eliminates an arbitrary additive function of ~. The above expressions for velocity field, stream function, and load, describe the behavior of the pertinent dependent variables for this problem. § 3. Results
Velocity distributions. The behavior of the functions p (5o), q (~0), P'(~0) and q'(Zo) is shown in Figs. 2a-2d and Figs. 3a-3d. Note from equations (23) that P'(Zo) is proportional to the in-phase (sine) component of radial velocity (the radial velocity at the instant of time cot = n/2), and q'(~o) is proportional to the out-of-phase (cosine) component of radial velocity (radial velocity at the instant of time cot = 0). The terms "in-phase" or "out-of-phase" indicate phase relationships with respect to the upper plate velocity V. Similarly, P(Zo) is proportional to the in-phase (sine) component of the axial velocity distribution, and q(5o) is proportional to the out-ofphase (cosine) component of the axial velocity.
VISCOELASTIC SQUEEZE FILMS
51 1
or (d) 1.0 -.12
0.5
1.0 [-
W
osL
I.--0
, .06
_..--'~ :o,,,,o~._~"..p:' ,
L°F
00/"
=0
IV/
J I. (c) I,O-0.6
0.5 [
~ ' " ~ R e
0:J.o
I 0.5
I
|
J2
/~----R~:O,t,,o- - - - - . &
I//"-
0 V" 0
0
~b)t.o -.o3~
// Re:IOO~
/ / < / ~ R e =0,I,I0
,
0.5 p(2'o)
o
.os~
J~7/-R~=o ( ( ~R~-, ~ k ~
,°:°,
1.0 -.035 (a)
n~-
0 q(Zo)
.055
Fig. 2. Behavior of the functions P(Zo) and q(Z7o).
F o r Newtonian fluids, the relaxation times 0~ are zero, and consequently the Deborah number 0 is zero. In this case, from equations (29) and (30), a = b = x / R e ~ , and the functions P(Zo) and q(io) become identical to expressions obtained by Terrill [4]. These results are displayed in Figs. 2a and 3a. F o r the non-inertial case, the Reynolds number R e is zero. In the limit as R e ~ O, the equations (26) reduce to p(%) = 3%z - 2~o3,
q(%) = o,
(36) P'(Zo) = 6(z0 - 2~),
q'(Zo) = 0,
512
J O H N A. TICHY I.O
~r ~__.__.~...__~ ~ i ~o / / ~
-LO
0
1.0
2.0
~ Re=1o
-.5
0
-2
o
.5
(d)
, -,.o
o-~
o ~ , o
,o
20
(e)
1.0~~ ~o
-,.0
t
~ R e
/J
0 I..,~I"- t 0
i
l.O
2.0
§ =o.~
=t(30
r
(b) -.2
0
' "O~------x.,
I
~°
4
Re=O,l,IO
O15 L.-
2
.2
~:~e
N R~:°''''°
O.Sh I ) IRe=tOO y
=,OO
/~,4
_
( ~p--Re=10 ""......~Re =O
0=0 ,
-I.0
ol..
0
.,
1.0
P'(Zo)
.
2.0
(a)
,
-.2
....
0
.2
q'( Zo)
Fig. 3. Behavior of the functions P'(Zo) and q'(~o), proportional to the sine and cosine components of radial velocity profile, respectively.
regardless of the values of/7. This result, the familiar parabolic radial velocity profile, is the same as the prediction of lubrication theory. The cosine components of velocity are identically zero. Kramer [5] also obtained this result for Lodge's "rubber-like liquid" by neglecting fluid inertia forces in his original formulation. These resulting velocity distributions for Re = 0 are shown on Figs. 2 and 3 for the various values of Deborah number 0. Inspection of Figs. 2 and 3 reveals that both Re and 0 affect p'(5o), the inphase component of radial velocity, in the same way. Below certain values of Re and 0, the parabolic velocity field is changed very little, then suddenly an unusual flow reversal in the centerline area occurs. The higher the
513
VISCOELASTIC SQUEEZE FILMS
Deborah number, the lower the value of Reynolds number at which this reversal occurs. The out-of-phase radial velocity component always exhibits this flow reversal phenomena unless Re = 0, in which case the q'(~o) is identically zero. The higher the values of Oand Re, the greater the magnitude of q' (~o). For low Deborah numbers, say O = 0.1 as shown in Figs. 2b and 3b, there is little qualitative difference in velocity field than that of Newtonian fluids, O--0. This would likely be the case for many applications involving slightly viscoelastic fluids, such as petroleum-polymer lubricant oils. To provide a pictorial representation of typical velocity fields, Fig. 4 has been prepared. Streamlines, as evaluated from (35), are also displayed. For the case Re = 10, O = 1, the cosine component of velocity field is shown in Fig. 4a (the velocity field at the instant of time cot = 0) and the negative in the sine component of the velocity field is shown in Fig. 4b (cot = 3~/2). Fig. 4a is especially interesting. Note that although the plates are motionless at this instant, the fluid is undergoing a peculiar circulating motion with the aforementioned flow reversal. Loads. In the case of loads, it is more instructive to consider the amplitude and phase relationships rather than the sine and cosine components F s and F c. Dimensionless load amplitude F* and phase angle
Re--Jo, ~--,
Fig, 4. Typical velocityprofilesand streamlines.
(a)
(b)
514
JOHN A. TICHY
6o t-
50-0=0.1
,L 3 o -
20--
I00,~--~]--'0
20
I
I
I
40
60
80
I00
Re
Fig. 5. Load amplitude versus Reynolds number for various Deborah numbers.
20.
~
0 0,0 I
II -~o.I0
I 20
40I
60I
80I
Re
Fig. 6. Phase angle versus Reynolds number for various Deborah numbers.
I I00
VISCOELASTICSQUEEZEFILMS
515
are defined by F * = (F 2 + F2~ ~
(37)
and
tan
==Fc Fs
(38)
such that F(t-) = F* sin(t+ ~)
(39)
The phase angle q) represents a shift between the load F and plate velocity V, equation (1). The amplitude and phase angle of load are evaluated from (32), (37) and (38) and displayed in Figs. 5 and 6. For the limit R e ~ O, if*
-
6 x/1 + 02
(40)
and tan q) = - 0
(41)
which can be shown to be equivalent to expressions obtained by Kramer. Inspection of Figs. 5 and 6 shows a smooth behavior of the functions for load amplitude and phase angle with Reynolds number, at low Deborah number. Then a wildly oscillating behavior at 0 = 10 is exhibited as repeated load resonances occur at various Reynolds numbers. Note from Fig. 5, that apart from the resonance phenomena a reduction in load amplitude with higher Deborah numbers is noted. An increase in load amplitude with increasing Reynolds number is evidenced due to fluid inertia forces adding to the resistance of the fluid viscous forces. § 4. Discussion
A solution to the equations of motion, including inertia forces, for a hypothetical viscoelastic fluid - Lodge's "rubber-like liquid" has been presented. The solution is exact for any Reynolds number provided the amplitude of oscillation is very small. The problem was readily formulated by using Kramer's convected coordinate ~ystem for the oscillatory squeeze film flow.
516
JOHNA. TICHY
The results of this study shed more light on several areas of application. This simple squeeze film flow is representative of the mode of lubrication in high-speed bearings undergoing oscillatory loading. Because modern multi-grade libricants are viscoelastic due to the addition of polymer additives for viscosity-index improvement, the inclusion of a viscoelastic fluid model in the analysis of a highly unsteady flow is pertinent. In particular, lubrication engineers have desired to know what effect, if any, viscoelastic lubricant properties have on load capacity. It is generally assumed, from an application standpoint, that a bearing running on a viscoelastic lubricant is operating in the low Deborah number regime. For this case, it is seen from Fig. 5 that there is very little change in load capacity due to viscoelastic effects. If however, a high-speed bearing were operating at high Deborah number, a variety of anomalous resonance effects would apparently occur, complicating matters considerably. In fact, the existing controversy in lubrication on the effect of viscoelasticity on bearing load capacity, and the discrepancies between experimental and analytical results, may be, in part, explained by such resonance phenomena. (See, for example, the discussions of Appeldoorn [8] or Harnoy [9].) It should perhaps be added that Kramer's prediction using Lodge's equation for the case on an applied step load on a squeeze film bearing, does not describe the data of certain other workers. Kramer predicts (for the step load case) a more rapid squeezing rate for viscoelastic liquids than for Newtonian fluids, while the data of Leider [10] and Brindley et al. E11] indicate a slower squeezing rate. The results of this study may prove useful to rheometry. Many measurements of rheological properties are limited in that, at high frequencies, fluid inertia effects may intervene and limit the range of measurement. In this study, there is no limit on Reynolds number. The only limitation is that the oscillation amplitude must be small. It would be extremely interesting, and relatively simple, to attempt to experimentally verify the predicted resonance behavior of the viscoelastic fluids. A study of this type would provide further information on the utility of Lodge's viscoelastic model, and on the actual behavior of viscoelastic fluids in a relatively simple, but non-viscometric flow. Received15 September1977
VISCOELASTIC SQUEEZE FILMS
517
REFERENCES [1] S~VAN,J., "Versuche fiber die scheinbare Adhesion," K. Akad. Wissenschaften, MathNaturwissenschaftliche Klasse, Wien. Sitzungsberichte, ~. 69, pp. 713, 1874. [2] LE~OZR,P. J. and R. B. BIRD, Ind. Eng. Chem. Fundam. 13, No. 4 (1974) 336. [3] Ktrr~, K. C. and C. C. YATES,ASME Trans. 7 (1964) 299. [4] TER~LL, R. M., Jour. Lub. Tech., Trans. ASME, Series F, 91 (1969) 126. [5] KRa~mR, J. M., Appl. Sci. Res. 30 (1974) 1. [6] LODGE,A. S., Elastic Liquids, Academic Press, London and New York, 1964. [7] GRIraM,R. J., Appl. Sci. Res. 32 (1976) 149. [8] APPELDOORN,J. K., J. Lub. Tech., Trans. ASME, Ser. F 90, n. 3 (1968) 526. [9] HARNOY,A. and W. PUILIVVOVF,ASLE Trans. 19, n. 4 (1976) 301. [10] LEIDER,P. J., Ind. Eng. Chem. Fundam. 13, n. 4 (1974) 336. [11] B~NDLEY,G., J. M. DAVmSand K. WALTERS,Jour. Non-Newtonian Fluid Mech. 1, n. 1 (1976).