Publicly Available Published by De Gruyter August 7, 2015

# A Note on Kirchhoff’s Junction Rule for Power-Law Fluids

• Eduard Marušić-Paloka and Igor Pažanin

## Abstract

We study a multiple pipe system filled with non-Newtonian fluid being in stationary regime and obeying the power law. Introducing the small parameter ε, being the ratio between the cross-section area and the length of the pipes, we find the effective flow via singular perturbation as ε tends to zero. The flow in each pipe of the system remains in the power-law Poiseuille form governed by the pressure drop between pipe’s end and the junction. Kirchhoff’s law providing the value of pressure at the junction point is derived, and its solvability is discussed. In the vicinity of junction, an interior layer appears, so we correct the asymptotic solution by introducing a non-Newtonian Leray’s problem. The error analysis is also provided.

## 1 Introduction

Non-Newtonian fluid flows through tubes that are either thin or long have many industrial and biological applications. It is mainly due to the fact that numerous real fluids (such as polymer melts and solutions, blood, pulps, cement, and mud) cannot be successfully described by the linear Newton’s law because its viscosity significantly changes with the shear rate (see, e.g. [1]). The simplest way to capture those phenomena is to introduce the shear-dependent viscosity. In this article, we study the case where the viscosity is given by the power law (or Ostwald–de Waele law):

(1)η(u)=μ|D(u)|r2. (1)

Here, u denotes the velocity, while D(u)=1/2[∇u + (∇u)T] is the strain rate tensor, i.e. the symmetric component of the velocity gradient. Two positive constants μ and r in (1) are called the flow consistency and the flow behaviour index, respectively. Depending on the value of the index r, we have shear-thinning (or pseudo-plastic) behaviour for 1<r<2 or shear-thickening (or dilatant) behaviour for r>2. Note that the case r=2 corresponds to the standard Newtonian behaviour. The power-law model has been extensively used in the engineering community because of its simplicity and ability to fit the experimental viscosity curves for many real fluids.

Our goal is to study a multiple pipe system filled with power-law fluid obeying (1). Inspired by numerous applications where several pipes are interconnected forming a branching structure, here we focus at the individual section of the network containing one junction point and m straight pipes with (possibly) different profiles (see Fig. 1). The same arguments of the presented analysis can be employed in the case of a general network with multiple junctions as well. Frequently, different scales naturally appear in these systems, e.g. if the extension in one space direction is much larger than in the other two directions. Thus, we introduce a small parameter ε representing the ratio between the cross-section area and the pipes’ length and find the effective behaviour of the flow via asymptotic analysis with respect to ε.

Figure 1:

Junction of thin pipes (for m=3).

The flow of classical Newtonian fluid through a multiple pipe system has so far been addressed in several papers; we refer the reader to [2–5]. The case of power-law flow through a single pipe has been extensively studied for various regimes of the flow, see, e.g. [6–9]. In particular, starting from the Stokes system with a shear-dependent viscosity, the power-law Poiseuille flow (PLPF) has been rigorously derived in [10] via two-scale convergence. Here, we formally confirm that in each pipe of our multiple pipe system, the flow remains in the form of PLPF, at least far from the junction (see Section 3). However, to determine the complete form of PLPF, one needs to deduce the pressure drop along the pipe, i.e. the value of the pressure at the junction point.1 Using the incompressibility of the flow, we obtain the equation satisfied by the junction pressure p0, see Section 4. In analogy with electrical circuits, we call it Kirchhoff’s junction rule or simply Kirchhoff’s law. In a case of a simple system consisting of two pipe connected at a junction, it turns out that Kirchhoff’s law can be explicitly solved yielding the explicit formula for p0. In general situation (the system of m pipes, m>2), we cannot derive the explicit solution of the obtained equation due to its nonlinearity. Nevertheless, we can prove that it generates the unique solution p0 and that represents our main contribution. Moreover, in the vicinity of junction, the interior layer has been observed in which the fluid exhibits behaviour different than far from the junction. For that reason, in Section 5, we correct the asymptotic approximation for dilatant flows by introducing non-Newtonian Leray’s problem and performing the appropriate matching procedure. Last but not least, the error analysis is provided in Section 6, indicating the order of accuracy of our method. To conclude, we believe that the presented result could be instrumental for creating more efficient numerical algorithms acknowledging the non-Newtonian effects on the pipe networks flow.

## 2 Description of the Problem

Let SεiR2 be a smooth-bounded domain denoting the cross section of the ith pipe, i=1, …, m, mN. The ith pipe of the system with length ℓi in (possibly) different coordinate systems (eji)j=1,2,3 is then defined as

Ωεi={x=(x1i,xi)R3:0<x1i<i,xi=(x2i,x3i)Sεi},i=1,m.

Motivated by the applications, we want to study the system in which the pipes are either very thin or very long. The usual way to describe that is to introduce the small (dimensionless) parameter ε as the ratio between pipe’s diameter and its length and put Sεi=εSi. Our thin (or long) pipes are connected at junction point O where we choose to put the origin. In view of that, we introduce Ωε0=εΩ0 as a small reservoir where the pipes meet, with Ω0 being a smooth-bounded domain containing O. The domain under consideration is then defined as (see Fig. 1)

(2)Ωε=i=0mΩεi. (2)

The end of the ith pipe is denoted by

Σεi={(i,xi)R3:xi=(x2i,x3i)εSi},   i=1,,m.

Finally, the lateral boundary of Ωε is denoted by Γε=Ωε\(i=1mΣεi).

We suppose that the pipe system Ωε is filled with non-Newtonian fluid with viscosity η obeying the power law:

η(u)=μ|D(u)|r2,  1<r<+,  μ=const.>0.

Here, D(u)=1/2[∇u + (∇u)T] is the symmetrised gradient of the velocity, namely

D(u)jk=12(ujxk+ukxj),  |D(u)|=(j,k=13|D(u)jk|2)1/2.

The equations to be considered have the following form:

(3)μdiv{|D(uε)|r2D(uε)}+pε=0 in  Ωε, (3)
(4)div  uε=0  in  Ωε. (4)

Remark 1As you can see, we assume that the flow is sufficiently slow so that the inertial term can be neglected. Such assumption is needed to have the well-posed problem2and not for the sake of our analysis. Indeed, if we write (3) in non-dimensional form and employ the generalised Reynolds number Regderived for cylindrical ducts of arbitrary cross section (see [11]), we get Reg= O (εγ) with γ being such that the inertial term does not contribute to the macroscopic model.

To complete the problem, suitable boundary conditions have to be added. We impose a standard no-slip condition for the fluid velocity on the lateral boundary and prescribe the values of pressures pi at the end of each pipe. In view of that, the boundary conditions read as follows:

(5)uε=0  on  Γε, (5)
(6)uε×e1i=0,   pε=pi  on  Σεi,  i=1,,m. (6)

Our goal is to investigate the effective behaviour of the flow described by (3–6), as ε → 0.

## 3 Power-Law Poiseuille Flow

Let r′ be the conjugate exponent of r, i.e. r′=r/(r – 1). It is reasonable to expect that in each pipe Ωεi of our system, the flow remains in the so-called PLPF, at least far from the junction. Indeed, by substituting

(7)uε=εrUi(x1i,xiε),  pε=Pi(x1i,xiε)+εQi(x1i,xiε) (7)

in the governing equations, we arrive at

(8)ε1:    yiPi=0      Pi=Pi(x1i), (8)
(9)1:    μdivyi{|Dyi(Ui)|r2Dyi(Ui)}+yiQi=Pix1ie1i, (9)
(10)εr1:  divyiUi=0. (10)

The scaling in (7) is suggested by the result from [10] and will be rigorously justified in Section 6 by proving the error estimates. Here, we introduce the fast variable yi=xi/ε, so the above equations are posed in ε-independent domains Ωi={(x1i,yi)R3:0<x1i<i,yi=(y2i,y3i)Si}. We also employ the following notations for the formal partial differential operators:

yiΦi=Φiy2ie2i+Φiy3ie3i,   divyiVi=V2iy2i+V3iy3i,  Vji=Vieji,j=1,2,3,Dyi(Vi)jk=12(Vjiyki+Vkiyji),j,k=2,3,Dyi(Vi)1j=Dy*i(Vi)j1=0,j=1,2,3.

Taking into account the no-slip boundary condition for the velocity and the values of the prescribed pressures at pipes ends, the system (8–10) can be solved by taking

Qi=0 and

(11)Ui(yi)=wi(yi)|pip0i|r2(p0pii)e1i,yiSi, (11)
(12)Pi(x1i)=p0+(pip0i)x1i,  x1i(0,i). (12)

Here, p0 denotes the (constant) value of the pressure at the junction point O, while wi is the solution of the auxiliary problem posed on the cross section of the ith pipe, namely:

(13)μ div(|ωi|r2ωi)=2r/2 in Si,ωi=0 on  Si. (13)

By a simple integration by parts, we deduce

(14)ωi:=Siωidyi=μ2r/2Si|ωi|rdyi>0. (14)

Remark 2If the pipes are assumed to have circular cross sections, namely Si= {y ∈ R2: |y|<Ri} (Ri=const.>0), then wican be explicitly computed as

(15)ωi(yi)=12r/21μr1r(Rir|yi|r). (15)

Then, in view of (11), the non-zero component of the velocity in the ithpipe has the form

Uxi(yi)=12r/21μr1r(Rir|yi|r)|pip0i|r2(p0pii),yiSi.

The obtained result is most conveniently illustrated using the dimensionless quantities defined by

U^xi=UxiVmaxi,   y^i=yiRi,

whereVmaxiis the maximum velocity (obtained fory2i=y3i=0). The dimensionless velocity fieldU^xi=1|y^i|r/(r1)is plotted in Figure 2 for different values of the flow behaviour index r. For r=2 (red line) we observe the standard parabolic profile of Newtonian fluid. For r<2 (i.e. for pseudo-plastic fluids), we see that the velocity profile becomes more flattened in the centre as r decreases. The observed phenomena, namely the flattening in the vicinity of the pipe’s centre line and relegation of the majority of shear toward the pipe’s walls could be indicative for practical implications, for instance, in polymer processing flows (see, e.g. [12]).

Figure 2:

Dimensionless velocity distribution U^xi as a function of |y^i| for different values of r.

## 4 Kirchhoff’s Junction Rule

In order to know the complete form of the PLPF derived in the previous section, note that we still need to determine the value of the junction pressure p0. Denoting Uεi=εrUi(xi/ε), we can use the incompressibility of our flow and (11) to obtain

(16)0=i=1mΣεiUεie1idxi=εr+2i=1mωi|pip0i|r2(pip0i) (16)

Let us introduce

(17)Hi:=|ωi|r2ωi>0 (17)

It can be easily verified that

|Hi|r2Hi=ωi

so that from (16) it follows:

(18)i=1m|Hii(pip0)|r2Hii(pip0)=0. (18)

If we consider a pipe system with only two pipes connected at the junction (i.e. m=2), then we can explicitly compute the junction pressure p0. Indeed, from (18), we deduce

|H11(p1p0)|r2H11(p1p0)=|H22(p2p0)|r2H22(p2p0).

The above equality is going to be fulfilled if and only if

H11(p1p0)=H22(p2p0)

implying

p0=H11p1+H22p2H11+H22.

Taking into account (17), we get the explicit formula for p0

(19)p0=i=12pii|ωi|r2ωii=121i|ωi|r2ωi (19)

In the case of Newtonian fluid (r=2), (19) reduces to standard Kirchhoff’s law derived in [2].

Remark 3Let us compute the formula (19) in the case of two circular pipes, i.e. for Si={y ∈ R2: |y|<Ri}, i=1, 2…. Taking into account (15), we get

ωi=π2r/21μr1(r+2)Rir+2

(20)p0=p11R13r2+p22R23r211R13r2+12R23r2 (20)

To illustrate the result, we investigate how the value of p0changes with respect to the flow behaviour index r. In view of that, the junction pressure p0is plotted in Figure 3 as a function of the index r. A wide range for r is taken (1<r<3), and results for different magnitudes of prescribed pressures p1and p2are presented.3We observe an interesting phenomena: if characteristic dimension of one pipe is much larger than the other one (green and yellow lines), the junction pressure significantly changes for 1<r<2 and then stabilises for r>2. It means that, for dilatant flows, the junction pressure remains (almost) the same no matter what r we choose. On the other hand, if we have two pipes of similar dimensions, no such phenomena are observed.

Figure 3:

The junction pressure p0 as a function of the index r.

In general case (for arbitrary, m>2), the junction pressure p0 is given by (18). Evidently, it is not reasonable to expect that we can explicitly solve such nonlinear algebraic equation, as we managed in the case m=2. The best we can do is to prove that there exists a unique p0 such that

(21)F(p1,p2,,pm,p0)=0,F(p1,p2,,pm,p0)=i=1m|Hii(pip0)|r2Hii(pip0). (21)

For given pressures p1, p2,…, pm, we consider a function tG(t) where

(22)G(t)=F(p1,p2,,pm,t). (22)

Note that

(23)limt+G(t)=,limtG(t)=+, (23)

so to prove the existence, it is sufficient to show that tG(t) is a strictly decreasing function. To accomplish that, we employ the following inequalities:

(24)for r>2:C1(|b|r2+|a|r2)|ba|2(|b|r2b|a|r2a)(ba), (24)
(25)for 1<r<2:C2|ba|2|b|2r+|a|2r(|b|r2b|a|r2a)(ba). (25)

The above inequalities hold true for any a, bRn, where C1 and C2 denote positive constants independent of a and b. For detailed proof of (24) and (25), we refer the reader to [13], Lemma 3.1. Depending on the value of the exponent r from the power law, we deduce as follows:

• Case 1<r<2 (i.e. r′>2):

(ts)[G(t)G(s)]==i=1m{|Hii(pit)|r2Hii(pit)|Hii(pis)|r2Hii(pis)}(st)C1i=1m{|Hii(pit)|r2+|Hii(pis)|r2}Hii|st|2=C1i=1mωii{|piti|r2+|pisi|r2}|st|20,for any t,sR.

Obviously,

(ts)[G(t)G(s)]=0  t=s,

so we conclude that G is a strictly decreasing function. Therefore, due to (23), there exists tR such that G (t)=0. That proves the existence of the junction pressure p0 satisfying (18). To prove the uniqueness, we assume the opposite. Then, there exist p0 and q0 such that G(p0)=G(q0)=0. But then (p0q0)[G(p0) – G(q0)]=0 implying p0=q0.

• Case r>2 (i.e. 1<r′<2):

Instead of (24), here we use the inequality (25) to obtain

(ts)[G(t)G(s)]=C2i=1m1|Hii(pit)|2r+|Hii(pis)|2rHii|st|2=C2i=1mωii1|piti|2r+|pisi|2r|st|20, for any t,sR.

Using the same arguments as above, we deduce the existence and uniqueness of p0 in this case as well.

## 5 Interior-Layer Correction

From the physical point of view, one should expect that the fluid behaviour is somewhat different in the vicinity of junction than far from it. Thus, in this section, we construct the interior-layer corrector capturing those effects and, consequently, improving the order of accuracy of our approximation.

We start by introducing the interior layer by

Gε=Ωε0(i=1mGεi),

where

Gεi={(x1i,xi)Ωεi:0<x1i<dεlnε}.

Parameter d>0 describes the thickness of the interior layer and, as such, it will play an important role later on in the error analysis (see Sec. 6). If we pass to the fast variable y=x/ε, namely if we substitute uε(x)=εrV(x/ε),pε(x)=εΠ(x/ε) into (3) and (4), we obtain

(26)μ div{|D(V)|r2D(V)}+Π=0 in G, (26)
(27)divV=0 in G. (27)

The above equations are posed in the domain G=Ω0i=1mΩi consisting of a bounded reservoir Ω0 and m infinitely long pipes Ωi={y=(y1i,yi)R3:y1i>0,yi=(y2i,y3i)Si}. Since far from the reservoir, the flow in pipes remains in the form of PLPF (see Sec. 3), we endow (26) and (27) with the conditions:

(28)V=0 on G, (28)
(29)|y|+limV=Ui in Ωi,i=1,,m, (29)

where Ui is given by (11). The Stokes problem (corresponding to the value r=2 in our case) in an infinite domain with asymptotic Poiseuille velocities is called Leray’s problem, and an overview of the related mathematical theory can be found in [14]. Here, the system (26–29) forms its non-Newtonian generalisation and has been theoretically analyzed in [15] by the first author of the present article. In the case of dilatant flow (r>2), the existence and uniqueness of the solution have been established as well as the exponential decay estimate for the solution. To our knowledge, the similar result for pseudo-plastic flows (1<r<2) is still missing and cannot be found in the existing literature. For that reason, we are in position to define the correction in the interior layer Gε only for r>2. Namely, it has the form

(30)Vε(x)=εrV(xε),Πε(x)=εΠ(xε), (30)

with (V, Π) being the solution of the non-Newtonian Leray’s problem (26–29).

Note we still need to perform some matching with the outer flow. Indeed, the PLPF (Uεi,Pi) representing the solution in

ωεi={(x1i,xi)Ωεi:dεlnε<x1i<i}

and the interior-layer corrector (Vεi,Πεi) in Gε are only matched asymptotically (for ε → 0), due to (29). Thus, on each interface,

γεi={(x1i,xi)Ωεi:x1i=dεlnε}

between Gε and ωεi, our asymptotic approximation is not continuous but it holds

Vε=Uεi+Eεi,  Eεi(x*i)=εr[V(dlnε,xiε)Ui(x˜iε)].

Thanks to the exponential decay result of the Leray solution (see [15], theorem 5.1), the residual Eεi can be controlled. Moreover, by a simple change of variables from (29), it follows

γεiEεie1idxi=εr+2Si[V(dlnε,yi)Ui(yi)]dyi=0.

That enables us to introduce hεi as the solution of the standard divergence problem (see, e.g. [14], chapter III.3):

(31)divhεi=0 in ωεi,hεi=εrEεi on γεi,hεi=0 on ωεi\γεi. (31)

Finally, we define our (matched) asymptotic solution as

(32)Uεapprox=εrUεapprox, where Uεapprox={Ui(xiε)+hεi in ωεiV(xε) in Gε, (32)
(33)Pεapprox={Pi in ωεiΠε in Gε (33)

## 6 Error Analysis

We restrict ourselves only on the case r>2. The error estimate can be derived in the similar way as in [2], using, in addition, the monotonicity of the differential operator given by the power law. Indeed, it can be verified by direct computation that, for some C, c>0 (independent on ε) it holds

(34)cΩε|D(Rε)|rμΩε[|D(uε)|r2D(uε)|D(Uεapprox)|r2D(Uεapprox)]D(Rε)==i=1mγεiJεiRεC|εlogε|11/ri=1m|Jεi|Lr(γεi)|D(Rε)|Lr(Ωε), (34)

where Rε=uεUεapprox and Jεi, in ith pipe, is the jump of the normal component of the stress tensor over the surface γεi. Essentially, we have

Jεi=ε[Π(d|lnε|,xiε)d|lnε|pip0i+μ|D(V(d|lnε|,xiε))|r2D(V(d|lnε|,xiε))μ|D(Ui(xiε))|r2D(Ui(xiε))]ei=O(εd+1).

Taking into account that the measure of the surface |γεi|=O(ε2) gives

|D(Rε)|Lr(Ωε)Cε(4+d3r)/(r1)|lnε|1/r.

Thus, following the similar procedure as in [2], we obtain for r′=r/(r – 1):

(35)1|Ωε|Ωε|D(uε)D(Uεapprox|rCεβ|lnε|,β=r(r1)(4+dr)3r12. (35)
(36)1|Ωε|Ωε|uεUεapprox|rCεβ+r|lnε| (36)
(37)1|Ωε|Ωε|pε𝒫εapprox|rCεσ|lnε|1/(r1)2,σ=r(4+d)3(r1)3r2. (37)

For that velocity estimate to make sense, we need to have β + r>0, which is true for any r>2 and d>0. For the pressure estimate, we need σ>0, leading to the more restrictive condition

d>1r[(r1)2(3r2)+3]4.

Thus, the thickness of the interior layer is depending on the flow index r and increases with it. The magnitude of the interior-layer corrector for the velocity is of the order ε1/r (in Lr norm) and for the pressure ε1+1r, in Lr′ norm. So, the correctors improve our approximation for the velocity only in the case d>r – 1. On the contrary, for the pressure, the corrector makes the difference only if d>5r2 – 7r + 1 – 3/r, which is very large for r>2. Thus, it is reasonable to keep the velocity corrector if we want better approximation, while the pressure corrector is just a technical aid needed for deriving the satisfactory error estimate.

## 7 Conclusion

The study of non-Newtonian fluid flow has gained much attention in recent years because of its numerous technological applications. Several models of non-Newtonian fluids have been proposed and among all, the simplest and the most common one is the power-law model. On the other hand, engineering practice requires extensive knowledge of the pipe network flow. Therefore, the aim of this article was to consider the steady flow of a power-law fluid through a junction of thin (or long) pipes with possibly different profiles. No serious assumptions are imposed in order to simplify the original 3D problem describing the real physical situation. Using asymptotic analysis with respect to the pipes thickness, the effective behaviour of the flow has been found. In particular, Kirchhoff’s law (18) has been derived yielding the value of the junction pressure. The flow in the vicinity of the junction has been studied in the case of dilatant flow, and the corrected asymptotic approximation (32) and (33) has been built. The accuracy of the obtained approximation has been discussed by providing the corresponding error estimates [see (36) and (37)]. Due to the lack of the existence and decay result concerning the pseudo-plastic flow in unbounded domain, we are not in position to correct the PLPF in the vicinity of the junction for 1<r<2. This is the subject of our current investigation. Nevertheless, we believe that the obtained result is important from both theoretical and practical points of view and, hopefully, could improve the known engineering practice.

Corresponding author: Igor Pažanin, Faculty of Science, Department of Mathematics, University of Zagreb, Bijenička cesta 30, 10000 Zagreb, Croatia, E-mail:

## Acknowledgments

The authors have been supported by the Croatian Science Foundation (project 3955: Mathematical modeling and numerical simulations of processes in thin or porous domains). The authors would like to thank the referees for their helpful comments and suggestions that helped to improve the paper.

## References

[1] S. L. Rosen, Fundamental Principles of Polymeric Materials, Wiley, New York, 1993.Search in Google Scholar

[2] E. Marušić-Paloka, Asymptotic Anal. 33, 51 (2003).Search in Google Scholar

[3] M. Lenzinger, Asymptotic Anal. 75, 1 (2011).Search in Google Scholar

[4] F. Blanc, O. Gipouloux, G. Panasenko, and A. M. Zine, Math. Mod. Meth. Appl. Sci. 9, 1351 (1999).Search in Google Scholar

[5] G. Panasenko, Appl. Anal. 76, 363 (2000).Search in Google Scholar

[6] R. T. Balmer and M. A. Florina, J. Non-Newtonian Fluid Mech. 7, 189 (1980).Search in Google Scholar

[7] L. C. Zheng, X. X. Zhang, and L. X. Ma, Chin. Phys. Lett. 25, 195 (2008).Search in Google Scholar

[8] T. G. Myers and J. Low, Int. J. Therm. Sci. 70, 127 (2013).Search in Google Scholar

[9] L. Leftona, D. Weib, and Y. Liu, Int. J. Comput. Fluid D. 28, 351 (2014).Search in Google Scholar

[10] S. Marušić and E. Marušić-Paloka, Asymptotic Anal. 23, 23 (2000).Search in Google Scholar

[11] F. Delplace and J. C Leuilet, Chem. Engng. J. 56, 33 (1995).Search in Google Scholar

[12] Z. Tadmor and C. G. Gogos, Principles of Polymer Processing, Wiley, Hoboken, New Jersey, 2006.Search in Google Scholar

[13] D. Sandri, Math. Anal. Numer. 27, 131 (1993).Search in Google Scholar

[14] G. Galdi, An Introduction to the Mathematical Theory of the Navier-Stokes Equations, Vol. I, Springer-Verlag, New York 1994.10.1007/978-1-4612-5364-8Search in Google Scholar

[15] E. Marušić-Paloka, Math. Mod. Meth. Appl. Sci. 10, 1425 (2000).Search in Google Scholar