Adiabatic shearing of one-dimensional thermoviscoelastic flows caused by boundary and inertial forces

Adiabatic shearing of one-dimensional thermoviscoelastic flows caused by boundary and inertial forces N. C. C H A R A L A b l B A K I S and E. N. H O U S T I S

Aristotle University of Thessaloniki, Thessaloniki, Greece We study the analytical and numerical behaviour of adiabatic shearing flows of viscoelastic fluids with temperature dependent viscosity, under several types of boundary conditions and an 'oscillatory' inertial force. We consider the shearing adiabatic flow of an incompressible Newtonian fluid with temperature dependent viscosity, caused by two types of shearing forces and a time dependent 'oscillatory' body force. Our first objective is to study the analytical behaviour of the solution of this problem. We show that the velocity and the stress are well behaved, provided the rate of change of viscosity with temperature satisfy appropriate bounds, that correspond to the type of the fluid. Our second objective is to devise a numerical solution of the problem and study its behaviour. The numerical results obtained indicate an agreement between the analytical and numerical behaviour. Key Words: incompressible Newtonian fluid, temperature dependent viscosity, adiabatic shearing,

a priori estimates, system of nonlinear partial differential equations, finite element collocation program

1. INTRODUCTION

In this paper, we consider the flow problem of adiabatic shearing of viscoelastic liquids or gases, due to combinations of time periodic or steady boundary shearing forces and an 'oscillatory' inertial force. If we normalise units so that the density of the fluid and its specific heat are unity, then from the conservation laws of momentum and energy we have for (x, t) E [0, 1]x [0, oo): V t ( x , t) = Ox(X , t ) + f ( t )

(1.1)

Ot(x, t) = o(x, t) Vx(X, t)

(1.2)

where v(x, t) is the velocity in the direction of the flow, f(t) the inertial force, O(x, t) the temperature and o(x, t) the stress for Newtonian fluid defined by:

o(x, t) = lz(O(x, t)) Vx(X, t)

(1.3)

with #(0) being the temperature dependent viscosity function. The inertial force is assumed to be oscillatory in the sense: t

f

f(r)

dr

0

is bounded for 0 ~< t < oo. The subscript indicates differentiation with respect to the variable indicated. In this paper we assume ta(0) is twice continuously differentiable and: Accepted May 1985. Discussion closes February 1986. 0264-682X/85/040205--6 $2.00 © 1985 Computational Mechanics Publications

(cl) the (v, 0) solution is subject to the following compatible initial and boundary conditions:

v(x, O) = re(x), O(x, O) = Oo(x)

for 0 ~
o(1, t) = 0(t) and o(0, t) = 0

for 0 ~< t < oo (1.5)

where re, Oo and O are given functions such that: Vo E W2'2(0, 1), 0oE W1'2(0, 1),

?l(1)=la(Oo(1))Vox(1)

0EL~oc(0, oo),

and Vex(0)=0

(c2) the viscosity function/1(0) satisfies one of the following relations for 0 < 0 < oo: (H1) /a(0)>0, # ' ( 0 ) > 0 ,

1<

u(O) u"(O) u'(0) 2

~
(H2) #(0) > 0, #'(0) ~< 0, --/a'(0) ~< ~u2(0), 0 ~< • < oo where assumption (HI) corresponds to the viscosity of typical gases and (H2) to typical liquids. The first objective of this paper is to study the analytical behaviour of the solution (v, 0) of (1.1)-(1.5) under the conditions (c 1) and (c2) and various choices of (O(t), f(t), #(0)). The discussion is presented in Section 2. The paper's second objective is to devise a numerical solution of (1.1)(1.5) and study its behaviour. The discussion is presented in Section 3.

Engineering Analysis, 1985, Vol. 2, No. 4

205

Adiabatic shearing of one-dimensional thermoviscoelastic flows." N. C Charalambakis and E. N. Houstis The problem of adiabatic shearing of an incompressible Newtonian fluid with temperature dependent viscosity was investigated by Dafermos and Hsiao t in which an asymptotic stability was established for the case of steady boundary velocities. The case of frictionless shearing caused by time dependent inertial force was studied (reference 2), where it was proved that the flow converges to a 'rigid' motion, in which the velocity depends only on the time and the temperature depends only on x.

1

- - ~ 1" [/a'(0(x, 0) 2 + I~(O(x, t)) ll"(O(x, t))] 0

x g"(O(x, t)) v6(x, t) dx (2.7)

= costox(1, t) On account of hypothesis (H1): U'2 ~< -- A(//2 + #"/a)

(2.8)

where # = ta(0) and 2. ANALYTICAL BEHAVIOUR OF "DIE SOLUTION

A>-Hence, using the inequality: 1

1

f ,##4 v x dx <<.-2i fe 0

if

l~'21av6 x dx + 2e 0

0 1

f , (/a2#+/~"/l)V6xdx+~eI f v2xdx 0

Assuming that g(0) satisfies hypothesis (H1) and

O(t) = sin t

(2.1)

mxin(u(0o(x)) ) > A '/2

(2.2)

0

where Vx =-Vx(X, t), we obtain from (2.7): 1

1

dtalf~2

where A > 1 / v - 1 and

ox dx +

0

t

1

¼f

II'l~v4xdx +

0

(2.3)

+

0

Io(x,t)l<<.K, Iv(x,t)l<~K,

Ivx(x,t)l<~K

f,

l~#VxdX--A

f

#V2xdx<~costox(1, t)

0

(2.4)

1

tox(1,t)l lf G(x,,) dx+K

O(x,O (2.5)

and 1

1

o2(x, t)<~ f o2x(x, t)dx<<.lfo2xx(x, t ) d x + K 0

0

1

#(O(x, t)) v2x(x, t) dx <. min (Oo(x)) f o2(x, t) dx

Proof. We start by showing that: 0

xE[0,1l

1

(2.13)

(2.6)

Using (A.3), (A.4) (see Appendix) and (2.1):

ldfo2x(x,t) dx+f#(O(x,t))v~t(x,t)dx 2 dt

--

dt

l tf

1

eXt) = ~ #'(O(x, t)) g(O(x, t)) v4x(x, t) dx

0

Engineering Analysis, 1985, Vol. 2, No. 4

+ K~t) <~K

(2.14)

where

0 1

206

de(t)

1

+4

d 0

(2.11), (2.12), (2.2), (H1) and (1.1), we obtain from (2.10) the differential inequality:

0

0

(2.12)

Using

Oo(x).

1

(2.11)

0

In the sequel, K will stand for a generic constant that depends on v or X and can be estimated from above in terms of lower bounds for min Oo(x) and upper bounds for x the W2'2(0, 1) norm of Vo(X) and the W1,2(0, 1) norm of

f O2x(X,t) dx <~K

(2.10)

where Ox - Ox(X, t), Vxt- Vxt(X, t), #' = #'(0). On account of (A.6), (A.7), (A.1) and (2.1):

for (x, t) E IR and

f la(~) d~ <~Kt Oo(x)

laV2xtdx

1

0

for arbitrary t > 0, then, for the solution (v(x, t), O(x, t)) of (1.1)-(1.5) in [lq = [0, 1] x [0, ~):

f 0

1

f f(r) dr <~constant

Pv2xdx

1

< 21e ( - - A )

Proposition 2.1

(2.9)

v-1

2.1. Shearing caused by time harmonic boundary force We consider the flow of a viscoelastic gas, caused by an oscillatory inertial force and a time periodic boundary shearing force. The analysis, which is summarised in Proposition 2.1, indicates that the stress, the velocity, the velocity gradient and the temperature remain bounded, provided the initial viscosity is sufficiently large.

1

1

o~x, t) dx + -2 u'(O(x, t)) u(O(x, t)) 0

x v~(x, t) dx

0

(2.15)

Adiabatic shearing of one-dimensional thermoviscoelastic flows: N. C Charalambakis and E. N. Houstis t 1

Hence ~t) < K

lf f

(2.16)

8

and recalling (2.15), we obtain (2.6). The first of (2.4) follows directly from (2.6) using (A.6). Next, using (1.2) and (1.3), we obtain Ivx(X, t)l ~
, dx d,

O0

1

if

(2.21)

dx

2

0

def

co(t) = minv2(x, t)

whence

x

we have:

1

v=(x, t) <. 26o(0 + 2 f v2x(x, t) dx

0

Hence, using (A.6), (1.3), (HI), (2.10) and (A.8), we obtain (2.19). Using (A.5), we obtain (2.20).

0 1

(2.22)

f o~(x, t) dx <~K

1

2

Proposition 2. 3 Assuming that #(0) satisfies (H2) and

0

0

whence, on account of (A.8), (2.1) and (2.3), we obtain Iv(x, t)l < K .

#(0) -x ~< constant x

#(~) d/j

p > 1 (2.23)

00

2.2. Shearing caused by steady boundary force We now consider the flow of a viscoelastic liquid or gas, caused by a steady boundary shearing force and an 'oscillatory' inertial force. The analytical behaviour of the solution (v, 0) of (1.1)-(1.5) is summarised in propositions 2.2 and 2.3 and indicates that the stress and the velocity gradient of the gas remain bounded independently of the initial data. On the other hand, for given time the velocity of the liquid remains bounded, provided the viscosity function satisfies appropriate bounds.

Proposition 2.2

1

-- t ~< v(1, t) <. Kt K

f la(~)d~ <~Kt Oo(x) Iv(x, t)l ~
(2.17)

t

(2.18) 1

f u(o(x, t)) v~(x, t) dx <<.Kv(1, t)

then, for the solution (v(x, t), O(x, t)) of (1.1)-(1.5) in Ilq =- [0, 1] x [0, 0o): Io(x, t)l <.K, IVx(X, t)l ~
0

0 (x, t)

u(O(x, t))-' <<.c[ f u(~)d$]V"<~Kv(1, t) 1/"

O(x,0 f la(~)d~ <<.Kt Oo(x) parts: 1

t 1

21 f Ox(X,2 t) d x +

ffu(O(x,r))~r(X,r)dxdr O0

for every x. Hence, on account of (2.28) and (2.29): 1

f o

v2x(x, t) dx <.

v(1, t) <. Kv(1, 0 (1/°)+1

min #(O(x, t)) x~[O, ll

(2.30) Using (A.8), (2.17), (2.18) and x

1

0

(2.29)

Oo(x)

(2.20)

Proof. Integrating (A.3) over t and after integration by

if

(2.28)

By (A.6), (A.5), (2.23) and (2.27):

for (x, t) E IR, and

+4

(2.27)

0 0

0

0

(2.26)

1

f f o (x,r)dx f(r) dr < constant

(2.25)

Using (A.2), (1.5), (2.17), (2.18), (H2)and

t

f

(2.24)

0(x, t)

Proof. (A.6):

Assuming that/a(0) satisfies hypothesis (H1) and:

O(t) = 1

and that (2.9) and (2.10) hold, then, for the solution (v(x, t), O(x, t)) of (1.1)-(1.5) in ~ = [0, 1] x [0, oo):

la(O(x, t)) p'(O(x, t)) v4(x, t) dx

v(x, t) = v(l, t) + f Vx(~, t) d~ 1

Engineering Analysis, 1985, Vol. 2, No. 4

207

Adiabatic shearing o f one-dimensional thermoviscoelastic flows: N. C Charalambakis and E. N. Houstis (a) Approximate space. For each t, each component of the solution vector is approximated by a cubic piecewise polynomial with C~([xl,xr]) smoothness. Assuming a partition A =- {xx = ~1
we obtain: 1

v(1, t) <~2 f [Vx(X, t)[ dx + Kt o

whence, using Schwarz's inequality, (2.30) and p > 1, we obtain the right inequality (2.24). Combining v(1, t)<~Kt and (2.27), we obtain (2.24). Inequality (2.25) follows directly from (A.6), (A.5) and (2.27). Finally, (2.26) follows from (2.30) and (2.24).

3. NUMERICAL SIMULATION In this section we consider the numerical solution of the thermomechanical model described by the relations (1.1)(1.5). The objective is to verify the analytical behaviour of the model. The relations (1.1) and (1.2) can be written in the form:

vt = #'(0) O~vx + u(O) vxx + f(t)

(3.1a)

0 t = #(0) V2x

(3.1b)

The boundary conditions (1.5) indicate that (v, 0) satisfies the conditions:

Vx(X,t)=Oatx=O

and #(0) V x = O ( t ) a t x = l

(3.2)

In the case of a null boundary condition as for the equation (3.1b), we assume that it is satisfied at x = 0 and x = 1 and that (v, 0) is subject to the initial conditions (1.4). For the numerical solution of the initial/boundary value problem defined by the equations (3.1) and subject to boundary conditions (3.2), (3.1b) at x = 0 and x = 1 and to initial conditions (1.2), we employ a method on lines based on a finite element semi-discretisation procedure.

2k+2

ui(x, t) ~- Ui(x, t) = Z cti,J(t) q)j(x) i = 1. . . . . N j=l where (~j(x) are the basis functions of the chosen approximation space. For the definition of ¢~'s see reference 3. (b) Semi-discretisation. We set U = (U1,..., UN). For each t, we force /) to satisfy the partial differential equation (3.3) at certain interior points (:ca < r2 < . . . < r2k+l < Xr), called collocation points:

dai,/ ~. dPf(rm)--~-t =fi(t, rm, U(rm, t), (]x('rm, t), (]xx(Tm, t))

2k+2 j=l

(3.6) for i = 1,2 . . . . . N and m = 2 , . . . , 2 k + 1. The choice of collocation points is a consequence of the optimum convergence of the method in reference 3. In this case, the ri's are chosen as the roots of a quadratic Legendre polynomial with respect to each subinterval of A. (c) Semi-discretisation o f the BC. For each t, we force 0 to satisfy the equations:

~ [ab k Obk "~ {--Dtu/+--Dtu/xl=Dtzk j = 1 \ ~Ui

~Ujx

t?(x, to) = ¢(x)

The numerical method used is designed to solve a general system of N nonlinear partial differential equations of the form:

D ~ = F ( t , x , u , ux, Uxx) (x,t)~[x~,x,] × [to, °°) (3.3) where

O-x= (Dxul . . . . . DxUN)

P = (f,,f2 . . . . . fN)

gtxx = (D2xUl, "'" , D2xUN)

The components of the vector solution ~ and its derivatives are functions of x and t. The right side of (3.3) is a vector whose components may depend on the indicated variables. The boundary conditions are imposed at x~ and/or xr (or not imposed in th absence of conditions) and must have the form: b(a, ax) = Z(t)

(3.4)

where b, g are some N-vector valued functions. We assume that t~ satisfies some initial conditions:

us(x, to) = t.bs(x)

i = 1,2 . . . . . N

Engineering Analysis, 1985, Vol. 2~ No. 4

(d) Software implementation. The method was realised by the appropriate selection of parameters of a general finite element collocation program PDECOL. 4 These parameters are K O R D - t h e order of the piecewise polynomials, NCC - the number of continuity condition at the knots of the piecewise polynomial, M F - the method for solving the system of ordinary differential equations and various accuracy specifications. We have considered the Hermite cubic piecewise polynomial approximation of u in space x (KORD = 4, NCC = 2) and Gear's stiff method to solve the ordinary differential equations (3.6) and (3.7), with the needed Jacobian matrix explicitly defined (ME = 21).

3.2. Numerical results We considered four instances of the thermomechanical model to study its numerical behaviour. These cases are characterised by the selection of the viscosity function ta(0) and the form of the stress O(t) at Xr. Throughout these experiments, it is assumed that f(t) = sin t.

(3.5)

that satisfy (3.4) at t = to. The model of the thermomechanical processes considered in this paper belongs to the above class of initialboundary value problems. The numerical method chosen to solve the above class of problems consists of the following components:

208

(3.7)

at x t and x r. It is worth noticing that the unknown coefficients as, j depend only on the time t and satisfy the ordinary differential equations (3.6) and (3.7) subject to initial conditions:

3.1. Numerical method used

t~ = (Ut . . . . . UN)

k = 1. . . . N

]

Case (i)." Is(O) = 07, 7 > 0 and 6(t) = sin t The method was applied for 0 < 3' < ~ that satisfy (H1) and initial conditions vo = 0 (and 1), 0o = [3'(1 -- 7)] -v2~. The stress function found to be periodic in time and its behaviour for each time interval (27r(i--1), 2rri] with i = 1 , 2 , . . . is summarised in Fig. 1. In this case, the

Adiabatic shearing o f one-dimensional thermoviscoelastic flows: N. C. Charalambakis and E. N. Houstis numerical simulation indicates that the stress, velocity and its gradient are bounded in [0, 1] x [0, ~). This is in agreement with Proposition 2.1.

1.000

0.600

Case (ii): la(O) = 07, 7 > 0 and ?fit) = 1 The method was applied for the initial conditions 00 = 1 and Vo=X2/2, chosen to agree with the boundary conditions ~, 0.200 (1.5) at t = 0. The parameter 7 was chosen to satisfy ( H I ) . . ~ , The numerical simulation has shown that the stress is independent of t and that the velocity and the temperature are ~, -0.200 increasing functions of t. Notice that the input parameters chosen satisfy assumptions of Proposition 2.2 and the numerical results obtained are in agreement with its assertion. Figures 2 and 3 depict stress and velocity behaviour. -0.600

J

t=7.85498 t=6,28419 t=4,71339 t=3,14259 t=1,57180 t=O

-1.000 0

1.000

t=12.5674 t=10.9966 t=9.42578

!

!

0.2oo

0400

1.00o X

Figure 3. Numerical behaviour o f velocity for 0(1, t) = 1, ?t = 0.3, 0o = 1, Vo =x2/2, #(0) = 07, and f ( t ) = sin(t)

0.600"

0.200" ~

Case (iii): #(0) = 0-7, 7 > 0 and O(t) = sin t The method was applied for the same initial and boundary conditions as in the case (i). The parameter 3' was chosen to satisfy (H2) for 0 = 0o. The numerical behaviour of the solution coincides with the one observed in case (i).

t=lr

-0.200

-0.600'

t : 3 ~ r / 2 ~

-1.000

.

0

0.200

,

,

,

0.400

0.600

0.800

1.000

Case (iv): #(0) = 0 -7, 7 > 0 and ~(t) = 1 The method was applied for the same input parameters as in the case (ii). The parameter 7 was chosen to satisfy (H2) for 0 = 0o. The numerical results indicate the same type of behaviour as in the case (ii), which agrees with the analytical behaviour described in Proposition 2.3.

x

Figure 1. Numerical behaviour o f stress in the case o f periodic shearing force with 7 = 0.3, Vo = O, 0o = 13.48, #(0) = 07 and f(t) = sin(t)

1.000

0.600

1 Dafermos, C. M. and Hsiao, L. Adiabatic shearing of incompressible fluids with temperature dependent viscosity. Quart. Appl. Math. 1983, 41, 45 2 Charalambakis, N. Adiabatic shearing flow caused by time dependent inertial force. Quart. Appl. Math. 1983, 41,275 3 Houstis, E. N. Application of method of collocation on lines for solving hyperbolic problems. Math. Comp. 1977, 31 (138), 443 4 Madsen, N. K. and Sincovec, R. F. PDECOL: General collocation software for PDE, ACN. Trans. Math. Softw. 1979, 5 (3), 3 2 6

APPENDIX

In this appendix we show a number of identities and inequalities that were used in Section 2. Throughout we assume that (v(x, t), O(x, t)) is a f'Lxed classical solution of (1.1)-(1.5) on [0, l]x [ 0 , ~ ) such that v(., t), Vx(., t), vt(., t), Vxx(., t), 0(., t), Ox(., t) are all in C°((0,~); L2(0, 1)) while Vxt(.,t) is in C°((0,~); L~(0, 1)) and vtt(., t) is in L~oc((0, ~); L2(0, 1)). We multiply (1.1) by Ox(X, t), integrate over x, integrate by parts and use (1.5) to obtain:

0.200 -

-0.200

REFERENCES

-

-0.600 -

1

-1.000

a(1, t)[V t(1, t) -- f(t)l -- | o(x, t) Vxt(X, t) dx 0

0.200

0.400

0.600

0.800

1.000

o

X

Figure 2. Numerical behaviour o f stress with steady boundary shearing force, 7 = 0.3, 0o = 1, Vo = X2/2, la(O) = OL and f(t) = sin(t)

= f o2x(x, t) dx J

Engineering Analysis, 1985, Vol. 2, No. 4

(A. 1)

209

Adiabatic shearing of one-dimensional thermoviscoelastic flows." N. C. Charalambakis and E. N. Houstis Integrating (A.1) over (0, t), integrating by parts and using (1.1)-(1.5) we arrive at: tl

tl

ff

'If

O2x(X,r) dx dr --

/~'(O(x, r)) g(O(x, r))

2

O0

la'(O(x, t)) #(O(x, t)) v3x(x, t) Vxt(X, t) la 1 - 4 ~t [/i(0(x, t))#(0(x, t)V4x(X, t)] - - 4 [~(0(x, t))) 2

+ u(O(x, t))u"(O(x, t))l u(O(x, t)) v6(x, t)

O0

Next, we use (1.2) and (1.3) to deduce:

1

x v4x(x, r) dx dr +

O(x,O

t

f u(Od~=fo2(x,r)dr

g(O(x, t)) V2x(X,t) dx o

Oo(x) t

=

f

(A.5)

0

1

o(1, r) Ox(1, r) dr +

o

I~(Oo(x))V~x(X) dx

(A.2)

On account of o(0, t) = 0 and Schwarz's inequality we have:

o

1

We multiply (1.1) by Oxt(X, t), integrate with respect to x over (0, 1), integrate by parts and use again (1.2), (1.3) and the boundary condition a(O, t) = 0 to obtain: 1

2 dt

f

(A.6)

o:(x, t) <<.f O2x(X,t) dx o

Finally, using o(0, t) = 0:

1

a2x(X,t) dx + o

U(O(x,t)) ~t(x, t) dx

1

o(1, t) = ox(l,

0

o

We note that, on account of (1.2) and (1.3):

EngineeringAnalysis, 1985, Vol. 2,No. 4

(A.7)

and integrating (1.1) over x and t:

+ f #'(O(x, t)) la(O(x, t)) Vax(X,t) Vxt(X, t) dx = [Vt(1, t ) - - f ( t ) ] ot(1, t)

t ) - f x o ~ ( x , t)dx o

1

210

(A.4)

1

t

t

1

(A.8)

(A.3) 0

0

0

0