Introduction to nonNewtonian ﬂuid mechanics
Josef Málek
October 2012
In this lecture we will try to answer the following questions.
Questions
Q1.What is mechanics?
Q2.What is a ﬂuid?
Q3.What is a Newtonian ﬂuid?
Q4.What are nonNewtonian ﬂuids?
Q5.Why mechanics of nonNewtonian ﬂuids is important?
Q6.What are the eﬀects that can not be described by classical models?
Q7.What will be the content of this course?
Q8.Is J.M.appropriate lecturer?
Q9.What are the unique features of the course?
Q10.Are there other approaches to describe such ”strange” features of materials?
Classical mechanics
Classical mechanics describes the evolution of mass particles which follow three Newton laws:
Law of inertia If the external force F = 0 then object in a state of uniform motion tends
to remain in that state of motion.
Law of force If there is an external force,it is equal to the time change of the momentum
mv
d
dt
(mv) = F:
Law of action and reaction For every action there is an equal and opposite reaction.
We show four examples of onedimensional linear spring.
1.Linear spring The displacement of the spring y satisﬁes ODE
m
d
2
y
dt
2
= F
E
+F
S
;
where F
S
is the spring force and F
E
is the external force.If the spring force F
S
= ky
and F
E
is gravitational force,then
m
d
2
y
dt
2
= mg ky:
The initialvalue problem has to be solved
d
2
y
dt
2
+
k
m
y = 0
y(0) = y
0
dy
dt
(0) = y
1
:
1
2.Spring with damping force We can add the damping force for the spring F
d
= b
dy
dt
d
2
y
dt
2
+
b
m
dy
dt
+
k
m
y = 0:
3.Generalized nonlinear spring Spring force is a nonlinear function of y:F
S
= F
S
(y),
for example Duﬃng spring F
S
(y) = y +y
3
.The wellknown Duﬃng oscillator with
the dumping satisﬁes
d
2
y
dt
2
+
b
m
dy
dt
1
m
(y +y
3
) = 0:
4.Further generalizations – implicit relation Finally,we can have a fully implicit re
lation for the between the spring force F
S
and the displacement y,and between the
damping force F
D
and the velocity
dy
dt
~
f (F
S
;y) = 0
~g
F
D
;
dy
dt
= 0:
Continuum mechanics
Instead of describing only a mass point we would like to describe whole body.To do this we
describe it as a continuum.Let B is a body,K
0
(B) a conﬁguration in the initial state and
K
t
(B) the conﬁguration at time t (see Figure).
The we deﬁne a motion which maps from K
0
(B) to K
t
(B):x = (X;t).Further we deﬁne
the following kinematical quantities
Velocity
v(X;t) =
@(X;t)
@t
=)v(x;t):= v(
1
(x;t);t)
Deformation gradient
F =
@(X;t)
@X
=
@x
@X
Left CauchyGreen tensor
B = FF
T
Right CauchyGreen tensor
C = F
T
F:
2
Balance equations
Balance of mass This balance is derived in the following way:(% density,v ﬂuid velocity)
m(P
t
) =
Z
P
t
% dx
d
dt
m(P
t
) = 0 8P
t
K
t
(B)
0 =
d
dt
m(P
t
) =
Z
P
t
%(x;t) dx =
Z
P
0
d
dt
%((X;t);t) det r
X
(X;t) dX
=
Z
P
0
@%
@t
+
@%
@
k
@
k
@t
+%div v
det r
X
(X;t) dX
=
Z
P
t
@%
@t
+r
x
%v +%div v
dx =
Z
P
t
@%
@t
+div (v%) dx
%
;t
+div (%v) = 0
Balance of linear momentum
(%v)
;t
+div (%v
v) div T = %b;
b volume force,T Cauchy stress tensor.
Balance of angular momentum
T = T
T
Balance of energy
@(%E)
@t
+div (%Ev) +div q div (Tv) = %b v +%r;
E energy,q thermal ﬂux,r thermal source.
Together with these balance equations we need to add some information to be able to
solve the problem:
1.Constitutive equations – characterize the response of class of materials at certain class
of processes.So far we have 8 equations and 14 unknowns still necessary to add some
constitutive equations.
2.Boundary conditions
3.Initial conditions
Comments:
Balance equations are formulated on averages (similar to calculus of variation) it sug
gests to work with weak formulation.
The last formulation is not equivalent in a weak sense (only classical sense)
3
What is a Newtonian ﬂuid?
We have to distinguish between compressible and incompressible ﬂuid.
Constitutive equation for the Cauchy stress tensor T for compressible Newtonian ﬂuid
T = p(%)I +2(%)D(v) +(%)(div v) I (1)
= p(%)I +2(%)D
d
(v) +
2(%) +3(%)
3
(div v) I
where A
d
:= A
1
3
(tr A) I.
Incompressibility
Deﬁnition 1.Volume of any chosen subset (at initial time t = 0) remains constant during
the motion.
0 =
d
dt
jP
t
j =
d
dt
Z
P
t
dx =
Z
P
t
div v dx =)div v = 0
in solid mechanics we can write det F = 1
Z
P
0
det F dX =
Z
P
t
dx
Consequences of incompressibility,balance of mass + incompressibility
@%
@t
+v r% +%div v = 0
div v = 0
=)
@%
@t
+v r% = 0
div v = 0
Nonhomogeneous incompressible Newtonian ﬂuid
T = mI +2(%)D
d
(v);(%) > 0:(2)
Homogeneous incompressible Newtonian ﬂuid
T = mI +2
D
d
(v);
> 0:(3)
The mean normal stress m= (tr T)=3 is the unknown.The pressure p in the compressible
model is in a similar situation as m in the incompressible model,but in the compressible
model the pressure has to be constitutively deﬁned,and it is not the unknown.
What is nonNewtonian ﬂuid?
A ﬂuid is nonNewtonian if its response to mechanical loading can not be describe either (1),
(2) or (3).
Concept of viscosity
Newton (1687):
"The resistance arising from the want of lubricity in parts of the ﬂuid,other
things being equal,is proportional to the velocity with which the parts are
separated from one another."
4
Deﬁnition 2 (Viscosity).Viscosity is the coeﬃcient proportionality between the shear stress
and the shearrate
Generalized viscosity
T
xy
=
F
z
A
v(h +y) v(y)
h
v
0
(y) = 2D
xy
g
():=
T
xy
()
;where = v
0
Simple shear ﬂow:v(x;y;z) =
0
@
v(y)
0
0
1
A
D=
1
2
0
@
0 v
0
0
v
0
0 0
0 0 0
1
A
:
bT
xy
= v
0
(y):
Experimental data show that the viscosity depends on the shearrate,pressure,etc.
One can understand the word ’proportional’ in the sense of arbitrary denpendence.Then
we can understand Newton’s statement in the sense of implicit relation
G(T
xy
;D
xy
;:::) = 0:
Physical dimension of the dynamic viscosity
[] =
N
m
2
m
m s
1
=
kg m s
2
m
2
s
= kg m
1
s
1
=:1000 cP;
and kinematic viscosity
%
= m
2
s
1
:
Kinematic viscosity of some ﬂuids
Fluid Viscosity [cP]
Air (at 18
C) 0.02638
Benzene 0.5
Water (at 18) 1
Olive oil (at 20) 84
Motor oil SAE 50 540
Honey 2000–3000
Ketchup 50000–70000
Peanut butter 150000–250000
Tar 3 10
10
Earth lower mantle 3 10
25
5
NonNewtonian phenomena and model for nonNewtonian
ﬂuids
Deﬁnition 3.The ﬂuid is nonNewtonian if and only if its behavior can not be described by
the model of NavierStokes ﬂuid
compressible T = p(%;)I +2(%;)D+(%;)(div v)I
incompressible T = pI +2(%;)D:
This deﬁnition is not useful.We describe the nonNewtonian phenomena.
NonNewtonian phenomena
(1) Shear thinning,shear thickening property
(2) Pressure thickening property
(3) Presence of activation/deactivation criteria (connected with the stress or with the shear
rate)
(4) Presence of nonzero normal stress diﬀerences in a simple shear ﬂow
(5) Stress relaxation
(6) Nonlinear creep
Consider a ﬂuid velocity v = (u(y);0;0) with the stress tensor in the form
T = pI +S:
We are interested in the shear rate jD
12
j and the shear stress jS
12
j.The shear rate jD
12
j is
in the literature denoted for example as neboor _ ,shear stress jS
12
j is denoted as or .
Ad (1) Shear thinning,shear thickening property
We are interested in the dependence of jS
12
j on jD
12
j.Let us deﬁne a generealized
viscosity
g
as
g
(jD
12
j):=
jS
12
(jD
12
j)j
jD
12
j
:
Newtonian ﬂuid
For the Newtonian ﬂuid it holds S
12
= D
12
and D
12
= u
0
(y),so jS
12
j is a linear
function of jD
12
j with the proportionality coeﬃcient (see Figure 1).
In the case of Newtonian ﬂuid
g
= = konst:,the graph
g
vs._ is in the Figure 2.
NonNewtonian ﬂuid
Shear thickening property (dilatant ﬂuids) means that jS
12
j is a superlinear function
of jD
12
j (Figure 3).The generalized viscosity
g
is a increasing function,in the Figure
4 is the example where the generalized viscosity is degenerate at the beginning,usually
the generealized visosity is positive in zero.
Shear thinning property (pseudoplastic ﬂuids) means that jS
12
j is a sublinear function
of jD
12
j (Figure 5).The generalized viscosity
g
is a decreasing function,in the Figure
6 is the example where the generalized viscosity is singular at the beginning,usually
the generealized visosity is ﬁnite in zero.
Generally we can say that there is no reason that jS
12
j is monotone function of jD
12
j
1
.
1
It does not have to be even a function.
6
Figure 1:Newtonian ﬂuid,
(shear stress/shear rate)
Figure 2:Newtonian ﬂuid,
(generalized viscosity/shear rate)
Figure 3:Shear thickening,
(shear stress/shear rate)
Figure 4:Shear thickening,
(generealized viscosity/shear rate)
Figure 5:Shear thinning,
(shear stress/shear rate)
Figure 6:Shear thinning,
(generalized viscosity/shear rate)
The best known models are powerlaw ﬂuids where the stress tensor look like
T = pI + 2
jD(v)j
r2

{z
}
g
(jDj)=~
g
(jDj
2
)
D(v):
For r = 2 (3) we speak about the Newtonian ﬂuid,for r > 2 (1) shear thickening ﬂuid
and for 1 < r < 2 (2) shear thinning ﬂuid,see Figures 7 and 8.
Shear thinning ﬂuid is for example the blood.Shear thickenning ﬂuid are for example
some painting colors,sauces (dressings),ketchups.Viscosity is declared as a quality
parameter for some industrial products.
7
Figure 7:Shear stress/shear rate.
Figure 8:Generalized viscosity/shear rate.
Ad (2) Pressure thickening property The generalized viscosity is not constant as in the
case of Newtonian ﬂuids
2
(see Figure 9) but it is a function of pressure p (see Figure
10).Experimental data say that the viscosity is an incresing function of the pressure,
not decreasing.
Figure 9:Newtonian ﬂuid
Figure 10:Viscosity dependent on pressure
Around the year 1890 Barus proposed a model
(p) =
0
exp(p); > 0:
The Nobel prize was given for the experimental validation of the pressure dependence.
These models are used for example in journal bearings where the pressure diﬀerences
are very high.
Compressible vs.incompressible ﬂuid For the compressible ﬂuid the pressure
p = ~p(%) is a function of the density %.Let us suppose that we can invert this relation
% = ~%(p),then the formulation
T = ~p(%)I +2(%)D+(%)(div v)I
is equivalent with
T = pI +2~(p)D+
~
(p)(div v)I;
which is the model with the viscosity dependent on the pressure.
Een if the material is compressible (but not extremely),the incompressible pressure
dependence model is used.Such models are diﬃcult to treat both analytically and
numerically.
2
Already Stokes in his work wrote about the dependence of the viscosity on the pressure.For the simplicity
he further supposed the viscosity to be constant.
8
Ad (3) Presence of activation/deactivation criteria
The ﬂuid starts to ﬂow when it reaches the critical value of the stress – called the yield
stress.in the case that the following dependence of the stress on the shear rate is linear,
respectively nonlinear,we call the ﬂuid Bingham ﬂuid,respectively HerschelBulkley
ﬂuid
3
.
The standard formulation is the following:
D= 0,jSj ;
T = pI +2
D
jDj
+ ~
g
(jDj
2
)D,jSj > :
If ~
g
is constant,it is the Bingham ﬂuid,if not,it is HerschelBulkley ﬂuid (see Figure
11).
Further example is the ﬂuid where the response is connected with the chemical pro
cesses.At some value jD
12
j the ﬂuid locks and does not ﬂow,see Figure 12.
Figure 11:HerschelBulkley
Figure 12:Locking
Except these explicit relations between the shear stress jS
12
j and shear rate jD
12
j we
can consider a general implicit relation in the form
g(jS
12
j;jD
12
j) = 0:
We do not have to just compute the tensor S from the knowledge of D and insert it
into the balance of linear momentum,it is also possible to deﬁne a tensor function G
G(D;S) = 2(jDj
2
)
+(jSj
)
+
D(jSj
)
+
S
that maps
R
dd
sym
R
dd
sym
!R
dd
sym
and this equation add to the system of balance equations.The advantage is that G is
continuous and no special tools has to be used.The disadvantage is that the system
we are solving is bigger.
Ad (4) Presence of nonzero normal stress diﬀerences in a simple shear ﬂow
In the three dimensional space we deﬁne three viscometric functions:
viscosity
N
1
:= T
11
T
22
ﬁrst normal stress diﬀerence
N
2
:= T
22
T
33
second normal stress diﬀerence
In the simple shear ﬂow we suppose a ﬂuid velocity in the form v = (u(y);0;0).
3
We call it ﬂuid although the deﬁnition of the ﬂuid says that the ﬂuid can not sustain the shear stres.
9
Newtonian ﬂuid Let us compute the Cauchy stress T
D=
1
2
0
@
0 u
0
(y) 0
u
0
(y) 0 0
0 0 0
1
A
;T = pI +2D=
0
@
p u
0
(y) 0
u
0
(y) p 0
0 0 p
1
A
:
We can see that N
1
= N
2
= 0,and so there are no normal stress diﬀerences in the
Newtonian.
NonNewtonian ﬂuid NonNewtonian ﬂuids like for example the amyl dissolved in
the water,exhibit the presence of nonzero stress diﬀerences in the simple shear ﬂow.
If we press the ﬂuid in one direction,it reacts in other direction,usually perpendicular.
These eﬀects are connected with this property:
Rod climbing
Die swell
Delayed die swell
Let us denote D the diameter of the pipe and D
E
the biggest diameter of the ﬂuid
after ﬂowing out from the pipe.For the Newtonian ﬂuids it was experimentally
found that
D
E
D
= 1:13:
For nonNewtonian ﬂuids is this ratio equal even to four and it holds
D
E
D
=
"
0:13 +
1 +
1
2
N
1
2S
w
2
#
1=6
;
where S is the shear stress and a subscript w means that it is a value on the wall.
Flow through the sloping channel
In the case of Newtonian ﬂuid the surface of the ﬂuid is smooth (free boundary),
in the case of nonNewtonian ﬂuid the surface is not smooth,there are bumps.
Inverted secondary ﬂow
Consider a container with a ﬂuid,the ﬂuid is covered with no gap between the
ﬂuid and the cover.Now we start to rotate with the cover.The secondary ﬂow
creates in the ﬂuid,this secondary ﬂow can be seen in the crosssection.In the
case of nonNewtonian ﬂuid is the direction of the secondary ﬂor oposite than in
the case of Newtonian ﬂuid and the direction depends on the size of N
2
.
Ad (5) Stress relaxation
Stress relaxation test is the following (see Figure 13):At time t = 0 we deform the
material with constant relative prolongation"and then we study the shear stress .
Consider two basic materials:linear elastic spring and linear dashpot (two coaxial
cylinders with almost same radii,between these cylinders is a Newtonian ﬂuid with the
viscosity ).
Linear spring is described with Hooke law = E",where E is the elastic modulus.
Linear dashpot is described by the relation = _",where is the viscosity.
The result of the stress relaxation test is in the case of linear spring in the Figure 14,in
the case of linear dashpot in the Figure 15.In the case of dashpot is the stress singular
in zero.
Most of materials are the combination of these two basic materials,their response is
depicted in the Figure 16
10
)
Figure 13:Stress relaxation test
Figure 14:Linear spring
Figure 15:Linear dashpot
Figure 16:Viscoelastic solid Viscoelastic ﬂuid
Ad (6) Nonlinear creep
Nonlinear creep test is the following (see Figure 17):At time t = 0 we deform the
material with constant stress ,at time t = t
we turn oﬀ the stress and study the
relative prolongation".
The result of the creep test is in the case of linear spring in Figure 18,in the case of
linear dashpot in Figure 19.
Again,most of materials are the combination of two basic materials,in Figure 20 is
the resultfor the material similar to viscoelastic solid a viscoelastic ﬂuid.
Except the algebraic models there exist also ratetype models,integral models and
stochastic models.A detailed overview of nonNewtonian models is in the following
section.
11
)
Figure 17:Nonlinear creep test
Figure 18:Linear spring
Figure 19:Linear dashpot
Figure 20:Viscoelastic solid Viscoelastic ﬂuid
NonNewtonian models
Models with variable viscosity
General form:
T = pI +2(D;T)D

{z
}
S
(4)
Particular models mainly developed by chemical engineers
.
Shear dependent viscosity
Ostwald–de Waele power law [24],[12]
(D) =
0
jDj
n1
(5)
12
Fits experimental data for:ball point pen ink,molten chocolate,aqueous dispersion of
polymer latex spheres
Carreau Carreau–Yasuda [8],[36]
(D) =
1
+
0
1
(1 +jDj
2
)
n
2
(6)
(D) =
1
+(
0
1
) (1 +jDj
a
)
n1
a
(7)
Fits experimental data for:molten polystyrene
Eyring [14],[29]
(D) =
1
+(
0
1
)
arcsinh(jDj)
jDj
(8)
(D) =
0
+
1
arcsinh(
1
jDj)
1
jDj
+
2
arcsinh(
2
jDj)
2
jDj
(9)
Fits experimental data for:napalm(coprecipitated aluminumsalts of naphthenic and palmitic
acids;jellied gasoline),1% nitrocelulose in 99% butyl acetate
Cross [11]
(D) =
1
+
0
1
1 +jDj
n
(10)
Fits experimental data for:aqueous polyvinyl acetate dispersion,aqueous limestone suspen
sion
Sisko [35]
(D) =
1
+jDj
n1
(11)
Fits experimental data for:lubricating greases
Models with pressure dependent viscosity
Barus [1]
(T) =
ref
exp (p p
ref
) (12)
Fits experimental data for:mineral oils
4
,organic liquids
5
Models with stress dependent viscosity
Ellis [21]
(T) =
0
1 +jT
j
n1
(13)
Fits experimental data for:0.6% w/w carboxymethyl cellulose (CMC) solution in water,
poly(vynil chloride)
6
4
[20]
5
[6]
6
[33]
13
Glen [16]
(T) = jT
j
n1
(14)
Fits experimental data for:ice
Seely [34]
(T) =
1
+(
0
1
) exp
jT
j
0
(15)
Fits experimental data for:polybutadiene solutions
Blatter [25],[5]
(T) =
A
(jT
j
2
+
2
0
)
n1
2
(16)
Fits experimental data for:ice
Models with discontinuous rheology
Bingham Herschel–Bulkley [4],[17]
jT
j >
if and only if T
=
D
jDj
+2(jDj)D
jT
j
if and only if D= 0
(17)
Fits experimental data for:paints,toothpaste,mango jam
7
Diﬀerential type models
Rivlin–Ericksen ﬂuids
Rivlin–Ericksen [31],[32]
General form:
T = pI +f(A
1
A
2
A
3
:::) (18)
where
A
1
= 2D (19a)
A
n
=
dA
n1
dt
+A
n1
L+L
T
A
n1
(19b)
where
d
dt
denotes the usual Lagrangean time derivative and L is the velocity gradient.
Criminale–Ericksen–Filbey [10]
T = pI +
1
A
1
+
2
A
2
+
3
A
2
1
(20)
Fits experimental data for:polymer melts (explains mormal stress diﬀerences)
Reiner–Rivlin [30]
T = pI +2D+
1
D
2
(21)
Fits experimental data for:N/A
7
[2]
14
Rate type models
Maxwell,Oldroyd,Burgers
Maxwell [22]
T = pI +S (22a)
S +
1
O
S = 2D (22b)
O
M:=
dM
dt
LMML
T
(23)
Fits experimental data for:N/A
OldroydB [23]
T = I +S (24a)
S +
O
S =
1
A
1
+
2
O
A
1
(24b)
Fits experimental data for:N/A
Oldroyd 8constants
[23]
T = I +S (25a)
S +
1
O
S +
3
2
(DS +SD) +
5
2
(tr S) D+
6
2
(S:D) I (25b)
=
D+
2
O
D+
4
D
2
+
7
2
(D:D) I
(25c)
Fits experimental data for:N/A
Burgers [7]
T = I +S (26a)
S +
1
O
S +
2
OO
S =
1
A
1
+
2
O
A
1
(26b)
Fits experimental data for:N/A
Giesekus [15]
T = I +S (27a)
S +
O
S
2
S
2
= D (27b)
Fits experimental data for:N/A
15
PhanThien–Tanner [26],[27]
T = I +S (28a)
Y S +
O
S +
2
(DS +SD) = D (28b)
Y = exp
"
tr S
(28c)
Fits experimental data for:N/A
Johnson–Segalman [19]
T = pI +S (29a)
S = 2D+S
0
(29b)
S
0
+
dS
0
dt
+S
0
(DaD) +(DaD)
T
S
0
= 2D (29c)
Fits experimental data for:spurt
Johnson–Tevaarwerk [18]
T = pI +S (30a)
O
S +sinh
S
s
0
= 2D (30b)
Fits experimental data for:lubricants
Integral type models
Kaye–Bernstein–Kearsley–Zapas [3],[9]
T =
Z
t
=1
@W
@I
C+
@W
@II
C
1
d (31)
Fits experimental data for:polyisobutylene,vulcanised rubber
16
Entropy and Newtonian ﬂuids
Deﬁnition 4.There exists a speciﬁc density of entropy (further just entropy) which is a
function of state variables,
= ~(y
0
;y
1
;:::);
where y
0
= e is a density of internal energy (further just energy),and it holds:
(i)
@~
@e
> 0 ) e = ~e(;y
1
;y
2
;:::),denote a temperature :=
@~e
@
(ii) !0
+
,if !0
+
(iii) S(t) =
Z
(%)(t;x) dx!S
max
for t!+1,if
is mechanically and thermally (ener
getically) isolated
For Newtonian ﬂuids it holds
= ~(e;%),e = ~e(;%):
Balance equations:
_% = %div v;
% _v = div T+%b;
%
_
E = div(Tv +q) +%b v:
Diﬀerentiate e,multiply by %,
%_e = %
@~e
@
_ +%
@~e
@%
_%
and insert it into the balance of linear momentum,
%b v +div(Tv +q) (div T) v %b v = %
@~e
@
_ %
2
@~e
@%
div v:
Now denote
=
@~e
@
;p = %
2
@~e
@%
and obtain
% _ = T rv +div q +pdiv v
symetrie T
= T D+div q +pdiv v:
Using the deviatoric parts we get
% _ = T
d
D
d
+
1
3
(tr T)
1
3
(div v) I I
{z}
3
+div q +pdiv v
= T
d
D
d
+
1
3
(tr T) +p
div v +div q:
Dividing by we get
% _ div
q
=
1
T
d
D
d
+
1
3
(tr T) +p
div v +
q r

{z
}
:
17
is called rate of entropy production.The second law of thermodynamics says that 0,
and so
:=
= % _ div
q
0:
If we denote m= (tr T)=3 the mean normal stress,then it holds
= T
d
D
d
+
1
3
(tr T) +p
div v+
q r
=
T
d
;m+p;q

{z
}
termodynamical ﬂuxes
D
d
;div v;
r

{z
}
termodynamical aﬁnities
(32)
What is for NavierStokes equations?
Compressible NavierStokesFourier equations
T = p(%;)I +2(%;)D+(%;)(div v)I;
q = K(%;)r:
It holds
m+p =
2 +3
3
div v
and
T
d
= 2D
d
:
Insert it into (32).
Option 1
= 2jD
d
j
2
+
2 +3
3
(div v)
2
+K
jrj
2
:
If 0,2 +3 0 a K 0,then 0.
Option 2
=
1
2
jT
d
j
2
+
3
2 +3
(m+p)
2
+
1
K
jqj
2
:
If > 0,2 +3 > 0 a K > 0,then > 0.
We ask:When is = 0?
Option 1 We do not want any restriction for the velocity or temperature,that is why
= 2 +3 = K = 0.
Option 2 It is possible only if T
d
= 0;m= p;q = 0.
Incompressible NavierStokesFourier equations
Option 1
= 2jD
d
j
2
+K
jrj
2
:
If 0 a K 0,then 0.
Option 2
=
1
2
jT
d
j
2
+
1
K
jqj
2
:
If > 0 a K > 0,then > 0.
18
Remark
=
@e
@
;p = %
2
@e
@%
:
Gibbs equation
d = de +pd
1
%
;
from this equation the relation for the pressure can be derived.This equation holds for
compressible gas,there is no reason why it should hold in general for more complicated
materials.
Question:Is there a way how from the knowledge for constitutive relations for and
get the constitutive relations for T and q?
Principle of maximimzation of the rate of entropy production
Our aim is to derive NavierStokes equations.It holds
= T
d
D
d
+(m+p) div v +
q r
;(33)
we choose a constitutive relation for
~
=
~
(D
d
;div v;r) = 2jD
d
j
2
+
2 +3
3
(div v)
2
+K
jrj
2
:
There is many possibilities when =
~
,we choose the one when
~
is maximal.That is why
we maximize it over D
d
;div v;r such that (33) holds
max
(33)+D
d
;div v;r
~
(D
d
;div v;r):
This heuristic principle is called principle of maximimzation of the rate of entropy production.
The maximization is done using Lagrangeo multipiers.We deﬁne the Lagrange function
L =
~
(D
d
;div v;r) +
~
(D
d
;div v;r)
T
d
D
d
+(m+p) div v +
q r
;
and we are looking for its extreme,i.e.
@L
@D
d
= 0;
@L
@ div v
= 0;
@L
@r
= 0:
19
We obtain
(1 +)
@
~
@D
d
= T
d
;
(1 +)
@
~
@ div v
= (m+p);
(1 +)
@
~
@r
=
q
:
(34)
Now we eliminate in this way.We multiply the ﬁrst equation by D
d
,the second by div v
and the third by r,
1 +
=
~
@
~
@D
d
D
d
+
@
~
@ div v
div v +
@
~
@r
r
=
1
2
:
If we insert this result into (34) and do the diﬀerentiation,we obtain NavierStokes equations
T
d
= 2D
d
;
m+p =
2 +3
3
div v;
q
= K
r
:
Maximization
~
for incompressible NavierStokes equations
For isothermal incompressible materials it holds
tr D= div v = 0; =
= const:
Constraint for is
= T D:(35)
First,let
~
is given by
~
= 2jDj
2
:
We use the principle of maximization of the rate of entropy production
max
(35)+D+tr D=0
~
(D):
We deﬁne the Lagrange function
L =
~
(D) +
1
~
(D) (T D)
+
2
tr D
and maximize it with respect to D,
(1 +
1
)
@
~
@D
=
1
T
2
I:
Multiply this equation by D and get
1 +
1
1
=
T D
@
~
@D
D
=
1
2
;
further take the trace of this equations and get
3
2
=
1
tr T )
2
1
=
1
3
tr T:
20
We have
T = 2D+
1
3
(tr T)I:
Incompressible NavierStokes equations can be also derived in other way.Let
~
=
1
2
jT
d
j
2
(q = 0) (36)
and we are maximizing
~
max
(36)+T
d
~
(T
d
):
We obtain
D
d
= D=
1
2
T
d
=
1
2
T
1
3
(tr T)I
;
which is less usual,but equivalent formulation of incompressible NavierStokes equations.
Integral form of entropy production
Integrate the balance of energy (b = 0,r = 0,v n = 0) over
d
dt
Z
%
e +
jvj
2
2
dx =
Z
div (Tv +q) dx =
Z
@
(Tv +q) n dS = 0
In order to have conservation of energy it suﬃces to assume pointwise
q n = Tv n = Tn v = ((Tn)
+(Tn)
n
) (v
+(v n)n) = (Tn)
v
For the total entropy S =
R
% dx we also have
d
dt
S(t) =
Z
!
1
T
:D
+:::
Z
@
q n
d
dt
S(t) =
T
;D
L
2
(
)
+(m+p;div v)
L
2
(
)
+
q;
r
L
2
(
)
+((Tn)
;v
)
L
2
(
)
The simplest proportional relation
T
d
= 2D
d
;
m+p =
2 +3
3
div v;
Q
= K
r
(Tn)
= v
with
0;2 +3 0;K 0; 0
guarantee the second law of thermodynamics.
In general we have
G
1
T
d
;D
d
= 0;
G
2
(m+p;div v) = 0;
G
3
Q
;
r
= 0;
G
4
((Tn)
;v
) = 0:
If D
= 0 or div v = 0 or r = 0 or v
= 0 we get constrained material (rigid),resp.
incompressible,resp.isothermal process,resp.noslip boundary condition.
No dissipation if T
= 0 and q = 0 and m= p and (Tn)
= 0.
21
Remark:If you have dissipation you have energy estimates
Cauchy model for ﬁnite elasticity
Kinematical quantities
F
K
0
=
@
@X
;
_
F
K
0
= LF
K
0
) B
K
0
= F
K
0
F
T
K
0
;
where B
K
0
is the left CauchyGreen tensor,we can compute its time derivative
_
B
K
0
= LB
K
0
+B
K
0
L
T
a tr
_
B
K
0
= 2D B
K
0
:
Let us have the entropy
= ~(e;B
K
0
),e = ~e(;B
K
0
);
For example the neoHookean material fullﬁls
~e = e
0
() +
2
(tr B
K
0
3):
Insert e into the balance of energy
% _ = div(Tv +q) div T v B
K
0
D
= T DB
K
0
D+div q
= (TB
K
0
) D+div q:
We get
% _ +div
q
=
1
(TB
K
0
) D+
q r

{z
}
produkce entropie
:
Elastic material does not disipate,that is why the rate of entropy production has to be zero,
and so
T = B
K
0
;
and for the incompressible material we have
T = I +B
K
0
:
Korteweg ﬂuids
We will use the principle of maximization of the rate of entropy production to derive the
model derived in 1901 by Dutch mathematician Korteweg.This model is used for describing
capilar eﬀects,multiphase and granulated materials.The model is formulated by the following
equations
_% = %div v;
% _v = div T;
T = p(%)I +2(%)D(v) +(%)(div v)I +K;
K=
2
((%
2
) jr%j
2
)I (r%
r%):
22
In the paper [13] Dunn and Serrin showed that this model is thermodynamically compatible.
We know that if we suppose this constitutive relation
= ~(e;%);
then for compressible NavierStokes equations we obtain
% _ +div
q
=
1
T
d
D
d
+(m+p) div v +
q r
=
1
T
d
dis
D
d
+t
dis
div v +
q
dis
r
:
We introduce disipative terms (with subscript dis) and write the constitutive relation for the
rate of entropy production
=
~
(T
d
dis
;t
dis
;q
dis
) =
1
2
jT
d
dis
j
2
+
2
2 +3
t
2
dis
+
1
K
jq
dis
j
2
and using the principle of maximization of the rate of entropy production we get Navier
StokesFourier system.
For the Korteweg model we consider the entropy dependent on the gradient of density %
= ~(e;%;r%) $e = ~e(;%;r%):
Diﬀerentiate
%_e = %
@~e
@
_ +%
@~e
@%
_% +%
@~e
@r%
_
r%;
we know that _% = %div v,and need to compute
_
r%,
_
r% = r(%div v) (rv)r%:
Let us denote
=
@~e
@
;p = %
2
@~e
@%
;
therefore we have
% _ = %
_
E % _v v p
_%
%
%
@~e
@r%
_
r% =
= div(Tv +q) div T v +pdiv v +%
@~e
@r%
r(%div v) +%
@~e
@r%
r%
rv =
= Trv+
%
@~e
@r%
r%
rv+div q+pdiv v+div
%
2
div v
@~e
@r%
%(div v) div
%
@~e
@r%
:
The principle of objectivity says that e is a function of the size of r%,
e = ~e(;%;jr%j);
this implies that
@~e
@r%
r%
is a symmetric tensor.We obtain
% _ =
T+%
@~e
@r%
r%
D+
p %div
%
@~e
@r%
div v+div
q +%
2
(div v)
@~e
@r%
=
=
T
d
+%
@~e
@r%
r%
d
!

{z
}
T
d
dis
D
d
+
p %div
%
@~e
@r%
+m+
%
3
@~e
@r%
r%

{z
}
t
dis
div v+div
q +%
2
(div v)
@~e
@r%

{z
}
q
dis
:
23
We get the equation
% _ div
q
dis
=
1
T
d
dis
D
d
+t
dis
div v +
q
dis
r
:
= 0,pokud T
d
dis
= 0;t
dis
= 0;q
dis
= 0,pak
T = T
d
+mI = %
@~e
@r%
r%
d
p(:::)I
%
3
@~e
@r%
r%
I +%div
%
@~e
@r%
I
= p(e;%;r%)I +%div
%
@~e
@r%
I %
@~e
@r%
r%
:
The special choice of internal energy
e = ~e(;%;r%) = e
0
(;%) +
2%
jr%j
2
) %
@~e
@r%
= r%
gives the Korteweg model
T = p(:::)I +% div %I (r%
r%)
= pI +
2
(%
2
) jr%j
2
I (r%
r%):
From the choice of constitutive equations for the rate of entropy production
=
~
(T
d
dis
;t
dis
;q
dis
) =
1
2
jT
d
dis
j
2
+
3
2 +3
t
2
dis
+
1
K
jq
dis
j
2
and the choice
8
T
d
dis
= 2D
d
;
t
dis
=
2 +3
3
div v;
q
dis
= Kr
we see that 0 and the second law of thermodynamics is satisﬁed.Then we get the
Korteweg model
q = Kr %
2
div v
@~e
@r%
= Kr %(div v)r%;
T = pI +2D+(div v)I +
2
(%
2
) jr%j
2
I (r%
r%):
Viscoelastic materials
Our aimis to model materials capable of describing nonNewtonian phenomena like nonlinear
stress relaxation,nonlinear creep or the normal stress diﬀerences in simple shear ﬂow.Such
materials are for example geomaterials,biological materials,polymers and other chemical
products like food for example.
Repeat what is the stress relaxation test and creep test.
Stress relaxation test Stress relaxation test is the following (see Figure 21):At time
t = 0 we deform the material with constant relative prolongation"and then we study the
shear stress .For linear spring it holds
= E";
24
)
Figure 21:Stress relaxation test
Figure 22:Linear spring
Figure 23:Linear dashpot
for the linear dashpot it holds
= _":
Let us deﬁne the stress relaxation function
G(t) =
(t)
"
0
:
Real materials exhibit the following bahavior
Figure 24:Viscoelastic solid Viscoelastic ﬂuid
Creep test Creep test is the following (see Figure 25):At time t = 0 we deform the
material with constant stress ,at time t = t
we turn oﬀ the stress and study the relative
prolongation".
8
It is possible to use the principle of maximization of the rate of entropy production,but for the simplicity
we can guess it,because the rate of entropy production is quadratic.
25
)
Figure 25:Creep test
Figure 26:Linear spring
Figure 27:Linear dashpot
Let us deﬁne the creep test function
J(t) =
"(t)
0
:
Real materials exhibit the following bahavior
Figure 28:Viscoelastic solid Viscoelastic ﬂuid
The result of stress relaxation test is in the case of linear spring in Figure 26,in the case
of linear dashpot in Figure 27.
Creation of simple models by combining linear springs and linear
dashpots
Maxwell element We get the Maxwell element if we connect the linear spring and linear
dashpot in series.Denote F resp.F
S
resp.F
D
the force in the whole element resp.in the
26
Figure 29:Maxwell element
spring resp.in the dashpot,and resp.
S
resp.
D
o the prolongation of the whole
element resp.of the spring resp.of the dashpot.Then for the linear spring holds
F
S
= E
S
;
S
= E"
S
;
for the linear dashpot
F
D
=
_
D
;
D
= _"
D
:
For the Maxwell element we do the computation in details.First we derive the constitutive
relation.In the series connection is the stress in the spring and dashpot the same,i.e.
D
=
S
,the whole prolongation is the sum of the prolongation of the spring and the
dashpot,i.e. =
S
+
D
.We substitute from the constitutive relation for the linear
spring and the linear dashpot
_
=
_
S
+
_
D
=
_
F
S
E
+
F
D
;
and get the constitutive equation fot the Maxwell element
E_"= _ +
E
,p
1
_ +p
0
= q
1
_":
From the initial condition we obtain
p
1
(0+) = q
1
"(0+):
Now copmute how the Maxwell element behaves in the creep test.Let (t) =
0
,then
"(t) =
0
q
1
(p
1
+p
0
t)
and the creep function is equal to
J(t) =
"(t)
0
=
1
q
1
(p
1
+p
0
t):
The response is linear it is not a satisfacory result.
Further compute the stress relaxation test for Maxwell element.Let"(t) ="
0
H(t),where
H(t) is the Heavisideova function
(t) =
q
1
p
1
"
0
e
p
0
p
1
t
and the stress relaxation function is equal to
G(t) =
(t)
"
0
=
q
1
p
1
e
p
0
p
1
t
:
This response looks reasonably.
27
Now compute how the stress depends on th deformation"
_
(t)e
p
0
p
1
t
=
q
1
p
1
_"e
p
0
p
1
t
(t) = (0+)e
p
0
p
1
t
+
q
1
p
1
Z
t
0
_"(t)e
p
0
p
1
(t)
d
initial cond.
=
q
1
p
1
"(0+)e
p
0
p
1
t
+
q
1
p
1
Z
t
0
_"(t)e
p
0
p
1
(t)
d
="(0+)G(t) +
q
1
p
1
Z
t
0
_"(t)G(t ) d:
In the higher dimension this result corresponds to the fact that we obtained the integral
relation between T and D.
Now compute how"depends on .
q
1
"(t) = q
1
"(0+) +p
1
Z
t
0
_() d +p
0
Z
t
0
() d
= q
1
"(0+) +
Z
t
0
_()(p
1
+p
0
(t )) d +p
0
(0+)t:
Using initial condition we obtain
"(t) =
p
1
q
1
+
p
0
q
1
t
(0+) +
Z
t
0
_()
p
1
q
1
+
p
0
q
1
(t )
d
= J(t)(0+) +
Z
t
0
_()J(t ) d:
We obtained the ratetype model and two integral type models,all equivalent.All models
hold in one dimensional space,but it is not clear howto generalize theminto three dimensional
space.
1.Linear dashpot is not an approximative theory,but the linear spring is.
2.What to do if the dashpot or the spring depend on the deformation in nonlinear way?
3.Partial time derivative is neither material time derivative nor objective time derivative,
the generalizations are not unique.
4.Many nonlinear 3D models can reduce in 1D to the same equation.So many models
can be used to capture 1D experimental data.
KelvinVoigt element We get the KelvinVoigt element if we paralelly connect the spring
and the dashpot.
The prolongations are the same for the spring and the dashpot.The total force is equal
to the sum of the forces in spring and dashpot,so
F = F
S
+F
D
; =
S
=
D
:
Using constitutive relations for the linear spring and the linear dashpot we get
p
0
= q
0
"+q
1
_";
"(0+) = 0:
The stress relaxation function is equal to
G(t) =
q
0
p
0
+
q
1
p
0
(t);
28
Figure 30:KelvinVoigtů element
this is disadvantage for KelvinVoigt.The creep function is
J(t) =
p
0
q
0
1 e
q
0
q
1
t
;
which is the advantage of KelvinVoigt element
Both derived models were onedimensional,we show teh modiﬁcation into the higher
dimension.Suppose that the ﬂuid is incompressible,tr D= 0.Consider the stress tensor in
the form
T = pI +S
and say that S satisﬁes the same equation as derived earlier,let us try
S +
_
S = 2D:
Suppose we have two coordinate systems – starred and unstarred,the starred one is
related to the unstarred by
x
= Q(x x
0
) +c;
where Q is orthogonal.We want the tensor S to transform as
S
= QSQ
T
:
In the starred system the material satisﬁes the equation
S
+
_
S
= 2D
;
which is equivalent to
QSQ
T
+
d
dt
QSQ
T
= 2QDQ
T
:
If it holded
QSQ
T
+
d
dt
QSQ
T
= Q
S +
_
S
Q
T
;
everything would be OK,but this does not hold because
d
dt
QSQ
T
=
_
QSQ
T
+Q
_
SQ
T
+QS
_
Q
T
and the sum of the ﬁrst and last term can never be equal zero.So,our try is bad –
_
S is not
objective.Let us denote the objective derivative of S by
S
satisfying
S
= Q
S
Q
T
:
29
Deﬁne now
S
=
dS
dt
LS SL
T
+(tr L)S
and verify that this derivative is objective.We want
S
= Q
S
Q
T
= Q
dS
dt
LS SL
T
+(tr L)S
Q
T
:
Compute.
dS
dt
L
S
S
(L
T
)
+(tr L
)S
=
QSQ
T
=
=
d
dt
QSQ
T
L
QSQ
T
QSQ
T
(L
T
)
+(tr L)QSQ
T
:
We need to know what is L
.After some algebraic manipulation it can be computed that
L
=
_
QQ
T
+QLQ
T
:
Inserting this into the previous calculation we get
d
dt
QSQ
T
L
QSQ
T
QSQ
T
(L
T
)
+(tr L)QSQ
T
=
=
_
QSQ
T
+Q
_
SQ
T
+QS
_
Q
T
_
QQ
T
+QLQ
T
QSQ
T
QSQ
T
_
QQ
T
+QLQ
T
+tr
_
QQ
T
+QLQ
T
QSQ
T
=
= Q
_
SQ
T
QLSQ
T
QSL
T
Q
T
+tr
_
QQ
T
+QLQ
T
QSQ
T
=
= Q
_
S LS SL
T
Q
T
+Q(tr L)SQ
T
;
because the trace is cyclic and
0 =
_
QQ
T
+Q
_
Q
T
:
The derivative is objective even if there is not the last term with the trace
O
S
=
dS
dt
LS SL
T
:
This derivative is called Oldroyd upper convected derivative.If we use it in the Maxwell
element we obtain the Maxwell model
T = pI +S;S +
O
S= 2D:
Incompressible ratetype ﬂuid models – derivation of ther
modynamically compatible models
Motivaation and the aim Combining two linear dashpots and one linear spring we obtain
the relation between the shear stress and shear rate
p
1
_
T
xy
+p
0
T
xy
= q
1
_
D
xy
+q
0
D
xy
:
30
In the previous section we explained why it is not trivial to generalize this model into three
dimensional space.A similar model is a Yeleswerapu model for blood:
div v = 0;
% _v = div T+%b;
T = pI +S;
S +
1
_
S LS SL
T
= (jDj)D+
2
(
_
DLDDL
T
);
(jDj) =
1
+(
0
1
)
1 +ln(1 +jDj)
1 +jDj
:
But is this model thermodynamically compatible satisfying the second law of thermodynam
ics?
First we derive three dimensional KelvinVoigt model.
KelvinVoigt model
Suppose the entropy in the form = ~(e;%;B
k
R
),e = ~e(;%;B
k
R
) = ^e(;%;tr B
k
R
) and
deﬁne kinematical quantities.First the deformation gradient
F
k
R
=
@
k
R
@X
left and right CauchyGreen tensor
B
k
R
= F
k
R
F
T
k
R
;C
k
R
= F
T
k
R
F
k
R
:
Compute
%
_
E % _v v = %_e = %
@^e
@
_ +%
@^e
@%
_% +%
@^e
@ (tr B
k
R
)
I
_
B
k
R
:
We know that
_
F
k
R
= LF
k
R
which implies
_
B
k
R
= LB
k
R
+B
k
R
L
T
) I
_
B
k
R
= 2B
k
R
D:
Continue computing
% _ div q = T D+pdiv v 2%
@^e
@ (tr B
k
R
)
B
k
R
D=
=
T2%
@^e
@ (tr B
k
R
)
B
k
R
d
D
d
+
m
2
3
%
@^e
@ (tr B
k
R
)
tr B
k
R
+p
div v:
Divide by temperature and get
% _ div
q
=
1
T
d
dis
D
d
+t
dis
div v +
q r

{z
}
:
The incompressibility and constant temperature gives
= T
d
dis
D
d
=
T2%
@^e
@ (tr B
k
R
)
B
k
R
d
D
d
:
Further for a general =
~
(T
dis
;D) choose the constitutive relation
~
=
1
2
jT
dis
j
2
;
31
or
~
= 2jDj
2
:
In the ﬁrst case we maximize with respect to T
dis
with no other constrains,in the second
case we maximize with respect to Dwith the incompressibility condition tr D= 0.We obtain
T = mI +2D+2%
@^e
@ (tr B
k
R
)
B
d
k
R
:
If we use neoHookean free energy
% =
2
(tr B
k
R
3);
the it holds
@^e
@ (tr B
k
R
)
=
2
I
and we obtain KelvinVoigt model
T = mI +2D+B
d
k
R
:
If the rate of entropy production is zero,then
T
d
dis
= 0;
p +m+
1
3
tr B
k
R
= 0;
q = 0
and we obtain
T = T
d
+mI = T
d
dis
+B
d
k
R
+mI = pI +B
k
R
:
If we do not require the incompressibility and use the relations as in NavierStokesFourier
T
d
dis
= 2D
d
;
p +m+
1
3
tr B
k
R
=
2 +3
3
div v;
q = kr;
we get
T = pI +2D+(div v)I +B
k
R
;
q = kr:
Viscoelastic model (K.R.Rajagopal,A.R.Srinivasa [28])
The body made from the viscoelastic material is at time t in conﬁguration (t),this con
ﬁguration is called the current conﬁguration.We will compare this conﬁguration with the
conﬁguration (r) which is a conﬁguration of the system in rest at the beginning,so called
the reference conﬁguration.
We deﬁne now the natural conﬁguration (p(t)).It is a conﬁguration of the system cor
responding to the current conﬁguration (t) the body would go to on the sudden removal
of external stimuli.Then the only diﬀerence between the natural and current conﬁgura
tion is purely elastic deformation.By deﬁning the natural conﬁguration we split the whole
deformation into purely elastic part and the dissipative part.
We motivate the notion of the natural conﬁguration using a simple example.Suppose we
have a simple linear viscoelastic material where in every material point is Maxwell element.
The spring represent the elastic part and the dashpot the viscous part In Figure 31a) jthe
material is at the beginning in rest,relaxed.The spring and the dashpot are not stretched.
32
Figure 31:Motivational example for the natural conﬁguration
The system is in the reference conﬁguration (r).Then we stretch it (Figure 31b)) and
the system goes to the current conﬁguration (t).Now we let the system relax,the spring
shrinks to its reference length but the dashpot remains stretched.(Figure 31c)).The system
is in the natural conﬁguration (p(t)).Because during the relaxation the deformation of the
dashpot does not change,we can think that the whole free energy of the system is hidden
only in the deformation of the natural conﬁguration.konﬁgurace.
Figure 32:Přirozená konﬁgurace
We deﬁne some kinematical quantities.First we deﬁne the deformation gradient (in
Figure 32 denoted by F)
F
k
R
= r
that maps inﬁnitesimal element from (r) to (t).We deﬁne standard left and right Cauchy
Green tensors.
B
k
R
= F
k
R
F
T
k
R
;C
k
R
= F
T
k
R
F
k
R
;
Further we deﬁne F
k
p(t)
as a mapping of the inﬁnitesimal element from (p(t)) to (t) (in
Figure.32 denoted by F
p(t)
) and Gthat maps from (r) to (p(t)).These two new variables
are related through the deformation gradient
G= F
1
k
p(t)
F
k
R
:
33
Deﬁne leftand right CauchyGreen tensor corresponding to (p(t))
B
k
p(t)
= F
k
p(t)
F
T
k
p(t)
;C
k
p(t)
= F
T
k
p(t)
F
k
p(t)
and the rate of the deformation of (p(t)) (in the analogy with the usual relation between
the velocity gradient and deformation gradient),resp.its symmetric part
L
k
p(t)
=
_
GG
1
;D
k
p(t)
=
L
k
p(t)
+L
T
k
p(t)
2
:
The material time derivative of B
k
p(t)
is equal to
_
B
k
p(t)
=
_
F
k
p(t)
F
T
k
p(t)
+F
k
p(t)
_
F
T
k
p(t)
=
_
F
k
R
G
1
F
T
k
p(t)
+F
k
R
_
G
1
F
T
k
p(t)
+F
k
p(t)
_
G
T
F
T
k
R
+F
k
p(t)
G
T
_
F
T
k
R
=
LB
k
p(t)
+B
k
p(t)
L
T
2F
k
p(t)
D
k
p(t)
F
T
k
p(t)
;
where we used
_
A
1
= A
1
_
AA
1
:
And we naturally obtained the Oldroyd upper convected time derivative
O
B
k
(p(t))
=
_
B
k
p(t)
LB
k
p(t)
B
k
p(t)
L
T
= 2F
k
p(t)
D
k
p(t)
F
T
k
p(t)
:
Note that
O
I
= 2D.It holds that
tr
_
B
k
p(t)
= 2B
k
p(t)
D2F
k
p(t)
D
k
p(t)
F
k
p(t)
:
Choose the entropy in the form = ~(e;%;B
k
p(t)
),we can invert it and use the internal
energy e = ~e(;%;B
k
p(t)
),concretely choose
e = e
0
(;%) +
2%
tr B
k
p(t)
3
and as in the case for KelvinVoigt we obtain
= (TB
k
p(t)
) D+(p + ~m) div v +
q r
+F
k
p(t)
D
k
p(t)
F
k
p(t)
:
For simplicity q = 0 and div v = 0 and like in [28],we suppose that the elastic response
is also incompressible,i.e.det B
k
p(t)
= 1.From this and incompressibility F
k
R
results that
det G= 1.Diﬀerentiating det G gives tr D
k
p(t)
= 0.
For the rate of entropy production holds (we the zero traces of D and D
k
p(t)
)
= (TB
k
p(t)
)
d
D+C
d
k
p(t)
D
k
p(t)
= T
d
dis
D+C
d
k
p(t)
D
k
p(t)
:
We choose the following constitutive relation for the rate of entropy production
=
~
(D;D
k(p(t)
;C
k
p(t)
) = 2
0
jDj
2
+2
1
D
k
p(t)
C
k
p(t)
D
k
p(t)
:
We observe that D
k
p(t)
B
k
p(t)
D
k
p(t)
= jD
k
p(t)
F
k
p(t)
j
2
0,and the second law of thermody
namics is satisﬁed.We deﬁne the Lagrange function
L(D;D
k
p(t)
) =
~
(:::) +
1
(
~
(:::) T
d
dis
DC
d
k
p(t)
D
k
p(t)
) +
2
tr D+
3
tr D
k
p(t)
:
and maximize.The necessary conditions are
@L
@D
= 0;
@L
@D
k
p(t)
= 0:
34
We substite for L
1 +
1
1
@
~
@D
= T
d
dis
+
2
1
I;(37)
1 +
1
1
@
~
@D
k
p(t)
= C
d
k
p(t)
+
3
1
I:(38)
First we do the diﬀerentiation
@
~
@D
= 4
0
D;
@
~
@D
k
p(t)
= 4
1
C
k
p(t)
D
k
p(t)
:
and eliminate the Lagrange multipliers.To do this we add (37) D+(38) D
k
p(t)
,and use
the constrains of incompressibility
1 +
1
1
=
T
d
dis
D+C
d
k
p(t)
D
k
p(t)
@
~
@D
D+
@
~
@D
k
p(t)
D
k
p(t)
=
~
2
~
=
1
2
:
We take the trace of (37) and obtain
0 = tr T
d
dis
3
2
1
)
2
1
= 0:
If we take the trace of (38),we obtain
3
1
=
2
3
1
D
k
p(t)
C
k
p(t)
:
We get
T
d
dis
= 2
0
D;(39)
2
1
C
k
p(t)
D
k
p(t)
1
3
(D
k
p(t)
C
k
p(t)
)I
= C
d
k
p(t)
:(40)
From the equations (39) it implies that
T = pI +2
0
D+B
d
k
p(t)
:
Before the further calculation we remind that
_
B
k
p(t)
= LB
k
p(t)
+B
k
p(t)
L
T
2F
k
p(t)
D
k
p(t)
F
T
k
p(t)
;
and compute the Oldroyd upper convected time derivative of the deviatoric part of B
k
p(t)
.
O
B
d
k
p(t)
=
O
B
k
p(t)
1
3
(tr B
k
p(t)
)I
=
O
B
k
p(t)
1
3
(tr
_
B
k
p(t)
)I
1
3
(tr B
k
p(t)
)
O
I
=
= 2F
k
p(t)
D
k
p(t)
F
T
k
p(t)
1
3
tr(LB
k
p(t)
+L
T
B
k
p(t)
)I+
+
2
3
tr(F
T
k
p(t)
F
k
p(t)
D
k
p(t)
)I +
2
3
(tr B
k
p(t)
)D=
= 2F
k
p(t)
D
k
p(t)
F
T
k
p(t)
2
3
(D B
k
p(t)
)I +
2
3
(C
k
p(t)
D
k
p(t)
)I +
2
3
(tr B
k
p(t)
)D:
Multiply (40) from the left by F
T
k
p(t)
and from the right by F
T
k
p(t)
,we obtain
2
1
F
k
p(t)
D
k
p(t)
F
T
k
p(t)
1
3
(C
k
p(t)
D
k
p(t)
)I
= C
d
k
p(t)
:
35
From the relations for the Oldroyd derivative of the deviatoric part of B
k
p(t)
we get
2
1
1
2
O
B
d
k
p(t)
1
3
(D B
d
k
p(t)
)I +
1
3
(tr B
k
p(t)
)D
!
= B
d
k
p(t)
:
Note that the trace both sides of the equation is zero,in the three dimensional space it is ﬁve
scalar equations for six unknowns.We complete the system of equations with the equation
for incompressibility of the elastic response
det
B
k
p(t)
+
1
3
(tr B
k
p(t)
)I
= 1:
We obtain the set of equations for nonlinear viscoelastic ratetype ﬂuid model
div v = 0;
% _v = div T;
T = pI +2
0
D+B
d
k
p(t)
;
2
1
1
2
O
B
d
k
p(t)
1
3
(D B
d
k
p(t)
)I +
1
3
(tr B
k
p(t)
)D
!
= B
d
k
p(t)
;(41)
det
B
k
p(t)
+
1
3
(tr B
k
p(t)
)I
= 1:(42)
Linearization of the viscoelastic model
By linearization of the elastic response we show that this model reduces to OldroydB model.
Suppose that the left CauchyGreen tensor is in the form B
k
p(t)
= I + A,where kAk =
";0 <" 1.By linearization we understand that after substitution into the nonlinear
viscoelastic model we will neglect all term of the order O("
2
) and higher.We use the last
equation for incompressibility (42) and do the Taylor expansion of determinant.We obtain
1 = det
B
k
p(t)
= det(I +A) = 1 +tr(A) +
1
2
(tr A)
2
tr
A
2
+O("
3
) )
tr(A) =
1
2
tr
A
2
(tr A)
2
+O("
3
)

{z
}
O("
2
)
:
From this it implies that
tr B
k
p(t)
= 3 +
1
2
tr
A
2
(tr A)
2
+O("
3
);
we obtain the relation for B
d
k
p(t)
,
B
d
k
p(t)
= B
k
p(t)
1
3
(tr B
k
p(t)
)I = A
1
6
tr
A
2
(tr A)
2
I +O("
3
)I:
Further we need to compute
O
B
d
k
p(t)
.
O
B
d
k
p(t)
=
_
B
d
k
p(t)
LB
d
k
p(t)
B
d
k
p(t)
L
T
=
O
A
1
3
_
A A(tr A)

{z
}
O("
2
)
(tr
_
A)
I +O("
2
):
Inserting all into (41) we get
A+
1
O
A
= 2
1
D
1
3
1
h
2D
_
A
A
i
I +O("
2
):(43)
36
This model is ver similar to OldroydB model,the diﬀerence is in the last term.We show
that this term is of the order O("
2
),and so we can neglect it.Take the scala product of (43)
with tensor A,we get
jAj
2

{z
}
O("
2
)
+
1
_
A A
1
(LA A+AL
T
A)

{z
}
O("
2
)
= 2
1
(D A)
1
3
1
h
2D
_
A
A
i
tr A

{z
}
O("
2
)
:
From this it implies that
9
2D
_
A
A= O("
2
):
We showed that the nonlinear viscoelatic model linearizes to OldroydB model.
div v = 0;
% _v = div T;
T = pI +2
0
D+A;
A+
1
O
A
= 2
1
D;
and so derived thermodynamically compatible nonlinear viscoelastic model is in some sense
consistent with the wellknown viscoelastic model.
Nonlinear twocomponent viscoelastic model
Let consider the material consisting of two components,to each component corresponds
one natural conﬁguration k
p
1
(t)
a k
p
2
(t)
(see Figure (33)).We choose the entropy =
Figure 33:Two natural conﬁgurations
~(e;%;B
k
p
1
(t)
;B
k
p
2
(t)
),we can invert it and use the internal energy e = ~e(;%;B
k
p
1
(t)
;B
k
p
2
(t)
)
and concretely we choose
e = e
0
(;%) +
1
2%
tr B
k
p
1
(t)
3
+
2
2%
tr B
k
p
2
(t)
3
;
9
Note that from this we can read
_
A= 2D+O(");
which is expected.
37
again we consider only isothermal processes and the incompressible material tr D= tr D
k
p
1
(t)
=
tr D
k
p
2
(t)
= 0.After some calculation we obtain a relation for the rate of entropy production
= (T
1
B
k
p
1
(t)
2
B
k
p
2
(t)
) D+
1
C
d
k
p
1
(t)
D
k
p
1
(t)
+
2
C
d
k
p
2
(t)
D
k
p
2
(t)
= T
d
dis
D+
1
C
d
k
p
1
(t)
D
k
p
1
(t)
+
2
C
d
k
p
2
(t)
D
k
p
2
(t)
:
We choose the following constitutive relation for the rate of entropy production
=
~
(D;D
k
p
1
(t)
;D
k
p
2
(t)
;B
k
p
1
(t)
;B
k
p
2
(t)
) = 2
0
jDj
2
+2
1
D
k
p
1
(t)
B
k
p
1
(t)
D
k
p
1
(t)
+
2
2
D
k
p
2
(t)
B
k
p
2
(t)
D
k
p
2
(t)
:
Again we observe that the rate of entropy production is nonnegative and so the second law
of thermodynamics is satiﬁed.We deﬁne the Lagrange function
L(D;D
k
p
1
(t)
;D
k
p
2
(t)
) =
~
(:::)+
1
(
~
(:::)T
d
dis
D
1
C
d
k
p
1
(t)
D
k
p
1
(t)
2
C
d
k
p
2
(t)
D
k
p
2
(t)
)+
2
tr D+
3
tr D
k
p
1
(t)
+
4
tr D
k
p
2
(t)
:
Necessary conditions for the maximum are
@L
@D
= 0;
@L
@D
k
p
1
(t)
= 0;
@L
@D
k
p
2
(t)
= 0:
We substitute for L
1 +
1
1
@
~
@D
= T
d
dis
+
2
1
I;(44)
1 +
1
1
@
~
@D
k
p
1
(t)
=
1
C
d
k
p
1
(t)
+
3
1
I;(45)
1 +
1
1
@
~
@D
k
p
2
(t)
=
1
C
d
k
p
2
(t)
+
4
1
I:(46)
Do the diﬀerentiation
@
~
@D
= 4
0
D;
@
~
@D
k
p
1
(t)
= 4
1
C
k
p
1
(t)
D
k
p
1
(t)
;
@
~
@D
k
p
2
(t)
= 4
2
C
k
p
2
(t)
D
k
p
2
(t)
:
and eliminate the Lagrange multipliers.To do this we add (44)D+(45)D
k
p
1
(t)
+(46)D
k
p
2
(t)
,
use the constrains of incompressibility and get
1 +
1
1
=
T
d
dis
D+
1
C
d
k
p
1
(t)
D
k
p
1
(t)
+
2
C
d
k
p
2
(t)
D
k
p
2
(t)
@
~
@D
D+
@
~
@D
k
p
1
(t)
D
k
p
1
(t)
+
@
~
@D
k
p
2
(t)
D
k
p
2
(t)
=
~
2
~
=
1
2
:
Now take the trace of (44) and obtain
0 = tr T
d
dis
3
2
1
)
2
1
= 0:
By taking the traces of (45) and (46),we obtain
3
1
=
2
3
1
D
k
p
1
(t)
C
k
p
1
(t)
;
4
1
=
2
3
2
D
k
p
2
(t)
C
k
p
2
(t)
:
38
We get
T
d
dis
= 2
0
D;(47)
2
1
C
k
p
1
(t)
D
k
p
1
(t)
1
3
(D
k
p
1
(t)
C
k
p
1
(t)
)I
=
1
C
d
k
p
1
(t)
;(48)
2
2
C
k
p
2
(t)
D
k
p
2
(t)
1
3
(D
k
p
2
(t)
C
k
p
2
(t)
)I
=
2
C
d
k
p
2
(t)
:(49)
From (47) it implies that
T = pI +2
0
D+
1
B
d
k
p
1
(t)
+
2
B
d
k
p
2
(t)
:
As in the case of model with one natural conﬁguration it holds
O
B
d
k
p
1
(t)
= 2F
k
p
1
(t)
D
k
p
1
(t)
F
T
k
p
1
(t)
2
3
(D B
k
p
1
(t)
)I +
2
3
(C
k
p
1
(t)
D
k
p
1
(t)
)I +
2
3
(tr B
k
p
1
(t)
)D;
O
B
d
k
p
2
(t)
= 2F
k
p
2
(t)
D
k
p
2
(t)
F
T
k
p
2
(t)
2
3
(D B
k
p
2
(t)
)I +
2
3
(C
k
p
2
(t)
D
k
p
2
(t)
)I +
2
3
(tr B
k
p
2
(t)
)D
Multiplying (48) from the left by F
T
k
p
1
(t)
and from the right by F
T
k
p
1
(t)
,and (49) from the left
by F
T
k
p
2
(t)
and from the right by F
T
k
p
2
(t)
,we get
2
1
F
k
p
1
(t)
D
k
p
1
(t)
F
T
k
p
1
(t)
1
3
(B
k
p
1
(t)
D
k
p
1
(t)
)I
=
1
B
d
k
p
1
(t)
;
2
2
F
k
p
2
(t)
D
k
p
2
(t)
F
T
k
p
2
(t)
1
3
(B
k
p
2
(t)
D
k
p
2
(t)
)I
=
2
B
d
k
p
2
(t)
:
From the relations for Oldroyd derivatives of B
d
k
p
1
(t)
a B
d
k
p
2
(t)
we obtain
2
1
1
2
O
B
d
k
p
1
(t)
1
3
(D B
d
k
p
1
(t)
)I +
1
3
(tr B
k
p
1
(t)
)D
!
=
1
B
d
k
p
1
(t)
;
2
2
1
2
O
B
d
k
p
2
(t)
1
3
(D B
d
k
p
2
(t)
)I +
1
3
(tr B
k
p
2
(t)
)D
!
=
2
B
d
k
p
2
(t)
:
The trace of both equations is zero,we have 2 5 scalar equations for 2 6 unknowns.
We complete the system of equations with two equations for incompressibility of two elastic
responses
det
B
k
p
1
(t)
+
1
3
(tr B
k
p
1
(t)
)I
= 1;
det
B
k
p
2
(t)
+
1
3
(tr B
k
p
2
(t)
)I
= 1:
39
We obtain the governing equations for the nonlinear twocomponent viscoelastic ratetype
ﬂuid model
div v = 0;
% _v = div T;
T = pI +2
0
D+
1
B
d
k
p
1
(t)
+
2
B
d
k
p
2
(t)
;
2
1
1
2
O
B
d
k
p
1
(t)
1
3
(D B
d
k
p
1
(t)
)I +
1
3
(tr B
k
p
1
(t)
)D
!
=
1
B
d
k
p
1
(t)
;
2
2
1
2
O
B
d
k
p
2
(t)
1
3
(D B
d
k
p
2
(t)
)I +
1
3
(tr B
k
p
2
(t)
)D
!
=
2
B
d
k
p
2
(t)
;
det
B
k
p
1
(t)
+
1
3
(tr B
k
p
1
(t)
)I
= 1;
det
B
k
p
2
(t)
+
1
3
(tr B
k
p
2
(t)
)I
= 1:
40
Homework
This is not full list of homework.It is just to complete the lecture notes.
Homework 1.Under the assumption
= T
d
D
d
+(m+p) div v +
q r
;(50)
derive NavierStokes equations by principle of maximization of entropy production
max
(50)+T
d
;m;q
~
(T
d
;m;q);
where
~
is given by
~
=
1
2
jT
d
j
2
+
3
2 +3
(m+p)
2
+
1
K
jqj
2
:
Homework 2.Consider the following model
T = pI +
1
+
0
1
(1 +
2
jDj
2
)
1n
2
D;
where
0
;
1
; > 0.Explain which nonNewtonian behaviour exhibits ﬂuid described by that
model for n 2 R.
Homework 3.Show that constitutive relation for Bingham ﬂuid
jSj
,D= 0;
jSj >
,S =
D
jDj
+2(jDj
2
)D
is equivalent to
2(jDj
2
)
+(jSj
)
+
D= (jSj
)
+
S;
where
> 0,():R
+
0
!R
+
0
and x
+
= maxf0;xg.
Homework 4.Use the principle of maximization of entropy production =
~
(T
d
;D
d
) for
T
d
2 A,where
A:=
T
d
2 R
33
sym
; = T
d
:D
d
:
In two cases for
~
given by:
1.
~
= jT
d
j
r=(r1)
;
2.
~
= jD
d
j
2r
jT
d
j
2
:
Homework 5.Stokes ﬂuid is given by constitutive relation
T = pI +2D+D
2
;
where and are constant.Explain which nonNewtonian behaviour exhibits Stokes ﬂuid.
Homework 6.Consider Korteweg’s model in unidirectional velocity ﬁeld v = (u(y);0;0).
Show that model exhibits normal stress diﬀerences
T
11
T
22
6= 0;
T
22
T
33
6= 0:
41
Bibliography
[1] C.Barus.Isotherms,isopiestics and isometrics relative to viscosity.Amer.J.Sci.,
45:87–96,1893.
[2] Santanu Basu and U.S.Shivhare.Rheological,textural,microstructural and sensory
properties of mango jam.J.Food Eng.,100(2):357–365,2010.
[3] B.Bernstein,E.A.Kearsley,and L.J.Zapas.A study of stress relaxation with ﬁnite
strain.Trans.Soc.Rheol.,7(1):391–410,1963.
[4] C.E.Bingham.Fluidity and plasticity.McGraw–Hill,New York,1922.
[5] H Blatter.Velocity and stressﬁelds in grounded glaciers – a simple algorithm for in
cluding deviatoric stress gradients.J.Glaciol.,41(138):333–344,1995.
[6] P.W.Bridgman.The eﬀect of pressure on the viscosity of fortyfour pure liquids.Proc.
Am.Acad.Art.Sci.,61(3/12):57–99,FEBNOV 1926.
[7] J.M.Burgers.Mechanical considerations – model systems – phenomenological theories
of relaxation and viscosity.In First report on viscosity and plasticity,chapter 1,pages
5–67.Nordemann Publishing,New York,1939.
[8] Pierre J.Carreau.Rheological equations from molecular network theories.J.Rheol.,
16(1):99–127,1972.
[9] IJen Chen and D.C.Bogue.Timedependent stress in polymer melts and review of
viscoelastic theory.Trans.Soc.Rheol.,16(1):59–78,1972.
[10] William O.Criminale,J.L.Ericksen,and G.L.Filbey.Steady shear ﬂow of non
Newtonian ﬂuids.Arch.Rat.Mech.Anal.,1:410–417,1957.
[11] Malcolm M.Cross.Rheology of nonnewtonian ﬂuids:A new ﬂow equation for pseudo
plastic systems.J.Colloid Sci.,20(5):417–437,1965.
[12] A.de Waele.Viscometry and plastometry.J.Oil Colour Chem.Assoc.,6:33–69,1923.
[13] J.E.Dunn and J.Serrin.On the thermomechanics of interstitial working.Archive for
Rational Mechanics and Analysis,88:95–133,1985.
[14] Henry Eyring.Viscosity,plasticity,and diﬀusion as examples of absolute reaction rates.
J.Chem.Phys.,4(4):283–291,1936.
[15] H.Giesekus.A simple constitutive equation for polymer ﬂuids based on the concept of
deformationdependent tensorial mobility.J.NonNewton.Fluid Mech.,11(12):69–109,
1982.
[16] J.W.Glen.The creep of polycrystalline ice.Proc.R.Soc.AMath.Phys.Eng.Sci.,
228(1175):519–538,1955.
42
[17] Winslow H.Herschel and Ronald Bulkley.Konsistenzmessungen von Gummi
Benzollösungen.Colloid Polym.Sci.,39(4):291–300,August 1926.
[18] K.L.Johnson and J.L.Tevaarwerk.Shear behaviour of elastohydrodynamic oil ﬁlms.
Proc.R.Soc.AMath.Phys.Eng.Sci.,356(1685):215–236,1977.
[19] M.W.Johnson and D.Segalman.A model for viscoelastic ﬂuid behavior which allows
nonaﬃne deformation.J.NonNewton.Fluid Mech.,2(3):255–270,1977.
[20] Michael M.Khonsari and E.Richard Booser.Applied Tribology:Bearing Design and
Lubrication.John Wiley & Sons Ltd,Chichester,second edition,2008.
[21] Seikichi Matsuhisa and R.Byron Bird.Analytical and numerical solutions for laminar
ﬂow of the nonNewtonian Ellis ﬂuid.AIChE J.,11(4):588–595,1965.
[22] J.Clerk Maxwell.On the dynamical theory of gases.Philos.Trans.R.Soc.,157:49–88,
1867.
[23] J.G.Oldroyd.On the formulation of rheological equations of state.Proc.R.Soc.
AMath.Phys.Eng.Sci.,200(1063):523–541,1950.
[24] Wolfgang Ostwald.Über die Geschwindigkeitsfunktion der Viskosität disperser Systeme.
I.Colloid Polym.Sci.,36:99–117,1925.
[25] Erin C.Pettit and Edwin D.Waddington.Ice ﬂow at low deviatoric stress.J.Glaciol.,
49(166):359–369,2003.
[26] N.Phan Thien.Nonlinear network viscoelastic model.J.Rheol.,22(3):259–283,1978.
[27] N.Phan Thien and Roger I.Tanner.A new constitutive equation derived from network
theory.J.NonNewton.Fluid Mech.,2(4):353–365,1977.
[28] K.R.Rajagopal and A.R.Srinivasa.A thermodynamic frame work for rate type ﬂuid
models.J.NonNewton.Fluid Mech.,88(3):207–227,2000.
[29] Francis Ree,Taikyue Ree,and Henry Eyring.Relaxation theory of transport problems
in condensed systems.Ind.Eng.Chem.,50(7):1036–1040,1958.
[30] M.Reiner.A mathematical theory of dilatancy.Am.J.Math.,67(3):350–362,1945.
[31] R.S.Rivlin and J.L.Ericksen.Stressdeformation relations for isotropic materials.J.
Ration.Mech.Anal.,4:323–425,1955.
[32] R.S.Rivlin and K.N.Sawyers.Nonlinear continuum mechanics of viscoelastic ﬂuids.
Annu.Rev.Fluid Mech.,3:117–146,1971.
[33] T.A.Savvas,N.C.Markatos,and C.D.Papaspyrides.On the ﬂow of nonnewtonian
polymer solutions.Appl.Math.Modelling,18(1):14–22,1994.
[34] Gilbert R.Seely.Nonnewtonian viscosity of polybutadiene solutions.AIChE J.,
10(1):56–60,1964.
[35] A.W.Sisko.The ﬂow of lubricating greases.Ind.Eng.Chem.,50(12):1789–1792,1958.
[36] Kenji Yasuda.Investigation of the analogies between viscometric and linear viscoelastic
properties of polystyrene ﬂuids.PhDthesis,Massachusetts Institute of Technology.Dept.
of Chemical Engineering.,1979.
43
Σχόλια 0
Συνδεθείτε για να κοινοποιήσετε σχόλιο