Microfluidic determination of minimum miscibility pressure (MMP) in dynamic CO₂/n-decane flow

Abstract

Carbon dioxide enhanced oil recovery (CO₂-EOR) has been recognized as a viable pathway for carbon capture, utilization, and storage (CCUS). Among its variants, miscible CO₂-EOR offers a considerable additional oil recovery of approximately 5–20%, making the determination of minimum miscibility pressure (MMP) a critical design consideration. In this study, we employ a high-pressure microfluidic platform to investigate the miscibility transition between CO₂ and n-decane at temperatures (T) of 40, 50, 70, and 90 °C. At T = 40 °C, with increasing pressure (P), microfluidic visualization reveals a series of distinct flow regimes: dripping, quasi-steady jetting, unsteady jetting, transitional, and ultimately diffusive regimes. In the diffusive regime, miscibility is achieved through intensive mixing, leading to the disappearance of the fluid–fluid interface. Based on these microfluidic observations, we propose a new criterion for MMP determination: the minimum pressure required to reach the diffusive regime for the dynamic CO₂–oil flow. The experimentally determined MMP values show good agreement with previous microfluidic studies and predictions from the Peng–Robinson equation of state (PR-EOS). Furthermore, the MMP increases linearly with temperature from 40 to 90 °C, consistent with the reduced solubility of CO₂ in n-decane at higher temperatures. This microfluidic method provides a rapid and visual approach to assess miscibility transitions in CO₂-EOR applications.

---

1. Introduction

Carbon dioxide enhanced oil recovery (CO₂-EOR) has been widely implemented as a large-scale and cost-effective method under carbon capture, utilization, and storage (CCUS) strategies to mitigate rising greenhouse gas emissions. In CO₂-EOR, captured CO₂ is injected into existing petroleum reservoirs to enhance oil recovery—a process that not only enables low-cost production but also reduces the net carbon footprint of the extracted oil.^1–4 Since the 1990s, the number of CO₂-EOR projects has steadily increased and is projected to continue growing through 2040.⁵ Various CO₂-EOR approaches exist, including immiscible injection, miscible injection, CO₂ foam-assisted EOR, and the huff-n-puff method.^4,6–8 Compared with immiscible CO₂-EOR, miscible and near-miscible operations have shown particular promise,^9–12 as they are more effective in improving sweep efficiency and the overall recovery factor.¹¹ Miscibility refers to the ability of two fluids to mix uniformly and form a single-phase mixture.¹³ In CO₂-EOR, achieving miscibility between injected CO₂ and in situ oil can significantly improve oil recovery through mechanisms such as oil swelling and viscosity reduction,^9,14,15 vaporizing and condensing gas drive processes,^16,17 and the mitigation of unfavorable flow phenomena including gravity override^18–21 and viscous fingering,^22–28 which commonly reduce recovery efficiency in conventional gas injection methods.²⁹

In miscible CO₂-EOR, the minimum miscibility pressure (MMP) is one of the most critical parameters, defining the pressure above which CO₂ and oil become fully miscible.¹³ MMP is influenced by several factors, such as temperature, oil composition, and injected gas composition.³⁰ Two commonly used definitions for miscibility are first-contact and multi-contact miscibility.^13,30,31 First-contact miscibility refers to a thermodynamic condition in which two fluids form a single-phase mixture upon contact, regardless of their mixing ratio.¹³ In contrast, multi-contact miscibility is more commonly encountered in practical CO₂-EOR operations, as reservoir oils typically contain multiple components, including light, intermediate, and heavy fractions.^32,33 Multi-contact miscibility is achieved through repeated contact and component exchange between the injected gas and reservoir oil, involving both vaporization and condensation mechanisms.^13,31

Several long-established, macroscopic experimental techniques have been developed since the 1980s to determine MMP, including the vanishing interfacial tension (VIT) method,^34,35 rising bubble apparatus (RBA),³⁶ and slim tube test (STT).^32,37 The VIT method involves measuring interfacial tension, typically using a pendant drop technique, to construct an IFT–pressure curve and extrapolate the pressure at which IFT approaches zero, which is taken as the MMP.^34,35 However, this method is sensitive to fluid composition and may result in deviations in MMP estimation.^31,38 The RBA involves injecting a gas bubble into a column of oil and observing its shape and rising behavior, which offers qualitative insight into miscibility conditions.³⁶ While relatively simple and fast to perform, RBA may not reliably capture the vaporizing and condensing mechanisms involved in multi-contact miscibility, due to the short interaction time and lack of continuous phase exchange between gas and oil.^31,39 The STT uses a coiled tube packed with porous media, such as sand or glass beads, to simulate reservoir flow conditions.^32,37 CO₂ is injected at varying pressures, and oil recovery is measured; MMP is operationally identified as the pressure at which the recovery factor transitions from a rapid increase to a slower rate of change, which typically occurs at a recovery of ≈90–95%.^13,32,37 The STT is widely considered the most accurate laboratory method for representing multi-contact miscibility in CO₂-EOR.^30,31,33 However, it is time-consuming (often requiring 2–6 weeks to complete), lacks standardized design in tube dimensions and packing material, and can be costly to perform.^30,31,33 In addition to the above-mentioned experimental techniques, MMP can also be estimated using empirical correlations^40–42 and equation-of-state- (EOS-) based methods.^43,44

In contrast, more recent microfluidic approaches offer complementary capabilities, including reduced sample volumes, faster measurements, and direct visualization at the pore scale. Microfluidic platforms, hence, have recently emerged as versatile tools for studying fluid displacement,^{12,25,45–54} mass transport,^55–66 and phase behavior^67–72 at the pore scale. When operated under high-pressure conditions, these systems can mimic the temperature and pressure environments of subsurface reservoirs, making them particularly useful for investigating CO₂–oil MMP.^73–82

Microfluidic systems, when combined with fluorescence microscopy, allow direct visualization of phase behavior at the CO₂–oil interface, with the oil phase typically labeled using a fluorescent dye. For example, Nguyen et al.⁷³ were among the first to examine MMP using a microfluidic T-junction, where CO₂ and oil were co-injected, and miscibility was evaluated by monitoring the fluorescence intensity at a fixed detection region.⁷³ Below MMP, the presence of segmented CO₂ bubbles led to periodic drops in fluorescence intensity; above MMP, a more uniform signal indicated complete mixing.⁷³ Microfluidic analogs of the vanishing interface technique (VIT) have also been developed, where oil is trapped in a dead-end microchannel and the average fluorescence intensity is monitored to determine whether a homogeneous mixture is formed.^74–78 Similarly, the “slim-tube-on-a-chip” design mimics the concept of conventional slim-tube tests by injecting CO₂ into oil-saturated microfluidic porous media, where multi-contact miscibility is evaluated based on the recovery factor.^79–82 Unlike conventional slim tubes, microfluidic platforms enable faster testing (within hours instead of weeks) and allow precise control over pore size, geometry, and permeability, offering more consistent and tunable experimental conditions.

Compared to conventional laboratory methods, microfluidic approaches offer several remarkable advantages for studying miscibility. These include low fluid consumption, real-time optical access, rapid testing, and precise control over the porous medium geometry and flow conditions.^83–87 Microfluidic platforms are particularly well suited for studying miscibility because their characteristic length scale, typically ranging from 100 nm to 100 μm,⁸⁸ is comparable to the pore scale of subsurface porous media.⁸⁴ At this scale, the reduced channel dimensions result in a high surface-to-volume ratio, which enhances interfacial contact and promotes mass transfer,⁸⁹ thereby making microfluidic systems more representative of multi-contact miscibility development in porous reservoirs. Moreover, recent microfluidic studies have extended into the nanofluidic regime, where channel depths on the order of 10 to 100 nm introduce strong molecular–wall interactions and confinement effects.^75,78 The measurement of MMP under such nanoconfinement is particularly relevant for CO₂-EOR in tight formations such as shale. For example, parachor-EOS modeling based on interfacial tension predictions suggests that reducing the channel depth from 1 μm to 10 nm can lower the MMP of CO₂–octane from 15.10 to 10.10 MPa at 160 °C.⁷⁸

Building on the advantages of microfluidics, this work introduces a new, dynamic, and visually intuitive method for determining MMP based on flow regime transitions in continuous flow. While previous microfluidic approaches using static geometries such as dead-end channels have provided valuable insights into equilibrium phase behavior,^74–76,78 our work introduces an original approach and a unique perspective that captures the transient evolution of miscibility under dynamic, flowing conditions. We develop a high-pressure microfluidic platform capable of operating under reservoir-like conditions (15 MPa, 90 °C) and use it to visualize interfacial dynamics in CO₂–n-decane systems across a range of temperatures and pressures. By identifying the transition from dripping, jetting, and transitional flows to a diffusive regime, we define MMP through direct visualization of miscibility with high precision. This approach enables real-time assessment of phase behavior and provides new insight into the dynamics of interfacial deformation and mixing near the miscibility threshold. The experimental MMP values are compared with predictions from the Peng–Robinson equation of state (PR-EOS) and previous microfluidic experimental data, confirming the method's accuracy, reproducibility, and potential for broader applications in energy and chemical systems.

2. Experiments

This study investigates the transition to miscibility between CO₂ and n-decane from 40 to 90 °C, where miscibility typically occurs at pressures on the order of 1 to 10¹ MPa. A robust microfluidic setup was developed to operate with current gas supply up to 15 MPa and 90 °C, enabling operation under reservoir-relevant conditions.

2.1. Materials

The materials used included n-decane (≥99%, CAS No. 124-18-5, Thermo Fisher Scientific), Nile red dye (CAS No. 7385-67-3, Sigma-Aldrich), and CO₂ (99.99% purity, 4.0 anaerobic grade, Linde). For microfabrication, prime-grade silicon (Si) wafers and round borofloat glass substrates, each with a diameter of 100 mm, were used.

2.2. Microfluidic chip fabrication

The microfluidic chip was fabricated using silicon and glass, materials commonly employed in high-pressure microfluidics.^90,91 On the front side, the silicon wafer was first cleaned with a piranha solution and primed with hexamethyldisilazane (HMDS) vapor to enhance photoresist adhesion (YES Oven). A layer of photoresist (AZ 1512) was spin-coated and patterned using direct-write photolithography with a maskless aligner (Heidelberg MLA150). Deep reactive ion etching (DRIE) was then performed using a high-rate single-stage Bosch recipe (Oxford Estrelas) to etch the microchannels to a depth of 50 μm.

Following the two-sided lithography method adapted from Qi et al.,⁹¹ the fluidic ports were etched through the remaining 450 μm-thick wafer from the backside using a bilayer of AZ4330 photoresist. This double-sided approach was adopted to eliminate the need for mechanical drilling, which previously caused silicon cracking under elevated pressure.

Residual photoresist was removed with a photoresist asher (PVA TePla Ion 10Q) followed by piranha cleaning. The patterned silicon wafer was then anodically bonded to a borofloat glass wafer at 370 °C with 600 V using an AML wafer bonder, forming a sealed chip. Finally, the bonded wafer stack was diced into individual chips using a dicing saw (Disco 3240). An image of the final silicon–glass microfluidic chip used in the experiments is shown in Fig. 1A, with overall chip dimensions of 25 × 15 mm. A microscopy view of the microfluidic design is shown in Fig. 1B, illustrating the CO₂–decane flow within the serpentine microchannel.

Fig. 1 (A) Photograph of the fabricated silicon–glass microfluidic chip. (B) Microscopy image showing the CO₂–decane flow in a serpentine microchannel. Scale bar: 200 μm. (C) Schematics of the microfluidic experimental setup. CO₂ and n-decane were injected using two separate high-pressure syringe pumps. Center: 3D model of the stainless steel microfluidic chip holder with a sapphire viewing window and the microfluidic chip positioned beneath it. Backpressure was controlled using a backpressure regulator with an electronic pilot controller.

2.3. Experimental setup

A stainless steel microfluidic chip holder was designed to securely house the microfluidic chip, as shown in Fig. 1C. The chip holder was fabricated using a CNC mill (Tormach PCNC 1100) and a milling–drilling machine (Optimum BF 30 Vario). Internal fluidic channels were precision-drilled to the interface with external 1/16″ stainless steel tubing. The chip was compressed and secured with four bolts, which also clamped a sapphire viewing window between the chip and the holder for real-time microscopy observations. The modular holder design enables replacement of the microfluidic chip with alternative layouts for other high-pressure experiments.

As illustrated in Fig. 1C, a high-pressure syringe pump (Teledyne ISCO 100DX) was used for CO₂ injection and pressure control. A pressure transducer (OMEGA PX309) was installed upstream of the CO₂ microfluidic inlet to provide accurate monitoring of the CO₂ inlet pressure.

For n-decane injection, a high-pressure syringe pump (Chemyx Fusion 6000X) was used with a 6 mL stainless steel syringe and an inline pressure gauge to verify that the target pressure was reached. Downstream pressure was regulated using a back pressure regulator (Equilibar ZF Series) with an electronic pilot controller (Equilibar EPR-3000), which maintained a stable back-pressure setpoint.

Temperature control was achieved using a temperature controller (OMEGA CN740). An aluminum heating plate was mounted on the stainless steel chip holder and secured with Kapton tape. A K-type thermocouple was inserted into a hole drilled 1 mm away from the microfluidic channel to provide accurate temperature readings. Thermal paste was applied beneath the heating plate and inside the thermocouple hole to ensure efficient thermal conductivity.

2.4. Experimental procedures

We investigated the miscibility transition for CO₂ and n-decane at temperatures of T = 40, 50, 70, and 90 °C. For experiments at 40 °C, an initial coarse pressure ramp was performed from 1 to 10 MPa in 1 MPa increments to observe flow behavior across a broad pressure range. Subsequently, a finer pressure ramp was performed from 7.5 to 8.5 MPa in 0.1 MPa increments to determine MMP more precisely. For temperatures of 50, 70, and 90 °C, only the fine pressure ramp (0.1 MPa increments) was performed to precisely determine the MMP values. For each pressure and temperature condition, the experiment was repeated at least three times to ensure reproducibility.

To enhance contrast between the oil and CO₂ phases, n-decane was mixed with Nile red at a concentration of 10 g m⁻³. The solution was stirred with a magnetic stirrer until fully homogenized, and then degassed in a vacuum desiccator for 45 minutes until no visible gas bubbles remained.

Prior to each experiment, the target temperature and back pressure were set using the temperature controller and back pressure regulator, respectively. CO₂ was gradually pressurized to the desired testing pressure using a high-pressure syringe pump and introduced into the microfluidic system. Simultaneously, the oil phase (n-decane) was also pressurized to the same target pressure with its valve kept closed. Once both pressure and temperature reached their targets, the decane valve was opened, and both CO₂ and decane were co-injected into the microfluidic chip. Both injection pumps were set to maintain a constant flow rate of 50 μL min⁻¹ throughout the experiment. The system was allowed to stabilize for at least 5 minutes before recording began.

After the experiments were completed at a given pressure, both CO₂ and n-decane lines were closed. The back pressure was then adjusted for the next setpoint, and each fluid was re-pressurized separately. Once stabilized at the new pressure, the fluids were reintroduced at the same flow rate of Q̇ = 50 μL min⁻¹.

2.5. Imaging acquisition and data analysis

The microfluidic chip holder was placed on an inverted microscope with 5× or 10× objectives. A mercury lamp (HXP 200C) was used for illumination. A high-speed camera was attached to the microscope side port to capture image sequences at 7500 frames per second (fps).

Initial image processing was conducted in ImageJ,⁹² including cropping the region of interest (ROI), applying a bandpass filter, and subtracting the background to improve illumination uniformity and reduce noise. Data processing, including equation of state calculations and figure generation, was carried out using MATLAB.

3. Results and discussion

3.1. Microfluidic flow regimes and visualization

We performed experiments by simultaneously injecting CO₂ and n-decane at a constant flow rate of 50 μL min⁻¹ into a microfluidic flow-focusing device, where CO₂ was the dispersed phase and n-decane was the continuous phase. A schematic of the microchannel layout is shown in Fig. 2A.

Fig. 2 (A) Microchannel schematic showing section 1 and section 2 as the regions of interest under a microscope with a 10× objective. CO₂ (dispersed phase) and n-decane (continuous phase) were co-injected at a constant flow rate of 50 μL min⁻¹ at the flow-focusing junction. Section 1 captures the junction where the CO₂–decane interface forms, while section 2 is located downstream (8.6 mm from the junction). (B) Phase diagram summarizing all the experimental conditions, showing flow regime transitions at T = 40, 50, 70, and 90 °C. (C) Representative snapshots at T = 40 °C for each flow regime, captured from section 1 (left) and section 2 (right) as indicated in (A). Experimental conditions are noted; scale bars: 200 μm.

At T = 40 °C, increasing the pressure from 1 to 10 MPa resulted in a series of distinct flow regime transitions, as illustrated by the data points in Fig. 2B and the corresponding snapshots in Fig. 2C. The observed flow regimes include dripping, quasi-steady jetting, unsteady jetting, transitional, and diffusive regimes. A detailed discussion of each regime is provided in sections 3.1.1 to 3.1.4.

For T = 50, 70, and 90 °C, finer pressure steps of 0.1 MPa were applied over a narrower pressure range to capture the miscibility transition. As shown in Fig. 2B, similar flow regime transitions were observed, from unsteady jetting to transitional, and ultimately to diffusive regimes.

3.1.1. Dripping regime. At T = 40 °C and P = 1–4 MPa, a dripping flow regime was observed, in which CO₂ droplets broke up immediately after the flow-focusing junction,⁹⁵ as shown in Fig. 2C. The two fluids were immiscible with a sharp and well-defined interface. The droplets were initially elongated and formed at regular intervals with uniform size.⁹⁶ As they traveled downstream into section 2 of the channel, the droplet size gradually decreased due to CO₂ dissolution into the surrounding decane.

3.1.2. Quasi-steady jetting regime. At T = 40 °C and P = 5–6 MPa, the flow exhibited a quasi-steady jetting regime. Although droplet formation remained periodic, breakup occurred further downstream and at higher frequency than in the dripping regime. The resulting droplets were smaller and more closely spaced, as shown in Fig. 2C (40 °C, 5.1 MPa).

The transition from dripping to jetting can be triggered when either inertial forces from the dispersed phase or viscous shear from the continuous phase become large enough to overcome interfacial tension.^95,97 Each mechanism is represented by a corresponding dimensionless number:

(1) Inertia-surface tension balance: when the Weber number of the dispersed phase (We_CO2) exceeds a critical value, the inertial force of the moving CO₂ thread overcomes the restoring force of surface tension, leading to the transition from dripping to jetting.

where ρ_CO2 is the CO₂ density extracted from the NIST database,⁹³ u = Q̇/A is the superficial velocity (i.e., the volumetric flow rate Q̇ divided by the microchannel cross-sectional area A), D_h = 2HW/(H + W) is the hydraulic diameter, and γ is the interfacial tension (IFT) between CO₂ and n-decane, estimated using the empirical model developed by Pan and Trusler.⁹⁴

(2) Viscous-surface tension balance: when the capillary number of the continuous phase (Ca_n-decane) exceeds a critical value, the viscous shear exerted by the flowing n-decane balances and overcomes surface tension, also triggering the transition to jetting.

where μ_n-decane is the viscosity of decane from the NIST database.⁹³

Both We_CO2 and Ca_n-decane increase with pressure (Fig. 3A and B). The increase in We is particularly pronounced during the transition from dripping to jetting, reflecting the growing inertia of the inner CO₂ stream as its density increases with pressure due to compression.

Fig. 3 Variation of dimensionless numbers and fluid properties with pressure at T = 40 °C. (A) Weber number of the CO₂ phase, calculated using eqn (1), showing its pressure dependence. (B) Capillary number of the n-decane phase, defined by eqn (2), as a function of pressure. (C) Reynolds number of the CO₂ phase, defined by eqn (3), as a function of pressure. (D) Pressure-dependent variation of fluid properties: density ratio ρ_CO2/ρ_n-decane (blue □, left axis),⁹³ viscosity ratio μ_CO2/μ_n-decane (red ○, right axis),⁹³ and IFT (inset, calculated using the empirical correlation from ref. 94 based on pendant drop measurements). The shaded region marks the pressure range where IFT approaches zero; We and Ca are not computed in this regime.

The change in flow behavior can be further explained by the pressure-dependent variation in fluid properties, illustrated in Fig. 3D. Between 1 and 5 MPa, the CO₂ density increases steadily, while the interfacial tension with decane decreases, as seen in the inset of Fig. 3D. This combined effect of higher density and lower interfacial tension causes the Weber number to increase, due to increasing inertial force and weaker surface tension. The steady rise in the Weber number promotes the transition from dripping to quasi-steady jetting.

3.1.3. Unsteady jetting regime. At T = 40 °C and P = 7–7.8 MPa, an unsteady jetting regime was observed. This regime was characterized by a large variation in droplet size and a breakup point located further downstream compared to the quasi-steady jetting regime, as shown in Fig. 2C (40 °C, 7.1 MPa).

In this pressure range, the Weber number continues to rise due to increasing CO₂ density and decreasing interfacial tension. The capillary number also increases more significantly as interfacial tension drops further, as seen in Fig. 3B. The increasing trend in We and Ca indicates an increasing influence of inertial and viscous forces relative to surface tension, leading to delayed droplet pinch-off. The Reynolds number of the CO₂ stream, defined as

Re_CO2 = ρ_CO2uD_h/μ_CO2, (3)

where μ_CO2 is the CO₂ viscosity,⁹³ increases nearly linearly with pressure (Fig. 3C). The rising Re reflects a shift toward inertia-dominated flow, contributing to jet instability and variation in droplet size.

To better visualize the jet-tip dynamics at elevated pressures, time-lapse snapshots were captured at 40 °C and 7.1 MPa to highlight the instability at the jet tip. As shown in Fig. 4A (left), frames acquired every 0.4 ms reveal a weak and unstable jet tip, with inefficient pinch-off and variation in droplet size.

Fig. 4 (A) A series of time-lapse snapshots at T = 40 °C, showing the evolution of the CO₂ jet tip at 7.1 MPa (unsteady jetting), 8.1 MPa (transitional), and 8.6 MPa (diffusive). Images were taken every 0.4 ms from t = 0 to t = 1.6 ms. (B) Representative images of the meander region, highlighting the jet front where droplet formation and mixing occur. In the transitional regime (8.1 MPa), the jet front appears wavy and unstable, producing many satellite droplets. In contrast, the diffusive regime (8.6 MPa) exhibits fewer droplets and enhanced mixing near the bend, indicative of miscible flow behavior. Scale bars: 200 μm.

Similar unsteady jetting behavior was consistently observed at 50, 70, and 90 °C (Fig. 5–7 (left)). At elevated temperatures, the onset of unsteady jetting occurred at correspondingly higher pressures, consistent with reduced interfacial tension and changes in fluid properties.

Fig. 5 Time-lapse snapshots at T = 50 °C showing the evolution of the CO₂ jet tip at 9.1 MPa (unsteady jetting), 9.6 MPa (transitional), and 10.0 MPa (diffusive). Scale bars: 200 μm.

Fig. 6 Time-lapse snapshots at T = 70 °C showing the evolution of the CO₂ jet tip at 11.6 MPa (unsteady jetting), 12.5 MPa (transitional), and 12.8 MPa (diffusive). Scale bars: 200 μm.

Fig. 7 Time-lapse snapshots at T = 90 °C showing the evolution of the CO₂ jet tip at 13.7 MPa (unsteady jetting), 14.7 MPa (transitional), and 15.1 MPa (diffusive). Scale bars: 200 μm.

3.1.4. Transitional regime. At T = 40 °C and P = 7.9–8.4 MPa, the flow entered a transitional regime. As seen in Fig. 4A (middle), the CO₂ jet developed wavy perturbations and frequently generated small satellite droplets near the jet tip. Tails occasionally formed at the rear of the droplets, but over time the droplets became rounded due to the surface tension effect minimizing surface area. A representative image near the first meander is shown in Fig. 4B (left). The wavy interface extended through the meander and eventually broke into droplets downstream. Although the CO₂–decane interface became less distinct in the meander section, the jet retained its overall structure and continued to form droplets after exiting the curved region.

In the pressure range corresponding to transitional flow, the inset in Fig. 3D shows IFT values of approximately 1 to 2 mN m⁻¹, a range commonly regarded as near-miscible in EOR applications.⁹⁸ Based on both the IFT values and observed flow characteristics, this regime is interpreted as near-miscible.

Transitional flow was also observed at elevated temperatures of 50, 70, and 90 °C. Compared to the transitional regime at 40 °C, droplets at higher temperatures (50–90 °C) generally appear smaller, although no consistent size trend is evident due to the unstable nature of the flow. As T increased, the CO₂ jet also appeared thinner and more laminar, with breakup occurring further downstream – sometimes even near the channel exit. Although the flow pattern appears different from that at 40 °C, it is still considered transitional, since a visible interface was present and discrete droplets eventually formed, consistent with low-IFT conditions.

These changes with increasing temperature can be explained by variations in fluid properties. Under the pressure and temperature conditions corresponding to the snapshots in Fig. 4A (middle) and 7 (middle), the viscosity of decane decreased by 37.5%, while the density of CO₂ increased by 24.1%. The reduction in viscous drag from the outer decane stream, combined with the increased inertial force from the inner CO₂ stream, led to a more extended CO₂ jet and a delayed breakup point further downstream.

3.1.5. Diffusive regime. At pressures exceeding 8.5 MPa, a diffusive regime was observed. In Fig. 2C (40 °C, 8.6 MPa), near the flow-focusing junction, the interface between CO₂ and n-decane appeared blurred. By the time the flow reached section 2 of the microchannel, the interface was no longer distinguishable, indicating complete mixing and the onset of full miscibility. Section 2 was selected as a practical observation point far enough downstream for the intense interface mixing to develop and complete, yet upstream of the regime where complete drop dissolution could be misinterpreted for miscibility.

The time-lapse images at 40 °C and 8.5 MPa captured the mixing evolution, as shown in Fig. 4B (right). As the CO₂ jet passed through the meandering region, streak-like patterns emerged and the jet appeared stretched and filamented. This pattern suggested intensive mixing with only a few satellite droplets. Further downstream, the interface disappeared entirely, suggesting that the two fluids transitioned into a homogeneous phase.

The miscibility was not achieved immediately upon contact near the junction. This observation is consistent with the findings of Cubaud and Notaro,⁹⁹ where miscible threads exhibited different flow behaviors depending on the Péclet number (Pe), defined as:

where Q̇₁ and Q̇₂ are the flow rates of the two fluids, D is the diffusivity, and h is the channel height.⁹⁹ The Péclet number represents the ratio of advective to diffusive transport. In our study, Pe ≈ 6 × 10³ at 40 °C, which exceeds the ultra-diffusive threshold of Pe_D = 1.5 × 10³ reported by Cubaud and Notaro.⁹⁹ The higher Pe suggests that the diffusive rate was insufficient to produce an ultra-diffusive interface immediately at the junction. Instead, we observed a gradual development of miscibility, characterized by a transient mixing zone that evolved as the fluids flowed downstream.

At 50 and 70 °C, a similar diffusive flow was observed. At 90 °C, the mixing dynamics became even more pronounced. As shown in Fig. 7 (right), chaotic advection led to tangled streak patterns and rapid interface breakdown.

In this study, the onset of the diffusive regime is used to define the minimum miscibility pressure (MMP). Microfluidic visualization shows that mixing does not begin immediately after the flow-focusing junction due to limited diffusion. A diffusive thread gradually forms downstream, with a blurred interface and intensified mixing enhanced by recirculation in the meandering channel. The appearance of this flow structure marks the transition to a single phase and provides a clear visual indicator of miscibility.

In addition, we conducted a supplementary grayscale intensity analysis at section 2 to compare the unsteady, transitional, and diffusive regimes (Fig. 8). Fig. 8B presents the grayscale profiles extracted from the corresponding experimental snapshots in Fig. 8A. Grayscale values were normalized to the range of 0–1, averaged along the vertical direction, and plotted against distance along the flow direction. As shown in Fig. 8B, sharp drops in intensity indicate the locations of the dark CO₂–oil interfaces—pronounced in the unsteady jetting regime, less frequent in the transitional regime, and absent in the diffusive regime, where the grayscale profile becomes flat, resulting in a more uniform gray color. Fig. 8C shows the histogram analysis of the same ROI. Compared with the unsteady jetting and transitional regimes, the probability of darker intensities (approximately 0.2–0.4) drops by orders of magnitude in the diffusive regime. This change reflects the disappearance of distinct interfaces as the pressure exceeds the MMP. Both methods—direct visualization of the diffusive regimes and grayscale intensity analysis—yield consistent MMP results. While MMP values in this work can be conveniently identified through direct visualization, the quantitative grayscale intensity analysis (e.g., Fig. 8B and C) also serves as an alternative approach.

Fig. 8 (A) Experimental snapshots at section 2 for unsteady jetting (40 °C, 7.1 MPa), transitional (40 °C, 8.4 MPa) and diffusive regimes (40 °C, 8.5 MPa). (B) Grayscale profiles extracted from the corresponding snapshots in (A). Grayscale values were normalized to the range of 0–1. A smaller grayscale value corresponds to a darker region in the snapshot; for example, CO₂–oil interfaces appear as valleys in the profile. For each profile, the grayscale values are averaged along the vertical direction of the region of interest (ROI) and plotted against distance along the flow direction. (C) Histogram calculated from the same ROI. The probability for each grayscale value was obtained by dividing the pixel count at that value by the total count within the ROI. The histogram therefore represents the probability distribution of normalized grayscale values within the ROI.

To evaluate whether the visual criterion for miscibility is influenced by flow conditions, we conducted additional experiments at 40 °C using different combinations of CO₂ and n-decane flow rates, ranging from 25 to 100 μL min⁻¹, at pressures between 8 and 9 MPa. While varying the flow rate altered specific aspects of the two-phase flow morphology—including the droplet size, breakup location, and thread width—the pressure at which the flow transitioned from transitional to fully diffusive regimes remained unchanged. For example, in the diffusive regime, reducing the CO₂ flow rate from 100 μL min⁻¹ to 25 μL min⁻¹ decreased the initial thread width from 60.5 μm to 36.0 μm and changed the interface curvature near the inlet. These results are consistent with previous observations of Re and Pe dependence of the viscous thread morphology of miscible flows in a flow-focusing microchannel.⁹⁹

The agreement of the regime transition pressure for different flow rates confirms that full miscibility is a thermodynamic effect, driven by a vanishing interfacial tension at the MMP, rather than by hydrodynamic effects. Therefore, the fixed flow rate of 50 μL min⁻¹ used in the main experiments provides a valid basis for determining MMP, while further studies could explore flow-rate effects, particularly changes in the flow morphology, under other regimes such as dripping or jetting.

3.2. Minimum miscibility pressure (MMP)

We compared the minimum miscibility pressure (MMP) determined from our microfluidic experiments with predictions from a thermodynamic model based on the Peng–Robinson equation of state (PR-EOS).¹⁰⁰ The PR-EOS is widely used in the petroleum industry for modeling the phase behavior and thermophysical properties of CO₂–hydrocarbon mixtures, offering a practical balance between simplicity and accuracy.^33,101,102 The following section provides a brief overview and computational workflow of the PR-EOS approach, while detailed equations and implementation steps are presented in Appendix A.

3.2.1. Thermodynamic modeling using PR-EOS. We developed a custom MATLAB implementation of the Peng–Robinson EOS for analyzing the thermodynamic phase behavior of the CO₂–decane system at experimental temperatures (T = 40–90 °C), for which published data are limited.^94,103,104 The model is generalizable and can be adapted to other CO₂–hydrocarbon mixtures.

PR-EOS is a thermodynamic model using a cubic equation of state that describes the pressure–volume–temperature (P–V–T) relationships for real fluids, accounting for intermolecular forces and molecular volume.^100,102,105 It defines pressure (P) as a function of molar volume (v) and temperature (T):¹⁰⁰

where R = 8.314 J mol⁻¹ K⁻¹ is the gas constant. The coefficients a (attraction parameter) and b (van der Waals covolume) are functions of the critical properties.¹⁰⁰ For multi-component systems, the parameters a and b are determined using the classical van der Waals mixing rules as outlined in Appendix A. Eqn (5) can be reformulated into a cubic equation in compressibility factor Z = Pv/RT (details in Appendix A).

In our study, PR-EOS was used to compute the pressure–composition (P–x) phase envelope for the CO₂–decane system. The envelope is bounded by two curves: the bubble point curve (open ○ in Fig. 9A), constructed by varying the CO₂ mole fraction in the liquid phase (x₁), which marks the pressure at which the first bubble of vapor forms from a liquid mixture, and the dew point curve (filled ◆ in Fig. 9A), constructed by varying the CO₂ mole fraction in the vapor phase (y₁), which marks the pressure at which the first drop of liquid condenses from a vapor mixture.¹³ These two curves are computed by satisfying vapor–liquid equilibrium,¹³ which ensures that there is no net mass transfer between the phases, i.e., thermodynamic equilibrium is achieved.¹⁰⁶

Fig. 9 Prediction of the minimum miscibility pressure (MMP) using the Peng–Robinson equation of state (PR-EOS). (A) Phase envelopes of CO₂/n-decane mixtures at T = 40, 50, 70, and 90 °C, computed using PR-EOS. The phase envelope consists of the bubble point curve (open ○), which plots the bubble point pressure P_bubble against the CO₂ mole fraction in the liquid phase (x₁), and the dew point curve (filled ◆), which plots the dew point pressure P_dew against the CO₂ mole fraction in the vapor phase (y₁). MMP is determined at the highest point of P on the phase envelope, above which the two components form a fully miscible, single-phase mixture at all proportions. (B) Flow chart illustrating the bubble point calculation procedure. CO₂ is treated as component 1 and n-decane as component 2. For each specific CO₂ liquid mole fraction x₁, the corresponding bubble point pressure P_bubble and vapor mole fraction y₁ are obtained using a Newton–Raphson method based on PR-EOS with the van der Waals mixing rule. The iteration converges when the fugacity residual between liquid and vapor phases, ‖F = f_L − f_V‖, falls below a tolerance of 10⁻⁶, satisfying vapor–liquid equilibrium (VLE). By the same token, dew point pressures (P_dew) are calculated analogously by fixing y₁ and solving for x₁ and P_dew. Each converged solution contributes to the construction of the phase envelopes in Fig. 9A.

Together, the bubble point and dew point curves form the phase envelope. In Fig. 9A, they appear as a single continuous boundary (liquid or vapor) enclosing the two-phase region, where vapor and liquid coexist. Outside this region, only a single phase exists. The peak of the envelope corresponds to the pressure above which the mixture becomes single-phase at all compositions, also known as the first-contact minimum miscibility pressure (FC-MMP).³¹

The workflow in Fig. 9B outlines the bubble point calculation procedure. Following common convention, we define x₁ and x₂ as the mole fractions of CO₂ (component 1) and n-decane (component 2) in the liquid phase, respectively, such that x₁ + x₂ = 1. Similarly, y₁ and y₂ denote the mole fractions of CO₂ and n-decane in the vapor phase. The subscript L refers to the liquid phase and V to the vapor phase. For a binary mixture of CO₂ and n-decane under isothermal conditions, we fix the temperature T and x₁, and solve for two unknowns: the bubble point pressure P_bubble and y₁. These parameters must satisfy vapor–liquid equilibrium (VLE), where the fugacity of each component is equal in both phases: f_1,L = f_1,V and f_2,L = f_2,V.

Under the PR-EOS formulation, fugacity depends nonlinearly on pressure and composition through the compressibility factor Z and other mixture-specific EOS parameters.¹⁰⁰ As a result, the VLE conditions form a coupled system of two nonlinear equations in two unknowns: P_bubble and y₁. To enforce these conditions numerically, we define a residual vector F = f_L − f_V, and solve the resulting nonlinear system using the Newton–Raphson method. Convergence is reached when the residual norm ‖F‖ falls below 10⁻⁶, indicating that thermodynamic equilibrium has been achieved within an acceptable tolerance. Initial guesses for P_bubble and y₁ are estimated using the Wilson correlation¹⁰⁷ for better convergence. Additional computational details are provided in Appendix A. Each converged solution yields a point on the phase envelope shown in Fig. 9A.

By repeating this calculation over a range of compositions, we traced the full phase envelopes at each temperature. The resulting envelopes for T = 40, 50, 70, and 90 °C are shown in Fig. 9A. As the temperature increases, the envelopes shift toward higher pressures, with the corresponding MMP increasing from 8.26 MPa at 40 °C to 9.73 MPa at 50 °C, 12.88 MPa at 70 °C, and 15.43 MPa at 90 °C. The CO₂–decane phase envelopes determined by PR-EOS in this study are benchmarked against the experimental data from Jiménez-Gallegos et al.¹⁰³ and EOS models reported by Nascimento et al.¹⁰⁴ and Pan and Trusler⁹⁴ at 40 and 70 °C, showing good agreement as presented in Fig. 11 in Appendix B.

3.2.2. Temperature dependence of MMP. Following the discussion in section 3.1, we define MMP as the pressure at which the diffusive regime first appears. The experimentally determined MMP values at four temperatures are plotted in Fig. 10. The experimental results are strongly reproducible, with small error bars calculated from three independent runs using 0.5 × (MAX − MIN) to quantify variability. Literature data from Bao et al.,⁷⁵ Shi et al.,⁷⁶ Pan et al.,⁷⁷ and Tao et al.⁷⁸ for the CO₂–decane system, shown as gray symbols in Fig. 10, were obtained using microfluidic dead-end channel geometries and are also consistent, agreeing well with both our experiments and PR-EOS predictions.

Fig. 10 Comparison of minimum miscibility pressure (MMP) determined from our microfluidic experiments using dynamic flow () and Peng–Robinson EOS predictions () at T = 40, 50, 70, and 90 °C. Literature data using dead-end microfluidic channels to determine CO₂–decane MMP from Bao et al.⁷⁵ (), Shi et al.⁷⁶ (), Pan et al.⁷⁷ () and Tao et al.⁷⁸ (). Error bars represent the scatter across three independent experimental datasets, calculated as 0.5 × (MAX − MIN). A linear regression of the experimental data is shown as a dashed line, given by MMP (MPa) = 0.136 × T (°C) + 2.976, with a coefficient of determination R² = 0.997.

Fig. 11 Benchmarking of PR-EOS-predicted phase envelopes in this work for the CO₂–n-decane system against experimental data from Jiménez-Gallegos et al.¹⁰³ at 71.6 °C, and PR-EOS models from Pan and Trusler⁹⁴ and Nascimento et al.¹⁰⁴ at 40 and 70.9 °C. The present model shows good agreement across the full two-phase region.

Within the tested temperature range (T = 40–90 °C), MMP increases nearly linearly with temperature. This strong correlation from microfluidic measurement aligns with previous studies identifying temperature as a dominant factor affecting MMP, based on sensitivity analysis.⁴²

The observed increase in MMP with temperature reflects the reduced solubility of CO₂ in n-decane at elevated temperatures.¹⁰⁸ Since CO₂ dissolution is an exothermic process, adding thermal energy shifts the equilibrium against solubility, thereby requiring higher pressure to achieve single-phase conditions.^6,108 Our observed slope of MMP increasing with temperature (ΔMMP/ΔT = 0.136 MPa °C⁻¹) is comparable to that reported for light crude oils in VIT-based MMP studies, where ΔMMP/ΔT = 0.1448 MPa °C⁻¹.¹⁰⁹ As noted by Lashkarbolooki et al.,¹⁰⁹ this slope may increase for a heavier crude oil.

A recent micro/nanofluidic study by Tao et al.⁷⁸ explored a wider temperature range (50–180 °C) and reported a non-monotonic trend in MMP. Based on thermodynamic modeling, the authors found that for a 1 μm deep channel, the MMP increases with temperature up to 140 °C, and then begins to decrease beyond this transitional point.⁷⁸ This behavior was attributed to changes in the CO₂–oil phase equilibrium: at lower T, MMP increases due to decreasing CO₂ solubility, while at higher T, the oil phase exhibits gas-like behavior that enhances mass transfer, reduces interfacial tension, and lowers the pressure required for miscibility.⁷⁸ Additionally, the transitional temperature decreased from 140 °C to 110 °C as the channel depth was reduced from 1 μm to 10 nm, highlighting the pronounced effect of nanoconfinement.⁷⁸

3.2.3. Comparison of EOS and microfluidic MMP. The MMP values of CO₂ and n-decane determined from microfluidic experiments are 8.33 ± 0.08, 9.77 ± 0.11, 12.71 ± 0.12, and 15.08 ± 0.06 MPa at 40, 50, 70, and 90 °C, respectively. As shown in Fig. 9A, the PR-EOS predicts the corresponding MMP values of 8.26, 9.73, 12.88, and 15.43 MPa. The microfluidic MMP measurements agree well with PR-EOS predictions across the full temperature range tested (Fig. 10), with a maximum deviation of 2.25% at 90 °C. The good agreement indicates that microfluidic visualization of flow regime transitions offers a reliable method for determining MMP, even under dynamic conditions where thermodynamic equilibrium is not strictly maintained.

The close agreement between microfluidic MMP measurements and PR-EOS predictions may be attributed to enhanced mixing and mass transfer at the fluid interface. The small length scale of the microfluidic platform, along with the meandering channel geometry, promotes passive mixing and increases interfacial contact between CO₂ and decane. The resulting rapid and continuous component exchange likely allows the local composition to approach near-equilibrium, leading to flow regime transitions that closely match PR-EOS predictions, despite the system being dynamically operated.

Nguyen et al.⁷³ pioneeringly introduced a microfluidic flow-based method for quantifying MMP in microfluidics, using fluorescence intensity profiles in a T-junction geometry to monitor composition changes.⁷³ In contrast, our approach determines MMP by identifying hydrodynamic flow regime transitions—from transitional to diffusive flow—in a flow-focusing geometry. This strategy eliminates the need for fluorescent tracers, enables real-time assessment with standard optical imaging, and directly captures miscibility-driven morphological changes. Furthermore, our method spans a broader temperature range (40–90 °C) and is validated against Peng–Robinson EOS predictions, demonstrating applicability under diverse thermodynamic conditions.

While the applicability of this approach to more complex fluid systems, such as CO₂ with multi-component crude oil, remains to be investigated, the microfluidic method provides valuable and unique insight into the transient evolution of fluid interfaces, including interface formation, droplet breakup, and mixing behavior, which cannot be captured by thermodynamic EOS models. In contrast to EOS-based predictions, which require careful model selection and parameter tuning,^31,102 microfluidic visualization provides a direct visual indicator for miscibility transitions under reservoir-relevant conditions.

4. Conclusions

In this study, we developed and utilized a high-pressure microfluidic platform with a flow-focusing junction to investigate the miscibility transition between CO₂ and a model light oil, n-decane, under dynamic co-flow conditions at temperatures of 40, 50, 70, and 90 °C. At 40 °C, five distinct flow regimes were observed as the pressure increased from 1 to 10 MPa: dripping, quasi-steady jetting, unsteady jetting, transitional, and diffusive regimes. At higher temperatures (50, 70, and 90 °C), we focused on the latter three regimes, where similar transitions were observed with increasing pressure. The flow regime transitions were driven by the increasing dominance of inertial and viscous forces over surface tension, arising from pressure-dependent changes in CO₂ properties. Specifically, increasing pressure led to higher density and viscosity of CO₂, along with a reduction in interfacial tension between CO₂ and decane.

The minimum pressure at which the diffusive regime first appeared was used as the criterion for determining the MMP experimentally. The diffusive regime was characterized by intense mixing and the disappearance of a sharp interface downstream, indicating fully miscible conditions. Miscibility was not achieved immediately at the flow-focusing junction but developed progressively along the microchannel, enhanced by mixing induced in the meander region.

The experimentally determined MMP values showed good agreement with Peng–Robinson EOS predictions with deviations less than 2.25%. Over the tested range of 40 to 90 °C, MMP increased linearly with temperature, consistent with the decreasing solubility of CO₂ in decane at higher temperatures.

These findings highlight the potential of microfluidic flow visualization as a reliable and rapid method for identifying miscibility transitions. In addition to MMP determination, the transitional regime observed just below the diffusive regime provides further insight into near-miscible conditions. Within this pressure range, the low IFT values are likely sufficient to suppress capillary trapping and enhance sweep efficiency. In CO₂-EOR applications where full miscibility cannot be achieved due to operational or reservoir constraints, operating near the transitional regime could still enable effective oil displacement.

A Peng–Robinson EOS-based prediction of minimum miscibility pressure (MMP)

This appendix provides the thermodynamic formulation and solution algorithm used to calculate the phase envelope of the CO₂–n-decane system using the Peng–Robinson equation of state (PR-EOS).¹⁰⁰

The overall solution algorithm follows the steps illustrated in Fig. 9B. For each temperature T, we construct the bubble point curve (open ○ in Fig. 9A) by incrementally varying the CO₂ mole fraction in the liquid phase (x₁). At each x₁, the goal is to determine the corresponding bubble point pressure (P_bubble) and CO₂ mole fraction in the vapor phase (y₁) that satisfy vapor–liquid equilibrium (VLE) conditions for both components. This requires solving the Peng–Robinson EOS to compute fugacities in both phases and enforcing f_i,L = f_i,V for i = 1 and 2, resulting in a nonlinear system with two equations and two unknowns. The system of nonlinear equations is solved using the Newton–Raphson method.

A1 Input properties and initial guess

We begin by specifying fluid properties for CO₂ (component 1) and n-decane (component 2). For each component, the critical temperature T_c,i, critical pressure P_c,i, and acentric factor ω_i are sourced from the DIPPR 801 database.¹¹⁰ The Wilson correlation¹⁰⁷ provides an initial estimate of the vapor mole fraction (y₁) and bubble point pressure (P_bubble) for improved convergence in the Newton–Raphson iteration.

A2 Peng–Robinson EOS parameters

For implementation, the cubic form of the PR-EOS is used, where is the compressibility factor, reflecting how much a real fluid deviates from ideal gas behavior:¹⁰⁰

Z³ − (1 − B)Z² + (A − 3B² − 2B)Z − (AB − B² − B³) = 0. (6)

Eqn (6) is the same as eqn (5) but expressed in terms of compressibility factor Z, where

A3 Pure component properties

To compute A and B, the mixture parameters a and b must first be calculated from pure-component properties. For each component i in the mixture, the pure-component parameters a_i and b_i are defined by:¹⁰⁰andwhere κ_i is a constant related to the acentric factor ω_i by

κ_i = 0.37464 + 1.54226ω_i − 0.26992ω_i². (10)

A4 van der Waals mixing rule

The mixture parameters a and b are calculated based on the pure component parameters (a_i, b_i), the component mole fractions (x_i, x_j), and the intermolecular interactions (a_ij). These are evaluated using the classical van der Waals mixing rules:¹⁰⁰where

a_ij = (1 − k_ij)a_i^1/2a_j^1/2, (12)

andIn eqn (12), k_ij is a binary interaction parameter typically determined empirically. In this work, k_ij is obtained by interpolating the data from Li et al.¹⁰¹

With the mixture parameters a and b determined, the constants A and B are computed using eqn (7). These values are then used to solve the cubic equation in eqn (6) for the compressibility factor Z. In the two-phase region, eqn (6) yields three real roots: the smallest root corresponds to the liquid phase compressibility Z_L, and the largest to the vapor phase compressibility Z_V.¹⁰⁰

A5 Fugacity and vapor–liquid equilibrium

The liquid and vapor compressibility factors (Z_L and Z_V) are used to determine the fugacity of each component (k = 1, 2) in the liquid and vapor phases, respectively.¹⁰⁰ Fugacity is the effective partial pressure of a component in a real mixture, given by:¹⁰⁰Here, f_k is the fugacity of component k in the mixture, x_k is its mole fraction, P is the pressure, and Z is the compressibility factor for the phase in which the fugacity is being evaluated. Eqn (14) is applied using Z_L and Z_V to compute f_k,L and f_k,V for the liquid and vapor phases, respectively. Here, k identifies the component being evaluated, while i serves as the summation index.

At vapor–liquid equilibrium, the fugacity of each component must be equal in both phases.¹⁰⁰ To maintain consistency with the rest of the appendix, we use i to represent component indices in the following discussion:¹⁰⁰

f_i,L = f_i,V for component i = 1, 2. (15)

We define the fugacity residual vector:

which quantifies deviation from equilibrium and is used in the Newton–Raphson solver.

A6 Newton–Raphson method

The Newton–Raphson method is used to solve the nonlinear system of equations F = 0, where F is defined in eqn (16). The variable vector is X = [y₁, P_bubble]^T, and an initial guess is provided using the Wilson correlation¹⁰⁷ to improve convergence. The residual vector F is evaluated at each iteration, and the solution is updated iteratively using the Newton–Raphson scheme:¹¹¹

X_n+1 = X_n − J⁻¹F(X_n), (17)

where J is the Jacobian matrix of partial derivatives. Iteration continues until ‖F‖ < 10⁻⁶. Repeating this process for a range of x₁ at fixed T yields the full bubble point curve. The dew point curve is computed analogously by fixing y₁ and solving for x₁ and P_dew.

Finally, the MMP is identified as the highest pressure on the phase envelope, above which the mixture becomes single-phase at all compositions.

B Benchmark of PR-EOS phase envelopes

This appendix presents a benchmark of the calculated phase envelopes using the PR-EOS formulation described in Appendix A. Model predictions for the CO₂–n-decane system are compared against experimental data¹⁰³ and previously published results^94,104 at 40 °C and approximately 70 °C. Data at 70.9 °C and 71.6 °C were used as the closest available references.^94,103 This comparison confirms that the present implementation (shown in solid lines) captures the two-phase behavior with high accuracy across the full range of conditions. Following successful benchmarking, we extended the PR-EOS calculations to the temperature ranges of interest that have not been previously reported in the literature.

C Temperature-dependent variation of dimensionless parameters at MMP

Fig. 12 illustrates the variation of dimensionless parameters at the minimum miscibility pressure (MMP) across different temperatures, indicated by different colors. The vertical dashed lines mark the MMP for each temperature, and the corresponding horizontal lines denote the value of each dimensionless parameter at that transition point.

Fig. 12 Variation of dimensionless parameters at different temperatures and pressures. (A) The density ratio ρ_CO2/ρ_n-decane, (B) the viscosity ratio μ_CO2/μ_n-decane, and (C) the Reynolds number of CO₂.

In Fig. 12A, the density ratio ρ_CO2/ρ_n-decane at MMP shows no clear trend with temperature and remains approximately 0.5. In Fig. 12B, the viscosity ratio μ_CO2/μ_n-decane at MMP increases with temperature. In Fig. 12C, despite changes in pressure and temperature, the Reynolds number of the CO₂ stream at MMP remains close to 150 for different temperatures. Based on the definition of the Reynolds number and the fixed flow velocity u and characteristic length L in our setup, this nearly constant value of Re suggests that the kinematic viscosity of CO₂ (μ/ρ) at the miscibility transition changes very little with temperatures.

D Grayscale intensity analysis of flow regimes

We performed a grayscale intensity analysis for one set of experiments at 40 °C, covering three characteristic flow regimes: unsteady jetting (7.1 MPa), transitional (8.4 MPa), and diffusive (8.5 MPa) regimes. Grayscale profiles were extracted from both section 1 and section 2 of the microfluidic channel using ImageJ. The selected regions of interest (ROI) correspond directly to the snapshots shown in Fig. 13. At each x-location along the flow direction, grayscale values were averaged across the channel width (vertical direction) and then normalized between 0 and 1 for comparison. In these profiles, smaller grayscale values correspond to darker regions in the snapshot, such as CO₂–oil interfaces, which appear as valleys in the curve.

Fig. 13 (A) Microchannel schematic indicating section 1 and section 2 as the defined regions of interest (ROI). Section 1 is positioned at the junction where the CO₂–decane interface forms, while section 2 is located 8.6 mm downstream. Scale bars: 200 μm. (B) Experimental snapshots at section 1 (left) and section 2 (right) for unsteady jetting (40 °C, 7.1 MPa), transitional (40 °C, 8.4 MPa), and diffusive (40 °C, 8.5 MPa) regimes. (C) Grayscale intensity profiles (normalized to 0–1) along the flow direction for each ROI, averaged across the channel width from the corresponding snapshots in (B).

As shown in Fig. 13, the section 1 profiles close to the initial flow junction were similar for all three regimes. In section 2, however, the unsteady jetting case exhibited large oscillations and multiple local minima due to distinct CO₂ bubble boundaries appearing dark, the transitional regime showed fewer fluctuations from smaller, less frequent bubbles, and the diffusive regime displayed an almost flat profile, indicating the absence of visible interfaces.

Author contributions

Junyi Yang: investigation, methodology, software, formal analysis, writing – original draft, visualization. Peichun Amy Tsai: conceptualization, methodology, validation, software, formal analysis, writing – revision, visualization, supervision, funding acquisition.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information: Video S1. Experimental video clips showing the complete set of flow regimes at T = 40 °C. See DOI: https://doi.org/10.1039/D5LC00616C.

The data supporting this article have been included as part of the SI.

Acknowledgements

We gratefully acknowledge the support from the Canada First Research Excellence Fund (CFREF), the Future Energy System (FES T02-P05 CCUS projects) at the University of Alberta, and the Canada Foundation for Innovation (CFI 34546). J. Y. acknowledges the support of the Graduate Student Scholarship from Alberta Innovates (AI-GSS) and the Future Energy System Legacy Builders Program (FES T02-P05-LBP1). P. A. T. holds a Canada Research Chair (CRC) in Fluids and Interfaces and gratefully acknowledges funding from the Natural Sciences and Engineering Research Council of Canada (NSERC) and Alberta Innovates (AI) (ADVANCE 202102818), in particular the NSERC Canada Research Chairs Program (CRC 233147) and Discovery Grant (RGPIN-2020-05511). We thank the nanoFAB facility at the University of Alberta for support with microfabrication, particularly Aaron Hryciw for process consultation, and Breanna Cherkawski, Gustavo de Oliveira Luiz, and Scott Munro for training and technical assistance. We also acknowledge the Elko Engineering Garage at the University of Alberta for support with microfluidic chip holder fabrication, and thank Nicolas Marchand for training and assistance.

Notes and references

1. International Energy Agency (IEA), Can CO₂-EOR really provide carbon-negative oil?, 2019, https://www.iea.org/commentaries/can-co2-eor-really-provide-carbon-negative-oil, accessed: 2025-01-01.
2. J. H. Ching, P. Chen and P. A. Tsai, Phys. Rev. Fluids, 2017, 2, 1–9.
3. N. Moradpour, J. Yang and P. A. Tsai, Curr. Opin. Colloid Interface Sci., 2024, 101845.
4. N. Moradpour, R. Azadi and P. A. Tsai, Colloids Surf., A, 2025, 705, 135533.
5. International Energy Agency (IEA), Whatever happened to enhanced oil recovery?, 2018, https://www.iea.org/commentaries/whatever-happened-to-enhanced-oil-recovery, accessed: 2025-01-01.
6. S. K. Prasad, J. S. Sangwai and H. S. Byun, J. CO2 Util., 2023, 72, 102479.
7. A. Alinejad and H. Dehghanpour, Fuel, 2024, 359, 130387.
8. A. Alinejad, J. Behseresht and H. Dehghanpour, SPE J., 2025, 30, 251–271.
9. N. Mungan, SPE West. Reg. Meetin, Long Beach, California, 1991, pp. 113–122.
10. F. B. Thomas, N. Holowach, X. Zhou, D. B. Bennion and D. W. Bennion, Proc. - SPE Symp. Improv. Oil Recover, 1994, pp. 307–320.
11. N. Zhang, M. Yin, M. Wei and B. Bai, Fuel, 2019, 241, 459–467.
12. J. Sun, L. Sun and B. Bao, Fuel, 2025, 392, 134835.
13. C. Whitson and M. Brule, Phase Behavior, Society of Petroleum Engineers Inc., 2000, vol. 20.
14. A. Abedini and F. Torabi, Energy Fuels, 2014, 28, 774–784.
15. X. Zhang, L. Li, Y. Su, Q. Da, J. Fu, R. Wang and F. Chen, Appl. Energy, 2023, 348, 121518.
16. A. Zick, 61st Annu. Tech. Conf. Exhib. Soc. Pet. Eng., 1986, p. 15493.
17. F. I. Stalkup, 62nd Annu. Tech. Conf. Exhib. Soc. Pet. Eng., 1987, p. 16715.
18. J. Han, M. Lee, W. Lee, Y. Lee and W. Sung, Appl. Energy, 2016, 161, 85–91.
19. X. Chen, Y. Li, X. Tang, H. Qi, X. Sun and J. Luo, Energy, 2021, 226, 120294.
20. Y. Hao, Z. Li, Y. Su, C. Kong, H. Chen and Y. Meng, Energy, 2022, 254, 124349.
21. S. Yu, L. Yiqiang, Q. Huan, C. Jinxin and T. Xuechen, 16th Int. Conf. Appl. Energy, 2025.
22. P. Saffman and G. Taylor, Proc. R. Soc. London, Ser. A, 1958, 245, 312–329.
23. R. L. Chuoke, P. Von Meurs, J. Member and A. C. Von Der Poel, Transactions of the AIME, 1959, 216, 188–194.
24. T. T. Al-Housseiny, P. A. Tsai and H. A. Stone, Nat. Phys., 2012, 8, 747–750.
25. Y. Guo, F. Liu, J. Qiu, Z. Xu and B. Bao, Energy, 2022, 256, 124524.
26. J. Ren, W. Xiao, W. Pu, Y. Tang, Y. Bernabé, Q. Cheng and L. Zheng, Energy, 2024, 297, 131211.
27. G. Bongrand and P. A. Tsai, Phys. Rev. E, 2018, 97, 2–5.
28. A. Pouplard and P. A. Tsai, Sci. Rep., 2024, 14, 10–17.
29. N. Kumar, M. Augusto Sampaio, K. Ojha, H. Hoteit and A. Mandal, Fuel, 2022, 330, 125633.
30. M. Almobarak, Z. Wu, D. Zhou, K. Fan, Y. Liu and Q. Xie, Petroleum, 2021, 7, 245–253.
31. B. Dindoruk, R. Johns and F. M. Orr, SPE Reservoir Eval. Eng., 2021, 24, 367–389.
32. W. F. Yellig and R. S. Metcalfe, JPT, J. Pet. Technol., 1980, 32, 160–168.
33. G. Song, Y. Meng, C. Zhang, Z. Zhao and Q. Yang, ACS Omega, 2024, 9, 14747–14765.
34. D. N. Rao, Fluid Phase Equilib., 1997, 139, 311–324.
35. S. C. Ayirala and D. N. Rao, SPE Soc. Pet. Eng., 2006, 1, 349–363.
36. R. L. Christiansen and H. K. Haines, SPE Reservoir Eng., 1987, 2, 523–527.
37. K. Zhang, N. Jia, F. Zeng, S. Li and L. Liu, J. Pet. Sci. Eng., 2019, 183, 106366.
38. K. Jessen and F. M. Orr, SPE Reservoir Eval. Eng., 2008, 933–939.
39. D. Zhou and F. M. Orr, SPE J., 1998, 3, 19–25.
40. R. B. Alston, G. P. Kokolis and C. F. James, SPEJ, Soc. Pet. Eng. J., 1985, 25, 268–274.
41. O. S. Glaso, SPEJ, Soc. Pet. Eng. J., 1985, 25, 927–934.
42. H. Chen, B. Li, I. Duncan, M. Elkhider and X. Liu, Fuel, 2020, 278, 118272.
43. E. H. Benmekki and G. A. Mansoori, SPE Reservoir Eng., 1988, 3, 559–564.
44. S. N. Ashrafizadeh and A. Ameri Ghasrodashti, Chem. Eng. Res. Des., 2011, 89, 690–696.
45. K. Ma, R. Liontas, C. A. Conn, G. J. Hirasaki and S. L. Biswal, Soft Matter, 2012, 8, 10669–10675.
46. C. A. Conn, K. Ma, G. J. Hirasaki and S. L. Biswal, Lab Chip, 2014, 14, 3968–3977.
47. X. Zheng, N. Mahabadi, T. S. Yun and J. Jang, J. Geophys. Res.:Solid Earth, 2017, 122, 1634–1647.
48. Y. Li, F. Kazemifar, G. Blois and K. T. Christensen, Water Resour. Res., 2017, 5, 6178–6196.
49. M. Saadat, J. Yang, M. Dudek, G. Øye and P. A. Tsai, J. Pet. Sci. Eng., 2021, 203, 108647.
50. J. Yang and P. A. Tsai, Chem. Eng. Sci., 2024, 300, 120543.
51. S. Yang, H. Li, S. Suo and Z. Wu, Adv. Water Resour., 2024, 191, 104757.
52. L. Li, J. Zheng, Y. Shi, Y. Su, Y. Hao and Z. Chen, Energy Fuels, 2024, 38, 7997–8008.
53. J. Wang, J. Sun, J. Shi and B. Bao, Lab Chip, 2025, 25, 1981–1992.
54. N. Moradpour and P. A. Tsai, Microfluidic Investigation of CO₂ Foam Flow in a Heterogeneous Porous Medium, Lab Chip, 2025 10.1039/D5LC00544B , in press.
55. N. Liu, C. Aymonier, C. Lecoutre, Y. Garrabos and S. Marre, Chem. Phys. Lett., 2012, 551, 139–143.
56. M. Abolhasani, M. Singh, E. Kumacheva and A. Günther, Lab Chip, 2012, 12, 1611–1618.
57. A. Sell, H. Fadaei, M. Kim and D. Sinton, Environ. Sci. Technol., 2013, 47, 71–78.
58. N. Qin, J. Z. Wen and C. L. Ren, Phys. Rev. E, 2017, 95, 1–15.
59. N. Qin, J. Z. Wen, B. Chen and C. L. Ren, Appl. Phys. Lett., 2018, 113, 033703.
60. S. G. Lefortier, P. J. Hamersma, A. Bardow and M. T. Kreutzer, Lab Chip, 2012, 12, 3387–3391.
61. T. H. M. Ho, J. Yang and P. A. Tsai, Lab Chip, 2021, 21, 3942–3951.
62. T. H. M. Ho, D. Sameoto and P. A. Tsai, Chem. Eng. Res. Des., 2021, 174, 116–126.
63. J. Yang and P. A. Tsai, Biomicrofluidics, 2024, 18, 051301.
64. S. Yang, G. Kong, Z. Cao and Z. Wu, Chem. Eng. J. Adv., 2023, 16, 100518.
65. S. Virappane, R. Azadi, N. Mukherjee and P. A. Tsai, Phys. Fluids, 2024, 36, 102019.
66. Y. Guo, F. Li, S. Saber, M. Zargartalebi, S. S. Sun, Y. C. Xiao, B. Bao, Z. Xu and D. Sinton, Lab Chip, 2025, 25(12), 2918–2925,  10.1039/D4LC00772G.
67. F. Mostowfi, S. Molla and P. Tabeling, Lab Chip, 2012, 12, 4381–4387.
68. W. Song, H. Fadaei and D. Sinton, Environ. Sci. Technol., 2014, 48, 3567–3574.
69. B. Bao, J. Riordon, Y. Xu, H. Li and D. Sinton, Anal. Chem., 2016, 88, 6986–6989.
70. Y. Xu, J. Riordon, X. Cheng, B. Bao and D. Sinton, Angew. Chem., 2017, 129, 14150–14155.
71. B. Bao, J. Riordon, F. Mostowfi and D. Sinton, Lab Chip, 2017, 17, 2740–2759.
72. M. J. Fornerod, E. Amstad and S. Guldin, Mol. Syst. Des. Eng., 2020, 5, 358–365.
73. P. Nguyen, D. Mohaddes, J. Riordon, H. Fadaei, P. Lele and D. Sinton, Anal. Chem., 2015, 87, 3160–3164.
74. A. Sharbatian, A. Abedini, Z. Qi and D. Sinton, Anal. Chem., 2018, 90, 2461–2467.
75. B. Bao, J. Feng, J. Qiu and S. Zhao, ACS Omega, 2021, 6, 943–953.
76. J. Shi, L. Tao, Y. Guo, X. He, Y. Li and B. Bao, Fuel, 2024, 362, 130876.
77. X. Pan, L. Sun, F. Chen, X. Huo, Y. Wang, C. Feng, X. Zheng and Z. Zhang, Energy Fuels, 2024, 38, 10904–10913.
78. L. Tao, W. Liu, J. Shi, Y. Guo, W. Qin and B. Bao, Chem. Eng. Sci., 2025, 302, 120828.
79. H. Zou, A. C. Slim and A. Neild, Anal. Chem., 2019, 91, 13681–13687.
80. H. Zou, H. Kang, A. C. Slim and A. Neild, Lab Chip, 2020, 20, 3582–3590.
81. F. Ungar, S. Ahitan, S. Worthing, A. Abedini, K. Uleberg and T. Yang, J. Pet. Sci. Eng., 2022, 208, 109415.
82. D. Pereponov, M. Tarkhov, D. B. Dorhjie, A. Rykov, I. Filippov, E. Zenova, V. Krutko, A. Cheremisin and E. Shilov, Energies, 2023, 16, 4994.
83. D. Qin, Y. Xia, J. A. Rogers, R. J. Jackman and G. M. Whitesides, Microsyst. Technol. Chem. Life Sci., Springer, Berlin, Heidelberg, 1998, vol. 194, pp. 2–20.
84. V. A. Lifton, Lab Chip, 2016, 16, 1777–1796.
85. D. Sinton, Lab Chip, 2014, 14, 3127–3134.
86. A. Abedini, S. Ahitan, Z. Barikbin, V. Soni, J. Ratulowski and D. Sinton, Energy Fuels, 2022, 36, 8578–8590.
87. S. S. Datta, I. Battiato, M. A. Fernø, R. Juanes, S. Parsa, V. Prigiobbe, E. Santanach-Carreras, W. Song, S. L. Biswal and D. Sinton, Lab Chip, 2023, 23, 1358–1375.
88. N. Convery and N. Gadegaard, Micro Nano Eng., 2019, 2, 76–91.
89. S. Marre and K. F. Jensen, Chem. Soc. Rev., 2010, 39, 1183–1202.
90. S. Marre, A. Adamo, S. Basak, C. Aymonier and K. F. Jensen, Ind. Eng. Chem. Res., 2010, 49, 11310–11320.
91. Z. B. Qi, L. Xu, Y. Xu, J. Zhong, A. Abedini, X. Cheng and D. Sinton, Lab Chip, 2018, 18, 3872–3880.
92. C. A. Schneider, W. S. Rasband and K. W. Eliceiri, Nat. Methods, 2012, 9, 671–675.
93. National Institute of Standards and Technology (NIST), NIST Chemistry WebBook, SRD 69, 2025, https://webbook.nist.gov/chemistry/, accessed: 2025-01-01.
94. Z. Pan and J. P. M. Trusler, Fluid Phase Equilib., 2023, 566, 113700.
95. J. K. Nunes, S. S. Tsai, J. Wan and H. A. Stone, J. Phys. D: Appl. Phys., 2013, 46, 114002.
96. T. Cubaud and T. G. Mason, Phys. Fluids, 2008, 20, 053302.
97. A. S. Utada, A. Fernandez-Nieves, H. A. Stone and D. A. Weitz, Phys. Rev. Lett., 2007, 99, 1–4.
98. A. Skauge and K. Serbie, Soc. Pet. Eng. - SPE EOR Conf. Oil Gas West Asia 2014 Driv. Integr. Innov. EOR, 2014, pp. 891–905.
99. T. Cubaud and S. Notaro, Phys. Fluids, 2014, 26, 122005.
100. D.-Y. Peng and D. B. Robinson, Ind. Eng. Chem. Fundam., 1976, 15, 59–64.
101. X. Li, D. Yang, X. Zhang, G. Zhang and J. Gao, Fluid Phase Equilib., 2016, 417, 77–86.
102. J. S. Lopez-Echeverry, S. Reif-Acherman and E. Araujo-Lopez, Fluid Phase Equilib., 2017, 447, 39–71.
103. R. Jiménez-Gallegos, L. A. Galicia-Luna and O. Elizalde-Solis, J. Chem. Eng. Data, 2006, 51, 1624–1628.
104. F. P. Nascimento, M. L. Paredes, A. P. D. Bernardes and F. L. Pessoa, J. Supercrit. Fluids, 2019, 154, 1–9.
105. Y. S. Wei and R. J. Sadus, AIChE J., 2000, 46, 169–196.
106. İ. Tosun, The Thermodynamics of Phase and Reaction Equilibria, Elsevier, 2nd edn, 2021, pp. 365–478.
107. G. Wilson, AIChE Annual Meeting, Cleveland, Ohio, 1968.
108. M. Abdurrahman, A. K. Permadi, A. Arsad, A. F. Abdul Rahman, W. Bae, U. Z. Husna, A. L. Pang and R. Fauzi, ACS Omega, 2023, 8, 8703–8711.
109. M. Lashkarbolooki, M. J. Eftekhari, S. Najimi and S. Ayatollahi, J. Supercrit. Fluids, 2017, 127, 121–128.
110. Design Institute for Physical Property Research/AIChE, DIPPR's Project 801 Database, 2025, https://www.aiche.org/dippr, accessed: 2025-01-01.
111. J. M. Ortega, Iterative solution of nonlinear equations in several variables, Academic Press, 1970.

---