Sagyn Omirbekov

Hossein Davarzani
email: h.davarzani@brgm.fr

Stéfan Colombano

Azita Ahmadi-Senichault

Experimental and numerical upscaling of foam flow in highly permeable porous media

Keywords: bulk foam, porous media, rheology, non-Newtonian fluid, yield stress, upscaling

Foam in porous media has been studied as a tool for various applications. Recently, the technology has become relevant for contaminated-aquifer remediation, where porous media are highly permeable.

Therefore, the behavior of foam flow in high permeability porous media still raises numerous questions.

In particular, upscaling of the foam flow from pore to Darcy scale is still under debate. Since the behavior of bulk foam has been studied principally in the food and cosmetics industries, and foam flow in porous media has mainly been investigated in the oil industry, the link between bulk-foam behavior and foam flow in porous media is still missing. The upscaling of foam flow from the pore scale to the laboratory scale could give valuable insight for understanding foam flow in aquifers.

We studied the behavior of pre-generated foam with different foam qualities through the rheological characterization of bulk foam using a rheometer and also when flowing in a porous medium composed of 1 mm glass beads. Foam was formed by co-injecting surfactant solution and nitrogen gas through a porous column filled by fine sand. The homogenization method is used to study macroscopic foam flow properties in porous media by solving the non-linear boundary value problem. The rheology of bulk foam is then used as an input in the upscaling procedure for foam flow in different periodic model 2D and 3D unit cells.

From our experiments, we found that the bulk foam is a yield-stress fluid and that the yield-stress values increase with foam quality. Moreover, the rheology of bulk foam corresponds well to the yield stress (Herschel-Bulkley-Papanastasiou) model. We found that foam behaves as a continuous yield-stress fluid in highly permeable porous media. It was also shown that the apparent foam viscosity in porous media increases with the foam quality at the same total flow rate. The results obtained from the rheometer successfully match the outcomes of apparent foam viscosity obtained by flow in porous media by a shifting parameter for the same foam quality. The apparent foam viscosity found in 1 mm glass-bead packing was much higher than bulk foam viscosity.

Experimental results were compared to numerical results on simple unit cells. Although we observed considerable differences between the experimental and numerical results of upscaling, the general trend was identical. The differences can be explained by the complexity of the foam flow in porous media, especially foam compressibility. We found that foam flow at low capillary numbers is influenced by the trapping effect and at high pressure gradients by the compressibility. Compressibility was estimated for foam flow in 1 mm glass-bead packing. When foam compressibility is insignificant, the upscaling model can predict foam-flow behavior well at the Darcy scale.

Introduction

Foam is a two-phase system where gas bubbles are dispersed in the continuous liquid phase. The liquid phase in foam is generally an aqueous surfactant solution that plays a crucial role in stabilizing the liquid films between bubbles by increasing the surface tension.

Foam is used in a wide variety of industrial applications due to its mechanical behavior combined with its sizeable specific area relative to the density (Prud'homme & [START_REF] Prud'homme | FOAMS. Theory, Measurements, and Applications[END_REF]. For instance, foam in the oil industry can be used as a drilling fluid, as a mobility-control agent in enhanced oil recovery processes [START_REF] Chen | Switchable nonionic to cationic ethoxylated amine surfactants for CO2 enhanced oil recovery in high-temperature, high-salinity carbonate reservoirs[END_REF], in matrix-acidization treatments [START_REF] Rossen | Modeling foams for acid diversion[END_REF], and in gas-leakage prevention [START_REF] Bernard | Model study of foam as a sealant for leaks in gas storage reservoirs[END_REF]. In soil remediation processes, foam can be used as a sweeping pollutant fluid [START_REF] Hirasaki | Surfactant/Foam Process for Aquifer Remediation[END_REF], as fluid that transports soil-cleaning additives or gases [START_REF] Maire | Targeted delivery of hydrogen for the bioremediation of aquifers contaminated by dissolved chlorinated compounds[END_REF][START_REF] Choi | Enhancement of aerobic biodegradation in an oxygen-limiting environment using a saponin-based microbubble suspension[END_REF]. Moreover, in ordinary life foam can be found in many food and cosmetic products.

The mechanical behavior of bulk foam is complicated due to elastic, plastic, and viscous properties.

Therefore, their complexity is analyzed using rheology [START_REF] Coussot | Rheometry of Pastes, Suspensions, and Granular Materials: Applications in Industry and Environment[END_REF][START_REF] Larson | The Structure and Rheology of Complex Fluids. s.l[END_REF][START_REF] Macosko | Rheology: Principles, Measurements, and Applications[END_REF].

Many authors have experimentally studied the rheology of bulk foam in the last sixty years [START_REF] Cohen-Addad | Flow in Foams and Flowing Foams[END_REF][START_REF] Denkov | The role of surfactant type and bubble surface mobility in foam rheology[END_REF][START_REF] Hohler | Rheology of liquid foam[END_REF][START_REF] Katgert | The jamming perspective on wet foams[END_REF][START_REF] Kraynik | Foam Flows[END_REF][START_REF] Dollet | Rheology of aqueous foams[END_REF][START_REF] Khan | Foam Rheology: III. Measurement of Shear Flow Properties[END_REF][START_REF] Herzhaft | Rheology of Aqueous Foams: a Literature Review of some Experimental Works[END_REF][START_REF] Denkov | Wall slip and viscous dissipation in sheared foams: Effect of surface mobility[END_REF]. In most of the studies, bulk foam was identified as a yield-stress fluid and the majority of the results are commonly well fitted by the Herschel-Bulkley law. For instance, [START_REF] Khan | Foam Rheology: III. Measurement of Shear Flow Properties[END_REF] investigated bulk foam through a rheometer with parallel-plate geometry. The geometry was covered with sandpaper to eliminate the wall slip velocity.

The average diameter of foam bubbles was 65 μm, as measured by an optical method. They demonstrated the existence of a yield stress for foams with 92%, 95% and 97% of gas fractions (𝑓 𝑔 , foam quality), with the yield stress increasing with 𝑓 𝑔 .

A full understanding of foam flow in porous media is hindered due to the complex behavior of foam associated with the complex microstructure of porous media. There are also apparent contradictions in the conclusions of the studies on foam flow presented in the literature. For instance, some authors showed Newtonian [START_REF] Persoff | A Laboratory Investigation of Foam Flow in Sandstone at Elevated Pressure[END_REF] foam flow behavior in porous media, others shear-thinning [START_REF] Hirasaki | Mechanisms of Foam Flow in Porous Media: Apparent Viscosity in Smooth Capillaries[END_REF][START_REF] Falls | The apparent viscosity of foams in homogeneous bead packs[END_REF], a mixture of Newtonian and shear-thinning [START_REF] Vassenden | Experimental Foundation for Relative Permeability Modeling of Foam[END_REF] or related to foam quality [START_REF] Alvarez | Unified Model for Steady-State Foam Behavior at High and Low Foam Qualities[END_REF]. [START_REF] Hirasaki | Mechanisms of Foam Flow in Porous Media: Apparent Viscosity in Smooth Capillaries[END_REF] reported non-linear behavior in foam when flowing through a smooth capillary tube and pointed out the impact of foam texture (bubble size) on apparent foam viscosity. Shear-thinning behavior of foam flow in homogenous bead packs was shown by [START_REF] Falls | The apparent viscosity of foams in homogeneous bead packs[END_REF]. [START_REF] Rossen | Modeling foams for acid diversion[END_REF] modeled foam as a Bingham plastic flowing in sandstones. They assumed bubble size was fixed and roughly equal to pore size at lowquality regimes. Some authors also considered the presence of yield stress based on a threshold pressure gradient, which was linked with the types of gas and surfactant, surfactant concentration, and petrophysical properties of porous media [START_REF] Rossen | Percolation theory of creation and mobilization of foams in porous media[END_REF][START_REF] Dicksen | Conditions for Foam Generation in Homogeneous Porous Media[END_REF][START_REF] Ransohoff | Mechanisms of Foam Generation in Glass-Bead Packs[END_REF][START_REF] Gauglitz | Foam generation in homogeneous porous media[END_REF][START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF]. However, in most of the foam modeling investigations in porous media, foam was characterized as a pure power-law fluid without assuming yield stress [START_REF] Friedmann | Experimental and simulation study of hightemperature foam displacement in porous media[END_REF][START_REF] Kovscek | Mechanistic foam flow simulation in heterogeneous and multidimensional porous media[END_REF][START_REF] Bertin | Development of a Bubble-Population Correlation for Foam-Flow Modeling in Porous Media[END_REF][START_REF] Myers | Transient Foam Displacement in the Presence of Residual Oil: Experiment and Simulation Using a Population-Balance Model[END_REF][START_REF] Kam | Improved mechanistic foam simulation with foam catastrophe theory[END_REF][START_REF] Kovscek | Modeling Foam Displacement With the Local-Equilibrium Approximation: Theory and Experimental Verification[END_REF][START_REF] Ashoori | Multiple Foam States and Long-Distance Foam Propagation in EOR Displacements[END_REF]. Recently, [START_REF] Aranda | Experimental Study of Foam Flow in Highly Permeable Porous Media for Soil Remediation[END_REF] experimentally found the shear-thinning behavior of foam when studying the influence of permeability on the apparent foam viscosity in highly permeable porous media [START_REF] Aranda | Experimental Study of Foam Flow in Highly Permeable Porous Media for Soil Remediation[END_REF].

Several foam-modeling techniques have been developed since the 1980s, to understand and predict the complex behavior of foam flow in porous media [START_REF] Ma | Modeling Techniques for Foam Flow in Porous Media[END_REF]. [START_REF] Ma | Modeling Techniques for Foam Flow in Porous Media[END_REF] provided a complete overview of existing methods for modeling foam in porous media. They categorized the modeling techniques into three groups depending on implicit or dynamic foam textures: local-equilibrium models, population-balance models, and others. They also indicated the possibility of using all models for direct incorporation into reservoir simulators based on material balance and Darcy law. They highlighted some studies on the upscaling of foam from laboratory-scale to field-scale in which the upscaling problems were pointed out due to entrance and end effects in laboratory core experiments.

On the other hand, a link between pore or bubble scale (discontinuous) physics and Darcy-scale (continuum) models is missing. This link is often ascertained through a variety of upscaling approaches [START_REF] Berryman | Comparison of upscaling methods in poroelasticity and its generalizations[END_REF]. For instance, Darcy's law is obtained from the creeping incompressible flow of a Newtonian fluid via the volume-averaging method with a no-slip boundary condition at the liquid-solid interfaces [START_REF] Quintard | Transport in ordered and disordered porous media: volume-averaged equations, closure problems, and comparison with experiment[END_REF][START_REF] Quintard | Transport in ordered and disordered porous media III: Closure and comparison between theory and experiment[END_REF]. In this method, intrinsic permeability is calculated as a function of closure variables for a given geometry (porous media) [START_REF] Whitaker | The Method of Volume Averaging[END_REF]. However, upscaling is still a challenge for non-Newtonian fluids, which are described by non-linear partial differential equations. There are very few experimental, numerical, or analytical studies on the upscaling of the non-Newtonian fluids flow in porous media. Several authors theoretically/numerically investigated upscaling power-law fluids [START_REF] Idris | Microstructural effects on the flow law of power-law fluids through fibrous media[END_REF][START_REF] Wang | Derivation of the Darcy-scale filtration equation for power-law fluids with the volume averaging method[END_REF][START_REF] Woods | Creeping flows of power-law fluids through periodic arrays of elliptical cylinders[END_REF][START_REF] Orgéas | Modelling the flow of power-law fluids through anisotropic porous media at lowpore Reynolds number[END_REF], viscoelastic fluids, yield-stress fluids, or considered the case of generalized non-Newtonian fluids [START_REF] Orgéas | Upscaling the flow of generalised Newtonian fluids through anisotropic porous media[END_REF]. Therefore, the upscaling of foam flow in porous media is still questionable [START_REF] Ma | Modeling Techniques for Foam Flow in Porous Media[END_REF], especially in soil remediation, where the permeability of porous media is high. Some authors also empirically proposed a Darcy-scale equation for non-Newtonian fluids [START_REF] Larry | Fundamentals of Enhanced Oil Recovery[END_REF][START_REF] Chauveteau | Rodlike polymer solution flow through fine pores: influence of pore size on rheological behavior[END_REF][START_REF] Chevalier | Darcy's law for yield stress fluid flowing through a porous medium[END_REF][START_REF] Rodríguez De Castro | Extending Darcy's law to the flow of yield stress fluids in packed beds: Method and experiments[END_REF]. Put simply, the extension of the Darcy scale equation for yield-stress fluid in anisotropic porous media is very complex and still open.

To our knowledge, the modeling of non-Newtonian fluids, especially yield-stress fluid flow in porous media, has received less attention due to strong coupling between the porous medium and the fluid's rheology. Most of the previous studies on foam in porous media were related to the oil industry, where bubble size was considered to be roughly equal to pore size. However, in soil remediation cases, the porous medium is highly permeable and consists of large pores, and there the bubbles can be much smaller than the characteristic length in porous media.

Consequently, the main objective of this study was to show whether foam flow in highly permeable porous media can be deduced from the bulk foam behavior and if we can consider foam as a single-phase yield-stress fluid. We also aimed to characterize the rheological behavior of bulk foam experimentally and to scale-up analytically and numerically the foam flow from pore to laboratory-Darcy scale.

Theoretical considerations

One of the important characteristics of foam is its foam quality, which is the ratio of the injected gas volume to the total liquid and gas volume, and can be expressed as

𝑓 𝑔 = 𝑄 𝐺 𝑄 𝐺 + 𝑄 𝐿 (1)
in which 𝑄 𝐺 (mL/min) and 𝑄 𝐿 (mL/min) are the volumetric gas and liquid flow rates, respectively. The sum of 𝑄 𝐺 and 𝑄 𝐿 provides the value of total flow rate 𝑄. According to the value of gas fraction, bulk foam can be dry (𝑓 𝑔 >99%), wet (64%<𝑓 𝑔 <99%), or considered as a bubbly liquid (𝑓 𝑔 <64%) [START_REF] Langevin | Aqueous foams and foam films stabilised by surfactants. Gravity-free studies[END_REF]. Moreover, foam bubble shape can be transformed from spherical to polyhedral by increasing the gas fraction [START_REF] Kovscek | Fundamentals of foam transport in porous media[END_REF]. In porous media, [START_REF] Osterloh | Effects of gas and liquid velocity on steady-state foam flow at high temperature[END_REF] distinguished two particular foam flow regimes depending on foam quality in sand-pack experiments at a steady state. They found a low-quality regime (wet), in which the pressure gradient was constant regardless of the liquid flow rate, and a high-quality regime (dry), in which the pressure gradient was independent of the gas flow rate. A transition foam quality ( 𝑓 𝑔 * ) divides these two regimes, and depends on porous media characteristics, types of surfactants, and gas [START_REF] Alvarez | Unified Model for Steady-State Foam Behavior at High and Low Foam Qualities[END_REF].

Newtonian-fluid flow in porous media can be considered as models of flow through a bundle of tortuous capillary tubes which have uniform radii, 𝑅. By the theory of Kozeny-Carman (Kozeny, 1927;[START_REF] Carman | Fluid flow through granular beds[END_REF], tortuosity in capillaries is defined as the ratio between the effective length (𝐿 𝑒𝑓𝑓 ), i.e. length of the tortuous capillaries and the lengths of equivalent straight capillary tubes (𝐿). Hence, it was assumed that fluid entering into tortuous capillaries flows faster than the fluid flowing through the equivalent straight capillaries in order to reach to the end of the tube at the same time. Using this concept, flow in porous media can be presented as follows:

𝑢 = ∅𝑅 2 ∆𝑃 8𝜇𝐿 ( 𝐿 𝐿 𝑒𝑓𝑓 ) 2 (2)
𝑢 (m/s) and ∆𝑃 (Pa) are Darcy velocity and pressure drop in the direction parallel to length 𝐿 (m), 𝜇 (Pa.s)

is the dynamic fluid viscosity, ∅ is the porous media porosity. Eq. 2 is then identified with the Darcy velocity for the flow of Newtonian fluids in porous media [START_REF] Darcy | Les Fontaines publiques de la ville de Dijon[END_REF].

𝑢 = 𝑄 𝐴 = 𝐾 ∆𝑃 𝜇 𝐿 (3)
where 𝑄 (m 3 /s) and 𝐴 (m 2 ) are the flow rate and cross-section surface area of the porous media. The mean pore radius could be calculated as the radius of straight capillaries composing a bundle giving the same flow:

𝑅 𝑒𝑞 2 = (𝑅 𝐿 𝐿 𝑒𝑓𝑓 ) 2 = 8𝐾 ∅ (4)
In the preceding equation, 𝑅 𝑒𝑞 (m) is the average pore throat radius, 𝐾 (m 2 ) is the intrinsic permeability of the porous media. The term ( 𝐿 𝐿 𝑒𝑓𝑓 ) 2 represents the tortuosity in porous media. Hence, an equivalent shear rate (𝛾̇𝑒 𝑞 ) for flow in porous media can be expressed as [START_REF] Darby | Chemical engineering fluid mechanics[END_REF],

𝛾̇𝑒 𝑞 = 4α𝑢 ∅ ⁄ 𝑅 𝑒𝑞 (5)
in which α is an empirical shift parameter associated with the bulk rheology of the fluid and the tortuosity of porous media [START_REF] Chauveteau | Basic rheological behavior of xanthan polysaccharide solutions in porous media: effects of pore size and polymer concentration[END_REF]. The tortuosity of packed spheres was assumed to be 25/12 from the derivation of the Blake-Kozeny equation using the capillary model [START_REF] Christopher | Power-law flow through a packed tube[END_REF][START_REF] Hirasaki | Analysis of factors influencing mobility and adsorption in the flow of polymer solution through porous media[END_REF]. As a consequence, [START_REF] Hirasaki | Analysis of factors influencing mobility and adsorption in the flow of polymer solution through porous media[END_REF] proposed to use α ≈ 0.69 through 𝛼 = 1/ √ 𝐿/𝐿 𝑒𝑓𝑓 for Newtonian fluids. Recently, Rodriguez de Castro (2019) experimentally found 𝛼 equal to 0.68 in glass-bead packings.

If we consider Darcy's law as applicable for the flow of a non-Newtonian fluid in porous media, the effect of non-linear rheology of the fluid can be incorporated in the definition of an apparent viscosity given by

𝜇 𝑎𝑝𝑝 = 𝐾 ∆𝑃 𝑢 𝐿 (6)
This apparent viscosity 𝜇 𝑎𝑝𝑝 is not an intrinsic property of fluid and varies with the flow rate. It is a nonlinear function of the equivalent shear rate given in Eq. ( 5), which can be varied relying on the flow rate.

Moreover, the shear rate varies considerably as a function of the structure and the heterogeneity of porous media that change interstitial velocity.

The rheological behavior of yield-stress fluid for simple shear can be described through the Herschel-Bulkley model [START_REF] Herschel | Konsistenzmessungen von gummi-benzollösungen[END_REF]:

𝜏 < 𝜏 0 ⇒ 𝛾̇= 0 (solid regime); 𝜏 > 𝜏 0 ⇒ 𝜏 = 𝜏 0 + 𝑎𝛾̇𝑛 (liquid regime) (7)
where 𝜏 (Pa) is the shear stress, 𝛾̇ (1/s) is the shear rate, 𝜏 0 (Pa) is the yield stress, 𝑎 (Pa.s n ), and 𝑛 (-) are the consistency and flow indexes, respectively. The continuous, viscous model proposed by

Papanastasiou [START_REF] Papanastasiou | Flows of Materials with Yield[END_REF] is used to prevent the main difficulty caused by the Herschel-Bulkley model in numerical studies due to discontinuous behavior at insufficient shear rates that tend the apparent viscosity to infinity. The Herschel-Bulkley-Papanastasiou (H-B-P) model can be written as follows:

𝜇(𝛾) = 𝑎(𝛾) 𝑛-1 + 𝜏 0 𝛾̇ [1 -𝑒𝑥𝑝(-𝑚𝛾)] (8)
in which 𝑚 (s) is the exponent index. [START_REF] Denkov | Wall slip and viscous dissipation in sheared foams: Effect of surface mobility[END_REF] studied foam rheology and wall slip velocity with 90% gas fraction using a rheometer.

The foam used was generated using a syringe with a needle with id (inner diameter) of 2.5 mm. They found the power-law index equal to 0.25 and 0.42 for tangentially mobile and immobile bubble surfaces, respectively. [START_REF] Denkov | The role of surfactant type and bubble surface mobility in foam rheology[END_REF] studied the effects of surfactant type and bubble surfaces on bulk foam rheology with 𝑓 𝑔 ≥80%. They classified the rheological behavior of foam into two different classes depending on values of power-law index 𝑛, by taking into account the viscous friction between bubbles and also between bubbles and solid walls, qualitatively. The results with 𝑛≈0.5 referred to the first system, which was described as the dominant friction in foam films, and the second type 𝑛<0.5 (mostly between 0.2 and 0.25) was defined as the essential energy dissipation on the bubble surfaces. These studies show solid, plastic or viscoelastic behavior of foam below yield stress and non-Newtonian liquid regime above the yield stress. The transition from solid-like to liquid-like mechanical behavior is called yielding. The yielding can be explained through the complex shear modulus that appears at the transition from a mostly elastic regime where 𝐺 ′ (Storage modulus) > 𝐺 ′′ (Loss modulus) towards predominantly viscous behavior 𝐺 ′ < 𝐺 ′′ [START_REF] Larson | The Structure and Rheology of Complex Fluids. s.l[END_REF]. The behavior of inverse transition from liquid-like to solid-like is called jamming, which was studied by slightly decreasing the shear rates and measuring the stress [START_REF] Liu | Jamming is not just cool any more[END_REF]. The apparent viscosity of foam is several times larger than the viscosity of the continuous liquid phase, even at low shear rates. This can be explained from a description of foam at the molecular and bubble scale [START_REF] Hohler | Rheology of liquid foam[END_REF].

Experimental approach

Here we present the fluids and materials used. The experimental setups and procedures were provided to characterize the rheology of foam in bulk form and in porous media.

Fluids and materials

Surfactant solution

C14-16 Alpha-olefin sulfonate (AOS) based surfactant Rhodacal® LSS-40/AX (Solvay Novecare) and demineralized/degassed water were used to prepare a surfactant solution. The surfactant contained 40 wt% of active materials in the aqueous solution. The chemical surfactant was chosen on the basis of several factors like soil biodegradability [START_REF] Tuvell | AOS ---An anionic surfactant system: Its manufacture, composition, properties, and potential application[END_REF][START_REF] Talmage | Environmental and human safety of major surfactants: alcohol ethoxylates and alkylphenol ethoxylates[END_REF], market accessibility [START_REF] Cserháti | Biological activity and environmental impact of anionic surfactants[END_REF], and field tests for soil remediation purposes [START_REF] Svab | Soil flushing by surfactant solution: Pilotscale demonstration of complete technology[END_REF]. The critical micelle concentration (CMC) of the surfactant is 1.8±0.1 g/L, which is measured by the pendant-drop method [START_REF] Stauffer | The measurement of surface tension by the pendant drop technique[END_REF]) using a drop shape analyzer (DSA-100S, KRUSS). The surfactant solution was prepared at a concentration of four times CMC to avoid adsorption [START_REF] Paria | Surfactant-enhanced remediation of organic contaminated soil and water[END_REF][START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF][START_REF] Aranda | Experimental Study of Foam Flow in Highly Permeable Porous Media for Soil Remediation[END_REF] at low concentrations and respect the biodegradability at high concentrations.

Gas

The gas used to generate foam in all experiments was nitrogen, provided by Air Liquide, with a purity of 99.98%. N2 was selected after considering several studies presented in the literature, including foam stability [START_REF] Farajzadeh | Comparative Study of CO2 and N2 Foams in Porous Media at Low and High Pressure-Temperatures[END_REF][START_REF] Zeng | Role of gas type on foam transport in porous media[END_REF] and solubility [START_REF] Sander | Compilation of Henry's law constants for inorganic and organic species of potential importance in environmental chemistry[END_REF] of gases. Also, carbon dioxide (purity > 99.7%), provided by Air Liquide, is used for pre-saturation procedures to avoid trapping of gas bubbles.

Porous media

Unconsolidated porous media made from calibrated silica sand (BR-37) and 1 mm glass beads are used, which are provided by Sibelco™ and Sigma-Aldrich companies, respectively. The sand-grading characteristics were as follows: the uniformity coefficient (Cu) is 0.72, the curvature coefficient (Cc) is 0.98; the effective size (d10) and mean grain size (d50) are 0.180 and 0.135 mm, respectively. Fine-sand packing was used to generate foam, while foam flow was examined in 1 mm glass-bead (GB) packing.

The characteristics of the porous media are presented in Table 1.

Experimental setups

The setups used to conduct the experiments are shown in Fig. 1. To generate foam, N2 and the surfactant solution were co-injected into the foam generator column. Then the rheology of pre-generated foam was studied with different methods: A) in bulk form using a rheometer (bulk); B) in the main porous media column (p.m.).

The mass flow controller El-Flow Select F-201CV (Bronkhorst), providing a range of 0.16-10 mLn/min (±0.5 % reading, plus ±0.1 % full scale), was used to control the gas flow rate and ensure the flow stability from a gas bottle. A piston pump Model 1HM (Eldex), with a precision of ±0.3 %, was used to inject the surfactant solution at a constant flow rate. The foam generator column was 10 cm long, with 4 cm inside diameter, and made of transparent polyvinyl chloride (PVC). The sand-pack was retained using a metallic grid with 42 µm cell size. The following setups were used to investigate the pre-generated foam properties.

A)

A rotational rheometer Haake Mars 60, Thermofisher (min. torque rotation 0.02 µNm, torque resolution 0.1 nNm) with serrated parallel plate geometry (P60/Ti/SE and TMP60/SE) was used to examine the rheology of pre-generated foam in bulk form.

Fig. 1 Schematic of the experimental setups used to characterize the rheology of A) bulk foam and B) foam in porous media

B)

The main column (made of PVC material) with a length of 40 cm and an inner diameter of 4 cm packed by 1 mm glass beads was used to study foam from the foam generator column. Metallic grids with 150 µm cell size were used on each side of the column to contain the glass-bead packing. Pressure drop along the main column was measured through a differential pressure transmitter Rosemount 2051 (Emerson) with a 0-2500 mbar range (±7 mbar at the maximum value). An electronic balance model STX 6201 (OHAUS) with a precision of 0.1 g was used to measure the effluent weight. The maximum pressure limit of the experimental installation was up to 6 bars controlled by the pump pressure sensor.

Experimental procedures

First, the porous columns were carefully packed and checked for leaks. The columns were flushed by CO2 gas to remove any air. Then, degassed/demineralized water was injected into the columns, in a vertical position, with a 0.5 mL/min for at least five pore volumes (PV) to dissolve any CO2 and thus to saturate the column thoroughly. As CO2 is highly soluble in the water, this step was used to provide 100% of porous media saturation. To measure pore volume and porosity, the columns were weighed before and after the water saturation. The subsequent experimental steps were carried out with the column in a horizontal position. The permeability of the porous media was measured through the injection of degassed water with different flow rates while measuring the pressure difference and calculated via Darcy's law (Eq. 3). The measured properties of porous media were tabulated in Table 1. Next, the column was flushed with 3 PV of surfactant solution to complete adsorption. Finally, to generate foam, the surfactant solution and the nitrogen gas were simultaneously injected into the foam generator. 5 PV of fluids was co-injected to obtain a stable pre-generated foam. The pre-generated foam was further investigated A) using the rheometer and B) in glass-bead packing. Pre-generated foams with 60%, 85%, and 95% foam qualities, generated at the same total flow rates ( 𝑄 =2 mL/min, 𝑢 =2.65×10 -05 m/s), were investigated using the rheometer. The mechanical properties of foam were probed through the rheometer, considering previous studies [START_REF] Larson | The Structure and Rheology of Complex Fluids. s.l[END_REF][START_REF] Macosko | Rheology: Principles, Measurements, and Applications[END_REF]. The parallel-plate geometry with serrated surfaces was used to avoid wall-slip velocity [START_REF] Herzhaft | Rheology of Aqueous Foams: a Literature Review of some Experimental Works[END_REF][START_REF] Marze | Aqueous foam slip and shear regimes determined by rheometry and multiple light scattering[END_REF][START_REF] Princen | Rheology of foams and highly concentrated emulsions: IV. An experimental study of the shear viscosity and yield stress of concentrated emulsions[END_REF]. The shear gap of the geometry was 1 mm, which was chosen by taking into account previous studies [START_REF] Herzhaft | Rheology of Aqueous Foams: a Literature Review of some Experimental Works[END_REF]. Moreover, before starting the experiments, the dynamic stability of the foam at each quality was verified, based on the study of [START_REF] Khan | Foam Rheology: III. Measurement of Shear Flow Properties[END_REF], by examining samples at a constant shear rate. The idea was to observe how long the foam could withstand at constant torque. The decrease in torque values indicates the change in the foam structure. The volume of 3.1 mL of the foam sample was set in the gap, and the experiments were run in triplicate for each foam quality.

B)

Firstly, the pre-generated foam was connected to the main column to study the impact of 𝑓 𝑔 and to determine the transition foam quality (𝑓 𝑔 * ). Foam with a fixed total flow rate (2 mL/min) was studied for various foam qualities (40%-99%). The fractions of gas and liquid injected were controlled by the mass flow controller and the pump. The foam flow experiments were examined through the flow rate and pressure drop measurements along the column. Consequently, the apparent foam viscosity was derived from Darcy's law (Eq. 3). After determining 𝑓 𝑔 * , foam flow in 1 mm GB packing was investigated at lowquality regimes. Foam flow with 𝑓 𝑔 =60%, 85%, and 95% was studied for different total flow rates that were increased step by step from minimum ( 𝑄 =0.2 mL/min, 𝑢 =2.65×10 -06 m/s) to maximum (𝑄=2.2 mL/min, 𝑢=2.92×10 -05 m/s) technically possible values. Thus, 905 and 82 minutes were required to inject 1 PV of pre-generated foam through the main column at the minimum and maximum flow rates, respectively. Experiments were duplicated by at least a descending flow rate experience. At the same time, the electronic balance measured the liquid effluent to determine the change of surfactant solution saturation (Sw) inside the main column. All experiments were carried out at 20°C, controlled room temperature, and atmospheric pressure.

Theoretical approach: upscaling technique

Since the direct upscaling of foam flow from pore to laboratory scale is complex, we used upscaling technique to scale up foam flow by using the rheological data of bulk foam as an input. The objective was to examine whether foam-flow behavior in porous media can be predicted from its bulk behavior. For field application in aquifers, a second upscaling from laboratory to aquifer scale should be considered.

That was beyond the scope of this study.

In order to numerically study the macroscopic behavior of foam flow, it was first necessary to establish the geometry of representative elementary volume (REV) to correspond to the porous media permeability obtained from the laboratory experiments. To do so, we solved the closure problem for creeping flow of steady-state Newtonian fluids developed through the volume-averaging method [START_REF] Quintard | Transport in ordered and disordered porous media III: Closure and comparison between theory and experiment[END_REF][START_REF] Whitaker | The Method of Volume Averaging[END_REF] to define the permeability of the two-dimensional (2D) and three-dimensional (3D) unit cells (see Fig. 3). The size of solid inclusion in each kind of unit cell was adjusted to obtain the same permeability as in the 1D column experiment. Then the boundary value problem developed using the homogenization method [START_REF] Auriault | Heterogeneous medium. Is an equivalent macroscopic description possible?[END_REF] was solved to study macroscopic foam-flow behavior further.

The targets of the upscaling were to study the impact of 𝑓 𝑔 and the effect of different unit cell geometries on foam flow. Then the outcome of the 1D porous column experiments was used to correlate and validate the findings from upscaling.

Description of foam flow at the pore scale

The concept of porous media is presented in Fig. 2, in which we consider the foam flow through a periodic REV with 𝛺 being the total volume of REV. where 𝐯 is the fluid velocity field and 𝑃 is the pressure. 𝜇 is the dynamic fluid viscosity, which can be defined by the H-B-P model (Eq. 8) for our case. Hence, the set of equations is the pore-scale description of the foam flow problem in highly permeable porous media. The edge size of the periodic REV was selected to be 1 mm for all cases in order to be close to 1 mm GB size, and the permeability was set equal to 8.3 × 10 -10 m 2 . The permeability of each porous medium has been calculated by solving the closure problem for Newtonian fluid flow through the volume-averaging method [START_REF] Quintard | Transport in ordered and disordered porous media III: Closure and comparison between theory and experiment[END_REF][START_REF] Whitaker | The Method of Volume Averaging[END_REF]. It should be noted that ideally, the unit cells chosen must have the same porosity and permeability as the experimental sample. Moreover, they should feature similar microstructures in terms of disorder. Although all unit cells considered here are ordered geometries, the porosity values for cells SPC, SPS, SCP and BCP are relatively close to each other and also close to the value found for the glass-bead packing. The porosity for CPC and CPS are significantly higher for the same value of the permeability.

Characteristics of the microstructure

Upscaling

Separation of scales

In order to meet the requirements of the homogenization method [START_REF] Idris | Microstructural effects on the flow law of power-law fluids through fibrous media[END_REF][START_REF] Orgéas | Upscaling the flow of generalised Newtonian fluids through anisotropic porous media[END_REF][START_REF] Auriault | Heterogeneous medium. Is an equivalent macroscopic description possible?[END_REF][START_REF] Papanicolau | Asymptotic Analysis for Periodic Structures[END_REF][START_REF] Sanchez-Palencia | Non-Homogeneous Media and Vibration Theory Lecture Notes in Physics. s.l[END_REF], the condition of separation of the scales must be satisfied for the porous media considered. As Fig. 2 describes, the characteristic REV length 𝑙 is supposed to be small compared to the characteristic macroscopic length 𝐿, which can be expressed as

𝜖 = 𝑙 𝐿 ≪ 1 (12)
where 𝜖 represents the scale-separation parameter.

Boundary value problem

In this work, the homogenization procedure proposed by Orgeas (2007) has been used. Thus, a generalized Newtonian fluid through 2D and 3D elementary unit cells is studied using numerical simulations. To study foam flow in porous media, a boundary value problem concerning the first order velocity 𝐯 (0) also, to the second-order pressure ∇P (1) is presented in the following way:

𝛁 • (-𝜖𝑃 (1) + ( 𝜏 0 |𝛾| [1 -exp(-𝑚|𝛾|)] + 𝑎|𝛾| 𝑛-1 ) (𝛁𝐯 (0) + (𝛁𝐯 (0) ) T )) = 𝛁𝑃 x (0)
in 𝛺 𝛽 (13) 𝛁 • 𝐯 (0) = 0 in 𝛺 𝛽 ( 14)

𝐯 (0) = 0 on 𝐴 𝛽𝜎 (15)
Here, the unknowns 𝜖𝑃 (1) and 𝐯 (0) are periodic, and the macroscopic pressure gradient 𝛁𝑃 x (0) is given on the whole REV as a source term. Foam was considered as a yield-stress fluid following the rheological model H-B-P (Eq. 8), where the fluid viscosity is supposed to be a function of the magnitude of the shear rate |𝛾| given by

|𝛾| = √2𝐒: 𝐒 (16)
where the strain rate 𝐒 is defined as a function of the velocity field 𝐯 at the pore scale,

𝐒 = 1 2 (𝛁𝐯 (0) + (𝛁𝐯 (0) ) T ) (17)
The boundary value problem given above was solved in the 2D periodic unit-cell of circular cross-section (Fig. 3, a) for 𝑓 𝑔 =60%, 85% and 95%, where the rheological parameters of the continuous, viscous model were obtained from the experimental study. Foam with 𝑓 𝑔 =85% was investigated in all geometries presented in Fig. 3. For all simulations in the ordered isotropic porous media, the macroscopic pressure gradient imposed 𝛁𝑃 x (0) was in the direction of the x-axis. The non-linear boundary value problem (Eqs.

13-15) was solved using Comsol Multiphysics Finite Element with a mixed pressure-velocity (P1-P2)

formulation. The flow problems were solved through the Creeping flow module, which is designed for solving Stokes flow problems. The dynamic viscosity of the module was user-defined via the Equation View section. Solving this problem was similar to solving the direct flow problem in simple unit cells.

The finite-triangular elements were used to mesh the 2D geometries. An example of mesh and boundary conditions used are in Fig. 4.

Fig. 4 The mesh and boundary conditions used to run the simulations

Results presented in the literature for the same geometries [START_REF] Idris | Microstructural effects on the flow law of power-law fluids through fibrous media[END_REF][START_REF] Orgéas | Upscaling the flow of generalised Newtonian fluids through anisotropic porous media[END_REF] were used to verify the numerical procedure adopted here in the case of power-law and generalized Newtonian fluids.

Results and discussion

Rheology of bulk foam

Fig. 5 shows the rheology of pre-generated foam for different foam qualities, in which the shear stress 𝜏 is displayed as a function of the shear rate 𝛾. The dots represent the average value of experimental results, and the lines describe the results of the H-B-P model. Table 2 shows the fitting parameters with coefficients of determination. The exponent index 𝑚 was chosen to be 10000 in order to better represent the yield stress behavior in numerical studies. From these results, a non-Newtonian fluid behavior with yield stress is observed that was well fitted by the Herschel-Bulkley-Papanastasiou model. Moreover, we saw that yield stress increases with increasing foam quality thereby raising the bulk foam viscosity and the increasing yield-stress value with the foam quality, which agrees with the literature results [START_REF] Khan | Foam Rheology: III. Measurement of Shear Flow Properties[END_REF], where the storage modulus (𝐺') and 𝜏 0 as a function of foam quality [START_REF] Marze | Oscillatory rheology of aqueous foams: surfactant, liquid fraction, experimental protocol and aging effects[END_REF][START_REF] Saint-Jalmes | Vanishing elasticity for wet foams: Equivalence with emulsions and role of polydispersity[END_REF] were measured. The value of the yield stress increases by a factor 2 from 𝑓 𝑔 =85% to 𝑓 𝑔 =95%, which can be explained by the evolution of the bubble shape from spherical to polyhedralform, thereby increasing the friction between bubbles [START_REF] Denkov | The role of surfactant type and bubble surface mobility in foam rheology[END_REF]. It can therefore be concluded that in our experiments, bulk foam can be clearly considered as a yieldstress fluid. This is in agreement with previous studies in which the value of the flow index, 𝑛, obtained experimentally was indicated close to 0.5 [START_REF] Denkov | The role of surfactant type and bubble surface mobility in foam rheology[END_REF][START_REF] Marze | Aqueous foam slip and shear regimes determined by rheometry and multiple light scattering[END_REF][START_REF] Princen | Rheology of foams and highly concentrated emulsions: IV. An experimental study of the shear viscosity and yield stress of concentrated emulsions[END_REF][START_REF] Ovarlez | Wide-gap Couette flows of dense emulsions: Local concentration measurements, and comparison between macroscopic and local constitutive law measurements through magnetic resonance imaging[END_REF][START_REF] Ovarlez | Investigation of shear banding in three-dimensional foams[END_REF][START_REF] Tcholakova | Theoretical model of viscous friction inside steadily sheared foams and concentrated emulsions[END_REF].

It should be noted that the measurement of the rheology for the foam with the quality of 𝑓 𝑔 =60% was found to be complicated due to the drainage of foam liquid-phase from the rheometer geometry. Therefore, to avoid the liquid draining from the plates, we loaded the exact necessary amount of the foam sample (3.1 mL) by setting the gap to 3 mm, where the foam sample was placed between two plates. Then, the gap was decreased slowly until 1 mm to perform tests at different shear rates. The foam with the quality of 𝑓 𝑔 =60% was studied with great care at each shear rate by performing additional experiments to ensure accurate results.

These results are analyzed with the outcomes of the experimental study on porous media and further used as inputs for upscaling of foam flow in highly permeable porous media by considering foam as a yieldstress fluid.

Rheology of foam in porous media (Darcy scale)

Since foam was generated in the sand-pack, foam bubble size was expected to be smaller than the pore size of the glass-bead packing. Therefore, the bubble sizes were estimated according to previous studies.

For instance, [START_REF] Marsden | Use of foam in petroleum operations[END_REF] noted that some foams generated in porous media could flow with bubbles much smaller than pores, as well as a series of foam films of various sizes, shapes, and configurations depending on the medium. They studied a pre-generated foam through the series of four identical porous media consisting of sand with a porosity of 35% and a permeability of 5.2 d. The foam texture was studied through optical cells and using a microscope. They showed an increase in bubble size with a decrease in foam mobility (the ratio of effective permeability to apparent viscosity). However, in later work, the opposite was revealed [START_REF] Hirasaki | Mechanisms of Foam Flow in Porous Media: Apparent Viscosity in Smooth Capillaries[END_REF][START_REF] Falls | The apparent viscosity of foams in homogeneous bead packs[END_REF]. From the results of [START_REF] Marsden | Use of foam in petroleum operations[END_REF], we observed that the bubble size was smaller at high flow rates compared to the pore size, which was estimated to be five times smaller than the grain diameter.

Friedmann and Jensen (1986) found that the size of foam bubbles decreased with increasing velocities at constant foam quality by observing the texture of the foam exiting from several different porous media.

Ettinger and Radke (1992) studied the foam texture of pre-generated foam in Berea sandstone. They also observed a decrease in bubble size with an increase in flow rate.

However, according to studies by [START_REF] Falls | The apparent viscosity of foams in homogeneous bead packs[END_REF], the average size of the foam bubbles remained unchanged for the foam, which was pre-generated in 3 mm glass-bead packing and introduced into the second pack with the same grain size. The average size of the bubbles was estimated to be the same as the equivalent capillary radius of the packs. Consequently, in this study, we roughly assumed that the size of foam bubbles is equivalent to the pore size of the sand-pack (see Table 1).

The apparent foam viscosity calculated through Darcy's law (Eq. 6) applied to the experimental study of foam flow in the porous sample (1 mm glass-bead packing) as a function of foam quality is presented in Fig. 6. From these results, the transition foam quality (𝑓 𝑔 * ) delimiting low and high-quality regimes is determined.

Fig. 6 Apparent foam viscosity in porous media (1 mm GB) as a function of foam quality at 𝑄=2 mL/min

We observed that the transition foam quality was quite high and equal to 97% (dashed line). Hence, this result is in good agreement with the results of some authors [START_REF] Alvarez | Unified Model for Steady-State Foam Behavior at High and Low Foam Qualities[END_REF][START_REF] Briceno | Self-lubricated transport of aqueous foams in horizontal conduits[END_REF][START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF]. For instance, [START_REF] Alvarez | Unified Model for Steady-State Foam Behavior at High and Low Foam Qualities[END_REF] found a value of 𝑓 𝑔 * =97% for an experiment study of foam flow in a sand-pack with a permeability of 3.1 Darcy, where the foam was pregenerated through bronze wool [START_REF] Alvarez | Unified Model for Steady-State Foam Behavior at High and Low Foam Qualities[END_REF]. They used two types of surfactants, an octyl ethoxylated alcohol and AOS, to generate foams. However, the bubble size of the pre-generated foam was not mentioned in their work. [START_REF] Briceno | Self-lubricated transport of aqueous foams in horizontal conduits[END_REF] investigated pre-generated foam flow in a 1.2 m long channel (1 in. tall, 1/4 in. wide) and in a pipe (5/8 in. inner diameter). The foam was formed using an aqueous solution based on a mixture of AOS, polyacrylic acid polymer, cosurfactants (dodecanol), solvent (butanol), and water.

They found 𝑓 𝑔 * =97% through the pressure gradient vs. foam velocity plot where they discovered a shift in the results at 𝑓 𝑔 >97%. It was also visually shown that there was a change in the structure of foam at this transition zone, in which foam flow behavior changed from uniform to slug flow.

Recently, [START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF] also found that transition-foam quality is independent of the porousmedia grain size when the average foam-bubble size is smaller than the pore size of the porous medium.

The value of transition-foam quality, for the same foam composition and pre-generated foam used in this study, was 97% for 2 mm, 4 mm, and 8 mm glass-bead packing. Nevertheless, the foam's bubble size can be changed depending on time. This is discussed in Appendix A. With these results in mind, all our further studies were chosen to be conducted in low-quality regime (𝑓 𝑔 =85%) to avoid the instability around 𝑓 𝑔 * .

Fig. 7 shows the experimental results in terms of the pressure gradient along the porous column packed by 1 mm GB as a function of the foam flow velocity, for 𝑓 𝑔 =60%, 85%, and 95%. Non-linear behavior was observed where the pressure gradient of flow increased with foam quality. We were unable to obtain ∇𝑃 values above foam injection flowrate of 3-4 mL/min due to the maximum pressure limit of the experimental configuration for 𝑓 𝑔 =60% and 85% of foams.

Fig. 7 The pressure gradient along the column as a function of the superficial velocity in GB, for

𝑓 𝑔 =60%, 85%, 95%
To investigate the foam rheology, the apparent viscosity variations versus the equivalent shear rate (Eq. 5) were analyzed. The apparent foam viscosity as a function of the equivalent shear rate (Eq. 6) is plotted in Fig. 8, in which points represent experimental data, and the lines correspond to the fitting curves by the H-B-P model. Results show that foams in porous media behave as non-Newtonian yield-stress fluids, and the apparent viscosity of foam increases with foam quality. Therefore, we found that foam flow in porous media behaved as a bulk foam if the bubbles were smaller than the pore size and for a low-quality foam-flow regime. In addition, we observed that the yield-stress values for foams in porous media were about twice as high as those found for bulk foams (see Table 3). Moreover, the fluidity index 𝑛 increased when the gas fractions decreased and had a value of 0.65 for 𝑓 𝑔 =60%. In the literature, 𝑛>0.5 is found for low-quality bulk foams ( 𝑓 𝑔 <75%) (Otsubo & Prud'homme, 1994;Otsubo & Prud'homme, 1994;[START_REF] Larmignat | Rheology of colloidal gas aphrons (microfoams)[END_REF], so it agrees with our results. However, note that the power-law index 𝑛 was around 0.5 for all bulk foams. This may be a consequence of the surface roughness since our porous medium consisted of glass-bead packing, and the pore walls are smooth compared to the serrated geometry of the rheometer. Also, it may be a result of the enduring sample preparation time for the foam with the quality of 𝑓 𝑔 =60% mentioned in section 5.1. Consequently, the foams can be classified into two groups depending on foam quality and values of the power-law index, where we observed 𝑛≈0.5 for 𝑓 𝑔 =85% and 95% and 𝑛>0.5 for 𝑓 𝑔 =60% foams.

Fig. 8 Apparent viscosity as a function of the equivalent shear rate: 𝑓 𝑔 =60%, 85% and 95% In Fig. 9a, the apparent viscosity of the bulk foam defined as 𝜏/𝛾̇ is plotted as a function of shear rate for various foam qualities on a log-log scale. The apparent viscosity obtained from experiments in porous media is also presented in the same figure. The range of the shear rate was restricted in order to facilitate the comparison of the results. Fig. 9b shows the fitted results with the shift parameter 𝛼 for each foam quality (Eq. 5). In this figure, the colored dots represent the outcomes of foam flow in porous media, and star dots are the result of the bulk foam examined through the rheometer. The fitted shift parameters 𝛼 with the coefficient of determination are tabulated in Table 4 for each foam quality. Here, the shift parameter is found to allow a good superposition of the apparent viscosities of bulk and in-situ foam.

Moreover, 𝛼 becomes more significant with increasing foam quality. In previous studies related to yieldstress fluids, it was found that the shift parameter 𝛼>1 for sandstones [START_REF] Fletcher | Measurements of polysaccharide polymer properties in porous media[END_REF] and bead packings [START_REF] Rodríguez De Castro | Extending Darcy's law to the flow of yield stress fluids in packed beds: Method and experiments[END_REF][START_REF] Zitha | Permeability ̃ Dependent Propagation of Polyacrylamides Under Near-Wellbore Flow Conditions[END_REF]. Also, 𝛼 was mostly assumed to be constant and independent of flow rate. Recently, Rodriguez de Castro (2019) showed empirically, in glass-bead packings, the dependence of 𝛼 on flow rates where 𝛼>1 at slow flow rates and 𝛼<1 for high flow rates with a constant value at very high velocities (Rodríguez de Castro, 2019). However, note that all these studies were related to the flow of polymer solutions in porous media, and there is a lack of data for foams.

In this study, it was found that the values of 𝛼 for foams are even smaller than those for Newtonian fluids in glass-bead packings (Rodríguez de Castro, 2019). Also, we observed that apparent in situ foam viscosity is higher than bulk foam viscosity at the same shear rate. Since bulk foam was examined at ambient conditions, this phenomenon may occur due to the effect of foam compressibility.

Upscaling of foam flow

Influence of foam quality

The outcome of the numerical study on foam flow in the periodic 2D SPC unit cell is compared to the experimental results in terms of the velocity/pressure gradient relationship in Fig. 10, for different foam qualities. In order to compare model results and experimental data, Fig. 10 was plotted on a semilogarithmic scale. The dashed lines represent numerical results, in which the size of the solid cylinder was set in order to fit the experimental data in terms of permeability (𝐾=8.3×10 -10 m 2 ). The colored dots correspond to experimental results in 1 mm GB packing. Thus, the flow of foam with quality values of 60, 85, and 95% were examined to study the influence of foam quality on foam flow behavior. As experimental results, it is observed that the foam flow viscosity in porous media increases with 𝑓 𝑔 .

However, results from the numerical study show higher foam mobility even at lower pressure gradients compared to experimental results. The global trend of the curves and the behavior in terms of 𝑓 𝑔 follow the same tendency as those of the experimental results.

Fig. 10 Evaluation of the superficial velocity as a function of the macroscopic pressure gradient in 2D

periodic SPC geometry, for different foam qualities (𝑓 𝑔 =60, 85 and 95%)

Effect of porous media geometry

Fig. 11 shows the superficial velocity of foam as a function of the pressure gradient in different 2D and 3D geometries for 𝑓 𝑔 =85%. The dashed lines correspond to the numerical results where the unit-cell geometry was adjusted to match the 𝐾=8.3×10 -10 m 2 value in Newtonian case, and the colored dots represent the corresponding experimental results of GB packing for 𝑓 𝑔 =85%

As Fig. 11 shows, the flow is non-linear and depends on the geometry of porous media, which plays a crucial role in the flow of yield stress fluids. Moreover, flow in 2D-SPC and 2D-SPS geometries are close to each other but differ from the outcomes of other unit cells. In addition, the results for more complex 2D unit cells CPC and CPS are very close to those of 3D unit cells SCP and BCP. Since the upscaling results were plotted on the semi-logarithmic scale, the small difference between the results was not explicitly shown. The porosity of 2D-SPC and 2D-CPS cells respectively equal 0.49 and 0.63, which is higher than the porosity of other unit cells where the values are close to 0.4. Therefore, we conclude that although the foam's behavior depends on the porous media structure. The results for the complex 2D and 3D geometries are closer to the experimental data than the simple 2D geometries. The predicted superficial foam velocities are however still a few order of magnitudes higher than the experimental data. Fig. 11 The macroscopic pressure gradient vs. superficial velocity of foam flow in 2D SPC, SPS, CPC,

CPS and 3D SCP, BCP unit cells (𝑓 𝑔 =85%)

This discrepancy can be attributed to the effects of compressibility and gas trapping, which should be considered during the study of foam flow in porous media.

To estimate the trapping effect in foam injection experiments, we calculated the capillary number 𝑁 𝑐𝑎 using the following equation [START_REF] Satter | Reservoir engineering: the fundamentals, simulation, and management of conventional and unconventional recoveries[END_REF],

𝑁 𝑐𝑎 = 𝐾∇𝑃 𝛾 (18)
where 𝛾=0.036 (N/m) is gas-liquid interfacial tension. Consequently, the capillary numbers as a function of pressure gradient can be written as 𝑁 𝑐𝑎 = 2.31 × 10 -8 ∇𝑃. The minimum value of 𝑁 𝑐𝑎 is then 2.55×10 - 3 which corresponds to our minimum measured pressure gradient (110782 Pa/m) where foam quality was 60%. This capillary number is even greater than the critical capillary number of 10 -5 , above which the capillary forces become negligible according to several studies [START_REF] Chatzis | Correlation of capillary number relationships for sandstone[END_REF][START_REF] Ding | Capillary number correlations for gas-liquid systems[END_REF]. Therefore, gas trapping can be negligible in our experiments.

To better observe the effects of gas trapping and compressibility of foam with 𝑓 𝑔 =85% in the porous medium, we closed the outlet of the generator column and measured the pressure drop along the main column at the end of the experiment. From the result presented in Fig. 12, it is shown that the pressure in the column reduced gradually even after the foam supply was stopped. During the first two hours, the pressure gradient dropped sharply from 140,000 Pa/m to 37,000 Pa/m, after gradually decreasing, and the experiment lasted four days. After closing the foam-generator outlet, foam continues flowing through the column. This can be explained by the foam compressibility effect. However, the pressure gradient falls first very fast and then more slowly because of the gas trapping. From the values of the capillary number shown in Fig. 12, we can say that the gas trapping occurs below a critical capillary-pressure number, 𝑁 𝑐𝑎 =7×10 -4 . This shows that our foam flow experiments in 1 mm glass-bead packing (2.55×10 -3 <𝑁 𝑐𝑎 <5.42×10 -3 ) are not affected by the gas trapping.

We also observed the foam's compressibility and gas trapping effects by measuring the pressure gradient along the column and weighing the effluent. The pre-generated foam column was first connected to the GB column that was initially saturated by the surfactant solution. Fig. 13 shows a plot of the effluent mass and pressure gradient as a function of the main column PV. The PV of the GB packing was 180.9 mL, which is shown in Fig. 13 by the continuous black line. The total-injection rate into the foam generator was 2 mL/min, including 1.20 mL/min of N2 gas and 0.80 mL/min of surfactant solution. However, the gas breakthrough happened after producing of 2.94 PV of liquid, which corresponds to a slope change in the figure (dotted black line). Hence, 2.94 times more liquid was recovered than initially stated in the column. This indicates compressibility and trapping of gas in porous media that delayed the gas breakthrough time.

From Fig. 12, we discovered that the gas was trapped considerably below 37,000 Pa/m. Thus the liquid flowed within foam lamellae at the low pressure gradients. This phenomenon was also observed by [START_REF] Falls | The apparent viscosity of foams in homogeneous bead packs[END_REF]. They found the flow of foam lamellae through a single-flow path where stationary lamellae blocked most of the cross-sectional area of the bead-pack at sufficiently low pressure gradients. Since the gas was trapped in the porous medium, all liquids were transported through the lamellae, thus delaying the gas breakthrough. Therefore, we recovered 385 mL of liquid when the gas breakthrough occurred, that is, 2.1 times more than PV. Moreover, based on the trend equation of liquid production for the first half of the curve, the liquid recovery rate was approximately 128.8 mL per injected PV, which is 24% less than the injection volume (180.9 mL). If we consider that the liquid is an incompressible fluid, the volume of gas is reduced to 40%. Therefore, we can conclude that foam flow at low capillary numbers is influenced by the trapping effect and at high pressure gradients by the compressibility, which tends to decrease the foam's gas volume. After the breakthrough, when the foam had reached the column outlet, the change of effluent weight corresponded to the mass of the injected fluid. As a consequence of these results, foam compressibility must be taken into account during the study of foam rheology in confined media.

Since most of the foam volume is gas, compressibility is one of the essential parameters. The compressibility is often described in terms of the expansion ratio 𝜀 which is defined as the ratio between the density of liquid phase 𝜌 𝑙 and the foam density 𝜌 𝑓 (𝜀=𝜌 𝑙 /𝜌 𝑓 ) [START_REF] Valko | Volume equalized constitutive equations for foamed polymer solutions[END_REF]. The expansion pressure should be determined for each pressure drop. Using the volume-equalization method proposed by [START_REF] Valko | Volume equalized constitutive equations for foamed polymer solutions[END_REF], one should obtain a single master curve by plotting normalized-shear stress (𝜏/𝜀) versus normalized-shear rate (𝛾/𝜀) regardless of the pressure and foam quality. According to [START_REF] Valko | Volume equalized constitutive equations for foamed polymer solutions[END_REF], this method is relevant for low expansion foams (𝜀<4) with spherical bubbles. [START_REF] Gardiner | Prediction of pressure losses in pipe flow of aqueous foams[END_REF] highlighted an increase in apparent-foam viscosity with the expansion ratio by analyzing and fitting the foam-viscosity data from previous research. They also found that the volume-equalization method proposed by [START_REF] Valko | Volume equalized constitutive equations for foamed polymer solutions[END_REF] is applicable to polyhedral (𝜀>5) and bubbly-to-polyhedral (5≥𝜀≥4) foams [START_REF] Gardiner | Rheology of fire-fighting foams[END_REF]. Most of these studies were done using simple geometries such as pipes. Therefore, the determination of the expansion ratio in a porous medium is more difficult, if only the pressure drop along the column is measured.

To avoid the compressibility effect, we investigated the upscaling of foam flow in a more highly permeable porous medium based on the results of 4 mm glass-bead packing from [START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF].

In this case, foam was pre-generated using the same surfactant solution and sand-pack. The maximum pressure gradient was 23,210 Pa/m at the flow rate of 3 mL/min at 𝑓 𝑔 =85%. Hence, the permeability of the unit cells was adjusted to the permeability of 4 mm glass-bead packing, 𝐾=1.09×10 -8 m 2 . No gas trapping was expected as the capillary numbers were 4.5×10 -3 and 7×10 -3 at the flow rates of 0.22 and 3 mL/min, respectively. Fig. 14 shows the upscaling results through the different 2D and 3D geometries with the experimental outcomes from [START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF]. In this case, the upscaling results are close to the experimental data.

Because the compressibility effect is insignificant in 4 mm glass-bead packing (𝜀≌1), we estimated the expansion ratio for 1 mm glass-bead packing using the volume equalization method of [START_REF] Valko | Volume equalized constitutive equations for foamed polymer solutions[END_REF]. In Fig. 15, we plotted the shear stress versus the shear rate of 𝑓 𝑔 = 85% foam in 1 mm glass-bead packings without and with mean-expansion factor, based on the 4 mm GB packing results. We found a mean-expansion factor of 𝜀 = 2.1 for 1 mm glass-bead packing. Because the bubble shape of the foam in porous media is more spherical, the method of [START_REF] Valko | Volume equalized constitutive equations for foamed polymer solutions[END_REF] for low expansion foams is applicable. We note also the simplicity of the unit cells considered here where the pore-scale structure is highly ordered which may also influence the results and partly explain the discrepancies observed. Consequently, upscaling is an encouraging tool if the input parameters are settled by taking into account not only the real structure of the porous medium but also all phenomena involved for foam flow in confined media. In this paper, we only investigated the upscaling of foam flow from pore-scale to laboratory-scale. Therefore, the second upscaling from laboratory-scale to field-scale should be explored in the future. Then foam stability, including coalescence and gravity segregation, needs to be investigated. We know that in field applications, foam texture and bubble sizes may be different, as foam-residence time can be much longer.

Empirical Darcy-scale law

To further study foam flow in porous media using the bulk foam viscosity as input, we used the purely empirical approach of [START_REF] Chevalier | Darcy's law for yield stress fluid flowing through a porous medium[END_REF] for the flow of yield-stress fluids (YSF) in glass-bead packings with a bead diameter 𝑑, which is presented in the following form [START_REF] Chevalier | Darcy's law for yield stress fluid flowing through a porous medium[END_REF]:

∇𝑃 = 𝜒𝜏 0 𝑑 + 𝜔𝑎 ( 𝑢 𝑑 ) 𝑛 𝑑 (19)
They considered a polymer (Carbopol aqueous solution) and an emulsion (water-in-oil emulsion with a surfactant solution) as yield-stress fluids. In their study, the relationship between the pressure gradient (∇𝑃) and the Darcy velocity (𝑢) is of the same form as the constitutive law of the fluid, which contains the yielding term. In the preceding equation, terms 𝜒 and 𝜔 are dimensionless coefficients that depend on the distribution of shear rate intensity and on the power-law index [START_REF] Chevalier | Breaking of non-Newtonian character in flows through a porous medium[END_REF]. To determine the dimensionless coefficients of YSF Darcy's law, we used the H-B-P fitting parameters of bulk foam and deduced 𝜒 and 𝜔 parameters by adjusting to 𝑢(∇𝑃) results of foam flow in GB packing. Fig. 16 presents the experimental results with the fitted YSF Darcy law, in which the colored dots represent the results of experiments through the GB packing and dotted lines are the results of the extended Darcy's law for yield-stress fluids with fitted parameters. The corresponding dimensionless coefficients and coefficient of determination for each foam quality are presented in Table 5. The values of the dimensionless coefficients determined are found to be very high in comparison to previous results of [START_REF] Chevalier | Darcy's law for yield stress fluid flowing through a porous medium[END_REF][START_REF] Chevalier | Breaking of non-Newtonian character in flows through a porous medium[END_REF] [START_REF] Chevalier | Darcy's law for yield stress fluid flowing through a porous medium[END_REF][START_REF] Chevalier | Breaking of non-Newtonian character in flows through a porous medium[END_REF], in which they found 𝜒 equal to 12 and 5.5 for the polymer solution and for the water-in-oil emulsion, respectively.

The dimensionless coefficient 𝜔 was 85 for both types of fluids.

Our results showed that 𝜒 decreases with increasing foam quality, and is 2-6 times higher than the values found for the polymer and the emulsion. Therefore, we recommend considering 𝜒 as a function depending on foam quality. Moreover, it was observed that the coefficient 𝜔 was 870 for foams with 85% and 95%

foam qualities, but it was 535.6 for 𝑓 𝑔 =60%. These results can be an outcome of power-law index values, which was high for 𝑓 𝑔 =60% foams. Thus in our case, 𝜔 should be considered individually for bubbly liquids (𝑓 𝑔 <64%), low-quality (64%≤𝑓 𝑔 ≤97%) and high-quality foams (𝑓 𝑔 >97%).

Consequently, we can conclude that the Darcy-scale law depends on foam quality as foam rheology changes with the fraction of gas and liquid.

Conclusion

This work investigated the upscaling of foam flow for very high permeability porous media. The impact of foam quality on the apparent foam viscosity, with a fixed flow rate, was examined to distinguish the low and high-quality regimes. The rheology of pre-generated foam through the fine sand was studied as a bulk foam in a rheometer. Foams of different qualities were also injected into a porous column filled with 1 mm glass beads, from which the foam rheology was deduced. The results of both methods were fitted by a rheological model (Herschel-Bulkley-Papanastasiou model) and analyzed by comparison with each other. Foam flow's macroscopic behavior was also numerically studied in different periodic 2D and 3D unit cells by using the rheological data of bulk foam as an input. We also used the empirical Darcy scale law proposed by [START_REF] Chevalier | Darcy's law for yield stress fluid flowing through a porous medium[END_REF] to derive dimensionless coefficients from our experimental results. We draw the following conclusions:

 Bulk-foam was found to behave as a yield stress fluid and the yield-stress values increased with foam quality. The rheological behavior of foam fits well with the Herschel-Bulkley model.

 Foam in highly permeable porous media (𝐾=8.3×10 -10 m 2 ) was also found to behave as a yieldstress shear-thinning fluid regardless of the foam quality. It was also shown that the apparent foam viscosity in porous media increases with the foam quality at the same total flow rate.

 Apparent foam viscosity was seen to be 3-4 times higher than bulk foam viscosity at the same shear rate, which may be due to the complexity of the foam flow in porous media. However, this difference was smaller for very highly permeable porous media. Secondly, the values of dimensionless coefficients of empirical extended Darcy's law exhibited high values for foam compared to other yield-stress fluids and varied depending on foam quality.  Considerable differences were observed in the experimental and numerical results of upscaling.

The global trend of the curves and the behavior in terms of 𝑓 𝑔 followed the same tendency as those of the experimental results. The experimental results for a very highly permeable porous medium, where the compressibility effect is not significant, fitted well to the experimental data.

In the absence of the foam trapping (e.g., a high capillary number), the difference between upscaling results and experimental data can be related to the foam compressibility. Therefore, foam rheological behavior in porous media can be easily predicted from the bulk foam in the case where the foam compressibility is negligible. The calculated compressibility coefficients should be increased with increasing pressure gradients.

 Upscaling can be an encouraging tool to predict foam-flow behavior from its bulk behavior.

However, further investigations are necessary to incorporate the complex nature of the foam flow in porous media into the upscaling model. These insights are a promising point for studying foam flow in highly permeable porous media.

Moreover, from Fig. 2A, we can see that the shape of the bulk foam bubbles was changed from spherical to polyhedral over time, while the shape of the foam bubbles in porous media did not change, and it took the shape of pores to a greater degree. In Fig. 3A, we plot the characteristic ratio of the bubble area (𝐵 𝑐 ) as a function of time, where 𝐵 𝑐 is the ratio of MBApm to MBAbulk at the same time. In the beginning, the ratio is close to one. However, the values of 𝐵 𝑐 decrease with time showing that the variation of foam bubbles is much higher than for foam in the porous medium. Nevertheless, it must be noted that these experiments were conducted at ambient pressure in the static case. In addition, the height of the system exceeds 200 mm, which tends to increase bubble coalescence due to gravity drainage. In the 1D column experiment, we expected the bubble change to be less than that detected by the foam analyzer. Because the foam was examined horizontally, therefore, the effect of gravity drainage was assumed to be weaker. Moreover, the foam in the 1D column flowed under a certain pressure, and according to the study by [START_REF] Holt | Effects of pressure on foam stability; Implications for foam screening[END_REF], the foam strength (made of AOS surfactant) increases with increasing system pressure.

In our recent study [START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF], we compared the behavior of foams in 1 mm GB packing, where the foams were pre-generated in a 1 mm bead-pack and in a sand-pack used in this study. From the plot of apparent viscosity versus total flow rate, we noticed that the apparent foam viscosities vary greatly even at low flow rates, despite the long residence time of the foam in the column. Thus, it can be concluded that foam stability in confined environments is much longer than in bulk form. Nevertheless, we confirm that there was a small change in the size of the bubbles in the main column outlet. This should be studied with great care in future studies.

Fig. 2 A

 2 Fig. 2 A schematic illustration of a periodic porous medium with the associated REV

Foam

  Fig. 3 2D and 3D periodic unit-cells for an array of parallel cylinders of circular (a, c), square (b, d) cross-sections, (e, f) simple cubic and cubic body-centered packings

Fig. 5

 5 Fig. 5 Rheology of bulk foam with different foam quality

Fig. 9

 9 Fig. 9 Apparent viscosity vs shear rate a) experimental results b) fitted results with 𝛼 shift parameter

  Fig. 12 Pressure gradient and capillary number as a function of time after the foam flow has stopped (𝑓 𝑔 =85%)

  Fig.13Liquid mass of effluent as a function of PV for GB packing at 𝑄 𝑡 =2 mL/min for 𝑓 𝑔 =60%

Fig. 14

 14 Fig.14The superficial velocity of foam flow vs. the macroscopic pressure gradient in 2D and 3D unit cells and experimental results in 4 mm GB packing ofOmirbekov et al. 2020 (𝑓 𝑔 =85%)

Fig. 16

 16 Fig. 16 Superficial velocity as a function of the pressure gradient for column experiments and the fitting curve by YSF Darcy's law for different foam qualities (𝑓 𝑔 =60%, 85% and 95%)

Fig. 2A

 2A Fig. 2A Change of foam structure as a function of time for bulk foam and foam in porous media (1 mm GB packing): liquid and solid phase in black, gas phase in color

Fig

  Fig. 3A The characteristic ratio of bubble area versus time (𝐵 𝑐 = 𝑀𝐵𝐴 𝑝𝑚 𝑀𝐵𝐴 ⁄ 𝑏𝑢𝑙𝑘 )

Table 1

 1 Properties of porous media

	Porous	Mean grain size	Pore Volume,	Porosity,	Permeability, 𝐾	Mean pore
	media	diameter, d (mm)	PV (mL)	∅ (%)	(m 2 )	radius, 𝑅 𝑒𝑞 (µm)
	Sand BR37	0.135	51±2	38±1	7.2 (±1)×10 -12	11.5
	GB	1	181±2	36±1	8.30 (±0.1)×10 -10	133.5

Table 2

 2 Fitting parameters of H-B-P model

	𝒇 𝒈 (%)	60%	85%	95%
	𝒏 (-)	0.48	0.52	0.54
	𝝉 𝟎 (Pa)	2.03	2.87	4.49
	𝒂 (Pa.s n )	0.94	0.76	0.7
	𝒎 (s)	10000 10000 10000
	R 2	0.99	0.99	0.99

Table 3

 3

	𝒇 𝒈 (%)	60%	85%	95%
	𝒏 (-)	0.65	0.53	0.49
	𝝉 𝟎 (Pa)	6.05	6.86	8.06
	𝒂 (Pa.s n )	3.15	3.88	3.90
	𝒎 (s)	10000	10000 10000
	R 2	0.99	0.99	0.99

Fitting parameters for H-B-P model

Table 4 Shift

 4

	GB size		1 mm
	𝒇 𝒈 (%)	60%	85%	95%
	𝜶 (-)	0.34		0.44
	R 2	0.99	0.99	0.99

parameter for various foam qualities obtained from 1 mm GB packing

Table 5

 5 The dimensionless coefficients determined for foams (𝑓 𝑔 =60%, 85% and 95%)

	𝒇 𝒈 (%)	60%	85%	95%
	𝝌 (-)	38.2	31.8	26.4
	𝝎 (-)	535.6	870	870
	R 2	0.99	0.98	0.98

Acknowledgments

This study was performed as part of the "Famous" project. The authors would like to thank ADEME for co-funding the project under the "GESIPOL" program and Bolashaq International Scholarship from Kazakhstan government for providing the PhD grant for Sagyn Omirbekov. We gratefully acknowledge the financial support provided to the PIVOTS project by the "Région Centre -Val de Loire" and the European Regional Development Fund.

Appendix A

To study the variation of foam bubbles in bulk form (Fig. 1A a) and in 1 mm GB packing (Fig. 1A b), we used a dynamic foam analyzer (DFA-100, KRUSS). The foam analyzer consists of a glass column (250 mm long and 40 mm internal diameter), a camera, and a light source. The experiments were conducted as in studies by [START_REF] Omirbekov | Experimental Study of Non-Newtonian Behavior of Foam Flow in Highly Permeable Porous Media[END_REF], where N2 gas was injected into the surfactant solution (AOS, 4 × CMC) using a porous-glass disc. The porous-glass disc was chosen according to the foam generator, where the pore size varied from 10 to 16 µm. In the porous medium, the used volumetric amount of the surfactant solution and the gas injection time were reduced by 60% since the porosity of the bead-pack was about 40%. In both cases, the camera height was set to 100 mm to record the change in bubble size over time. 2A shows the change of structures in bulk foam and foam in 1 mm GB packing, where the mean bubble area (MBA) was compared after 60, 300, and 900 seconds. We observed that the bulk foam is less stable than foam in porous media. It is the consequence of the drainage process, which is greater for the bulk foam. For instance, the MBA of bulk foam doubled after 300 seconds, while the MBA of foam in 1 mm GB packing increased twice in 900 seconds, which was three times lower than the bulk foam. The bulk foam decays faster than the foam in porous media due to gravity drainage and gas diffusion, which tends to coalescence and coarsen average bubble size. While in porous media, the effect of gravity is small due to the short length scales within the bead-pack; therefore, the hydrostatic pressure differences are slight [START_REF] Jones | Surfactant screening for foam EOR: Correlation between bulk and core-flood experiments[END_REF]. [START_REF] Jones | Surfactant screening for foam EOR: Correlation between bulk and core-flood experiments[END_REF] also found a correlation between apparent foam viscosity in porous media and the stability of foam in bulk form.