Nonlinear Flow Behavior in Packed Beds of Natural and Variably Graded Granular Materials

Abstract

Under certain flow conditions, fluid flow through porous media starts to deviate from the linear relationship between flow rate and hydraulic gradient. At such flow conditions, Darcy’s law for laminar flow can no longer be assumed and nonlinear relationships are required to predict flow in the Forchheimer regime. To date, most of the nonlinear flow behavior data is obtained from flow experiments on packed beds of uniformly graded granular materials (C_u = d₆₀/d₁₀ < 2) with various average grain sizes, ranging from sands to cobbles. However, natural deposits of sand and gravel in the subsurface could have a wide variety of grain size distributions. Therefore, in the present study we investigated the impact of variable grain size distributions on the extent of nonlinear flow behavior through 18 different packed beds of natural sand and gravel deposits, as well as composite filter sand and gravel mixtures within the investigated range of uniformity (2.0 < C_u < 17.35) and porosity values (0.23 < n < 0.36). Increased flow resistance is observed for the sand and gravel with high C_u values and low porosity values. The present study shows that for granular material with wider grain size distributions (C_u > 2), the d₁₀ instead of the average grain size (d₅₀) as characteristic pore length should be used. Ergun constants A and B with values of 63.1 and 1.72, respectively, resulted in a reasonable prediction of the Forchheimer coefficients for the investigated granular materials.

1 Introduction

In general, flow processes in groundwater flow studies can be described by Darcy’s law (1856):

$$i = - \frac{1}{K}q$$

(1)

where i [m/m] is the hydraulic gradient, K [m/s] is the hydraulic conductivity, and q [m/s] is the specific discharge. At higher flow velocities in porous media, Darcy’s law does not hold and the specific discharge starts to deviate from the linear Darcy’s law after a certain threshold hydraulic gradient. Forchheimer (1901) proposed alternative equations to describe post-Darcian flow behavior, which can be categorized in the two-term, the three-term and the power law (Bear 1972). Theoretical base for the different Forchheimer forms in a particular porous medium for given flow velocity ranges are provided in the literature (e.g., Hassanizadeh and Gray 1987; Whitaker 1996; Skjetne and Auriault 1999; Chen et al. 2001; Adler et al. 2013). The second-order polynomial two-term Forchheimer equation is widely used for groundwater flow problems and describes the nonlinear behavior with adding a quadratic term to Darcy’s law:

$$i = - aq - bq^{2}$$

(2)

where a [s/m] is a parameter equal to the reciprocal of the hydraulic conductivity (i.e., a = 1/K) and b [s²/m²] is the empirical Forchheimer coefficient. The onset of nonlinear laminar flow behavior is ascribed to different mechanisms of energy dissipation (e.g., Hassanizadeh and Gray 1987; Ma and Ruth 1993; Lage 1998; Nield 2000; Chen et al. 2001; Panfilov and Fourar 2006).

Ergun (1952) extended the Kozeny–Carman relationship to relate the Forchheimer coefficients to the grain size and porosity of a porous medium by:

$$i = - A\frac{{\left( {1 - n} \right)^{2} v}}{{gn^{3} d^{2} }} q - B\frac{{\left( {1 - n} \right)}}{{gn^{3} d}} q^{2}$$

(3)

where g [m/s²] is the acceleration due to gravity, d [m] is the characteristic pore length by means of the particle diameter, v [m²/s] is the kinematic viscosity of the fluid, n [−] is the porosity, A is the Ergun constant that equals 150 and B is the Ergun constant that equals 1.75 (Ergun 1952). Many variations on the Ergun relation for different kinds of porous media and flow regimes are given in the literature (e.g., Engelund 1953; Irmay 1964; MacDonald et al. 1979; Du Plessis 1994; Sedghi-Asl and Rahimi 2011; Erdim et al. 2015; Guo et al. 2019).

The understanding of high-velocity flows over steep hydraulic gradients through hydrogeological formations is important for a broad range of large scale flow problems. Such flow problems are, for example, the flow through waste deposits (d₁₀ = 10 mm) up to velocities of 10 m/s (Greenly and Joy 1996), transient flow due to steep head gradients by sudden head rise or drawdown in aquifers and constructed wetlands (Moutsopoulos and Tsihrintzis 2005; Moutsopoulos 2007), as well as flow through landfill caps or leachate collection systems up to velocities of 0.6 m/s (Bordier and Zimmer 2000). Moreover, in many groundwater applications the efficiency of injection/abstraction wells is affected by nonlinear head losses since flow velocities increase due to the divergence/convergence of flow lines in the vicinity of wells (e.g., Basak 1978; Mathias and Todman 2010; Yeh and Chang 2013; Wen et al. 2013; Houben 2015; Mathias and Moutsopoulos 2016). In general, nonlinear head losses will occur in the gravel pack while flow in the aquifer can be considered as linear Darcian flow for common well designs (Houben, 2015). However, for naturally developed wells with a well screen in direct contact with the naturally developed aquifer material, steep hydraulic gradients and nonlinear head losses occur in the aquifer material itself. Houben et al. (2018) stated that reduction in porosity in the gravel pack by clogging could have a significant influence on nonlinear head losses during well operation. Besides well hydraulics, the influence of nonlinear flow behavior through natural granular materials is important for analyzing slug tests in highly permeable sandy aquifers of K > 50 m/d (McElwee and Zenner 1998), as well as for analyzing packer test data (e.g., Yamada et al. 2005).

To date, a broad range of experimental studies on nonlinear flow behavior in different kind of porous media is available in the literature. Experimental datasets on nonlinear flow behavior through various kinds of packed beds of granular material show that the Forchheimer equation (Eq. 2) can be used to predict flow over a wide range flow velocities (e.g., Moutsopoulos et al. 2009; Sedghi-Asl et al. 2014; Ferdos et al. 2015; Salahi et al. 2015; Li et al. 2017; van Lopik et al. 2017). Most literature provide data on flow through packed beds of coarser granular materials, such as gravels and cobbles with an uniform grain size distribution (e.g., Sidiropoulou et al. 2007; Moutsopoulos et al. 2009; Sedghi-Asl et al. 2014; Salahi et al. 2015). The experimental studies of van Lopik et al. (2017) and Li et al. (2019) focus on small variations in gradation by adding only a small amount of both finer and coarser material to a given uniformly graded reference sand.

In contrast, however, for natural sand and gravel deposits the grain size distribution is usually non-uniformly graded. Only a limited amount of studies provides information about both the grain size distribution of the material and the nonlinear flow behavior through non-uniformly graded granular materials (e.g., Yamada et al. 2005; Dan et al. 2016; Wang et al. 2019). Hence, the understanding of nonlinear flow behavior through packed beds of graded natural sand and gravel deposits requires a much broader experimental dataset than provided in the present literature.

The aim of this study is to provide insight in the differences in nonlinear flow behavior through packed beds of a wide range of graded sand and gravel types. Therefore, in the present study the effect of nonlinear flow behavior through packed beds of granular material with wider grain size distribution and lower corresponding porosities is investigated. This is done by experimental tests on packed beds of natural sands and gravels from unconsolidated aquifers in the Netherlands and Germany, as well as artificial composite mixtures of filter sands and gravels considering a broad range of grain size distributions and porosities at a given reference average grain size (d₅₀) and compaction grade. The dataset of the present study on nonlinear flow behavior through packed beds of variably graded granular material is compared to the existing datasets in the literature on uniformly graded granular material. This is done to obtain insight in the difference in nonlinear flow behavior between these types of granular materials.

2 Experimental Data from Previous Studies

A wide variety of experimental studies on nonlinear flow behavior through packed beds of all kinds of material is provided in the literature. Investigation of fluid/gas flow through packed beds of artificial compounds is done, such as cubes, cylinders, ellipsoids, fibers, cubic arrays of spheres, as well rounded spheres of glass or metal (e.g., MacDonald et al. 1979; Sidiropoulou et al. 2007; Allen et al. 2013; Huang et al. 2013; Bağci et al. 2014; Erdim et al. 2015). Besides artificial compounds, the literature provides a broad dataset on nonlinear flow behavior through all kinds of granular material, ranging from sand to cobbles. In the present study, the selection of the experimental data on nonlinear flow behavior through granular material in Table 1 is considered in order to compare the results between uniformly graded material and non-uniformly graded material.

Table 1 Selection experimental datasets from the literature on nonlinear flow behavior through packed beds of (uniformly graded) granular material

Most studies have focused on flow through packed beds with uniform-sized grains (e.g., Sedghi-Asl et al. 2014; Salahi et al. 2015; Li et al. 2017). For the experimental studies on packed beds of granular material with wider grain size distributions, the data provides information on gradation by means of the coefficient of uniformity C_u [−] and the coefficient of curvature C_c [−] (e.g., Yamada et al. 2005):

$$C_{U} = \frac{{d_{60} }}{{d_{10} }}$$

(4)

and

$$C_{C} = \frac{{\left( {d_{30} } \right)^{2} }}{{d_{10} *d_{60} }}$$

(5)

In the studies of van Lopik et al. (2017) and Li et al. (2019), mixtures of granular material were used to obtain slightly more well-graded composite sands. Wider grain size distributions (higher C_u values) for a material at a given average grain size and porosity resulted in increased flow resistance and, hence, increased Forchheimer coefficients a and b.

Besides gradation, the compaction grade (e.g., Dan et al. 2016; Banerjee et al. 2018) and the grain shape (e.g., Allen et al. 2013; Salahi et al. 2015) of the granular material has a significant effect on nonlinear flow behavior. Packed beds of angular grains resulted in increased flow resistance compared to rounded grains (Salahi et al. 2015). The influence of particle shape is also shown for randomly packed beds of cubes and cylinders instead of rounded spheres with a diameter of 2.7 mm (Allen et al. 2013). The flow resistance for such materials is drastically increased, resulting in the highest Forchheimer coefficients for the packed beds of cubes.

3 Methods

The experimental setup as shown in Fig. 1 is used for the one-dimensional packed column flow experiments. This setup is similar to the one used for the nonlinear flow experiments in the study of van Lopik et al. (2017) and consists of the following main components:

Fig. 1

Schematic overview of the experimental setup: (1) inlet compartment, (2) outlet compartment, (3) packed column with porous medium, (4) pressure sensor, (5) circulating tank with pump, (6) thermometer and (7) volumetric flow meter

• The inlet and outlet compartment (1 and 2). During the flow experiments the pressure at the inlet compartment (1) is measured, while at the outlet compartment (2) the hydraulic head is fixed at 0.505 m.

• The packed column (length is 0.507 m, diameter is 0.098 m) with the porous medium (3). Besides the inlet and outlet compartment, the pressure is monitored at 0.13, 0.26, and 0.38 m distance from the inlet with a BD Sensors pressure transmitter combined with the BD Sensors PA 430 display (4). The interpolation of the measured water heads at a given specific flow rate resulted in a linear hydraulic head gradient over the entire column length for all experiments.

• The circulation tank of 1 m³ with a submersible pump (5). Due to pumping, water was significantly heated during the experimental process. Therefore, the water temperature is measured in the outlet compartment (6) and cold water is added to the circulation tank during the flow experiment to ensure constant water temperature at 20 °C (± 0.5 °C) and constant dynamic viscosity of approximately 0.001 Pa s. The volumetric flow rate Q [m³/h] is measured using a SVM F2 Energy Meter (7). At a given hydraulic gradient, the flow rate is determined after steady-state flow conditions were obtained.

3.1 Characteristics of Investigated Granular Material

Two different sets of porous media are investigated in the present study:

• Natural sand and gravel deposits obtained from drilling samples from aquifers in the Netherlands (M.1–3) and Germany (W.1–2).

• Composite mixtures of sand and gravel obtained by mixing of selected portions of different uniformly graded filter sands and gravels (S.1–3).

3.1.1 Natural Sand and Gravel Deposit Samples

All natural sand and gravel samples are obtained from drilling cores in aquifers. The natural sands M.1–3 are obtained by cable tool drilling at a location near the Münsterlaan at the Utrecht University Campus (the Netherlands). The samples are from a sandy unconsolidated aquifer and are fluvial Pleistocene deposits classified as Kreftenheye formation. The samples are obtained from two different drilling locations at only 12 m distance of each other. Sand sample M.1 is obtained at a depth of 12.50–14.50 mbgl. At the other location, the sampling depth is the same for M.2 (12.50–14.50 mbgl), while the depth is 9.50–10.50 mbgl for M.3. The natural sand W.1 and gravel W.2 are reverse-circulation drilling disturbed soil samples at depths of, respectively, 20–21 mbgl, and 19–20 mbgl. The samples are Pleistocene river terrace sediments obtained from a sandy unconsolidated heterogeneous-layered aquifer near Wyenhütte in Rheindahlen, about 10 km southwest of Mönchengladbach (Germany). The grain size distributions for the natural sand and gravel samples (M.1–3 and W.1–2) are obtained by sieve analysis using sieve sizes 0.063, 0.125, 0.25, 0.355, 0.5, 0.71, 1.0, 2.0, 5.3, 8, 16.0, 22.6 and 32.0 mm.

3.1.2 Composite Mixtures of Sand and Gravel

Similar to the study of van Lopik et al. (2017), uniformly graded filter sands (with a C_u < 1.5) are mixed to obtain sand and gravel material with the desired gradation (see “Appendix 1” for the mass fraction mixing ratios of the filter sands). For the composite mixture samples with wider grain size distribution, the required gravel (grain sizes larger than 10 mm) is obtained by selection at minimum grain sizes of 10, 15, 20 and 30 mm. In the present study, the filter sands with an average grain size diameter of 1.0 (Sample 1), 1.5 (Sample 2) and 6.3 mm (Sample 3) are used as reference material. The reference samples 1, 2 and 3 are, respectively, samples 7, 9 and 11 in the dataset of van Lopik et al. (2017).

3.1.3 Analyzation of the Granular Material

The gradation for both kinds of materials is described by means of the coefficient of uniformity C_u and the coefficient of curvature C_c (Eqs. 4, 5). The porosity of each packed sample is determined by weighing the dry mass of granular sand or gravel in the column. The solid phase density of the sand and gravel used for the composite mixtures equals 2.66 kg/m³ (data by FILCOM, Papendrecht, the Netherlands). Moreover, the porosity of each material was determined by measuring the used granular material volume from the packed column under fully saturated conditions in a given water volume. Equal porosities are obtained for the composite mixtures of sand and gravel [Samples 1.(1–5), 2.(1–5) and 3.(1–3)] with above-mentioned methods. The characteristics of the natural as well as the composite mixtures of sand and gravel are given in Table 2.

Table 2 The grain size distribution characterization and porosity values for all sand and gravel samples

Additionally, the grain shape properties of both the natural sand and composite sands are obtained using the CAMSIZER^® P4 (Retsch, Germany) image analysis technique (see “Appendix 2”). For all samples, the grain sizes in the range of 0.063–2.0 mm were selected in order to determine the grain shape of the sand. The grain shape properties of the used filter sands for the artificially mixed composite sands as well as for the natural deposit sands are similar. All investigated granular material is subangular-subrounded.

3.2 Preparation Procedure for the Packed Columns

Similar to the procedure in the study of van Lopik et al. (2017), the consistent packing method of Rietdijk et al. (2010) is used to obtain reproducible and homogenous packed porous media. To ensure homogeneous packing over the entire length of the column, small sand/gravel portions with the required mixture of fine and coarse material are carefully poured into the fully saturated column. The column is prepared, packed and installed under fully saturated conditions in a large container filled with water. In order to obtain the highest compaction grade, the deposited sand in the column is manually compressed throughout the packing process.

The consistent use of the Rietdijk et al. (2010) method resulted in the highest grain size density for a given sample (van Lopik et al. 2017). In Fig. 2, the effect of gradation on the porosity is shown for the samples listed in Table 2. For sand and gravels with wide grain size distributions (high C_u values), the fines tend to fill the void spaces between the larger sand particles, resulting in lower porosity values (Das 2008). For the relatively uniformly graded material (C_u < 3), the porosity values are higher than 0.33. The porosity is reduced for granular material samples with wider grain size distributions (C_u > 6), resulting in porosity values lower than 0.27.

Fig. 2

Relationship between the coefficient of uniformity (C_u) and the porosity values for all prepared column samples with the Rietdijk et al. (2010) packing method

3.3 Ergun Friction Factor and Reynolds Number

The Reynolds number is used to indicate whether fluid flow is laminar or turbulent. The most common definition of the Reynolds number is defined as (Bear 1972):

$${\text{Re}} = \frac{\rho dq}{\mu } = \frac{\text{inertial forces}}{\text{viscous forces}}$$

(6)

where ρ [kg/m³] is the fluid density and µ [kg/m s] is the dynamic viscosity. In the present study, the onset of nonlinear flow behavior by means of a critical Reynolds number Re_cr is defined as a significant nonlinear deviation in specific discharge of (> 5%) from Darcy’s law.

For the Ergun equation (Eq. 3), the following form of the friction factor f and Reynolds number Re_p is obtained:

$$f = i\frac{gd}{{q^{2} }}\left( {\frac{{n^{3} }}{1 - n}} \right)\;{\text{and}} \;{\text{Re}}_{p} = \frac{\rho dq}{{\left( {1 - n} \right)\mu }} = \frac{\text{Re}}{1 - n}$$

(7)

Hence, the relation between friction factor f and Reynolds number Re_p is given by:

$$f = \frac{A}{{{\text{Re}}_{p} }} + B$$

(8)

where A and B are the Ergun constants (which are equal to 150 and 1.75 for the Ergun equation).

3.4 Empirical Correlation for the Nonlinear Flow Forchheimer Coefficient

The Ergun equation (Eq. 3) is the extended form of the Kozeny–Carman equation for nonlinear flow behavior and can be used to estimate both the Forchheimer coefficient a and b based on the average grain size diameter and porosity. In the literature, a wide variety of empirical relationships are provided which correlates the Forchheimer b coefficient to the intrinsic permeability (k), porosity and tortuosity (see Table 3).

Table 3 Different correlation forms for the Forchheimer coefficient β, where β = bg

Carman (1937) already recognized that porous media exists of complex pore structures, and streamlines are far from being completely straight and parallel to each other. This effect can be described by a dimensionless parameter $\tau$ called the hydraulic tortuosity (Bear 1972):

$$\tau = \frac{L}{{L_{\text{e}} }}$$

(9)

where L [m] is the straight line length that connects the two ends of the tortuous tube of length L_e [m]. Tortuosity is difficult to measure and hence, this parameter is less practical than porous media properties such as grain size and porosity. A broad range of empirical tortuosity equations exists for simpler pore geometries, while the available information in the literature on tortuosity for granular material that consists of mixed particle beds is scarce (e.g., Mota et al. 2001; Ghanbarian et al. 2013). Mota et al. (2001) proposed an empirical tortuosity–porosity power law for binary mixtures of spherical particles:

$$\tau = n^{ - 0.4}$$

(10)

In the present study, Eq. (10) is used to investigate if the correlation forms for β of Liu et al. (1995) and Thauvin and Mohanty (1998) with the tortuosity model of Mota et al. (2001) work for the investigated non-uniformly graded granular material.

4 Results

In the present study, a nonlinear relationship between the hydraulic gradient i and the specific discharge q is found for the natural sand and gravel samples [M.(1–3) and W.(1–2)], as well as for the composite mixture samples (S.1–3) (Figs. 3, 4, 5, 6, 7). Forchheimer coefficients are determined by fitting the second-order polynomial function (Eq. 3) over the specific discharges and corresponding hydraulic gradients. The total amount of realizations for the samples ranges between 25 and 33. Excellent Forchheimer fits were obtained for all samples (R² > 0.99). The Forchheimer coefficients a and b are listed in Table 4 for the natural sand and gravel samples, as well as for the composite mixture samples. In the present study, the determined critical Reynolds numbers range between 2.17 and 12.5. This is in line with the values of 1–15 suggested in the literature (Bear 1972; Hassanizadeh and Gray 1987; Ma and Ruth 1993; Comiti et al. 2000).

Fig. 3

Nonlinear relationships between the hydraulic gradient and the specific discharge and the grain size distributions for M.1–3

Fig. 4

Nonlinear relationships between the hydraulic gradient and the specific discharge and the grain size distributions for W.1–2

Fig. 5

Nonlinear relationships between the hydraulic gradient and the specific discharge and the grain size distributions for S.1.1–5

Fig. 6

Nonlinear relationships between the hydraulic gradient and the specific discharge and the grain size distributions for S2.1–5

Fig. 7

Nonlinear relationships between the hydraulic gradient and the specific discharge and the grain size distributions for S.3.1–3

Table 4 The Forchheimer coefficients a and b for the natural (M.1–3 and W.1–2), as well as for the composite mixtures of sand and gravel (S.1–3)

The measured ranges of Re_min–Re_max for the packed column experiments on the natural sand and gravel samples are shown in Table 4. For the sand samples from the Münsterlaan (M.1–3), multiple data points are in the Darcy regime before flow is increased above the critical Reynolds number Re_cr and hence, the onset of nonlinear flow is observed (respectively, 7, 5 and 3 datapoints for samples M.1, M.2 and M.3). For the composite mixtures, multiple datapoints are measured in the Darcy regime for sample S.1(1–5) and S.2(3–5). The measured datapoints are solely in the nonlinear flow regime for the remaining composite mixtures and the samples from Wyenhütte W.(1–2).

4.1 Nonlinear Flow in Natural Sands and Gravels

The same nonlinear flow behavior characteristics are observed for the samples M.1 and M.2 with approximately similar gradations, average grain sizes and porosity values. Compared to M.(1–2), the samples M.3 and W.1 have a wider grain size distribution (with C_u values of 3.86, and 4.84) and lower porosity values. The flow resistance was reduced for these samples with a higher average grain size, which resulted in lower Forchheimer coefficients a of 1334.9 (M.3) and 1326.1 s/m (W.1) compared to the values obtained for samples M.1–2 (Table 4). However, the Forchheimer coefficients b for these samples are, respectively, a factor 1.5 (sample M.3) and 2.5 (sample W.1) higher than for sample M.2. For these samples, reduced porosity values (up to 0.25) due to wider grain size distributions resulted in significantly higher coefficients b, while only a small difference is observed between the coefficients a.

4.2 Nonlinear Flow in the Composite Mixtures of Sand and Gravel

The results show that gradation of a granular material has a large impact on nonlinear flow behavior (Figs. 5, 6, 7). For all tested composite mixture samples, a wider grain size distribution (i.e., increased C_u value) resulted in increased flow resistance with respect to the uniformly graded reference sand. This results in a poor correlation between the average grain size d₅₀ and the Forchheimer coefficients a and b (Fig. 8a–b). Using the d₁₀ as characteristic poor length, the correlation is significantly improved (Fig. 8c–d).

Fig. 8

a Relationship between the median grain size (d₅₀) and Forchheimer coefficient a, as well as for b Forchheimer coefficient b, for the investigated sands. The increase in flow resistance for the sands with a higher C_u value is also indicated. c Relationship between the grain size (d₁₀) and Forchheimer coefficient a, as well as for d Forchheimer coefficient b, for the investigated sands. Note that the shown relationships between the grain size and Forchheimer coefficients are calculated based on d [units in m]

The d₃₀ is decreased by a factor of 0.80 and 0.64 with respect to the reference sample S.1 for, respectively, S.1.1 and S.1.2 (Table 2). Consequently, this results, respectively, in increased coefficients a by a factor of 1.64 and 2.04, as well as increased coefficients b by a factor of 1.47 and 1.84 (Table 4). The porosity was only slightly reduced to a value of 0.34, compared to the value of 0.35 for the reference sample S.1. For S.1.5, the d₃₀ is decreased by a factor of 0.44 with respect to reference sample S.1. This resulted in an increase by a factor up to 5.71 for coefficient a and up to 6.24 for coefficient b. In this case, the porosity was significantly lower than the reference sample with a value of 0.25. Similar trends are also observed for the composite mixtures dataset S.2.(1–5). For the composite gravel S.3.3, the relative increase in flow resistance, and hence Forchheimer coefficients a and b, compared to its reference sample S.3 was the highest of all samples with increases by a factor of 36.1 and 14.6, respectively.

Equal amounts of finer material were used for the samples S.1.3 and S.1.4 (Table 2 and Fig. 5). Composite mixture S.1.4 shows a relative increase in coarser material (i.e., tailing in the grain size distribution) with respect to S.1.3, which resulted in a porosity of 0.29 for S.1.4 instead of the 0.31 for S.1.3. Consequently, for S.1.4 the coefficient a is increased by a factor of 1.43 and the coefficient b is increased by a factor of 1.36 compared to the coefficients of S.1.3. A similar trend is observed for the composite mixtures S.2.3 and S.2.4, where also an equal amount of fine material was used (at d₃₀ of 0.7 mm). In this case, the porosity was reduced to 0.27 for S.2.4, compared to a value of 0.31 for S.2.3, and coefficients a and b were increased by a factor of 1.41 and 1.45, respectively.

4.3 Ergun Constants for Granular Material

In the present study, we have used the Ergun form of the friction factor (Eqs. 7, 8) to analyze the dataset from the literature (see Table 1), as well as the dataset from the present study on non-uniformly graded material (Fig. 9). For the characteristic pore length, the grain size d₅₀ was used to calculate the friction factor f and Re_p. Van Lopik et al. (2017) obtained the following relationship for the friction factor on their dataset of 20 different uniformly graded filter sands and gravels (C_u < 2):

Fig. 9

a Correlation of the friction factor f and the calculated Re_p using grain size d₅₀. The dashed correlation line is based on the dataset of 20 uniformly graded filter sands (van Lopik et al. 2017), which includes the reference samples 1–3 of the present study. The solid correlation line is based on the dataset of the present study. b The friction factor based on nonlinear flow through packed beds of different kinds of granular material using the literature dataset listed in Table 1. The studies that provide the ranges of Reynolds numbers over which the experiments were conducted are plotted. The envelope in both figures shows the upper and lower limits of the literature dataset on uniformly graded material

$$f = \frac{233.5}{{{\text{Re}}_{p} }} + 2.88; \quad {\text{for}}\; d_{50} \;{\text{and}}\; C_{u} < 2$$

(11)

Similarly, a relationship can be obtained from the results from present study:

$$f = \frac{498.7}{{{\text{Re}}_{p} }} + 5.03; \quad {\text{for}} \;d_{50} \;{\text{and}} \;Cu > 2$$

(12)

For the present study, a poor relationship between friction factor and Reynolds number is obtained for the entire dataset of the present study while using d₅₀ as characteristic pore length (see Fig. 9a). A wide range of Ergun constants of A are obtained, ranging from 254 for sample M.1 up to 6477 for sample W.2 (Table 5). Ergun constants of B are in the range of 3.5 (S1.1 and S2.1) and 19.3 (W.1). Especially for the coarser granular material samples (d₅₀ > 1.5 mm) with wider grain size distributions (C_u > 5), the obtained Ergun constants are significantly higher (A > 750 and B > 6) compared to uniformly graded materials.

Table 5 The resulting Ergun constants A and B from the present study based on the grain sizes d₁₀, d₃₀ and d₅₀

Using the dataset from the literature (Table 1) and the average grain size d₅₀ of the tested porous medium, the resulting Ergun constants A and B are obtained from the relationship between Reynolds number and friction factor (Fig. 9b). The literature provides Reynolds number ranges over which the packed beds were tested on nonlinear flow behavior. However, the definition of the Reynolds number differs for the given studies in Table 1, and therefore we recalculated these values to equivalent Re_p values (based on the d₅₀ and n, see Eq. 3). The individual analysis of the different packed beds with granular material provides a nice overview of the range of possible Ergun constants for different kind of uniformly graded granular material (see Fig. 10).

Fig. 10

Overview of dataset from the literature (see Table 1) of the resulting Ergun constants obtained from the provided grain sizes, porosity values and Forchheimer coefficients a and b in the literature for a Ergun constant A and b Ergun constant B

A wide variety of A and B values is obtained for granular material at given average grain size d₅₀. Comparing the dataset of uniformly graded granular material, the Ergun constant A seems to increase with increasing grain size (see Venkataraman and Rao 1998; Moutsopoulos et al. 2009; Salahi et al. 2015; Sedghi-Asl et al. 2014; Li et al. 2017; van Lopik et al. 2017). For example, estimating nonlinear flow behavior through packed beds of coarse gravel (d₅₀ > 20 mm) requires Ergun constants A that are much higher than for other granular material (values up to 4000, see Fig. 10a). The Ergun constant B is in general higher than the value of 1.75 (average value B of 2.8 for Ranganadha Rao and Suresh 1970; Moutsopoulos et al. 2009; Salahi et al. 2015; Li et al. 2017; van Lopik et al. 2017) (Fig. 10b). This is in line with the values suggested by Kovács (1981) (B = 2.4) for sediments, as well as the values from Engelund (1953), who stated that for smooth particles B is 1.8, while for angular particles B is 3.8.

The Ergun constants for granular materials with wider grain size distributions are notably higher compared to the dataset on uniformly graded granular material (Table 5; Fig. 10). Hence, the present dataset of graded granular material is also analyzed using the Ergun form of the friction factor (Eqs. 7, 8) while using the d₃₀ and d₁₀ as characteristic pore length (Fig. 11). For uniformly graded filter sands and gravels (C_u < 2) using d₃₀ and d₁₀, the friction factors are (van Lopik et al. 2017):

Fig. 11

a Correlation of the friction factor f and the calculated Re_p with a characteristic pore length of d₅₀, b as well as for a characteristic pore length of d₅₀. The dashed correlation line is based on the dataset of 20 uniformly graded filter sands (van Lopik et al. 2017), which include the reference samples 1–3 of the present study. The solid correlation line is based on the dataset of the present study

$$f = \frac{183.8}{{{\text{Re}}_{p30} }} + 2.53; \quad {\text{for}} \;d_{30} \;{\text{and}}\; C_{u} < 2$$

(13)

$$f = \frac{120.3}{{{\text{Re}}_{p10} }} + 2.03; \quad {\text{for}} \;d_{10} \;{\text{and}} \;C_{u} < 2$$

(14)

While for the dataset of the present study with granular material of wider grain size distributions (C_u > 2), the friction factor is:

$$f = \frac{168.0}{{{\text{Re}}_{p30} }} + 2.86;\quad {\text{for}}\; d_{30} \; {\text{and}} \;Cu > 2$$

(15)

$$f = \frac{63.1}{{{\text{Re}}_{p10} }} + 1.72;\quad {\text{for}} \;d_{10} \;{\text{and}} \;Cu > 2$$

(16)

Considering the entire dataset of the present study, the overall correlation between friction factor and Re_p with d₁₀ is good (Fig. 11b). Hence, we suggest the use of 63.1 and 1.72 for, respectively, Ergun constants A and B, while using the d₁₀ as characteristic pore length for sand and gravel with wider grain size distributions.

4.4 Correlation Forms for Forchheimer Coefficient b

The differences between the experimentally obtained Forchheimer coefficients b and the derived coefficients from empirical relationships listed in Table 3 are shown in Fig. 12. Using the correlation form of Geertsma (1974) for Forchheimer β, the estimation of the experimentally obtained values is poor (Fig. 12a). The derived Forchheimer coefficients b are highly overestimated, and the scattering around the linear regression slope is high with a R²-coefficient of 0.70. Due to the relative importance of the porosity term (n^−5.5) in the Geertsma (1974) equation in comparison with the other equations in Table 3, the low porosity values (0.25 < n < 0.33) of the samples in this study results in this significant overestimation. For Janicek and Katz (1955), the fit is reasonable with a R²-coefficient of, 0.76 (Fig. 12b).

Fig. 12

Linear regression analysis of the Forchheimer coefficient b using the correlation form of a Geertsma (1974), b Janicek and Katz (1955), c Liu et al. (1995) and Thauvin and Mohanty (1998) and d the modified correlation forms Eqs. (17)–(18)

The scattering around the linear regression slope is significantly reduced for the empirical relationships that included the tortuosity model of Mota et al. (2001) for mixed particle beds (Fig. 12c). Combining this tortuosity model with the correlation forms of Liu et al. (1995) and Thauvin and Mohanty (1998), the R²-coefficients are, respectively, 0.86 and 0.88. In order to obtain a better correlation with the experimental data, the coefficients in these correlation forms are modified (Fig. 12d). The modified correlation form of Liu et al. (1995) is:

$$\beta = 3.1*10^{6} \frac{\tau }{kn} = 3.1*10^{6} k^{ - 1} n^{ - 1.4} \quad {\text{with}}\; R^{2} = 0.86$$

(17)

where k [mD] is the intrinsic permeability.

The modified correlation form of Thauvin and Mohanty (1998) is:

$$\beta = 2.0*10^{3} k^{ - 0.98} n^{ - 0.29} \tau^{3.35} = 2*10^{3} k^{ - 0.98} n^{ - 1.63} \quad {\text{with}}\; R^{2} = 0.88$$

(18)

where k [D] is the intrinsic permeability.

5 Discussion

Nonlinear flow behavior in porous media is controlled by the microscale processes at the fluid–solid interfaces of the pore structures and pore geometry (e.g., Dybbs and Edward 1984; Hassanizadeh and Gray 1987; Ma and Ruth 1993; Panfilov and Fourar 2006). Moreover, these microscale processes differ over the different Reynold number ranges. Commonly, post-Darcian nonlinear flow is separated in two flow regimes. These are the nonlinear laminar flow regime (also known as Forchheimer regime or transition flow regime) and the fully turbulent flow regime (Andersen and Burcharth 1995; Lage 1998; Houben 2015). This is done by flow visualization studies in porous media that determined the onset of turbulent flow (e.g., Dybbs and Edwards 1984; Seguin et al. 1998), as well as analysis of macroscale flow experiments on packed beds (e.g., Venkataraman and Rao 1998; Sedghi-Asl and Rahimi 2011; Bağci et al. 2014). This transition between the different flow regimes in the post-Darcian flow regime is hard to define for natural irregular and graded granular materials such as used in this study (Andersen and Burcharth 1995).

In practice, the traditional Ergun-type formulations are still widely used to provide estimates for the Forchheimer coefficients a and b (Eq. 3), since these are based on the easy measurable parameters such as a particle diameter for the characteristic pore length and the porosity. Previous studies have shown that besides the particle diameter and porosity, also the particle shape and packing arrangement has a significant impact on nonlinear flow behavior (e.g., Engelund 1953; Allen et al. 2013; Salahi et al. 2015). Experimental studies on nonlinear flow behavior on cubic arrays of smooth glass spheres show that the flow resistance is drastically lower compared to granular material (Huang et al. 2013). This study showed that the Ergun constants provided by Irmay (1964) with values A and B of, respectively, 180 and 0.6 suffice to predict the experimental data. The cubic arrangement of the smooth spheres allowed for the use of low B values of 0.6. However, for randomly packed beds of smooth particles, resulting in non-cubic arrangements, the resulting Ergun constants are much higher. Even for randomly packed beds of smooth spheres a wide variety of Ergun constants is obtained in the literature due to the varying pore structures (e.g., MacDonald et al. 1979; Fand et al. 1987; Bagci et al. 2014; Erdim et al. 2015). For angular shaped particles, which usually vary in shape for each granular particle, it is even harder to obtain one single specific characteristic pore length to account for the entire pore geometry of the packed bed of granular material. Moreover, during packing of asymmetric elongated particles, no isotropic packing arrangement will be obtained, since the particles tend to align during packing in such way that the larger areas are parallel to the bottom of the outlet (Allen et al. 2013). Packing the grains in such way that they align with the flow direction, significant reduction in flow resistance is observed. In the present study slightly elongated subangular-subrounded sands were used (“Appendix 2”), and anisotropy in packed beds of the mixed granular material will have an impact on the flow resistance. However, the differences in grain shape between the used sands in this study are very small (“Appendix 2”) and hence, the variations in the flow resistance are mainly due to differences in the pore geometry by different grain sizes and packing arrangement of the non-uniformly graded granular material.

The present study shows that nonlinear flow behavior through packed beds of granular material with wide grain size distributions (C_u > 2) deviates drastically from flow through uniformly graded material (Figs. 9, 10). The packing arrangements of such granular material differ from uniformly graded material, which is reflected by the decreased porosities of the tested granular material (Fig. 2). The heterogeneity of the packed bed with a wide grain size distribution results in increased complexity of the pore structure, smaller pore spaces and hence, higher tortuosity values. The use of the correlation forms for the nonlinear Forchheimer coefficient β of Liu et al. (1995) and Thauvin and Mohanty (1998) in combination with the tortuosity model of Mota et al. (2001) for mixed spherical particles beds resulted in a reasonable prediction of the experimental dataset. However, the actual values for tortuosity in both Darcian and non-Darcian flow regimes are difficult to measure for a given pore structure of the packed bed. Therefore, one heavily relies on the different empirical relationships for tortuosity provided in the literature (e.g., Ghanbarian et al. 2013).

6 Conclusions

Nonlinear flow behavior through packed beds of 5 different types of natural sand and gravel from unconsolidated aquifers, as well as 13 different composite mixtures of uniformly graded filter sands at different grain size distributions and porosity values, is investigated. This broad set provides data on the hydraulic characteristics for sand and gravel types with C_u values of 2.00–17.5 and porosity values of 0.23–0.36. The comparison of the dataset of the present study on sand and gravel with wide grain size distributions (C_u > 2) with datasets in the literature on uniformly graded sand and gravel (C_u < 2), shows that the gradation has a significant effect on nonlinear flow behavior. For graded gravels (d₅₀ = 6.3 mm, C_u = 12.3), Forchheimer coefficients a and b are increased by a factor of 36 and 15, respectively, with respect to uniformly graded gravel (d₅₀ = 6.3 mm, C_u = 1.27).

Considering the Ergun relationship for empirical correlation of the experimental dataset while using the average grain size d₅₀, a wide range of Ergun constants A and B is obtained with values of, respectively, 233.5 and 2.88 for uniformly graded filter sands up to values of 2192 and 8.36 for poorly gap-graded gravel (C_u = 12.33). In order to estimate the Forchheimer coefficients of granular material with wide grain size distributions (C_u > 2), the d₁₀ should be used as characteristic pore length with the Ergun constants of A and B of 63.1 and 1.72, respectively.

Moreover, for granular material with wider grain size distributions (C_u values of 3.86–17.35) and low porosities (n < 0.30), the Forchheimer coefficient b for a given coefficient a is significantly higher than for more uniformly graded granular material. Hence, the quadratic term in the Forchheimer that accounts for increase in nonlinear flow resistance by additional energy dissipation is higher for such granular material.

Appendices

Appendix 1: Composition of the Artificially Mixed Filter Sands

See Table 6.

Table 6 Mass percentage of the uniformly graded filter sands of the artificially mixed filter sand samples in this study

Appendix 2: Grain Shape Analysis

The grain shape properties of the filter sands that are used for the artificially mixed filter sands (see “Appendix 1”), as well as the natural sands M.1 and W.1 are determined with the CAMSIZER^® P4 (Retsch, Germany). The use of the Camsizer imaging techniques enables to detect the various dimensions of irregularly shaped particles by means of the sphericity (S^Kr) and roundness (R^Kr) based on the Krumbein scale (Krumbein and Sloss 1963), sphericity (S^QP), symmetry (symm), convexity (C_x), and elongation (e) (e.g., Altuhafi et al. 2013):

$$S^{QP} = \frac{{\sqrt {4\pi A_{\text{s}} } }}{P}$$

(B1)

where A_s [m²] is the surface area, and P [m] is the perimeter of the sphere.

$${\text{symm}} = \frac{1}{2}\left[ {1 + \hbox{min} \frac{{r_{1} }}{{r_{2} }}} \right]$$

(B2)

where r₁ und r₂ [m] are the distances from the center of the grain to the outer borders in the measuring direction. For asymmetric particles Symm is < 1.

$$C_{x} = \sqrt {\frac{{A_{\text{real}} }}{{A_{\text{conv}} }}}$$

(B3)

where A_real [m²] is the real area of the particle projection and A_conv [m²] is the convex area of particle projection.

$$e = \frac{L}{W}$$

(B4)

where L [m] is the maximum grain diameter and the W [m] is minimum grain diameter.

The grain shape analysis in Table 7 shows that the grain characteristics of the filter sands and natural sands are similar and the granular material can be classified as subangular-subrounded sand.

Table 7 Grain shape properties of the filter sands with a grain size smaller than 2.0 mm and the natural sand M.1 and W.1