Population-balance modeling of CO2 foam for CCUS using nanoparticles

Foam implementation for carbon capture, utilization and storage (CCUS) can greatly improve CO2 mobility control, resulting in enhanced hydrocarbon production and carbon storage capacity. The use of nanoparticles (NP) to create robust foam structures has recently gained attention. Local foam generation and coalescence dynamics can be described by mathematical models. Here we address knowledge gaps for NP foam in porous media by tracking the bubble density (nf ) of NP foam data spatially and temporally using an established surfactant (SF) population-balance model. We suggest a reduced shear-thinning effect, compared to SF, to accurately model NP CO2 foam flow in both the highand low-quality regime. A NP foam rheological transition appeared at gas fraction fg 1⁄4 0.85. The nf parameter increased linearly with distance from inlet for NP foam, with reduced CO2 mobility and improved displacement efficiency compared to co-injections of water and CO2.


Introduction
Carbon geo-sequestration can benefit from improved CO 2 mobility control during fluid displacements and flow in porous media. Reduced gravity override and viscous fingering due to CO 2 foam generation leads to improved volumetric sweep in hydrocarbon reservoirs and carbon storage sites. Thus, a successful implementation of CO 2 foam flow enhances carbon storage capacity and hydrocarbon production compared to conventional low-viscous gas flooding. A key success criterion for CO 2 foam displacement is developing stable and strong foam structures. Surfactants (SF) are typically used to stabilize the lamellae separating individual CO 2 bubbles, demonstrating excellent properties as a foaming agent, albeit suffering from great adsorption loss and often becoming weakened by challenging reservoir conditions. Therefore, the emerging trend of using tailored nanoparticles (NP) to create more robust foam structures has received a lot of attention recently (Binks and Lumsdon, 2000;Dickson et al., 2004;Sun et al., 2017;Yekeen et al., 2018;Rognmo et al., 2019). While NP are usually less effective at generating foam compared with SF, the mechanical and thermal stability of NP make them more tolerant against high temperatures, pressures, shear and salinity (Bennetzen and Mogensen, 2014). NP at interfaces are associated with high surface adsorption energies that practically eliminate them from being desorbed (Binks, 2002), adding stability to CO 2 foam flow in porous media. Further, co-injecting silica NP solution with CO 2 lowers the pH of the mixture, reported to reduce the overall risk of NP aggregation . SF foam has been characterized extensively, and much is known about the foam behavior and rheology, particularly without the presence of oil (Bernard et al., 1980;Hirasaki and Lawson, 1985;Khatib et al., 1988). Foam is created by leave-behind, lamellae division and snap-off mechanisms, and the stability of foam in porous media is limited by capillary pressures (Khatib et al., 1988). A detailed description of foam creation and destruction mechanisms can be found elsewhere (Nguyen et al., 2000). Flowing SF foam is frequently reported to exhibit a shear-thinning behavior (Hirasaki and Lawson, 1985;Khatib et al., 1988;Marsden and Khan, 1966;Heller and Kuntamukkula, 1987;Falls et al., 2007;Fernø et al., 2016), explained by bubbles slipping on the pore wall and against each other. NP foam is less studied and the literature contain sparse information on NP flow properties, especially in porous media. Several authors report near-Newtonian NP foam behavior (Lotfollahi et al., 2016a;AlYousif et al., 2018), whereas others state a clear shear-thinning effect (Griffith et al., 2016;Worthen et al., 2015), in contrast to reported shear-thickening effects (Mo et al., 2012).
Predictive mathematical models are needed to describe foam flow in porous media accurately, currently underdeveloped for NP foam, for large-scale implementation of CO 2 foam in CCUS applications. Specifically, population-balance models have the unique potential to describe both transient and steady-state conditions for foam flow by spatially and temporally tracking the bubble density (n f ). The n f parameter is a function of generation and coalescence mechanisms, and is of great significance as bubble density/texture is directly coupled to flow resistance (Ettinger and Radke, 1992). This paper expands the application of an established population-balance model by including NP CO 2 foam data to improve characterization of dynamic bubble densities containing nanoparticles in the pore network. The paper is structured to first introduce the experimental data used for core-scale history matching. Then, the population-balance model and key parameters are presented. Further, we tuned the mechanistic model by adjusting specific kinetic terms introduced in section 1.2 to match the laboratory observations. Lastly, we characterized NP foam properties and directly compared the results with a baseline (no foaming agent) and a SF foam flood.

Experimental data
The experimental CO 2 foam flow data in this paper are based on earlier work from the same research group (Rognmo et al., 2017;Horjen, 2015;Rognmo, 2018). Co-injections of liquid CO 2 and a foaming solution (NP or SF dispersed in brine) were performed in two Bentheimer sandstone core samples (Table 1) to generate either NP CO 2 foam or SF CO 2 foam for a range of gas fractions and total injection rates. The pore pressure and temperature were 90 bar and 20 � C, respectively.
The sandstone core samples were fully saturated with foaming solutions (Table 2) prior to co-injections. The SF solution consisted of C 14-16 alpha olefin sulfonate (AOS, Stepan Company) dispersed in brine. To generate NP foam in porous media, particle concentration and shear rate should exceed a threshold value to maintain foam stability (Espinoza et al., 2010). Moreover, particle hydrophobicity can affect NP foam generation rate and bubble size/texture (Yu et al., 2014), and the silane-modified silica particles used herein (Levasil CC301, provided by Nouryon) obtained a hydrophilic surface (Rognmo et al., 2017). The small mean particle diameter (Table 2) allowed NP to flow unrestricted through the Bentheimer pore network, and NP foam appeared stable (low coalescence rate) at concentrations ranging from 0.15 to 0.5 wt% particles dispersed in brine at experimental conditions (Rognmo, 2018).
The effective (apparent) gas viscosity, μ f , was derived from experimental data by the expression: where k is permeability of the porous medium, Δp is the pressure gradient across the medium, and u l and u g is velocity of liquid and gas, respectively. Gas fraction (f g ) is commonly referred to as foam quality, and is the ratio of flowing gas rate (q g ) to total injection rate (q g þ q l ): Co-injections with NP solution and CO 2 used three different total volumetric injection rates (120, 180 and 240 ml/h) at f g ¼ 0.1, 0.2, 0.35, 0.5, 0.7, 0.85, 0.9 and 0.99 during drainage (increasing gas fraction). Co-injections with SF solution and CO 2 foam used two different total volumetric injection rates (120 and 180 ml/h), at f g ¼ 0.1, 0.5, 0.7, 0.9, 0.95 and 1.0 during drainage. The maximum effective viscosity for NP foam quality scan, expressed as effective viscosity as function of increasing f g , was observed at f g ¼ 0.7 (Fig. 2), compared with f g ¼ 0.9 for SF foam. Further, NP foam showed near-Newtonian flow (similar effective viscosity for all injection rates), whereas SF foam demonstrated the expected shear-thinning behavior (lower effective viscosity for higher injection rates). The effective viscosity was approximately two orders of magnitude higher compared than NP foam (Rognmo et al., 2017).
We propose that the reported experimental results of NP CO 2 foam flow is due to a reduced number of foam bubbles and reduced surface interaction compared with SF CO 2 foam. This is corroborated by comparing pore-scale bubble density between SF and NP CO 2 foam at similar experimental conditions using a high-pressure silicon micromodel, since reliable in situ bubble density in opaque core systems are very difficult to obtain. Direct comparison of foaming agents (Fig. 1) demonstrated quantitatively a substantial reduced bubble count for NP CO 2 foam (Benali, 2019) at equal PV CO 2 injected. Foam bubble density analysis (baseline ¼ 45 bubbles, SF foam ¼ 506 bubbles, NP foam ¼ 366 bubbles) revealed that both surfactants and nanoparticles were able to generate strong foam compared to the baseline without foaming agent (no stable bubbles). NP CO 2 foam bubbles were spherical in shape and heterogeneously distributed, whereas SF CO 2 foam generated a higher number of pore-spanning lamellae that were uniformly distributed in the pore network. Bubble size generally defines the foam texture (n f ) which significantly affects the foam flow properties. Finer foam texture (higher n f ) implies lower gas mobility (Kovscek and Radke, 1994).

Model description
In this study we use a well-established SF population-balance model to simulate CO 2 foam flow (Chen et al., 2010;Kovscek et al., 1995). The model keeps track of the number of foam bubbles in the porous media and adjust the effective gas viscosity accordingly. In a one-dimensional porous medium, the mass balance of foam bubbles can be expressed as: where φ is porosity, S f is flowing gas saturation, n f is foam texture/ bubbles per unit of flowing gas, S t is trapped gas saturation, n t is foam bubbles per unit of trapped gas, and u f is Darcy velocity. S g is the sum of flowing and trapped gas saturation (S g ¼ S f þ S t ), k 1 is foam generation rate, k À 1 is foam coalescence rate, v w is interstitial velocity of water, v f is interstitial velocity of gas, and Q b is a source/sink term for foam bubbles. Because no foam was pre-generated in the system, Q b is set to 0. Expressions for standard multiphase form of Darcy's law and relative permeability formulations with standard Corey exponent models can be found in Supplementary material and in the original model (Kovscek et al., 1995).
Foam generation rate is written as: where k 1 is the foam generation rate, k 0 1 is a foam generation constant, n f is foam texture, n * is limiting foam texture and ω determines the function type of foam generation. Because limiting foam texture depends on Table 1 Rock properties from (Rognmo et al., 2017) used as input in the population-balance model.  rock properties and we assume one bubble per pore (Ettinger and Radke, 1992;Bertin et al., 1998), n * was kept constant at 4000 mm À 3 for all simulations, both SF and NP CO 2 foam. Foam coalescence rate can be expressed as: where k À 1 is foam coalescence rate, k 0 À 1 is a scaling constant, P c is capillary pressure and P * c is limiting capillary pressure for foam coalescence. The limiting capillary pressure is dependent on foaming agent concentration and can be expressed as: where P * c;max is a limiting value for P * c , C s is the foaming agent concentration and C 0 s is a concentration threshold for strong foam generation. When the core is pre-saturated with a foaming solution, it is assumed that the concentration is uniform and rock adsorption is satisfied. In this instance, the foaming agent mass balance is automatically satisfied (Kovscek et al., 1995). The effective (apparent) gas viscosity can be expressed as: where μ f is effective gas viscosity (subscript f denotes flowing foam), μ g is gas viscosity, α is the viscosity proportionality constant, n f is foam bubbles/texture, v f is gas velocity and c is a velocity exponent which is theoretically close to 1 3 for shear thinning foams (Hirasaki and Lawson, 1985;Kolb and Cerro, 2002). This relationship in Eq. (7) means that increasing number of foam bubbles increases effective foam viscosity, and inversely, at insignificant bubble number and high gas velocity the initial gas viscosity is recovered (μ f ¼ μ g ). The fraction of foam trapped, X t , can be written as: where X t;max is maximum fraction of trapped foam, n t is foam bubbles per unit of trapped gas, and β is a trapping parameter. Capillary pressure  (Benali, 2019), showing a baseline (CO 2 and waterleft image) in sandstone pore network (grains ¼ black shapes) without foaming agent, SF CO 2 foam (middle) and NP CO 2 foam (right). Quantitative bubble density analysis (baseline ¼ 45 bubbles, SF foam ¼ 506 bubbles, NP foam ¼ 366 bubbles) revealed that both surfactants and nanoparticles were able to generate strong foam compared to the baseline (no foaming agent) where the CO 2 phase remained continuous in the pore system. Degree of match is displayed as mean squared error (MSE). Lowest (f g ¼0.1) and highest gas fraction (f g ¼1.0) has not been emphasized in the history matching due to experimental uncertainties. The highest effective viscosity in the foam quality scans was observed at f g ¼ 0.7, independent of total flow rate.
was modelled using the Leverett-J function (Leverett, 1941): where P c and σ are capillary pressure and interfacial tension respectively, between the phases CO 2 and water, k is permeability and φ is porosity. The aqueous solution is assumed to be incompressible, whereas the gas density scales linearly with pressure. Several different population-balance models exist, e.g. Chang et al., 1990;Kam and Rossen, 2003), but the key characteristics are similar. The differences between various models are typically the dominant foam generation mechanism , and the effect of oil on foam stability (Hematpur et al., 2018). The experimental data presented in this paper did not include an oil phase. All existing foam models, population-balance versions and implicit-texture versions, provide non-unique solutions (Lotfollahi et al., 2016b), however, the parameters chosen here are derived from physical observations.

Simulation results and discussion
A one dimensional (1D) model with 50 cells was created to replicate the core sample reported in the laboratory co-injections using rock and fluid parameters in Tables 1 and 2. The injection strategy reflected reported experimental conditions, described as co-injection of either 1) water and CO 2 (baseline; no foaming agent), 2) NP solution and CO 2 , or 3) SF solution and CO 2 . Limiting foam texture was set to 4000 mm À 3 based on pore size of 0.07 mm for pores with foam, and cubical packing of the foam bubbles. Cubical packing was used because the limiting foam texture is associated with wet foam where the shape tends to be spherical (Ettinger and Radke, 1992;Alvarez et al., 2001;Belyadi et al., 2017). The average pore size in the model was based on pore size distribution measurements in Bentheimer sandstone (Peksa et al., 2015), corroborating values reported elsewhere (Liaw et al., 1996;Barifcani et al., 2015). Overview of parameters used in the population-balance model (Kovscek et al., 1995) are presented in Table 3.

History-matching experimental data
The matching procedure for the population-balance model required an effective viscosity profile obtained at steady state (presented in Fig. 2), along with the accompanying subset data of pressure gradient vs volumetric gas rate at fixed water rate (presented in Fig. 3a), and the pressure gradient vs volumetric water rate at fixed gas rate (presented in Fig. 3b). Often, these data sets can be obtained within one experimental run (Kovscek et al., 1995). Steady-state conditions were not met for endpoint f g in the available laboratory co-injection raw data (Horjen, 2015), see Supplementary material for transient pressure data. Here, the differential pressure was still declining when f g was changed; hence, we expect the effective viscosity, μ f , to be lower than reported in (Rognmo et al., 2017) at f g ¼ 0.1 and f g ¼ 0.99. The experimentally measured data points without true steady-state conditions are indicated using smaller-sized points to reflect the lower confidence for the reported effective viscosities for these f g . NP foam effective viscosities simulated using the theoretical velocity exponent of c ¼ 1/3 matched the reported experimental data set moderately, with μ f gradually increasing with increasing gas fraction (Fig. 2 -broken lines). Maximum gas mobility reduction occurred at f g ¼ 0.7, followed by a gradual destruction of foam bubbles dictated by the capillary forces as foam quality approaches unity. The population-balance model reproduce the experimental trends with lesser degree of match at f g ¼ 0.1 and f g ¼ 0.99, however this is likely due to the aforementioned unfulfilled steady-state conditions. The experimentally measured effective viscosities (marked as data points) appear to be insensitive to total injection rate for each f g , whereas a flow rate dependency was observed in simulated data with decreasing viscosity with increasing total flow rate at fixed foam quality for most of the gas fraction interval.
The effective gas viscosity, μ f ; in the model is calculated using Eq.
(7), where the default value gives a rate dependency and non-Newtonian behavior of the flowing foam. As stated in the introduction section, shear-thinning behavior is observed for most SF foam where the pore wall is coated with surfactants causing the bubbles to slip at high shear. Shear thinning is also observed in other emulsion systems, both in bulk and porous media, for example oil-in-water emulsions (Zhang et al., 2010;Karambeigi et al., 2015;Zhou et al., 2017). In addition to bubble slipping, shear thinning has also been attributed to bubble size dependence on flow rate (Zhou et al., 2017). At low injection rate and large bubble size, high μ f was observed as bubbles blocked pore throats, in contrast to higher rates and reduced bubble size where the emulsions could more easily flow in the porous media. Recent literature data on flow rate dependency on NP foam is not coherent: shear-thinning properties in glass bead pack is reported (Xiao et al., 2017), but others found little evidence of shear-thinning viscosity with particle diameters above 12 nm (Kim et al., 2016). Low sensitivity of flow rates was reported for both shear-thinning and shear-thickening behavior for NP foam (AlYousif et al., 2018).
In accordance with our own experimental data and recent observations reported by others, the effect of flow rate on effective viscosity (expressed through the velocity exponent c) was reduced from the default value in the model to replicate reported NP CO 2 foam behavior. Decreasing the exponent c from 1/3 to 1/5.5 required an increase in the proportionality constant α from Eq. (7), to compensate for the overall decrease in μ f . By keeping the exponent above zero we retain some shear dependence, as observed in experiments at similar flow rates (Horjen, 2015), and the degree of match was improved (Fig. 2 -solid lines). The proportionality constant α was increased from 4.05 � 10 -17 for shear-thinning foam systems (Fig. 2 -broken lines) to 1.10 � 10 -16 for near-Newtonian NP foam (Fig. 2 -solid lines). The change in viscosity and subsequent increase of α, implies that when NP CO 2 foam becomes less shear thinning it also becomes stronger at equivalent flow rate and gas fraction. The rheology of NP CO 2 foam is further discussed in section 2.2. The experimental data from (Horjen, 2015) has also been simulated using a different population-balance approach, focusing on the gas fraction interval f g ¼ 0.5 to 1.0 (Ortiz et al., 2019). The authors achieved good agreement between experimental and simulated data, however the model relied on limiting water saturation to separate the high and low-quality regime, using different velocity exponents to separate the degree of non-Newtonian viscosity in the two regimes. We attempt to reproduce the entire foam quality scan, f g ¼ 0.1 to 1.0, by using a single velocity exponent (c¼1/5.5). The degree of match is shown in Fig. 3, displaying pressure gradient profiles as function of gas rate (3a), and function of water rate (3b), and a contour plot showing effective viscosity μ f as a function of both water and gas rates (Fig. 4). At a constant total flow rate, μ f increases with decreasing gas fraction in the high-quality regime, shifting to decreasing μ f with further decreasing gas fraction in the low-quality regime. As most of the data points correlates with the contour color, we have obtained a good match throughout the different rates and foam qualities.

NP CO 2 foam characterization by population-balance modeling
Previously, two distinct regimes have been identified for foam flow (Osterloh and Jante, 1992;Vassenden and Holt, 1998), a low-quality (wet foam) regime and a high-quality (dry foam) regime separated by a transition zone. Foam viscosity has been reported to be a function of these flow regimes (Alvarez et al., 2001) partially explaining the inconsistent rheological behavior reported in the literature. We find a moderate shear-thickening behavior at high NP foam quality (f g ¼0.9) and a shear-thinning behavior at low NP foam quality (f g ¼0.5) when plotting effective gas viscosity against total injection rate (Fig. 5), corroborating previous observations of SF foam (Alvarez et al., 2001). The transition from shear-thinning to shear-thickening NP foam flow emerges around gas fraction f g ¼ 0.85 (Fig. 2). Further, we observe that the effective viscosity μ f of NP foam stabilizes as the total flow rate increases (near-well conditions). It appears that the experimental data that we are matching, are from a less-sensitive range of the spectrum as illustrated in Fig. 5, where the experimental f g data (points) and simulation data f g (lines) are displayed in corresponding colors. More experimental data are needed for a wider range of flow rates to completely map the influence of foam quality and flow rate on NP CO 2 flow properties. Our preliminary results indicate that the viscosity of NP CO 2 foam depends on foam quality below a certain threshold rate (�200 ml/h for f g ¼ 0.9). Above this threshold, flowing NP foam appears near-Newtonian and less sensitive to injection rates and gas fractions. There have been reported near-Newtonian properties for SF foam (Ettinger and Radke, 1992;Persoff et al., 1991) due to increased snap-off frequency at higher liquid velocities. This effect is partially included in the simulation model, where the foam generation is dependent on gas and liquid velocities.
The effective viscosity μ f during NP foam flow was about an order of magnitude higher than water, and the residual saturation values during steady state NP CO 2 foam approached the connate water saturation for a range of gas fractions. At high f g (0.9), the baseline co-injection resulted in final water saturation of S w ¼ 0.55, compared with S w ¼ 0.30 for NP foam. Dynamic water profiles during the baseline fluid displacements (Fig. 6) were horizontal at all time steps during the co-injection. Hence, there was not a well-defined displacement front and the poor mobility ratio between CO 2 and water resulted in viscous fingering, leaving behind unswept zones of high water content in the pore system. In contrast, a distinct fluid displacement front developed during NP foam progression, and the front sharpened over time to become more pistonlike.  water and gas injection rates. Effective viscosity increases with decreasing gas fraction in the high-quality regime and decreases with further decreasing gas fraction in the low-quality regime. Degree of match is indicated by the color map.
Based on the dynamic water profiles in Fig. 6, it is clear that the NP foam generation is not instantaneous. A certain time is needed to generate strong foam to efficiently displacing the initial water. The development of bubbles (n f ) along the core length during NP foam flow provided additional important information on the foam generation process (Fig. 7). A shift in n f occurred between 200 s and 500 s, indicating increased foam generation kinetics and finer foam texture. This change in kinetics was directly linked to the development of a distinct displacement front in Fig. 6 occurring at the same time step. During the transient flow period, some n f profiles exceeded the steady-state texture, likely due to gas compressibility effects (Kovscek et al., 1995). Coalescence forces coarsened the foam texture over time and n f profiles at steady-state conditions (5000 s ¼ 2.2 PV and 99999 s ¼ 43 PV) retained similar shapes as time step 500 s. After foam breakthrough at the outlet (normalized length ¼ 1.0), a non-zero n f value is present at the inlet (normalized length ¼ 0.0), possibly due to a backward front movement of finer foam texture (Simjoo and Zitha, 2020; Apaydin and Kovscek, 2001;Almajid et al., 2019). At steady state, bubbles appeared across the entire core length, with a close-to-linear increase toward the outlet end. In proximity to the outlet of the core sample, the n f value approached the theoretical foam density limit n * at steady-state conditions. Steady state NP and SF effective viscosity profiles highlighted that at current conditions SF CO 2 foam shifted maximum μ f toward higher foam quality (f g ¼0.9) and significantly enhanced the effective gas viscosity (see Fig S2 in Supplementary material). A previous study successfully modelled combined NP and SF foam flow at near-well conditions using a SF foam mechanistic framework, with lamella division as the dominant foam generation mechanism (Prigiobbe et al., 2016). They reported a synergetic effect of enhanced gas mobility reduction from combining the two foaming agents. We assume the main differences between the two systems (NP CO 2 foam versus SF CO 2 foam) in our paper are the average foam texture (Fig. 8) and the limiting capillary pressure, P * c . The average foam texture, n f , obtained at steady state is consistently lower (factor Effective viscosity, μ f , versus injection rate for c ¼ 1/5.5 at two different gas fractions, f g ¼ 0.5 (blue) and f g ¼ 0.9 (orange), indicating shear-thinning trend at low NP foam quality and shear-thickening trend at high NP foam quality. The sensitivity decreases with increasing total flow rate. Fig. 6. Simulated water saturation as a function of normalized core length at constant f g ¼ 0.9, for NP foam (blue line) and for baseline (orange) without foaming agent. NP foam developed a distinct displacement front eventually (at 1000 s ¼ 0.43 PVpore volume), and improved the volumetric sweep compared to coinjection of CO 2 and water (no foaming agent).
2.1) in the low-quality regime for NP foam. Further, the foam texture profile for NP foam collapses earlier than for SF foam due to a lower P * c value represented here. Because P * c strongly depends on foaming agent concentrations (Apaydin and Kovscek, 2001;Aronson et al., 1994), this effect was considerable for the current experimental data. Instead of a direct performance comparison between NP and SF foam, this study emphasizes dynamic foam textures to accurately describe NP CO 2 foam behavior and mobility control potential at low particle concentration (0.15 wt%). Modeling NP foam instead of SF foam impacted the kinetic expressions k 1 (foam creation) and k À 1 (foam destruction), as well as the proportionality constant α and the velocity exponent c (detailed in Table 3).

Conclusions
Dynamic properties of CO 2 foam stabilized with nanoparticles (NP) were predicted using an established surfactant (SF) population-balance model, and led to the following key observations: � The population-balance model matched experimental data of NP foam effective viscosity over a range of gas fractions, by reducing the shear-thinning effect commonly used for SF foam. This was supported by observations of near-Newtonian NP foam behavior in the laboratory at specific conditions. Foam viscosity appeared sensitive to quality regimes at low total flow rates (shear thinning for gas fraction f g ¼ 0.5, and shear thickening for f g ¼ 0.9). � Simulations revealed that high-quality NP foam developed a distinct displacement front and improved local and overall sweep compared to co-injection of CO 2 and water (baseline; contained no foaming agents). Characteristics of the flowing foam were increased effective viscosity and decreased residual water saturation along the core length.

Fig. 7.
Foam texture (n f ) development as function of normalized core length for NP foam at constant f g ¼ 0.9. The first time steps (20, 50, 100 and 200s) follow a similar foam profile with a gradual increase in number of lamellae (� 10 12 m -3 ). During the transient flow period, the density profiles temporarily exceeded the steady-state texture (500 s < t < 3000 s; not shown in Fig. 7), due to gas compressibility effects. Fig. 8. Average foam texture, n f [m -3 ], simulated at steady-state conditions as a function of gas fractions for NP foam (blue) and SF foam (orange). The foam texture is consistently higher by a factor of 2.1 in favor of SF foam. Further, SF foam behavedmore stable than NP foam in the high-quality regime. In both cases, the average foam texture appeared constant over a wide range of gas fractions.
� At steady-state conditions, NP foam bubbles appeared across the entire core length, and the profiles showed increasing density of lamellae with distance within the sample. This led to reduced CO 2 mobility and improved displacement efficiency of low-concentration NP foam compared to baseline co-injections. Compared to highconcentration SF foam, the NP foam texture appeared coarser with lower effective viscosity, mainly ascribed to reduced bubble density and reduced limiting capillary pressure.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.