ORIGINAL ARTICLE: H. Binous, Introducing non-Newtonian fluid mechanics computations with Mathematica® in the undergraduate curriculum, Chemical Engineering Education, 41 (1), pp . 59 - 64, 2007. https://journals.flvc.org/cee/article/view/122478
Article Abstract: We study four non-Newtonian fluid mechanics problems using Mathematica®. Constitutive equations describing the behavior of power-law, Bingham and Carreau models are recalled. The velocity profile is obtained for the horizontal flow of power-law fluids in pipes and annuli. For the vertical laminar film flow of a Bingham fluid we determine the velocity profile. Both problems involve the use of the shooting techniques because they have split boundary conditions. Since Mathematica permits symbolic computations, we determine analytical expressions of volumetric flow rates for pipe flow of the Bingham and power-law fluids. We use the built-in optimization command of Mathematica, FindMinimum, in order to find the non-Newtonian fluid model from representative data of flow rates measured under different applied pressure gradients in a horizontal pipe. These pedagogic problems are used to introduce the field of non-Newtonian fluid mechanics to students at the National Institute of Applied Sciences in Tunis. The Mathematica® notebooks are available from the corresponding author upon request or at Wolfram Research.
We propose to find the velocity profiles for the laminar flow of Newtonian, power-law and Carreau fluids in an annulus, shown in below. We use the following values for the pressure difference \[Laplacian]P, the exponent n, the Newtonian fluid viscosity η, the consistency index k,
the infinite-shear viscosity
η
∞
, the zero-shear viscosity
η
0
, the relaxation parameter λ, the pipe length L, and radius R:
​
\[Laplacian]P =100; L = 50; R1=0.02 and R2 = 0.05
​
Newtonian fluid:
η=8.9×
-4
10
​
​
Dilatant fluid:
n=1.2
and
κ=4.7×
-4
10
​
​
Pseudo-plastic fluid:
n=0.5
and
κ=4.5×
-3
10
​
​
Carreau fluid:
n=0.5
,
λ=0.2
,
η
0
=2.04×
-3
10
, and
η
∞
=0
Flow of dilatant, pseudoplastic, Carreau and Newtonian fluids in an annulus.
Pseudoplastic Fluid
In[]:=
K=4.510^-3;n=0.5;ΔP=100;L=50;R1=0.02;R2=0.05;
In[]:=
system[Ω_]={D[r
τ
rx
[r],{r,1}]ΔP/Lr,D[
v
x
[r],{r,1}]If[
τ
rx
[r]≥0,-(
τ
rx
[r]/K)^(1/n),(-
τ
rx
[r]/K)^(1/n)],
τ
rx
[R1]Ω,
v
x
[R1]0};​​myODEsoln[Ω_?NumericQ]:=NDSolve[system[Ω],{
v
x
,
τ
rx
},{r,R1,R2}]​​yend[Ω_?NumericQ]:=Flatten[(
v
x
[r]/.myODEsoln[Ω])/.rR2]​​bc=FindRoot[yend[Ω]0,{Ω,-2,2}][[1,2]]
Out[]=
-0.035363
In[]:=
sol1=myODEsoln[bc]
Out[]=

v
x
InterpolatingFunction
Domain: {{0.02,0.05}}
Output: scalar
,
τ
rx
InterpolatingFunction
Domain: {{0.02,0.05}}
Output: scalar

In[]:=
plt1=Plot[
v
x
[r]/.sol1,{r,R1,R2},PlotStyle{Thick,Darker@Green},Frame->True,AspectRatio->1,FrameLabel->{Style["position",14],Style["velocity",14]},GridLines->Automatic,ImageSize->350{1,1},PlotRange->{{R1,R2},All}]
Out[]=
Dilatant Fluid
In[]:=
K=4.710^-4;n=1.2;ΔP=100;L=50;R1=0.02;R2=0.05;
In[]:=
system[Ω_]={D[r
τ
rx
[r],{r,1}]ΔP/Lr,D[
v
x
[r],{r,1}]If[
τ
rx
[r]>0,-(
τ
rx
[r]/K)^(1/n),(-
τ
rx
[r]/K)^(1/n)],
τ
rx
[R1]Ω,
v
x
[R1]0};​​myODEsoln[Ω_?NumericQ]:=NDSolve[system[Ω],{
v
x
,
τ
rx
},{r,R1,R2},MaxSteps->Infinity]​​yend[Ω_?NumericQ]:=Re@Flatten[(
v
x
[r]/.myODEsoln[Ω])/.rR2]​​bc=FindRoot[yend[Ω]0,{Ω,-2,2}][[1,2]]
Out[]=
-0.0377572
In[]:=
sol2=myODEsoln[bc]
Out[]=

v
x
InterpolatingFunction
Domain: {{0.02,0.05}}
Output: scalar
,
τ
rx
InterpolatingFunction
Domain: {{0.02,0.05}}
Output: scalar

In[]:=
plt1=Plot[Re[
v
x
[r]/.sol2],{r,R1,R2},PlotStyle{Thick,Darker@Magenta},Frame->True,AspectRatio->1,FrameLabel->{Style["position",14],Style["velocity",14]},GridLines->Automatic,ImageSize->350{1,1},PlotRange->{{R1,R2},All}]
Out[]=
Newtonian Fluid
In[]:=
μ=8.910^-4;n=1;ΔP=100;L=50;R1=0.02;R2=0.05;
In[]:=
system[Ω_]={D[r
τ
rx
[r],{r,1}]ΔP/Lr,D[
v
x
[r],{r,1}]-(
τ
rx
[r]/μ),
τ
rx
[R1]Ω,
v
x
[R1]0};​​myODEsoln[Ω_]:=NDSolve[system[Ω],{
v
x
,
τ
rx
},{r,R1,R2}]​​yend[Ω_?NumericQ]:=Flatten[(
v
x
[r]/.myODEsoln[Ω])/.rR2]​​bc=FindRoot[yend[Ω]0,{Ω,0,0.5}][[1,2]]
Out[]=
-0.0372962
In[]:=
sol3=myODEsoln[bc]
Out[]=

v
x
InterpolatingFunction
Domain: {{0.02,0.05}}
Output: scalar
,
τ
rx
InterpolatingFunction
Domain: {{0.02,0.05}}
Output: scalar

In[]:=
plt4=Plot[
v
x
[r]/.sol3,{r,R1,R2},PlotStyle{Thick,Darker@Magenta},Frame->True,AspectRatio->1,FrameLabel->{Style["position",14],Style["velocity",14]},GridLines->Automatic,ImageSize->350{1,1}]
Out[]=
Analytical Solution is available only for Newtonian Fluids.​Comparing numerical result (solid line) with analytical solution (dashed curve) for Newtonian fluid​​
v
x
(r)=
ΔP
4Lμ
2
R2
-
2
r
+
(
2
R2
-
2
R1
)Log
r
R2

Log
R2
R1

​​
In[]:=
v
analytical
[r_]:=ΔP/(4μL)(R2^2-r^2+(R2^2-R1^2)/Log[R2/R1]Log[r/R2])
In[]:=
plt5=Plot[
v
analytical
[r],{r,R1,R2},PlotStyle{Thick,Darker@Cyan,Dashed},Frame->True,AspectRatio->1,FrameLabel->{Style["position",14],Style["velocity",14]},GridLines->Automatic,ImageSize->350{1,1}]
Out[]=
Comparing velocity profiles obtained for the four fluids: dilatant, pseudoplastic, Carreau and Newtonian. These profiles are obtained under equal volumetric flow conditions.
The velocity near the wall is higher for Carreau and pseudo-plastic fluids than for Newtonian and dilatant fluids. This results in higher heat transfer rates due to a higher convection.
​
Final calculation checks : we find that all four volumetric flow are equal !

Acknowledgement

Since in this laminar flow problem we have split boundary conditions. Problem’s resolution is based on the Shooting Method developed in Mathematica by Professor Brian Higgins (University of California, Davis).

References

[1] H . Binous, Introducing non-Newtonian fluid mechanics computations with Mathematica® in the undergraduate curriculum, Chemical Engineering Education, 41 (1), pp . 59 - 64, 2007.
[2] J. O. Wilkes, Fluid Mechanics for Chemical Engineers, Upper Saddle River, NJ: Prentice Hall, 1999.

CITE THIS NOTEBOOK

Introducing non-Newtonian fluid mechanics computations in the undergraduate curriculum​
by Housam Binous and Ahmed Bellagi ​
Wolfram Community, STAFF PICKS, June 14, 2023
​https://community.wolfram.com/groups/-/m/t/2937252