Mathematical Analysis of Non-Newtonian Fluid Flow through Porous Media

Mathematical Analysis of Non-Newtonian Fluid Flow through Porous Media

Gitte Vainkat Taterao

Department of Mathematics

N.H. College of Engineering, Parli Vaijnath Dist.Beed, Maharashtra -431515

Abstract – The movement of non-Newtonian fluids through porous surfaces regulates numerous natural and engineering tasks such as polymer-enhanced oil recovery, hydraulic fracturing, contamination of soil and aquifers, filtering, and biological fluids movement. Since such fluids exhibit a nonlinear behaviour, the classical linear Darcys law is no longer applicable, and new macroscopic transport equations need to be derived from the pore-scale rheology. In this paper, mathematical analysis of single-phase creeping flow of inelastic and weakly-viscoelastic fluids flowing through porous media is conducted. The mathematical treatment starts from Stokes-type pore-scale equations with shear-rate dependent viscous behaviour. Closed-form flow relations for power law, Bingham, and Herschel-Bulkley fluids in capillaries are obtained, after which upscaling with tortuous capillary bundle closure allows deriving generalized Darcy law with explicit expression for effective viscosity. In order to demonstrate the connection between the Darcy law and the phenomenon of a state in the context of the power law, we will show that the problem reduces to a p-Laplace problem with known theorems on the solvability and regularity based on the properties of monotone operators. Yield-stress fluids lead to variational inequalities and thresholds of the percolation type. We also perform numerical tests on a poly-disperse capillary network illustrating the Darcy law of the fluid as well as the transition process of power law fluid and also the threshold pressure gradient of yield-stress fluids. Finally, we discuss the new results from pore-scale studies as well as the experiments providing evidence on the elastic instabilities which limit the application of viscosity models.

Keywords: non-Newtonian fluids; porous media; generalised Darcy law; power-law fluid; Herschel Bulkley fluid; pore-network modelling; p-Laplacian; yield stress; elastic turbulence

1. INTRODUCTION

   The experiments conducted by Darcy on the fountains of Dijon led to the first macroscopic identification and description of the flow of fluids in a porous medium, resulting in the establishment of a relationship between the discharge and pressure gradient [1]. Later this relation was derived from the Stokes equations by using methods of volume averaging and homogenization [2]. The importance of this derivation may be realized from the fact that it forms the basis of hydrogeology and reservoir engineering. The derivation of the above law assumes that the fluid is a Newtonian one, possessing constant viscosity. However, numerous industrial fluids violate this assumption: the viscosity of polymer solutions, drilling muds, cement slurries, foams, emulsions, colloidal suspensions, and biological fluids depends on either the rate of deformation or the deformation history or both of them, in addition to being able to present yield stress and/or elastic energy storage capabilities [3].

   This topic is of great significance in the field of improved oil production, filtering polymeric solutions, and the working out of methods for the removal of liquid contaminants from soils [3]. In unconventional oil fields, the rheology of fracturing fluids is a determinant for the proppant transport and leak-off [4], while in geotechnical engineering, shear-thinning polymeric solutions serve as fluid support solutions while excavating and their penetration into the surrounding soil must be predicted as well [5].

   Savins [6] reviewed the preceding work in the area of classical literature, while the earliest theoretical advancements were derived from capillary-bundle models [7]. In Sochi, four main approaches to the modelling of single-phase flow have been identified continuum models, models using bundles of tubes, numerical methods and pores scale network modelling [3]. The network concept has been used by Sorbie et al. [8] in connection with pseudoplastic fluids and later extended to packed beds [9].

   There are still three questions that have only been partially addressed. One is how does bulk rheology relate to the viscosity at Darcy's scale (in situ)? The empirical shift factor that links these two quantities is said to vary by three orders of magnitude depending on the particular fluidmedium pair [10]. The second question is how is the yield-stress fluid able to invade disordered voids? New studies show that the flow exhibits continuous transition controlled by the applied pressure gradient obeying a nonlinear law of Darcy near the threshold [11, 12]. Third, what occurs when the fluid is viscoelastic? In case of a polymer solution, the macroscopic resistance increases abruptly when flow exceeds some critical rate in the porous media but does not change in bulk. This effect is due to elastic instability making the flow exhibit some turbulence-like oscillations at very low Reynolds number [13]. Answers to these questions have been accelerating in recent times with development of pore-scale simulations of polymer flow [14], a generalized Darcys law based on pore net validated on sphere packings [5] and quite a new pore network model for yield-stress fluid undergoing slip at the wall [15].

   The current work consolidates these advances within a unified mathematical framework and provides concrete derivations and computations. The main achievements consist of: (i) explicit capillary relationships for power-law and Herschel-Bulkley fluids and their upscaled versions, including the effective viscosity resulting in the cancellation of tortuosity; (ii) the well-posedness study of the macroscopic power-law problems and Bingham problem; and (iii) numerical examples of the flow regimes obtained and the flow-channeling associated with shear-thinning. Section 2 gives a summary of the rheological models, Section 3 summarizes the governing equations, Section 4 presents the analysis, Section 5 describes the computational framework, and Section 6 gives the conclusions and outlines the limitations of the obtained results.

2. RHEOLOGICAL MODELS

   For a generalised Newtonian fluid, the deviatoric stress is = 2()D, where D is the rate-of-strain tensor and its scalar magnitude; the classification into time-independent, time-dependent and viscoelastic fluids follows Sochi [3]. The models being discussed in the paper are shown in Table 1. The power-law (Ostwaldde Waele) model is the most basic explanation of shear thinning (n < 1) and shear thickening (n > 1) but it does not reproduce the Newtonian plateaus that polymer solutions show at low shear rates and the Ellis model predicts very low viscosity at high shear rates whilst the Carreau [16] and similar models give the complete S-shaped plot of the viscosity as a function of shear rate [10]. However, the Bingham and HerschelBulkley models [17] account for yield stress effects. Viscoelastic effects are introduced into the models by means of an additional parameter called relaxation time that can be included, for instance, in the Oldroyd-B model [18].

   Table 1. Constitutive models used in the analysis (: shear stress, : shear rate).

   Model	Constitutive relation	Parameters	Behaviour
   Newtonian	=		Linear
   Power-law	= K n1	K, n	n<1 thinning; n>1 thickening
   Bingham	= y + p ( > y); = 0 ( y)	y, p	Yield stress, linear post-yield
   Herschel Bulkley	= y + K n ( > y)	y, K, n	Yield stress, power- law post-yield
   Carreau	= + (0 )[1 + ()2](n1)/2	0, , , n	Two Newtonian plateaus
   Cross	= + (0 ) / [1 + ()m]	0, , , m	Two Newtonian plateaus
   Oldroyd-B	= s + p; p + p = 2pD	s, p,	Elastic, normal stresses

3. Governing Equations and Dimensionless Groups

   Consider the steady or slowly varying flow of an incompressible fluid in the pore space p of a rigid porous medium. Mass and momentum conservation read

   (1)

   (2)

   with u the velocity, p the pressure, the density and D = ½(u + uT). At the low Reynolds numbers typical of subsurface flow the inertial term is negligible and Eq. (2) reduces to the Stokes-type system

   (3)

   accompanied by slip conditions on the surface of solid (wall slip is discussed in Section 5). The procedure to solve Eq. (3) in a realistic three-dimensional pore shape is costly; hence the need for upscaling. In the case of Newtonian fluids, volume averaging or homogenization results in the Darcy law[2].

   (4)

   In this case, q relates to the Darcy velocity (also known as superficial velocity), k refers to the intrinsic permeability and denotes viscosity. In terms of non-Newtonian fluids, the filtration law in its homogenised form is a nonlinear one and does not allow separating the pore scale and macroscopic law into independent

   issues for any rheological behaviour. The common practice to deal with non-Newtonian filtration is to preserve the structure of Darcy law but to introduce effective velocity-dependent viscosity [3, 5, 10]:

   (5)

   In this case, q relates to the Darcy velocity (also known as superficial velocity), k refers to the intrinsic permeability and denotes viscosity. In terms of non-Newtonian fluids, the filtration law in its homogenised form is a nonlinear one and does not allow separating the pore scale and macroscopic law into independent issues for any rheological behaviour. The common practice to deal with non-Newtonian filtration is to preserve the structure of Darcy law but to introduce effective velocity-dependent viscosity:

   with the inertial coefficient. Creeping flow is retained provided the permeability-based Reynolds number

   (6)

   (7)

   is limited in scope; the Newtonian principle's validity in terms of Reynolds number still needs to be studied for non-Newtonian flow. The problem in regards to the discussed non-Newtonian flow is represented by the following dimensionless groups: n (the flow index); Bn = y/[K(U/)n] (a dimensionless number that relates essential and viscous stress); and Wi = U/ (the ratio of the polymer relaxation time to the flow time scale at some characteristic velocity).

4. Mathematical Analysis

   1. Capillary flow of a power-law fluid

      Let a fluid with = Kn1 flow steadily through a straight capillary of radius R under a pressure gradient G

      = P/L. A momentum balance on a cylindrical fluid element gives the linear stress distribution rz = Gr/2, so that

      (8)

      Integrating with u(R) = 0 gives the velocity profile

      (9)

      and a second integration over the cross-section gives the volumetric flow rate and mean velocity

      (10)

      (As seen above, the numerical integration of formula (9) is consistent with formula (10) up to machine precision). Then, in order to normalize the obtained profile with respect to the average velocity, the shear rate at the boundary is calculated (taking into account the Rabinowitsch correction to the Newtonians' value 4/R),

      (11)

      We find ourselves confronted with two important consequences. First, when the values of n begin to decrease, it will be noted that the profile assumes a flatter shape: the ratio of the centerline velocity to the mean velocity, given by the expression (3n + 1)/(n + 1), will equal 1.46 for n = 0.3, 1.67 for n = 0.5, 2 for the Newtonian fluid, and 2.2 for n equal to 1.5. Second, the scaling of the volume flow rate is also changed in terms of the variable R: for a Newtonian fluid, it being represented as R(3n + 1)/n; for instance, R4 for a Newtonian fluid, but R5 for n = 0.5 and R6.3 for the case when n equals 0.3, at the same pressure gradient. Thus, shear-thinning fluids can be called more sensitive to the size of the pore spaces as compared to Newtonian fluids, which explains the channeling process that was mentioned in Section 6.

   2. Upscaling by a tortuous bundle of capillaries

      In its most basic form, the method of upscaling simply involves substituting the porous medium with a series of capillary tubes with a tortuosity of T, an average radius of R, and a porosity of . The pressure gradient in the capillaries is G/T, and the Darcy velocity is calculated as the interstitial speed multiplied by /T. In the case of a Newtonian fluid, the permeability can be calculated from the Hagen-Poiseuille lawk = R²/(8T²):

      (12)

      Inserting the mean velocity of Eq. (10) with G replaced by G/T gives, for a power-law fluid,

      (13)

      so that the Darcy velocity is a power of the pressure gradient with exponent 1/n. Defining the effective viscosity through eff = kG/q and eliminating R and T in favour of k and , all tortuosity dependence cancels and one obtains

      (14)

      (the algebra was verified numerically for several values of T). Equation (14) has the structure eff K qn1(k)(1n)/2 which is known from the classic analysis of packed-bed systems, where the numerical constant is a function of the package closure. Comparing it to equation (5), the bundle model instead gives the displacement factor

      (15)

      which is variable within narrow limits between 1.10 (n = 1) and 1.16 (n = 0.3). In practical cases, however, the does not remain constant; numerical simulation at the pore level shows that for permeable media typically

      increases, as the shear rate and local viscosity change very much in the area of low permeability due to the tortuosity effect, which is not accounted for in the capillary bundle model. This means that equation (14) should be used only as an approximation allowing a crude estimate of the real sysstems behavior in terms of order-of- magnitude eff qn1.

   3. Fluids with Newtonian plateaus: Darcy, transition and power-law regimes

      An ideal power-law fluid indicates zero viscosity when still and infinite viscosity in a state of high shear, which is unrealistic the he Carreau model

      (16)

      provides a more realistic description, as well as the very similar truncated power law, where the fluid has higher and lower viscosity levels [5]. The local flow curve slope can be defined as

      (17)

      All capillaries operate under high Newtonian viscosity regime with low pressure gradients (i.e., 0). The Darcy law is then applicable, with S = 1. At high pressure gradients, the majority of the cross-section is sheared under power-law conditions and S 1/n. The two regimes are separated by a transition between Poiseuille and two-layer profiles. Pore-network simulations on sixty sphere packings have illustrated just this behavior and enabled one to formulate a generalized Darcy law including the measurements of PHPA polymer solutions [5]. The polidispersed bundle simulation results represented in Fig. 2 also demonstrate the same stages of viscosity change: the value of the slope equals 1.00 at low values of G, reaches a maximum of .89 (close to 1/n = 2) in the power-law region, and decreases again (up to 1.43 at maximal gradient) on approaching the lower plateau.

   4. Yield-stress fluids and the threshold pressure gradient

      For a HerschelBulkley fluid, = y + Kn for > y, the fluid in a capillary can move only if the wall shear stress w = GR/2 exceeds y. The unyielded plug radius and the threshold gradient are

      (18)

      For G > Gc the sheared annulus rp < r < R satisfies du/dr = [(Gr/2 y)/K]1/n. Writing Q = r²(du/dr)dr and changing the variable to the local shear stress gives the closed form

      (19)

      which was checked against direct numerical quadrature for n = 0.6, 1 and 1.3. For n = 1 it reduces to the BuckinghamReiner equation

      (20)

      with p = K. In a bundle with an area-weighted pore-radius density fA(R), only capillaries wider than Rc = 2yT/G conduct, so that

      (21)

      with taken from Eq. (19). In the bundle the threshold is set by the widest capillary, Gc = 2y/Rmax, The flow curve begins with 0 and proceeds along with some concentration of pore types (Fig. 3). Thus, in the interconnected system such a case lacks completeness. The fluid flow needs to trace out a continuous succession of the pores that can bear the pressure, implying that the threshold is a property of the system as a whole, linked with the minimum-threshold path and percolation [11]. Flow in the two-dimensional disordered media undergoes continuous phase transition governed by the pressure gradient and exhibiting non-linear Darcy's law in the vicinity of threshold [11]. The types of flow can be connected to flow-field statistics [12].

      (22)

      with a non-trivial exponent . The concept of universal scale with respect to yielding in porous media has been there in Chaparians work [20]. Most recently, a pore-network model has been proposed which makes yielding as well as post yielding transport follow from the pressure-flow relation of a converging-diverging throat without fitting resistance parameters, with optional wall slip.

   5. Well-posedness of the macroscopic problems

      For a homogeneous medium with the power-law generalised Darcy law q = An|p|(1n)/np (the vector form of Eq. (13)) and ·q = 0, the pressure satisfies the p-Laplace equation

      (23)

      The equation is degenerate (the effective diffusivity |p|s2 vanishes where p = 0) for shear-thinning fluids, s > 2, and singular for shear-thickening fluids, 1 < s < 2.

      0

      Proposition 1. Let d be a bounded Lipschitz domain, n > 0, s = 1 + 1/n and pD W1, s(). Then the Dirichlet problem for Eq. (23) has a unique weak solution p pD + W 1, s(), which is the unique minimiser of

      (24)

      Moreover p C1,loc() for some (0,1).

      0

      Sketch of proof. The functional is convex, continuous and, by the Poincaré inequality, coercive on the reflexive space W 1, s() for 1 < s < ; The calculus of variations direct method provides a minimizer, while the strict convexity of |a|s leads to the uniqueness of the solution, and the Euler-Lagrange equation derived from is the weak form of Eq. (23). Moreover, the operator p ·(|p|s2p) is strictly monotone, for instance through

      (25)

      and coercive and hemi continuous, so that the MintyBrowder theorem applies [21]. The interior C1, regularity of p-harmonic functions is classical [22].

      A sharp result is obtained in the following sense: for shear-thinning fluids, the solution is usually not C2, as the degeneracy of the equation creates stagnant zones where the gradient of the pressure is low; this happens mathematically for the poorly swept region of a polymer flow. The significance of is that the infinite power- law flow minimizes the total power loss.

      In case of yield-stress fluids, the constitutive equation becomes set-valued for = 0 and the problem at the pore scale turns into a variational inequality: find u in the space of divergence-free fields that vanish on the solid boundary.

      (26)

      This holds good for all values of v that are acceptable, with f serving as the input and j being both convex and lacking any derivatives. The particular features of the function being minimized provide for the existence and uniqueness of u in this scenario (24). It is j's unresponsiveness that leads to the presence of the free boundary separating the yielded regions and the non-yielded ones and to the threshold behaviour described in Section 3.

   6. Viscoelastic effects

      Polymer solutions are often viscoelastic. In the Oldroyd-B model the polymer stress obeys

      (27)

      and the relevant dimensionless group is the Weissenberg number, together with the ratio of the measured to the purely viscous pressure drop:

      (28)

      Purely viscous models predict 1. Above a critical Weissenberg number Wic, however, increases sharply. Using direct three-dimensional visualisation methodology, researchers have established that the onset of an elastically unstable state causes strong spatio-temporal fluctuations comparable to the inertial turbulence at low Reynolds numbers with energy dissipated during such pore-scale fluctuations accounting quantitatively for the increase in flow resistance [13]. In later findings it was shown that stagnation point affects the chaotic fluctuations [24], while a review published recently outlines issues related to movement of complex fluids in

      porous media [25]. In this case, it seems evident that generalized Darcy law discussed in Sections 4.2-4.4 applies properly only in Wi < Wic.

5. Computational Framework

   1. Continuum discretisation

      The nonlinear problem of Eq. (23) can be solved by finite elements or finite volumes with a Picard (lagged- diffusivity) iteration. To avoid the degeneracy at p = 0 the diffusivity is regularised with a small parameter > 0,

      (29)

      Which converges monotonically to the minimization form (24). It is common to regularize yield stress fluids as well (in particular, by exponential viscosity regularization) or to treat them with variational inequality methods. On the pore scale, finite-volume solvers, finite elements, and lattice Boltzmann methods have been used for shear-thinning flows, and a newly developed open-source pore-scale simulator based on Open FOAM has analyzed the difference between the bulk and in situ viscosity of rock samples [14].

   2. Pore-network modelling

      A pore network includes void spaces represented as connected pores (nodes) by throats (links) [8, 9]. A system of equations representing conservation of mass at each node i is formulated based on the nonlinear throat relationship Qij(Pij), referred to from the capillary solutions provided in Section 4, and this system gets solved using the Newton-Raphson method [5].

      (30)

      where J is the Jacobian (conductance) matrix and R the residual vector. The procedure is:

      1. Extract or create the network, and set the radius, length, and shape of the throat.

      2. Decide on the rheology, and calculate the throat relation Qij(Pij) (Eq. (10), Eq. (19), or a numerical inversion method for Carreau fluids).

      3. Initialize the pressure field by using the Newtonian solution (linear).

      4. Construct the conductance matrix using gij = Qij/Pij and obtain the increment P.

      5. Repeat the above step until the residual and relative change of pressure fall below the acceptable limit.

      6. Calculate the value of q, eff = kG/q and the measured shift factor keeping the imposed gradient in a certain range in order to get the flow curve.

         A useful reference calculation is provided by Eq. (14). Table 2 lists the effective shear rate and viscsity obtained for a medium with k = 1010 m², = 0.35, a consistency index K = 1 Pa·sn and a Darcy velocity of 104 m/s.

         Table 2. Shift factor, effective shear rate and effective viscosity from the bundle closure (k = 1010 m², = 0.35,

         K = 1 Pa·sn, q = 104 m/s).

         n	n	eff (s1)	eff (Pa·s)
         0.3	1.161	19.6	0.124
         0.5	1.131	19.1	0.229
         0.7	1.115	18.9	0.414
         0.9	1.105	18.7	0.746
         1	1.101	18.6	1.000

6. Results and Discussion

   1. Velocity profile and shear rate

      The normalized profiles shown in Figure 1 are based on Equation (11). Shear thinning liquids have developed plug-like profiles where shear is mostly on the outer walls; at n = 0.3, the velocity at the center line is just at 46% of its average value whereas it equals 100% for Newtonian fluids. Shear thickening liquids exhibit the opposite behavior. The wall shear provided here exceeds its Newtonian equivalent (4/R) by a factor of (3n+1)/(4n) which is equal to 1.58 for n = 0.5. This equation demonstrates that shear rate derived from Darcy velocities cant be used for estimating shear at pores.

      Figure 1. Dimensionless velocity profiles of a power-law fluid in a capillary, Eq. (11), for different flow indices n.

   2. Flow curves of a polydisperse bundle

      Figure 2 shows the Darcy velocity of a bundle with a lognormal pore-radius distribution (median radius 40 m, logarithmic standard deviation of 0.35, = 0.35, T = 1) Newtonian fluid, power-law fluid of index n =

      0.5 and same index Carreau fluid. The Newtonian curve has a slope of 1 as dictated by Darcys law while the power-law curve has a slope of 1/n = 2. The Carreau fluid follows the Newtonian line at low gradients, coincides with the power-law line in the range in between, and diverges again at high gradients. The middle zone is limited because the two plateaus are two decades apart in terms of viscosity; greater separation between plateaus would lead to a wider range. This is what provides the basis for the three-regime concept of pore- network simulations [5] and explains why a single index cannot be used over the complete range of pressure gradients.

      Figure 2. Darcy velocity versus pressure gradient for a polydisperse capillary bundle: Newtonian, power-law and Carreau fluids (n = 0.5, 0 = 2 Pa·s, = 0.02 Pa·s, = 1 s). Computed from Eqs. (9)

      (10) and (16)(17).

   3. Yield-stress fluids

      Figure 3 compares Newtonian, Bingham and HerschelBulkley fluids of equal viscosity at a shear rate of 10 s1, with y = 1 Pa, in a bundle truncated at Rmax = 162 m. The threshold of Eq. (18) is Gc = 2y/Rmax = 1.23

      × 104 Pa/m. With the movement nothing flows below that point; from that point upward there are only couple of very big openings for flow to take place and the flow can only be recognized from about 2 to 3 units of Gc. After this moment the Bingham's curve remains somewhat straight with an offset in relation to what comes on the form of the HerschelBulkley's curve, which can show nonlinear characteristic. In case of a connected network, the threshold will be higher than in a case of a bundle of connected throats since the path of yielding has to exist [11]. But the value of the bundle is less than that of the network, so it cannot be an upper limit for the threshold.

      Figure 3. Flow curves of a bundle with truncated lognormal pore radii for Newtonian, Bingham (y = 1 Pa) and HerschelBulkley (n = 0.6, y = 1 Pa) fluids. The dotted line marks the bundle threshold Gc = 2y/Rmax. Computed from Eqs. (19)(21).

   4. Flow channelling by shear thinning

      In a system where capillary flow rate is related to the radii of pores, shear thinning favors broader pores. The figure shows the flow rate fraction carried by pores of different sizes while the pressure gradient in the system is fixed. Thus, the pores with the radius greater than 1.25 average size carry 78% for Newtonian fluids, while the flow rates for n = 0.7, 0.5, and 0.3 are 82%, 87%, and 94%, respectively. The poration process in these pores is characterized by the low values of effective viscosity, which enhances the phenomenon observed. Pore-network models support the claims made since the shear thinning leads to the changes in flow distribution, with high flows observed in areas of low viscosity values [5]. This phenomenon explains the reduction of the displacement efficiency of fluids with shear-thinning properties compared to Newtonian ones during the displacement with the lower efficiency than that predicted based on effective viscosity characteristics.

      Figure 4. Cumulative fraction of the total flow carried by pores of radius up to R at fixed pressure gradient, for power-law fluids with different n (lognormal pore-size distribution, median 40 m, = 0.35). The dotted line marks 1.25 times the median radius.

   5. Limitations

      The assessment has some shortcomings. The bundle closing neglects the pore connections and thus is not capable of portraying the percolation limit, the dependence of the tortuosity on the shift coefficient, or the dependence of the pore dimension on the throat position. [9,14] Interactions between the polymer and the wall

      such as adsorption, depletion zone, and mechanical destruction are not incorporated. The DarcyForchheimer term and the inertial regime are excluded from the mentioned numerical examples. [19] Above all, it is worth noting that the solely viscous approach is not appropriate when Wi is greater than Wic when the elastic instability leads to pressure drop that is significantly higher than the prediction of the generalized Newtonian approach. The experimental verification of the numerical trends should involve the use of rheometric data for the same liquid and for the range of the shear rates that are really available in the pores.[25]

7. CONCLUSIONS

This paper has analysed single-phase flow of non-Newtonian fluids in porous media from the pore scale to the Darcy scale. The main conclusions are:

For a power-law fluid the capillary flow rate scales as G1/nR(3n+1)/n; the tortuous-bundle upscaling gives the closed-form effective viscosity of Eq. (14), eff qn1(k)(1n)/2, in which tortuosity cancels, and an almost constant shift factor 1.10 1.16.

1. Fluids of Newtonian type have Darcy (S = 1), transitional and power-law (S 1/n) regimes at the same time, which agree with current results obtained with pore networks.

2. The macroscopic problem of the power-law type is the same as the p-Laplace equation with s = 1 + 1/n; it is well posed and has a unique variational solution together with regularity C1,.

3. Vinaptic fluids have variational characteristics in that they lead to the appearance of a variational inequality and a threshold value employed in the network which is a permeability characteristic.

4. Shear-thinning changes the flow due to the structure of the widest pores; for n = 0.3 pores wider than

   1.25 times median, i.e., passing 94% of the flow.

5. Instabilities of viscoelastic origin restrict the applicability limit of generalized Darcys laws (Wi < Wic).

Future research needs to integrate viscoelastic models of material behavior in conjunction with pore-network simulation methods, advance the yield-stress ethod to be applicable in three dimensions and for unsaturated soils, and devise new methods of determining the shift factor through pore-space statistics rather than use regression analysis for fitting.

REFERENCES

1. H. Darcy, Les Fontaines Publiques de la Ville de Dijon, Victor Dalmont, Paris (1856).

2. S. Whitaker, Flow in porous media I: a theoretical derivation of Darcy's law, Transp. Porous Media 1, 325 (1986).

3. T. Sochi, Flow of non-Newtonian fluids in porous media, J. Polym. Sci. B Polym. Phys. 48(23), 24372467 (2010). doi:10.1002/polb.22144

4. A.C. Barbati, J. Desroches, A. Robisson, G.H. McKinley, Complex fluids and hydraulic fracturing, Annu. Rev. Chem. Biomol. Eng.

   7, 415453 (2016).

5. S. Suo, S. Foroughi, M.J. Blunt, C. O'Sullivan, Pore-network modeling of polymer flow in porous media, Comput. Geotech. 182, 107142 (2025). doi:10.1016/j.compgeo.2025.107142

6. J.G. Savins, Non-newtonian flow through porous media, Ind. Eng. Chem. 61(10), 1847 (1969). doi:10.1021/ie50718a005

7. R.H. Christopher, S. Middleman, Power-law flow through a packed tube, Ind. Eng. Chem. Fundam. 4(4), 422426 (1965).

8. K.S. Sorbie, P.J. Clifford, E.R.W. Jones, The rheology of pseudoplastic fluids in porous media using network modeling, J. Colloid Interface Sci. 130(2), 508534 (1989).

9. M.T. Balhoff, K.E. Thompson, A macroscopic model for shear-thinning flow in packed beds based on network modeling, Chem. Eng. Sci. 61(2), 698719 (2006).

10. T. Shende, V.J. Niasar, M. Babaei, Upscaling non-Newtonian rheological fluid properties from pore-scale to Darcy's scale, Chem. Eng. Sci. 239, 116638 (2021).

11. C. Liu, A. De Luca, A. Rosso, L. Talon, Darcy's law for yield stress fluids, Phys. Rev. Lett. 122(24), 245502 (2019). doi:10.1103/PhysRevLett.122.245502

12. T. Chevalier, L. Talon, Generalization of Darcy's law for Bingham fluids in porous media: from flow-field statistics to the flow-rate regimes, Phys. Rev. E 91, 023011 (2015).

13. C.A. Browne, S.S. Datta, Elastic turbulence generates anomalous flow resistance in porous media, Sci. Adv. 7(45), eabj2619 (2021). doi:10.1126/sciadv.abj2619

14. M. Amiri, J. Qajar, A.Q. Raeini, A. Raoof, Assessing rheology effects and pore space complexity in polymer flow through porous media: a pore-scale simulation study, Water Resour. Res. 60(5), e2023WR036125 (2024). doi:10.1029/2023WR036125

15. C.P. Fonte, E. Sutton, K. Ohie, E. Doman, Y. Tasaka, A. Juel, Network modeling of yield-stress fluid flow in randomly disordered porous media, Appl. Phys. Lett. 128(24), 242703 (2026). doi:10.1063/5.0333161

16. P.J. Carreau, Rheological equations from molecular network theories, Trans. Soc. Rheol. 16(1), 99127 (1972).

17. W.H. Herschel, R. Bulkley, Konsistenzmessungen von Gummi-Benzollösungen, Kolloid-Z. 39, 291300 (1926).

18. R.B. Bird, R.C. Armstrong, O. Hassager, Dynamics of Polymeric Liquids, Vol. 1: Fluid Mechanics, 2nd ed., Wiley, New York (1987).

19. T. Tosco, D.L. Marchisio, F. Lince, R. Sethi, Extension of the DarcyForchheimer law for shear-thinning fluids and validation via pore-scale flow simulations, Transp. Porous Media 96(1), 120 (2013).

20. E. Chaparian, Yielding to percolation: a universal scale, J. Fluid Mech. 980, A14 (2024).

21. E. Zeidler, Nonlinear Functional Analysis and its Applications, Vol. II/B: Nonlinear Monotone Operators, Springer, New York (1990).

22. E. DiBenedetto, C1+ local regularity of weak solutions of degenerate elliptic equations, Nonlinear Anal. 7(8), 827850 (1983).

23. G. Duvaut, J.L. Lions, Inequalities in Mechanics and Physics, Springer, Berlin (1976).

24. S.J. Haward, C.C. Hopkins, A.Q. Shen, Stagnation points control chaotic fluctuations in viscoelastic porous media flow, Proc. Natl. Acad. Sci. U.S.A. 118(38), e2111651118 (2021).

25. M. Kumar, J.S. Guasto, A.M. Ardekani, Transport of complex and active fluids in porous media, J. Rheol. 66(2), 375397 (2022).