CO Convection in Hydrocarbon Under Flowing Conditions

Carbon-neutral oil production is one way to improve the sustainability of petroleum resources. The emissions from produced hydrocarbons can be offset by injecting capture CO 2 from a nearby point source into a saline aquifer for storage or a producing oil reservoir. The latter is referred to as enhanced oil recovery (EOR) and would enhance the economic viability of CO 2 sequestration. The injected CO 2 will interact with the oil and cause it to flow more freely within the reservoir. Consequently, the overall recovery of oil from the reservoir will increase. This enhanced oil recovery (EOR) technique is perceived as the most cost-effective method for disposing captured CO 2 emissions and has been performed for many decades with the focus on oil recovery. The interaction between existing oil and injected CO 2 needs to be fully understood to effectively manage CO 2 migration and storage efficiency. When CO 2 and oil mix in a fully miscible setting, the density can change non-linearly and cause density instabilities. These instabilities involve complex convectivediffusive processes, which are hard to model and simulate. The interactions occur at the sub-centimeter scale, and it is important to understand its implications for the field scale migration of CO 2 and oil. In this work, we simulate gravity effects, namely gravity override and convective mixing, during miscible displacement of CO 2 and oil. The flow behavior due to the competition between viscous and gravity effects is complex, and can only be accurately simulated with a very fine grid. We demonstrate that convection occurs rapidly, and has a strong effect on breakthrough of CO 2 at the outlet. This work for the first time quantifies these effects for a simple system under realistic conditions.


Introduction
Carbon dioxide (CO 2 ) injection in oil reservoirs to increase recovery and simultaneously store CO 2 is perceived as a highly cost-effective method for disposing captured CO 2 emissions (Heidug et al. 2015). This enhanced oil recovery (EOR) technique has been performed for many decades with the focus being on hydrocarbon recovery and not long-term CO 2 storage (Dooley et al. 2010). In CO 2 -EOR, CO 2 combines with existing hydrocarbons, increasing the efficiency of displacement and enabling greater recovery of the oil in place. CO 2 acts as a solvent and is highly miscible with oil at high pressure. This affects key fluid-flow processes in the reservoir during miscible flooding.
Recently, CO 2 -EOR has gained renewed attention with increasing focus on CO 2 storage during oil recovery (Elzinga et al. 2015). In any CO 2 -EOR flood, some amount of CO 2 is naturally lost through dissolution and spurious migration, which ultimately reduces the carbon footprint of the produced oil if the injected CO 2 is anthropogenically sourced. However, most injected CO 2 is recycled and any incidentally stored CO 2 is not enough to completely offset the emissions from the incremental oil produced in a conventional CO 2 -EOR operation (Nuñez-Lopez et al. 2019). Various techniques can be implemented to increase CO 2 stored during oil extraction and reduce the carbon footprint of oil production, socalled CO 2 -EOR+ technologies (Heidug et al. 2015). One such approach is to increase the efficiency of injected CO 2 for oil production, which can be achieved by increased reservoir sweep and reduced recycling.
Increasing the efficiency of CO 2 -EOR relies on improved understanding of CO 2 migration in the reservoir. Figure 1 gives a simplified sketch of a CO 2 flood, where CO 2 mixes with oil in the miscible zone, leading to increased oil recovery. This figure shows ideal conditions where CO 2 displaces oil in a piston-like manner. However, CO 2 is less dense and viscous than oil, resulting in an unfavorable mobility ratio that can lead to viscous fingering and a high density contrast that can lead to gravity override (Jenkins 1984). In the viscous setting, the concern is dominated by viscous fingering that destabilizes the miscible front (Stalkup 1983). The less viscous front can develop instabilities and finger through the reservoir, this phenomenon is widely studied, i.e., Bensimon et al. 1986;Homsy 1987. Gravity override, on the other hand, reduces reservoir sweep. Both viscous fingering and gravity override leads to bypassed oil and reduced production (Green and Willhite CO 2 Miscible zone Additional oil recovery Caprock Fig. 1 Injection of CO 2 in an oil reservoir in a CO 2 -EOR setting 1998). CO 2 breakthrough also occurs earlier, leading to more recycled CO 2 and less potential for stored CO 2 in the reservoir. Override can be minimized by increasing the viscous drive, usually by tighter well spacing. However, for offshore reservoirs, well placement is usually constrained by high drilling costs, and therefore tends to be more widely spaced. Greater distance between wells implies shifting flow regimes within the reservoir away from primarily viscous-driven flow, and as a result, density effects will become more important.
Of particular interest is the impact of CO 2 on oil density in the miscible zone. At high CO 2 concentrations, the oil mixture is less dense than the original oil and makes the oil swell. At lower CO 2 concentrations, the oil becomes more dense Both et al. 2015). This non-linearity in the density function of oil-CO 2 mixtures can lead to density inversions occurring as CO 2 mixes with oil in the miscible zone. Density inversion can trigger gravity fingers, resulting in density-driven convection in the miscible zone. This phenomenon can be observed in examples of field-scale injection of CO 2 into oil using both commercial (Ahmed et al. 2012;Both et al. 2015) and high-resolution simulations (Shahraeeni et al. 2015). In these studies, it is demonstrated that gravity fingers develop quickly in both homogeneous and heterogeneous settings.
Note that for certain high pressure and low temperature reservoirs, CO 2 can be denser than oil, which will lead to different type of instabilities. In , experimental work and numerical simulations report that unstable flow occurs when heavy CO 2 is injected above a lighter oil, or from below in combination with other phase-related phenomenon. It was shown that strong convection behavior can be expected, also with a different setting than described herein. Also, in a non-porous media setting, a large effect from mixing and convection was shown at a smaller scale in (Rongy et al. 2012).
Gravity fingering due to CO 2 dissolution is well-studied for CO 2 -brine systems. CO 2 is widely known to increase brine density, which again leads to convection. Convection in CO 2 -brine systems is widely studied in the literature (Elenius et al. 2014;Ennis-King and Paterson 2005;Farajzadeh et al. 2011;Kneafsey and Pruess 2010;Khosrokhavar et al. 2014;Lindeberg and Wessel-Berg 1997;Mykkeltvedt and Nordbotten 2012;Neufeld et al. 2010). In these systems CO 2 dissolution can be upscaled from fine-scale models into fieldscale relevant simulation (Gasda et al. 2011). At large spatial and time scales, convection has a significant impact on CO 2 plume and storage efficiency (Gasda et al. 2012). Most importantly, the act of convection impedes the migration of CO 2 and allows CO 2 to penetrate deeper into the reservoir due to negative buoyancy.
The known effects of gravity fingering in CO 2 -brine systems indicate that convection has the potential to improve sweep and reduce breakthrough in oil systems. Some evidence of this is demonstrated in the simulations by (Ahmed et al. 2012;Shahraeeni et al. 2015). Despite the similarities with CO 2 -brine, it is not evident that analysis and upscaling of CO 2 -hydrocarbon convective mixing can be simply extended from current knowledge. The CO 2 -hydrocarbon system differs from the CO 2 -brine system: The fluids are fully miscible, density is non-monotonic with CO 2 concentration, and significant density gradients occur over short distances (Elenius and Gasda 2021). In general, gravity fingering during miscible displacement of CO 2 and oil is much less studied compared to the viscous case. A recent study shows that non-monotonic density leads to complex convective behavior when gravity forces are isolated from viscous effects Elenius and Gasda 2021. This study characterized onset time and the speed of convection in dimensionless oil-CO 2 systems under idealized conditions. To our knowledge, no similar analysis has been carried out for convective mixing under flowing conditions, and a detailed analysis is required to quantify the balance between viscous and gravity forces in a miscible CO 2 flood. This work addresses this challenge through simulation of fine-scale coupled viscous and gravitational processes during miscible displacement. A series of numerical experiments are performed where the gravity-to-viscous ratio is varied for a simple system. The resulting effect on breakthrough of CO 2 at the outlet shows the importance of gravity even during typical flowing conditions. We also determine the necessary grid resolution to obtain a physical solution where we observe fingers.

Governing equations
We assume that only one phase is present, and thus all components n c exist in the same phase. The conservation equation for a single component i can be written as where denotes porosity (assumed constant), denotes density of the phase, and X i represents mass-fraction of component i in the phase. The flux i is the sum of the convective and diffusive flux, and is given by For simplicity, the diffusion coefficient D i is taken to be constant for all components. The volumetric flux is the given by Darcy's law for a single-phase: where k denotes permeability, represents viscosity, p is pressure, and is the gravity vector. Note that diffusion is often thought to be negligible, because the diffusion coefficients are small and the diffusive fluxes are generally slow compared to the convective flow.
The mass fraction X i of component i is defined using mole fractions x i as where M i is the molecular weight of component i. Note that the mole fractions x i for all the components i = 1, ⋯ , n c sum up to one, If the mole fractions x i sum up to unity, the mass fractions also sum up to unity. The density and viscosity of the mixture used in this work are dependent of pressure and mole fraction of the CO 2 component x CO 2 , In the well-studied CO 2 -brine system, the brine density increases approximately linearly with CO 2 concentration. For a CO 2 -hydrocarbon system thermodynamic models and available data indicate that the characteristic density curve is non-monotonic, i.e., a maximum mixture density occurs at an intermediate CO 2 concentration (Elenius and Gasda 2021). As a general rule, heavier oil components will exhibit greater non-monotonicity and (Elenius and Gasda 2021) find that, in a system where gravity effects are isolated, convection of CO 2 in oil is dependent on whether CO 2 originates from above or below the oil zone. From above, convection follows classic convective mixing but is accelerated by viscosity decrease with increasing CO 2 . From below, convection flows upward due to CO 2 buoyancy, but is countered by downward convection due to the heavier mixture density. In Elenius and Gasda (2021), the onset time is much earlier and the convection rates for the CO 2 -hydrocarbon systems studied there are 200-600 times greater than for the CO 2 -brine system under similar conditions. In this work, CO 2 is lighter than oil and is injected from the side and it is not known a priori how this flowing system relates to the gravity stable system in Elenius and Gasda (2021).
There is a scarcity of data for CO 2 -hydrocarbon systems at relevant temperatures and pressures. Also, different thermodynamic models will produce different density curves (Span et al. 2015). In  the density and viscosity measurement trends for oil samples from west Texas are reported. It shows that many different oils exhibit non-monotonic behavior when mixed with miscible CO 2 under miscible conditions. In this work, we have studied mixing of the binary system CO 2 and octane, and the non-monotonic density curve produced using the Trend library (Span et al. 2015) and the viscosity curve for CO 2 and octane remains mostly linear, see Figure 2.

Numerical experiments
To study convective mixing at the fine scale under flowing conditions, we have implemented necessary functionality in the existing reservoir simulator framework developed within the Open Porous Media (OPM) initiative (Rasmussen et al. 2020). This initiative (opm-project.org) aims to provide simulation tools for geological subsurface fluid flow problems. The module named opm-models is a numerical simulation framework focused on fully implicit finite volume solvers for the conservation equations of fluid flow in porous media (Lauser et al. 2018).
In this work, a new framework was created to study displacement influenced by convective mixing within a two-component single-phase system. This fluid system is intended to approximate the interaction between a hydrocarbon fluid and CO 2 within a fully miscible setting, i.e., within the mixing zone of a CO 2 flood where residual oil has been mobilized into an oil bank and miscibility has been fully established between the fluids. We, therefore, neglect the additional complexity induced by multiphase flow and multiple-contact miscibility. Given our desire to focus on a simple study of gravity effects in the competition between viscous and gravitational forces, we have limited the focus of this study, but the additional effects will be included in further work.
This section reports several numerical examples that show various settings where CO 2 is injected into a rectangular domain initially filled with octane. If not specified otherwise, the domain is 2.5 × 1 m large, with 2500 cartesian cells equidistant in both spatial dimensions. The rock properties used are porosity = 0.1 , permeability K = 100mD, diffusivity 10 −10 m 2 , and the fluid properties are given in Table 1 and Figure 2. To trigger the development of fingers at the boundary, the numerical diffusion is not sufficient. A small perturbation of the porosity is added, in addition to the diffusion coefficient mentioned. The flow will still be dominated by convection.
Furthermore, we consider reservoir conditions similar to the current Snøhvit CCS project in offshore Norway (approximately 200 bars and 80 • C) (Kaufmann and Skurtveit 2018). CO 2 is injected at the left (west) boundary, the right boundary (east) has free flow boundary conditions, and the top-and bottom-boundaries have no-flow conditions.

Sensitivity to Flow Regime
Under the above boundary conditions, a total of 28 test cases were performed to examine the sensitivity CO 2 migration and oil displacement to whether the regime is gravity-or viscous-dominated, or lies in some intermediate flow regime. We investigate changes in the impact of gravitational forces from zero gravity (angle 0 • with the xy-plane) to full gravitational effect (angle 90 • with the xy-plane).
The angle of the domain relative to the gravitational direction is shown in Figure 3. Increasing viscosity-dominated regime is imposed with increasing injection rate. The  Figs. 4 and 5, the CO 2 -plume is shown 4 and 8 days after injection started with increasing (impact of gravity) angle relative to the y-axis and increasing injection rate (impact of viscosity). The quantity shown in the figures is mole fraction of CO 2 in oil. Octane-rich areas are blue, and dark red means CO 2 rich. The four columns in the figures indicate increasing injection rates from left to right, and the seven rows represent increasing angle with the y-axis and thus affect of gravity from top to bottom. The top most row has zero impact of gravity, that means angle = 0 in Fig. 3. The simulations in rows 2-7 are performed with = 3, 6, 14, 30, 49, 90 • . Thus, the last row represents maximum gravity effect.
For the simulations with no gravitational forces (top panels in Figs. 4 and 5), the migration pattern represents classical viscous fingering throughout the domain. The increased injection rate only magnifies the viscous fingering, as expected. When the domain is tilted with increasing angle, the effect of the gravitational forces grow. For the small angles, we see a thin finger that migrates along the top of the domain, especially in the second row in Fig. 4 that represents early times. As the front advances over a short distance (approximately 1 m) more viscous fingers develop deeper into the domain, and the fingers advance in a manner that resembles the pure viscous case. The difference from the latter is that finger speeds are not uniform when gravity is stronger, but decreases from top to bottom. In Figure 5, we observe that fingers break through faster at the top of the domain than the bottom.
As the effect of the gravitational forces increases, we see that CO 2 convects downward from the top finger, and viscous fingers no longer develop deeper into the domain. Notable in these higher gravity cases is a bottom-most tongue forming at the inlet boundary from   Fig. 4 CO 2 migration 4 days after injection started for increasing injection rates (left to right), and increasing impact from gravitational forces (top to bottom). Dark red is CO 2 -filled, and dark blue is octane-filled early time. This tongue is a result of the dense CO 2 -oil mixture ( x CO 2 = 0.78 ) that sinks initially and migrates along the bottom of the domain, especially for the three largest angles at the latest time in Fig. 5. This dense gravity tongue flows at a slower rate than the light gravity tongue at the top of the domain due to viscosity differences. Eventually the bottom tongue merges with the gravity fingers convecting from above. Figures 4 and 5 show mole fraction of CO 2 in oil. In Fig. 6 the corresponding CO 2 plume is represented by density of the phase for the case with maximum impact of gravity and injection rate 8.12 kg/m 2 /day, 8 days after injection started. The density of the convecting fingers is nearly uniform at the maximum value of 740 kg/m 3 . The density   Fig. 6 CO 2 migration shown with density 8 days after injection started, using maximum impact from gravitational forces and injection rate 8.12 kg/m 2 /day contrast from the CO 2 -rich plume and the fingers is very sharp, pointing to the very strong convection occurring under the migrating CO 2 -rich finger.
To quantify the ratio of gravity to viscous forces, we define a gravity number N g .
where Q m is the CO 2 mass flow rate and is a scaling factor. Note that N g often is multiplied with the characteristic length so that the gravitational forces typically will dominate further away from the injection. Large gravity numbers ( N g >> 1 ) indicate that gravity forces dominate, and small gravity numbers ( N g << 1 ) indicate that viscous forces dominate.
As an intermediate state, where gravitational and viscous forces have similarly weighted effect on the system, we let = 15 • and the injection rate 8.12kg/m 2 / day. For this intermediate state, we scale such that N g = 1 . This gives a scaling factor = 3.7 × 10 −5 in Equation 7.
In Table 2, the gravity numbers N g are given for the same angles and injection rates that we used for the simulations in Fig. 5. These numbers confirm what the results in Fig. 5 told us, namely that N g >> 1 for large angles and low injection rates. Furthermore, N g << 1 for large injection rates and low angles . The values also indicate that an injection rate significantly larger than 12.18 kg/m 2 /day would be required to overcome the gravity forces for this system.
In Fig. 7, the effect of the gravitational forces on the plume migration are compared by the location of the tip of the plume at time 8 days after injection started for different rates and angles . The colour purple corresponds to the largest injection rate and the tip of the CO 2 has migrated the furthest. With increasing impact from the gravitational forces, we see that the location no longer increases, but rather stabilises around a migration distance below 250 cm. For the lowest injection rate shown in blue, we see that when the gravitational forces are more than 75%, the migration distance for the the tip of the plume decreases. This means that for some rates the location decreases when the angle is larger and the gravity dominates more. This trend is not unexpected, since we see from Figure 5 that for larger angels the CO 2 -rich fingers migrate more and more downward, increasing the areal sweep of the domain and changing the time of the breakthrough of CO 2 at the outlet. The time of CO 2 breakthrough at the distance x = 1.25 from the inlet is shown in Fig. 8 for different angels using the injection rate 8.12 kg/m 2 /day. The mass of CO 2 in the cells located at the column at x = 1.25 is summed and given at the y-axis. The case with the second largest angle (75% impact from gravity) has the earliest breakthrough, followed by the case with maximum impact from gravity at 90 • . This agrees well with the result from Figure 7 with the same injection rate, where the location of the tip of the plume is larger for 75% impact of gravity than for 100%. This observation is also visible and explained by the plume migration shown in Fig. 5 for the two bottom rows for the second largest injection rate. For the largest impact from gravity, more CO 2 will migrate downwards in the fingers, and less CO 2 will be available in the tip of the plume. Thus, the breakthrough time for the largest impact from gravity is smaller than for the second largest impact from gravity. Next is the lower angles, with the horizontal case having the latest breakthrough by far. We also  Fig. 8 Mass of CO 2 that has reached the cells located at x = 1.25 up to 16 days after injection started, for five different gravity scenarios with a fixed injection rate of 1.64 kg/m 2 /day observe that for the horizontal case the cells at x = 1.25 are filled completely with CO 2 after approximately 13 days, as expected given a piston-like front. For the cases with larger gravitational forces, we see that this process takes much longer, and that the mass of CO 2 at this location stabilises around some lower values. A significantly longer simulation time would be needed for the higher angle cases to fully sweep the domain.
In Fig. 9, we show a partition of the CO 2 migration around x = 1.25 after 6 and 16 days. The vertical magenta line indicates x = 1.25 , and the black vertical lines to the left and right, indicates x = 1 and x = 1.5 , respectively. The closer examination of this case shows the complexity of the fingers as they evolve and grow with increasing gravitational effects.

Sensitivity to Grid Refinement
Increasing grid refinements capture the mixing of CO 2 and oil differently. To illustrate this, we compare the CO 2 migration for six different grid refinements on 2 × 1 m rectangular domain. The injection rate used is 8.12kg/m 2 /day and the impact of gravity is at maximum ( =90 • ). In Fig. 10, the CO 2 migration is shown 4 days after injection started for increasing refinements of the grid.
We see that fewer number of grid cells the fingers grow more slowly, while the gravity tongue grows more quickly. We observe that 20000 cells is required to capture the fingering for this 2 × 1 m domain. The delay of convection for too coarse grid cells is consistent with previous results in the literature; however, the impact of grid refinement on resolving the dense gravity tongue is not previously quantified. This indicates that for field-scale domains with much larger grid-cells than the domain we are looking at here, this is a very challenging effect to capture. Even for meter-scales like this, we need a very fine grid resolution to capture the effect of the convective mixing correctly. To transfer this knowledge to a more realistic reservoir domain with much larger spatial and temporal scales will be hard. Also, connected with increased heterogeneity and other complicating factors, we can not say how dominating this effect will be. For this reason, this effect needs to be further studied, and the effects needs to be upscaled to a more realistic scale to be able to give information about the effect of this phenomena in a field scale case.

Dissolution Rate
To quantify the mixing, we study the mass of CO 2 convected into oil below some "interface." We let this interface be defined by some threshold of mole fraction of CO 2 in oil. For this study, we have chosen a threshold of 0.95, and stress that the results below are sensitive to this number. In this study, we have used a domain size of 2 × 1 m, and various numbers of cells in the discretisation.
In Fig. 11, the mass of dissolved CO 2 in oil is compared for increasing refinements of the grid. For the coarsest grid (5000 cells), we see that more CO 2 has convected into the oil below the 0.95 mass fraction threshold. The higher values for the coarser grids are due to the overestimation of the sinking gravity tongue at the inlet. Furthermore, for the other refinements, the mass of convected CO 2 is more similar as the time after injection started increases. The corresponding rates in Fig. 11 show that the rate corresponding to the coarsest grid oscillates more than the rates of the finer grids. In all cases, the rate is stable after 1 day.

Fig. 10
Migration of CO 2 shown with mole fraction of CO 2 in oil 4 days after injection started with maximum impact of gravity and injection rate 8.12kg/m 2 /day. From left to right and from the first to the second row, there is an increasing number of discretisation cells, meaning that the grid is more and more refined unfavourable mobility ratio between CO 2 and brine (Calabrese et al. 2005). When small scale heterogenities are considered, viscous and capillary forces dominate the flow and gravity forces are regarded as less important (Li and Benson 2014). Also, for these systems, we know that large-scale heterogeneity results in lower fluid-phase mobilities (Ott et al. 2015). Less is known for the CO 2 -hydrocarbon system, and the observations from the CO 2 -brine is not necessarily applicable here. The effect of heterogeneity in this setting is thus subject of future work. Also, the importance of diffusion in finger initiation and development is subject of future work. It is evident that diffusion is important especially in fractured and heterogeneous media .
The system is modeled in a single-phase, multi-component, framework. This is a simplification, since the hydrocarbon phase can divide in two phases (liquid and vapour) where multiple components can exist. In addition, if water is present, this will add another phase to the system. Fingering in multiphase compositional and compressible flow has received less attention in the literature, and high computational complexity is the reasoning for this. In multiphase flow, relative permeabilities affect the mobility contrast for a given viscosity ratio and phase behavior can also change local fluid properties, which can enhance or mitigate viscous and gravitational instabilities. This is subject of active research, where convective mixing will be exploited in a multiphase compositional framework within OPM.
Another simplification done here, is the reduced dimensions of the computations. Gravitational and viscous fingers are 3D naturally, and a 2D representation of these might lead to inaccurate results. In the literature, it is well known that extra dimensions result in higher numerical dispersion (Shahraeeni et al. 2015). For the 2D representation, it is assumed symmetry into the y-plane, and this assumption will not be applicable in many cases. The effect of the CO 2 plume will be different in 3D, and we expect that the front of the CO 2 plume will also be consumed by fingers in the third dimension, leading to increased retardation of CO 2 breakthrough at the outlet. To evaluate this, full 3D simulations need to be performed, and this is subject of future work.
We also note that the phenomena we are investigating in this work occurs on a very finescale. The results found in this work gives insights into gravitational and viscous fingers at the meter-scale, and the next step is to investigate the effect of convection at more realistic field scale. The works (Ahmed et al. 2012;Both et al. 2015;Shahraeeni et al. 2015) demonstrate that gravity fingers also develop quickly in both homogeneous and heterogeneous media under realistic reservoir injection conditions. The numerical representation of complex physical processes like viscous and gravitational fingering, where reducing numerical dispersion is important, has been discussed by several authors, i.e., (Shahraeeni et al. 2015). To capture these effects with lower order methods, very small cells-sizes should be used. The simple example in Section 3.2 with decreasing cell sizes emphasizes exactly this.

Summary and conclusions
In this work, we simulate gravity effects, namely gravity override and convective mixing, during miscible displacement of CO 2 and oil. The flow behavior due to the competition between viscous and gravity effects is complex due to non-monotonic density off the mixed fluids, and can only be accurately simulated with a very fine grid. We demonstrate and find that under a range of conditions: • convection occurs rapidly, and has a strong effect on breakthrough of CO 2 at the outlet, • CO 2 with convection can enhance the vertical sweep of the reservoir and lead to more efficient displacement process.
Especially, these effects are quantified for a simple system under realistic conditions. Increased efficiency is positive for increasing storage of CO 2 behind the oil front and maximising the ratio of CO 2 stored per oil produced. We emphasise that the knowledge about the effect of convective mixing on the field scale is limited. This is something the research community needs to focus more on, and the work done herein is a starting point to decrease this knowledge gap.