# Probing Current Sheet Instabilities from Flare Ribbon Dynamics

Ryan J. French  Mullard Space Science Laboratory, University College London,
Dorking, RH5 6NT, UK  Sarah A. Matthews  Mullard Space Science Laboratory,
University College London, Dorking, RH5 6NT, UK  I. Jonathan Rae  Department
of Mathematics, Physics and Electrical Engineering, Northumbria University,
Newcastle upon Tyne, UK  Andrew W. Smith  Mullard Space Science Laboratory,
University College London, Dorking, RH5 6NT, UK

(Accepted to ApJ, September 2021)

###### Abstract

The presence of current sheet instabilities, such as the tearing mode
instability, are needed to account for the observed rate of energy release in
solar flares. Insights into these current sheet dynamics can be revealed by
the behaviour of flare ribbon substructure, as magnetic reconnection
accelerates particles down newly reconnected field lines into the chromosphere
to mark the flare footpoints. Behaviour in the ribbons can therefore be used
to probe processes occurring in the current sheet.

In this study, we use high-cadence (1.7 s) IRIS Slit Jaw Imager observations
to probe for the growth and evolution of key spatial scales along the flare
ribbons - resulting from dynamics across the current sheet of a small solar
flare on December 6th 2016. Combining analysis of spatial scale growth with Si
IV non-thermal velocities, we piece together a timeline of flare onset for
this confined event, and provide evidence of the tearing-mode instability
triggering a cascade and inverse cascade towards a power spectrum consistent
with plasma turbulence.

sun — solar flares — magnetic reconnection

##  1 Introduction

In the standard eruptive flare model, magnetic reconnection originates in a
thin current sheet; created by the inflow of oppositely orientated magnetic
fields under a rising magnetic flux rope (Carmichael, 1964; Sturrock, 1968;
Hirayama, 1974; Kopp & Pneuman, 1976). The nature of magnetic reconnection in
the current sheet, and the processes leading to their formation, are not yet
fully understood. In order to explain the breakdown of ideal
magnetohydrodynamics (on small spatial scales) needed to enable fast magnetic
reconnection in the current sheet (leading to solar flare onset), recent
attention has turned to the role of the tearing-mode instability (Biskamp,
1986). In the tearing-mode instability, a classic Sweet-Parker current sheet
reaches an unstable aspect ratio, collapsing and reconnecting at multiple
points along the sheet to form magnetic islands, or plasmoids. These magnetic
islands form at a preferred spatial scale, and continue to collapse further to
produce islands at progressively smaller scales. Larger islands are also
produced via the reconnection and coalescence at island boundaries (Tenerani &
Velli, 2020). This critical aspect ratio and spatial scales are dependent on
the local magnetic conditions, in particular the Lunduist number. The cascade
to smaller island scales is predicted to produce plasma turbulence, with a
power-law similar to those observed in in-situ plasmas (Shaikh & Zank, 2010).
The exact relationship between the tearing-mode instability and turbulence
within a reconnecting current sheet is still a prominent area of ongoing
research.

Current sheets are notoriously difficult to observe directly, primarily due to
the low density and small width of the region. Heated regions of plasma around
current sheets (known as a plasma sheet) can occasionally be observed off-
limb, but the rarity of their observation means that spectroscopic and multi-
wavelength observations, needed to truly disentangle the energy release
process (and relationship with instabilities and plasma turbulence), are rare.
The most successful spectroscopic plasma sheet observation set to date
originated from the famous September 10th 2017 X-class flare (Warren et al.,
2018). This event was observed by multiple instruments across the entire
electromagnetic spectrum, facilitating insight into the nature of magnetic
reconnection in the event, and providing initial evidence for the presence of
turbulent or tearing-mode reconnection (Cheng et al., 2018; French et al.,
2019), persisting for the flare’s entire duration (French et al., 2020).
Although observations of this September 10th 2017 event provided many new
insights into the behaviour of current sheet dynamics, the uniqueness of such
a dataset requires alternative methods to be employed in order to make further
progress in understanding these processes. The field line connectivity between
flare ribbons that are routinely observed on disk and the reconnecting current
sheet, offers such a possibility. On-disk observations of flares are far more
common, as active regions can be observed for longer durations on-disk than
they can off-limb.

During a solar flare, high energy particles are accelerated from the
reconnection site down the flare loops to the chromosphere, depositing energy
at the magnetic footpoints to form flare ribbons (Cheng et al., 1983; Doschek
et al., 1983). Due to their magnetic connectivity, behaviour of flare ribbon
substructure must in some way reflect processes in the flaring current sheet
region (Forbes & Lin, 2000). Observational studies of flare ribbons have found
spectroscopic signatures suggestive of waves and/or turbulence (Brosius & Daw,
2015) \- consistent with the presence of either the tearing mode or Kelvin-
Helmholtz instability in the current sheet (Brannon et al., 2015). These
signatures of energy release in the chromosphere are detected on as fine a
scale as $\sim$ 150-300 km (Graham & Cauzzi, 2015; Jing et al., 2016).

The Earth’s magnetotail experiences similar cycles of energy storage and
release (Akasofu, 1964; McPherron, 1970). The energy builds in the Earth’s
magnetotail until it reaches a point of instability, at which point the energy
is explosively released through reconnection (Hones et al., 1984), causing
observable phenomena such as the rapid brightening of the aurora (Akasofu &
Lepping, 1977) (analogous to flare ribbon brightening on the Sun), and the
formation of auroral beads (Henderson, 1994, 2009). Fourier analysis of
auroral bead brightenings along the auroral arc has found that specific
spatial scales grow exponentially (Rae et al., 2010; Kalmoni et al., 2015).
Recently, high temporal and spatial resolution measurements of the aurora
enabled the comparison of observation and theory to show that non-local
effects, in the form of sheer Alfven waves, were responsible for plasma
instabilities in near-Earth space (Kalmoni et al., 2018).

Although fundamentally different in their temperature, density and collisional
timescales, there are distinct similarities between flare ribbons and
terrestrial aurorae driven by sub-storm activity. In this study, we take
inspiration from the magnetospheric work of Kalmoni et al. (2018), applying
the auroral beads methodology to flare ribbons of a small solar flare, to
determine whether similar evidence exists for the presence of a plasma
instability at flare onset.

##  2 Observations

We present IRIS (De Pontieu et al., 2014) observations of a confined B-class
flare in the core of AR 12615 on December 6th 2016. IRIS observed with a large
sit-and-stare Slit Jaw Imager (SJI) 1400 Å window, and a 1.7 s cadence. With
just the one slit-jaw channel, this is close to the fastest cadence IRIS can
observe. The SJI observations were 2x2 binned, resulting in a pixel width of
0.3327 ″ or 238 km. IRIS observed the flare continuously from pre-flare to
decay.

Figure 1 (top row) presents snapshots of the SJI 1400 Å evolution, with the
field of view cropped to center on the two flare ribbons. From hereon, the two
ribbons will be denoted as the east (left) and west (right) ribbon. Examining
the light curve of each region, we see the two ribbons brighten cotemporally,
and although the east ribbon is significantly larger than the west, maximum
intensities are comparable. Similar maximum intensities imply a similar level
of energy deposition in each ribbon.

Figure 1:  Top: Evolution of IRIS SJI 1400 Å. Second row: Evolution of AIA 131
Å, with IRIS ribbon contours to distinguish emission at high/low altitudes.
Third row: Evolution of HMI LoS B field, with IRIS ribbon contours. Bottom
Left: Ribbon light curve, of masked region in Fig 2. Vertical dashed lines
mark location of above panels. Bottom Right: Ribbon area over intensity
threshold, multiplied by HMI LoS field strength in each region. This is
analogous to magnetic flux (BxA) evolution within the bright ribbon region.
East ribbon flux is multiplied by -1.

The brightest ribbon regions are also observed in the 131 Å channel of the
Solar Dynamics Observatory Atmospheric Imaging Assembly (AIA, Lemen et al.,
2012). Following peak ribbon intensity, higher-altitude flare loops become
visible in this channel as they fill with hotter ablated chromospheric
material. The timing of these events is shown in Figure 1. Looking at the
photospheric Line of Sight (LOS) magnetic field measurements (B) in this
region from the Helioseismic and Magnetic Imager (HMI, Schou et al., 2012), we
see the flare ribbons (marked by red contour) appear either side of the
polarity inversion line, adjacent to the regions of highest field strength
(third row of Figure 1). At around 600 ″ from center disk, projection effects
have begun to affect the west of the HMI observations. Although the area
traced by the ribbons change throughout their evolution, the photospheric
field itself does not vary on the timescales studied. Examining the LoS
magnetic flux traced by the evolving (SJI 1400 Å) ribbons in the bottom-right
panel of Figure 1, we note that the flux is equal in each ribbon throughout
the flare. That is,

$B_{E}A_{E}=B_{W}A_{W}$ (1)

Matching flux in each ribbon implies connection of the two via a consistent
tube of flux, which does not interact significantly with any external field.
We can therefore imagine this flux tube passing from one ribbon to the other,
containing the reconnection regions within it. The flare is likely a result of
quadrapolar reconnection (Melrose, 1997; Pikelner & Livshits, 1977), with, due
to the small size of the event, footpoints of equal polarity merging together
to form one ribbon feature either side of the polarity inversion line. We
therefore expect the reconnecting current sheet to be somewhere between the
flare loops - perhaps even those seen in AIA 131 Å. The Active Region is
complex however, and other flaring configurations are possible. One other
possibility considered was the circular ribbon flare (Wang & Liu, 2012), which
would feature null point reconnection instead of the more-traditional laminar
current sheet. In this scenario, null point reconnection still contains layers
of current and thus would still be subject to instabilities (Pontin, 2011).

Spectral data of the flare ribbons were collected by the IRIS sit-and-stare
slit, centered across the larger (east) ribbon. The position of this vertical
slit can be seen as a slight enhancement in SJI image intensity. Jeffrey et
al. (2018) examine the spectral evolution of this flare ribbon, finding a
steep rise in Si IV 1402.77 Å non-thermal velocity (interpreted as a signature
of plasma turbulence) preceding the rise in Si IV intensity (marking plasma
heating flare onset).

In our study, we examine the spatial and temporal evolution of the two ribbons
in SJI 1400 Å observations, investigating the power spectra in the spatial
domain along the ribbon to provide insight into processes across the current
sheet, perpendicular to the magnetic field. We compare the power spectra at
multiple spatial scales with parameters predicted by theory, and the timings
of the Si IV non-thermal velocities introduced by Jeffrey et al. (2018). When
studying the ribbon spatial intensity evolution, we are not restricted by the
slit location and can therefore analyze both ribbons independently. Current
sheet processes may also affect structure across the ribbon, in addition to
the variation along the ribbon we probe in this study. Analyzing these
components together, along with the time-varying component, would provide the
full dispersion relation, but not possible with the resolution of this data
set. We therefore focus exclusively on searching for spatial scales along the
flare ribbon structures.

##  3 Ribbon Power Spectrum

###  3.1 Ribbon Tracking and Intensity Processing

Figure 2:  Top: Snapshots (at 10:37:37 UT) of the ribbon tracking centroid
position, along the masked data of the each ribbon. Second Row: Intensity
cross-sections along the centroids in the above panel.

In order to investigate intensity variation along the flare ribbons, we
developed a method to track a centroid along each flare ribbon as they evolved
in space and time. This method, and all subsequent analysis, is repeated
separately for both the east and west ribbon.

To begin, we produced a mask to remove all pixels never crossing an intensity
threshold (that never form part of the ribbon structure). Next, we take cross-
sections perpendicular to the main ribbon structure, and locate the position
of the brightest pixel at each of these cross-sections. We then take a moving
average of the pixel locations from this initial estimate, producing a smooth,
continuous centroid slit along the centroid of each ribbon. These centroids
change as the ribbons evolve, and a snapshot for each ribbon is shown in the
top row of Figure 2. Our main tracking limitation is the spatial resolution,
as each centroid position must still have an integer pixel position. The
tracking method works effectively however, capturing the change in curvature
of the ribbons as they evolve.

Figure 3:  Top: Ribbon intensity stack plots for both ribbons. Middle: Spatial
scale power stack plots. Bottom: Spatial scale power stack plots, normalized
at each spatial scale. All y-axes in units of Mm. (0.33 arcsec).

We measure the intensity cross-section along the ribbons, by taking the mean
intensity of a 3x3 pixel area around each pixel location in the centroid slit.
Measuring an average around each pixel provides a more reliable estimate of
intensity along the centroid (by not relying on a single pixel measurement),
and reduces the impact of higher intensities in the slit-jaw image at the sit-
and-stare slit location. Intensity measurements along this centroid are not
uniformly spaced, dependent on whether adjacent pixels in the centroid lay
horizontal/vertical or diagonally to one other. We interpolate along this
centroid to create an evenly spaced array - necessary to calculate a Fast
Fourier Transform (FFT). An example cross-section for each ribbon is shown in
the bottom panel of Figure 2. Here, we have converted position along the
centroid into units of Mm. The east ribbon is over twice as long as the west,
with centroid slit lengths of 17 and 8 Mm respectively. The intensity signal
displayed here is then processed by detrending the data with a moving average,
to avoid intensity jumps between the centroid slit ends dominating the power
spectrum. The moving average sizes were selected manually, to remove large
general trends in the data without inadvertently subtracting the smaller
variations we aim to detect. Finally, we apply a Hanning Window to the data
and introduce zero padding, to preserve the power and increase the spectral
resolution of the Fast Fourier Transform (FFT). The FFT returns the power at
spatial scales at intervals of L/2, L/3… L/n, where L is the length of the
array in real units and n the number of data points. Zero padding therefore
increases the spectral resolution of the FFT by calculating the power at
smaller intervals (n is larger), as well as minimizing the effects of
aliasing.

###  3.2 2D Power Spectrum

Intensity cross-sections along the ribbon centroid are measured and processed
for every 1.7 s time step, producing the intensity stack plots in the top row
of Figure 3. For each of these cross-sections (vertical slices along the stack
plot), we calculate the FFT of ribbon intensity, resulting in the power
spectrum in the spatial domain for each time step. This provides us with the
power of different spatial scales appearing along the flare ribbon. These
spatial scales (in Mm) are inverse to the dispersion relation wavenumber, $k$,
often used in the literature. We elect to work with units of spatial scale for
this study, as they are easier to picture in this instance. We combine each
spatial FFT into a single stack plot, presented in the middle row of Figure 3
for each ribbon. Maximum power varies significantly with spatial scale, so to
better see time evolution at each spatial scale we normalize the power
spectrum in time, to produce the bottom panel of the same figure. Between
10:37 and 10:38 UT, we see a steep rise in power across all spatial scales,
with growth rate and onset times varying across different scales. Before this
time, we see consistent signals at certain spatial scales - around 2.5 Mm in
the east ribbon (with a potential second, weaker peak higher around 5.5 Mm),
and 2 and 3.5 Mm in the west ribbon. These signals result from the spacing
between bright points in the pre-ribbon structure. This spacing can be
verified by studying the intensity stack plots in the same figure, with the
spacing between the bright horizontal features.

##  4 Spatial scales and Exponential Growth

Figure 4:  Top: Cross-section through power plot (Fig 3) at spatial scale of
1.75 Mm, showing exponential growth. Middle: Exponential growth rate at each
spatial scale. West scales are transformed to east scales. Bottom: Start/end
times of the fitted exponential period, for the east and scaled west ribbon
spatial scales.

###  4.1 Flux Conservation - Spatial Scaling

If the observed spatial scales in the flare ribbons result from processes
propagating down from the reconnection site, we can expect the length scales
to scale with the square root of the area of the flux tube connecting the two
ribbons. (Picture a wave with one mode, perpendicular to the magnetic field
direction. As the area of the flux tube changes, the wavelength scale
approximately changes by the square root of the area). Considering the flux
conservation shown in Figure 1, we can equivalently scale the length scales
through the ratios of magnetic field strength (equation 1).

As magnetic conditions vary along this tube of flux, we can apply this concept
to compare the spatial scales seen at each ribbon, by scaling the west ribbon
scales to the ’reference frame’ of the east ribbon through the square root of
the ribbon area ratios. The analysis below considers the original and scaled
spatial scales observed in the east and west ribbon respectively.

###  4.2 Growth with Spatial Scale

In MHD, exponential growth at multiple spatial scales is a classical signal of
plasma instability (Priest, 1985). We investigate the growth at a specific
spatial scale by taking a horizontal cross-section through the power spectrum
in Figure 3 \- tracking the evolution of power with time. Figure 4 (top panel)
presents an example cross-section for a scale of 1.75 Mm in the east ribbon.
The plot shows a noisy baseline power around $10^{0}$, which sharply increases
by nearly 5 orders of magnitude between 10:37:10 and 10:37:50 UT. We determine
the region of exponential growth in the data moving average, and fit an
exponential curve to the corresponding data points. This provides us with an
exponential growth rate and start/end times of the exponential phase.

By fitting an exponential curve for each spatial scale, we calculate the
variation in the start time, duration and rate of the exponential growth of
different spatial scales in each ribbon. The east ribbon growth rates are
shown as the cyan curve in the middle panel of Figure 4. The curve shows three
peaks around 1.2, 2.3 and 3.8 Mm, before tailing off at high spatial scales.
The error bars in this plot assume that the algorithm determining the
start/end point of exponential growth was correct within one data point. With
error bars considered, the plot is continuous with no major jumps in the data.
The growth rate itself varies by a factor of $\sim 2-3$, between $\sim
0.09-0.24$ $s^{-1}$. The scaled west ribbon scales provide a good alignment
with the east ribbon, and the location of the three peaks match reasonably
well.

We also plot the duration of exponential growth at each spatial scale for the
east and west ribbon in Figure 4 (bottom). Examining the east ribbon data, we
see exponential growth start initially at a single spatial scale (the 1.75 Mm
scale plotted previously), before beginning at all other scales up to 19
seconds later. This is suggestive of a process at a specific spatial scale
causing the growth at progressively shorter and longer scales through a
cascade and inverse cascade. The west ribbon produces noisier data (due to
fewer data points along the ribbon) - but still shows exponential growth
beginning first at the same spatial scale of 1.75 Mm. The end times of
exponential growth are not captured as clearly by our algorithm. Detecting key
spatial scales independently in both ribbons provides confidence that the
processes causing the features are likely linked - that is, they both
originate from instability processes at the reconnection site.

The east ribbon does not trace the strongest photospheric field of the active
region, resulting in a fairly low mean LoS field strength of 107.5 Gauss. We
do not expect the magnetic field of active region loops (Brooks et al., 2021;
Landi et al., 2021) to be significantly lower than this photospheric
measurement. Considering again flux conservation, and that the length scales
scale with the square root of the field strength ratios - the spatial scales
of these processes at the reconnection site (in the active regions loops) will
not be much larger than the scales observed in the flare ribbons (of order
1-10 Mm).

##  5 Discussion

Figure 5:  Left: Comparison between the non-thermal velocity timing, with the
start an end of exponential growth at different spatial scales (in the east
ribbon). Vertical ticks on the x axis mark the times shown in the right panel.
Right: East ribbon spatial power spectra at key times, compared with a power-
law distribution.

Our analysis of intensity variations along the flare ribbons give an insight
into the dynamics of the reconnecting current sheet. Our analysis has provided
us with observational constraints, allowing us to compare the observed
parameters with those predicted by reconnection instability theory. Our main
constraints are -

  * •

Exponential growth is seen across all wavelengths, with preferred scale at
1.75 Mm.

  * •

Exponential growth rate is in the order of $\approx 0.1-0.2$ $s^{-1}$.

  * •

Spatial scale growth times suggest a cascade and inverse cascade.

In (Tenerani & Velli, 2020), the authors use 2.5D simulations to study the
non-linear phase of tearing modes within a current sheet. The authors find
growth beginning at a specific wave mode, before leading to growth at lower
and higher modes (higher and lower spatial scales) through a simultaneous
cascade and inverse cascade. The combination of cascades to higher and lower
spatial scales simultaneously comes from the interplay of magnetic island
collapse and coalescence. The exponential growth rates in Tenerani & Velli
(2020), $\approx 0.45$ are also on the same order of magnitude as our
observations. The specific properties of the tearing mode instability within
this simulation are consistent with our observations. In particular, the
cascade/inverse cascade, which we would not expect from other instability
processes (such as Kelvin Helmholtz), suggests there is more at play than just
turbulence. We therefore conclude that the flare ribbon spatial scale growth
reflects the presence of the tearing mode instability in the current sheet.

With a higher spatial resolution, or larger flare with more data points along
the ribbon - it would be possible to measure the power at further larger and
smaller wavelengths. This would potentially allow for the measurement of a
more complete dispersion relation, including the time variation ($\omega$) and
second $k$ (variation across the ribbon) component, to allow further
comparison with tearing mode theory.

###  5.1 Flare timeline

Assuming that exponential growth at our key spatial scale is the signature of
tearing mode instability onset, we can compare this with the other observed
parameters for this same event. Jeffrey et al. (2018) find that the rise in
non-thermal velocity, believed to be a signature of turbulence at the flare
ribbon slit location, preceded the intensity increase associated with plasma
heating and flare onset. The authors note that there must be a driver for this
long-lived turbulence signal. High-cadence IRIS observations of other small
flares have found a similar pattern of enhanced non-thermal velocities
preceding rises intensity - consistent with signatures of turbulent
reconnection in the current sheet (Chitta & Lazarian, 2020). Comparing the
timing of the non-thermal velocities with the growth of spatial scales in
Figure 5 (left), we find growth at our key scale of 1.75 Mm begins some 15
seconds before the growth in non-thermal velocities. Similarly, the end of
exponential growth across all spatial scale occurs a similar time frame before
the drop of non-thermal velocity. This suggests that for this event, the
tearing mode instability is the driver of plasma turbulence in the region, as
a result of the cascade and inverse cascade to progressively smaller and
larger spatial scales.

Further evidence to support the breakdown to turbulence via a(n) (inverse)
cascade can be found by examining the power spectra for specific time frames
(taking a vertical cross-section through the middle row of Figure 3). Cross-
sections at important moments are shown in Figure 5 (right), and are as
follows:

  * •

10:36:50 - Large oscillations in the low-amplitude power spectra, equivalent
to the harmonics of the bright points in the pre-flare ribbon region.

  * •

10:37:05 - Key scale of 1.75 Mm appears for the first time.

  * •

10:37:13 - Exponential growth begins at the key scale (shown in left panel of
Figure 5).

  * •

10:37:31 - Exponential growth is well underway across all spatial scales, lead
by the key 1.75 Mm scale. The non-thermal velocity begins to rise at this
time.

  * •

10:37:51 - Power spectra have reached their peak, with all scales catching up
to the key scale to produce a power law of $\approx k^{-2.3}$.

  * •

10:38:00 onwards - Power law starts to break down, with a falling amplitude
and local maxima/minima forming. Non-thermal velocity is also returning to
background levels.

At the peak of the power spectra, the spatial scales form a power law of
approximately 2.3. A power law in this domain is a key prediction in plasma
turbulence, routinely seen in in-situ measurements. In the classic Kolmogorov-
like distribution of turbulence (Kolmogorov, 1941), the power law is predicted
to evolve to and end-state of 5/3 - slightly lower than our observed value of
$\approx 2.3$. There are two likely explanations for this disparity. Firstly,
it is possible that the tearing mode process did not persist for long enough
to drive the power law to the predicted Kolmogorov end-state of 5.3. Or, an
alternative turbulence model is a better fit for these conditions. Tenerani &
Velli (2020), for example, predict a power-law index of 2.2, when introducing
a model of intermittency (Frisch et al., 1978) to the collapsing magnetic
islands and beginning from a Kolmogorov state. Dong et al. (2018) also find an
index of 2.2, in their simulation of a plasmoid-mediated regime.

By comparing the power spectrum and non-thermal velocity evolution with the
timing of exponential growth onset, we have provided evidence for the tearing-
mode instability manifesting itself before plasma turbulence - similar to that
expected by simulation work, e.g. Dong et al. (2018); Tenerani & Velli (2020).
As shown in Jeffrey et al. (2018), the turbulence itself precedes the ribbon
intensity enhancements and flare onset. We thus have extended our knowledge of
the pre-flare timeline for this single event, and provided insight into the
nature of current-sheet collapse and onset of fast-reconnection via MHD
instabilities in a solar flare.

##  6 Conclusion

A major open question in solar flare physics is the sequence of events that
lead to the onset of the fast reconnection that is widely acknowledged to
cause the impulsive energy release. Central to this question is the formation
and disruption of the current sheet in which reconnection occurs. Current
sheet formation can occur through a variety of processes that include large
scale MHD instabilities that are a key element of eruptive flare models (e.g.
Green et al., 2018, and references therein), but also as the result of MHD
turbulence. The subsequent evolution of that turbulence can then lead to its
disruption through the onset of the tearing mode instability, which causes the
energy spectrum to steepen (Dong et al., 2018). The issue of the interplay and
feedback between reconnection, turbulence and current sheet disruption is
therefore complex, but critical to understanding flare onset.

For the event studied here we find evidence to support the development of the
tearing mode instability prior to the first detectable signatures of
turbulence, followed by the subsequent evolution of the energy spectrum to a
power-law dominated regime with slope $\approx 2.3$. While we cannot
completely rule out the existence of turbulence prior to the tearing mode
onset that we are just unable to detect with the available data, our analysis
strongly suggests current sheet formation by another mechanism, e.g. at the
interface between emerging and pre-existing magnetic fields. The subsequent
evolution of the instability and turbulent energy spectrum is consistent with
previous work by Dong et al. (2018) and Tenerani & Velli (2020) that supports
the presence of the combined effects of dynamic alignment of the magnetic
islands (or plasmoids), intermittency and recursive reconnection. We note that
recent work by French et al. (2019, 2020) provides an important new tool for
quantifying the internal sub-structure of current sheets that can be used to
further probe alignment and intermittency.

In conclusion, we find that for case study presented here the flare onset is
precipitated by current sheet formation, leading to the onset of the tearing
mode instability and subsequent development of a turbulent cascade.

IRIS is a NASA small explorer mission developed and operated by LMSAL with
mission operations executed at NASA Ames Research Center and major
contributions to downlink communications funded by ESA and the Norwegian Space
Centre. R.J.F. thanks the STFC for support via funding from the PhD
Studentship, ST/S50578X/1. S.A.M was supported by STFC Consolidated Grant (CG)
ST/S000240/1 and Hinode Operations (UKSA) ST/S006532/1. I.J.R was supported in
part by ST/V006320/1. AWS was supported by CG ST/S000240/1, and NERC
NE/P017150/1 and NE/V002724/1. The authors thank Gherardo Valori, Lidia van
Driel-Gesztelyi, Clare Watt and Philip Judge for their helpful discussions.

## References

  * Akasofu (1964) Akasofu, S. I. 1964, Planet. Space Sci., 12, 273
  * Akasofu & Lepping (1977) Akasofu, S. I., & Lepping, R. P. 1977, Planet. Space Sci., 25, 895
  * Biskamp (1986) Biskamp, D. 1986, Physics of Fluids, 29, 1520
  * Brannon et al. (2015) Brannon, S. R., Longcope, D. W., & Qiu, J. 2015, ApJ, 810, 4
  * Brooks et al. (2021) Brooks, D. H., Warren, H. P., & Landi, E. 2021, ApJ, 915, L24
  * Brosius & Daw (2015) Brosius, J. W., & Daw, A. N. 2015, ApJ, 810, 45
  * Carmichael (1964) Carmichael, H. 1964, NASA Special Publication, 50, 451
  * Cheng et al. (1983) Cheng, C. C., Oran, E. S., Doschek, G. A., Boris, J. P., & Mariska, J. T. 1983, ApJ, 265, 1090
  * Cheng et al. (2018) Cheng, X., Li, Y., Wan, L. F., et al. 2018, ApJ, 866, 64
  * Chitta & Lazarian (2020) Chitta, L. P., & Lazarian, A. 2020, ApJ, 890, L2
  * De Pontieu et al. (2014) De Pontieu, B., Title, A. M., Lemen, J. R., et al. 2014, Sol. Phys., 289, 2733
  * Dong et al. (2018) Dong, C., Wang, L., Huang, Y.-M., Comisso, L., & Bhattacharjee, A. 2018, Phys. Rev. Lett., 121, 165101
  * Doschek et al. (1983) Doschek, G. A., Cheng, C. C., Oran, E. S., Boris, J. P., & Mariska, J. T. 1983, ApJ, 265, 1103
  * Forbes & Lin (2000) Forbes, T. G., & Lin, J. 2000, Journal of Atmospheric and Solar-Terrestrial Physics, 62, 1499
  * French et al. (2019) French, R. J., Judge, P. G., Matthews, S. A., & van Driel-Gesztelyi, L. 2019, ApJ, 887, L34
  * French et al. (2020) French, R. J., Matthews, S. A., van Driel-Gesztelyi, L., Long, D. M., & Judge, P. G. 2020, ApJ, 900, 192
  * Frisch et al. (1978) Frisch, U., Sulem, P. L., & Nelkin, M. 1978, Journal of Fluid Mechanics, 87, 719
  * Graham & Cauzzi (2015) Graham, D. R., & Cauzzi, G. 2015, ApJ, 807, L22
  * Green et al. (2018) Green, L. M., Török, T., Vršnak, B., Manchester, W., & Veronig, A. 2018, Space Sci. Rev., 214, 46
  * Henderson (1994) Henderson, M. G. 1994, PhD thesis, UNIVERSITY OF CALGARY (CANADA).
  * Henderson (2009) —. 2009, Annales Geophysicae, 27, 2129
  * Hirayama (1974) Hirayama, T. 1974, Sol. Phys., 34, 323
  * Hones et al. (1984) Hones, E. W., J., Baker, D. N., Bame, S. J., et al. 1984, Geophys. Res. Lett., 11, 5
  * Jeffrey et al. (2018) Jeffrey, N. L. S., Fletcher, L., Labrosse, N., & Simões, P. J. A. 2018, Science Advances, 4, 2794
  * Jing et al. (2016) Jing, J., Xu, Y., Cao, W., et al. 2016, Scientific Reports, 6, 24319
  * Kalmoni et al. (2015) Kalmoni, N. M. E., Rae, I. J., Watt, C. E. J., et al. 2015, Journal of Geophysical Research (Space Physics), 120, 8503
  * Kalmoni et al. (2018) —. 2018, Nature Communications, 9, 4806
  * Kolmogorov (1941) Kolmogorov, A. 1941, Akademiia Nauk SSSR Doklady, 30, 301
  * Kopp & Pneuman (1976) Kopp, R. A., & Pneuman, G. W. 1976, Sol. Phys., 50, 85
  * Landi et al. (2021) Landi, E., Li, W., Brage, T., & Hutton, R. 2021, ApJ, 913, 1
  * Lemen et al. (2012) Lemen, J. R., Title, A. M., Akin, D. J., et al. 2012, Sol. Phys., 275, 17
  * McPherron (1970) McPherron, R. L. 1970, J. Geophys. Res., 75, 5592
  * Melrose (1997) Melrose, D. B. 1997, ApJ, 486, 521
  * Pikelner & Livshits (1977) Pikelner, S. B., & Livshits, M. A. 1977, Soviet Ast., 21, 601
  * Pontin (2011) Pontin, D. I. 2011, Advances in Space Research, 47, 1508
  * Priest (1985) Priest, E. R. 1985, Reports on Progress in Physics, 48, 955
  * Rae et al. (2010) Rae, I. J., Watt, C. E. J., Mann, I. R., et al. 2010, Journal of Geophysical Research (Space Physics), 115, A10222
  * Schou et al. (2012) Schou, J., Scherrer, P. H., Bush, R. I., et al. 2012, Sol. Phys., 275, 229
  * Shaikh & Zank (2010) Shaikh, D., & Zank, G. P. 2010, MNRAS, 402, 362
  * Sturrock (1968) Sturrock, P. A. 1968, in IAU Symposium, Vol. 35, Structure and Development of Solar Active Regions, ed. K. O. Kiepenheuer, 471
  * Tenerani & Velli (2020) Tenerani, A., & Velli, M. 2020, MNRAS, 491, 4267
  * Wang & Liu (2012) Wang, H., & Liu, C. 2012, ApJ, 760, 101
  * Warren et al. (2018) Warren, H. P., Brooks, D. H., Ugarte-Urra, I., et al. 2018, ApJ, 854, 122


# Pore-scale analysis of gas injection in gas-condensate reservoirs

Paula Kozlowski Pitombeira Reis  Marcio da Silveira Carvalho

###### Abstract

Condensate banking around wellbores can significantly shorten the production
from gas-condensate reservoirs. Different approaches to mitigate this issue
have been proposed in the literature, among which gas injection comes out with
promising results. With this method, pressure maintenance and condensate re-
vaporization can be achieved, lessening the flow blockage caused by liquid
dropout and accumulation in the porous medium. While gas injection in gas-
condensate reservoirs has been largely investigated at the meso and macro
scales, data regarding the method’s efficiency at the micro scale are scarce.
Therefore, the effects of local changes in gas and condensate properties
stemmed from the interaction between injected and reservoir fluids at the
pore-scale are poorly understood. In order to evaluate how these changes
affect the transport in porous media, a compositional pore-network model was
used to reproduce gas injection in a sandstone sample following condensate
accumulation. $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ and produced gas were tested
as the candidates for condensate banking remediation at different pressure
levels. After gas flooding, condensate saturation, heavy component recovery
and gas relative permeability were quantified to appraise the achieved gas
flow improvement. The results indicated that $C_{2}$ and $CO_{2}$ were the
most effective gases to clear the accumulated condensate and re-establish the
gas flow. Conversely, $C_{1}$ and $N_{2}$, especially mixed with the produced
fluids, displayed the least favorable results, and could even lead to gas flow
impairment.

###### keywords:

Pore-scale modeling , Gas-flooding , Condensate recovery

[inst1]organization=Department of Mechanical Engineering, Pontifícia
Universidade Católica do Rio de Janeiro,addressline=Rua Marques de Sao
Vicente, 225, city=Rio de Janeiro, postcode=22451-900, state=RJ,
country=Brazil

##  1 Introduction

Production from gas-condensate reservoirs can be significantly affected by
condensate banking. As the pressure around producing wells is lowered below
the reservoir mixture’s dew point, a liquid phase emerges in the porous medium
and partially blocks the gas flow. Several parameters affect the severity of
the damage, such as the reservoir permeability, fluid composition and
depletion level. Even lean gas mixtures can lead to significant flow blockage
once liquid dropout takes place [1]. The liquid and gas phases different
mobilities and molar contents induce a compositional shift in the two-phase
region, which enlarges the local heavier component fractions and can lead to
substantial liquid saturations [2].

In order to manage this production challenge, various production strategies
and EOR methods have been devised for gas-condensate reservoirs development
[3, 4, 5]. Among these methods, many involve the injection of gases in the
reservoir, aiming full or partial pressure maintenance, as well as re-
evaporation of accumulated condensate.

In full pressure maintenance, the reservoir pressure is kept above the dew
point pressure and liquid dropout is prevented. If feasible, this method can
maximize condensate recovery, as no liquid phase is deposited in the porous
medium and heavy components are produced in the gas. Luo et al. [6] conducted
long-core flooding experiments with lean gas injection to compare recoveries
above and below the dew point pressure. While injection below the dew point
managed to re-vaporize a fraction of both intermediate and heavy hydrocarbons
in the core, condensate recovery was considerably higher for the case above
the dew point. Implementing a full pressure maintenance scheme, however,
requires large amounts of injected gas and can be unpractical for cases in
which the initial reservoir pressure is close to the fluids dew point [6].
Also, special attention has to be directed to the injected gas composition.
Certain gases, such as nitrogen[6] and methane [7], have a tendency to
increase gas-condensate mixtures dew point and can lead to early condensation
in the reservoir.

For the case of partial pressure maintenance, gas can be injected following
liquid dropout in the reservoir. In this approach, the objectives are both
slowing down pressure depletion and re-vaporizing the accumulated condensate,
so that gas flowing paths are cleared and valuable heavy components are
recovered. Numerous studies have investigated the ability of different gas
compositions to support reservoir pressure and reduce liquid banking. Al-Anazi
et al. [8] evaluated the effectiveness of methane to recover condensate using
coreflooding experiments. They concluded that methane is able to re-vaporize
the liquid phase and restore the core’s original permeability, but that tens
to hundreds of injected pore volumes can be required for the method’s success.
High injection pressures and flow rates were pointed out as means to
accelerate the method’s outcomes. Gachuz-Muro et al. [9] performed core
flooding experiments to investigate the recovery obtained with gas injection
in naturally fractured gas-condensate reservoirs. A small gap was allowed
around the core to simulate a fracture and, after the core was depleted below
the dew point pressure, $N_{2}$, $CO_{2}$ and lean gas were injected. The
results indicated a condensate recovery of $51.7\%$ with the lean gas
injection, $34.78\%$ with $CO_{2}$ and only $18.7\%$ with the $N_{2}$
injection case. Al-Abri et al. [10] evaluated condensate recovery by the
injection of super-critical $CO_{2}$, $C_{1}$, and their mixtures. They
performed gas injection in a sandstone core following condensate flooding,
meaning that the liquid content was not established in the porous medium by
condensation. Unsteady-state relative permeability and recovery measurements
demonstrated that the injection of pure $CO_{2}$ delayed the gas breakthrough
and provided higher recovery of condensate when compared to the injection of
mixtures of $CO_{2}$ and pure $C_{1}$, at the same operating conditions.
Besides conventional gas flooding, huff-n-puff gas injection has also been
investigated experimentally, for application in shale gas-condensate
reservoirs. Due to shale ultra low permeabilities, this technique could be
beneficial to expedite both pressure boost and re-evaporation of condensate in
the producing wellbore vicinity [11]. Core flooding experiments in shale cores
displayed promising results with the huff-n-puff injection of produced gas
[12], $C_{1}$ [13] and $CO_{2}$ [14].

Additionally to core flooding tests, several reservoir-scale numerical studies
have been conducted to assess gas injection as an enhanced gas-condensate
recovery method. Marokane et al. [15] used a full-field compositional
reservoir simulation model to investigate the efficacy of a one-time lean gas
injection scheme to remove the condensate banking around producing wells.
Results indicated that, for a lean gas reservoir, maximum recovery was
obtained for injection at pressures below the maximum liquid dropout. For a
rich gas reservoir, however, starting gas injection before reaching the
maximum liquid dropout pressure produced better results. Linderman et al. [16]
used a similar approach to evaluate the feasibility of injecting nitrogen in a
gas-condensate reservoir. Their model predicted superior net gas production
with the injection of nitrogen than with the injection of produced gas, while
the condensate net recovery was virtually identical for the two scenarios.
They also compared the injection of $CO_{2}$ with $N_{2}$, finding no
significant differences in the achieved recoveries. Taheri et al. [17] modeled
enhanced condensate recovery in fractured gas-condensate reservoirs using
stock tank gas, $C_{1}$, $N_{2}$ and $CO_{2}$ injection. Results suggested
that $CO_{2}$ injection could generate considerably better results that the
other tested gases, and that no appreciable difference is obtained between the
injections of methane, nitrogen and stock-tank gas. It was also concluded that
injection at higher pressures boosted condensate recovery. Another numerical
study evaluating gas injection in fractured gas-condensate reservoirs was
conducted by Fath et al. [18]. After testing four different injected gas
compositions, they obtained the highest recovery with $CO_{2}$, followed by
$C_{1}$, produced gas and $N_{2}$. Their model also predicted that, for the
same injected volumes of the tested gases, nitrogen led to a significant
higher pressure average in the reservoir, while produced gas led to the
lowest. Recent numerical works also investigated the efficacy of $C_{2}$
injection for condensate recovery, with very favorable results. Sharma and
Sheng [19] evaluated the performance of huff-n-puff gas and solvent injection
in shale gas-condensate reservoirs. They concluded that ethane was the best
injection fluid on accounts of higher and much faster recovery, when compared
to methane, methanol and isopropanol. The positive results were attributed to
ethane’s ability to greatly reduce the gas-condensate mixture’s dew point
pressure, ensuring relatively low injecting volumes for condensate re-
vaporization. Zhang et al. [20] investigated the injection of $C_{1}$,
$C_{2}$, $CO_{2}$ and $N_{2}$ in a five-spot configuration for condensate
enhanced recovery. In their study, $C_{2}$ also led to the best condensate
recovery, of $56.8\%$ after the injection of 0.85 pore volumes, closely
followed by $CO_{2}$, with a recovery of $50.7\%$. Methane and nitrogen
displayed considerably worse performances, with $22.8\%$ and $21.3\%$ of
condensate recovery, respectively, for the same injected volume.

Although both experimental and numerical research indicate that gas injection
is an advantageous approach for condensate enhanced recovery, no micro-scale
analysis of the method has been reported in the literature. This represents a
significant gap in data for gas injection performance evaluation, especially
considering that pore-scale events are essential to understanding macro-scale
transport properties in porous media. During gas injection, local changes in
composition can alter significantly both bulk and interfacial properties of
gas and liquid phases, which in turn affect their flow, characterized by the
relative permeability curves. To address this gap, we used a compositional
pore-network model to evaluate gas injection in porous media after condensate
accumulation. The model has been validated against core flooding experiments
[21] and also used to evaluate wettability alteration as a gas-condensate
enhanced recovery method [22]. In the present study, first we replicated
condensate accumulation in the porous medium by flowing a representative gas-
condensate mixture through a sandstone-based network at different depletion
levels, starting from values just below the dew point, until pressures below
the maximum liquid dropout. Then $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ or
produced gas were injected in the condensate bearing network and the flow
improvement was evaluated. Final saturations, relative permeabilities and
recovery of heavy hydrocarbon components were quantified to compare the
efficacy of each injection scenario. The impact of the concentration of the
injected gases was also explored, by repeating the analysis with the injection
of mixtures containing $50\%$ in moles of $C_{1}$, $C_{2}$, $CO_{2}$ or
$N_{2}$, and $50\%$ the produced gas.

##  2 Pore-network model

The pore-network model used in the preset study was devised specifically to
represent gas-condensate flow in porous media. It was based on the fully
compositional pore-network model proposed by Santos and Carvalho [23] and
adapted to encompass realistic characteristics of in-situ condensate formation
and gas-condensate flow. In the following sections, the main aspects of gas-
condensate flow in porous media and their representation in the model are
described. For further details, we refer to Reis and Carvalho [21].

###  2.1 Pore-scale gas-condensate flow patterns

Gas-condensate mixtures undergo a transition from single-phase gas flow to
two-phase gas-condensate flow, as the dew point pressure is reached. Upon
condensation, it has been observed in micro model experiments [24, 25] that
the liquid phase wets completely the porous medium and tends to form a film
that flows adjacent to the walls, while the gas flows along the pore centers,
in a core-annular pattern. This flow configuration can be disrupted, however,
if the film becomes unstable and the snap-off of the gas happens. In this
case, a bridge of liquid is formed across the gas flowing channel and blocks
the flow. Once blocked, a balance between viscous and capillary forces
determines whether the channel is re-opened to gas flow. This intermittent
opening and closing of gas flowing channels is pointed out as the main
mechanism behind condensate blockage [26].

In our model, we represented this mechanism as illustrated in Figure 1. Figure
1(a) represents a flowing channel above the dew point pressure, where only gas
flows. Figure 1(b) depicts the flowing channel below the dew point pressure.
In this case, the phases flow in core-annular configuration, with condensate
at the walls and gas at the center. Figure 1(c) contains the phases
configurations once the snap-off of the gas happens and the bridge of liquid
is formed. Finally, Figure 1(d) illustrates the re-opening of the flowing
channel, with the mobilization of the condensate bridge.

Figure 1: Evolution of condensate configuration in a pore throat [21]

###  2.2 Pore Geometry

Networks of pore bodies connected by converging-diverging circular pore
throats were used in the model. Equation 1 expresses the radius of a pore
throat $j$, with length $L_{j}$, maximum radius $R_{max,j}$ and constricted
radius $R_{min,j}$. Control volumes for each pore $i$ were defined
encompassing half the volumes of the throats connected to them, as illustrated
in Figure 2. The constricted pore throat geometry was chosen so that an
accurate criterion for the snap-off of the gas phase could be implemented.
Apart from that, it allows the calculation of unique and continuously varying
capillary pressure as the the liquid bridge meniscii move in a pore throat
[27] (e.g. Figure 1(c) and (d) ). Consequently, the conditions for gas flow
blockage could be adequately represented.

$\displaystyle r_{j}(x)=\sqrt{a_{j}+b_{j}x^{2}}$ (1a) $\displaystyle
a_{j}=R_{min,j}^{2}$ (1b) $\displaystyle
b_{j}=\left(\frac{2}{L_{j}}\right)^{2}(R_{max,j}^{2}-R_{min,j}^{2})$ (1c)

Figure 2: Definition of a pore volume [22]

###  2.3 Fluid conductances and snap-off criteria

Regarding the two-phase flow calculations, it was considered that pore bodies
impose no restriction to the flow, so that conductances were devised only for
pore throats. For the annular flow developed after condensate dropout, the
conductances for the gas and liquid phases are presented in equations 2(a) and
2(b), respectively.

$\displaystyle g_{g}$
$\displaystyle=\frac{\pi}{8\mu_{g}L}(S_{g}R^{2})^{2}+\frac{\pi}{4\mu_{l}L}(R^{2}-S_{g}R^{2})S_{g}R^{2}$
(2a) $\displaystyle g_{l}$
$\displaystyle=\frac{\pi}{8\mu_{l}L}(R^{2}-S_{g}R^{2})^{2}$ (2b)

where, $\mu_{g}$ and $\mu_{l}$ are the viscosities of the gas and liquid
phases, $S_{g}$ is the gas saturation and $R$ and $L$ are, respectively, the
equivalent radius[28, 21] and length of the capillary.

The snap-off of the gas phase takes place in a pore throat once a critical
thickness of the liquid film is achieved [21] and the criterion proposed by
Beresnev et al. [29] is met. This criterion, given by equation 3, takes into
consideration the length, maximum and minimum radii of the throat, as well as
the condensate film thickness ($e$).

$L>2\pi\sqrt{(R_{min}-e)(R_{max}-e)}$ (3)

For the flow through a pore throat to be reestablished after the snap-off, the
pressure drop across the throat has to overcome the critical value given by
equation 4. It is a function of the interfacial tension, $\sigma$, between the
phases and the radii of the condensate bridge meniscii, $R_{1}$ and $R_{2}$.
The interfacial tension is calculated with the correlation proposed by Weinaug
and Katz [30], and the radii vary according to the throat geometry and the
condensate saturation.

$\Delta P_{crit}=2\sigma\left(\frac{1}{R_{1}}-\frac{1}{R_{2}}\right)$ (4)

###  2.4 Governing Equations

A system of non-linear equations relates the molar content and pressure of
each pore control volume $i$ in the model. It encompasses molar balance
equations (eq. 5), volume consistency equations (eq. 6), and equations to
enforce boundary conditions.

Equation 5(a) describes how the number of moles in a pore volume $i$ varies
due to the molar flow rate through its adjacent pore throats,
$\dot{n}^{k}_{j}$, and the flow at the network boundaries, $\dot{s}_{i}^{k}$.
In this equation, $c_{ij}$ are the entries of the incidence matrix $C$, which
maps which throats are connected to each pore volume in the network. The molar
flow through a pore throat, given by equation 5(b), converts the volumetric
flow, calculated with the conductances and the pressure drop, into the molar
flow, using the molar fraction of each component in the gas and liquid phases,
$y^{k}$ and $x^{k}$, and their molar densities, $\xi_{g}$ ans $\xi_{l}$.
Equation 5(c) represents the pressure difference that drives the flow through
the throats. The inclusion of the interface pressure difference ($\Delta
P_{j}^{int}$) in eq. 5(c) is controlled by the parameter $H^{int}$, which is
equal to $1$ when a liquid bridge is present and $0$, otherwise. In Equations
5, $i=1..n_{b}$, $j=1..n_{t}$, and $k=1..n_{c}$, represent, respectively, the
number of pore volumes, pore throats and fluid mixture components in a
network.

$\frac{\partial N^{k}_{i}}{\partial
t}=-\displaystyle\sum_{j=1}^{n_{t}}c_{ij}\dot{n}_{j}^{k}+\dot{s}_{i}^{k}$ (5a)

$\dot{n}_{j}^{k}=(y^{k}\xi_{g}g_{g}+x^{k}\xi_{l}g_{l})_{j}\Delta P_{j}\\\ $
(5b)

$\Delta P_{j}=\displaystyle\sum_{m=1}^{n_{b}}c_{mj}P_{m}-H_{j}^{int}\Delta
P_{j}^{int}$ (5c)

The volume consistency equations, given by eq. 6, are used to enforce
compatibility between the pore volumes and the volumes of the phases contained
in them. These equations are devised considering slightly compressible
networks, so that, for a given pressure $P_{i}$, the volume of a pore can be
approximated using the pore compressibility, $\nu_{i}$, along with reference
volume, $\overline{V_{i}}$, and pressure, $\overline{P_{i}}$ values. For the
calculation of the gas and liquid phases volumes, the pressure and temperature
in the pore, $P_{i}$ and $T$, are related with the fluid parameters
$\mathcal{L}_{i}$, $Z_{i}^{g}$, $Z_{i}^{l}$, $x_{i}^{k}$ and $y_{i}^{k}$
(respectively: the fraction of the $N_{i}$ moles in the liquid phase, the
compressibility factors and the molar fractions of each component $k$ in the
gas and liquid phases).

$N_{i}-\frac{\overline{V_{i}}[1+\overline{\nu_{i}}(P_{i}-\overline{P})]}{\mathcal{L}_{i}\left(\frac{Z_{i}^{l}RT}{P_{i}}-\sum_{k=1}^{n_{c}}v_{k}x_{i}^{k}\right)+(1-\mathcal{L}_{i})\left(\frac{Z_{i}^{g}RT}{P_{i}}-\sum_{k=1}^{n_{c}}v_{k}y_{i}^{k}\right)}=0$
(6)

Additionally to equations 5 and 6, simple equations are written to enforce the
boundary conditions for the system. In the presented model, the flow can be
driven either by imposing different pressure levels at the inlet and the
outlet, or by prescribing molar flow rate at the inlet and pressure at the
outlet. Other parameters that have to be imposed are the composition of the
fluid injected in the network and the temperature $T$.

The presented system of equations was based on the compositional formulation
proposed by Collins et al. [31], which decouples the flow equations from the
flash calculations. In this formulation, the system relating flow equations is
solved with the Newton-Raphson method, while thermodynamic equilibrium is
enforced at each newton iteration by performing phase equilibrium calculations
using the Peng-Robinson EoS and the trial solution variables. Details of the
implemented phase equilibrium calculations can be found in Santos and Carvalho
[23].

##  3 Results

###  3.1 Pore-network

The pore-network used in the present study, shown in Figure 3, represented a
$1.75mm^{3}$ cubic section of a sandstone with permeability of $169mD$ and
porosity of $17.1\%$. This network was first extracted from Micro-CT imaging
of a sandstone sample, using the network extraction algorithm proposed by Dong
[32], then adapted to meet the geometric criteria described in section 2.2.
Figure 4 represents the size distributions of the constricted and
unconstricted capillary radii of the adapted pore-network. More details of the
network adaptation process can be found in Reis and Carvalho [22].

Figure 3: Sandstone based network Figure 4: $R_{min}$ and $R_{max}$
distributions in the used pore-network

###  3.2 Fluids

A series of fluid mixtures was used in the present study in order to represent
condensate buildup in porous medium followed by gas injection and condensate
re-vaporization.

The composition of the fluid used during the condensate accumulation process
in the porous medium is presented in Table 1. This mixture represents a
typical fluid found in gas-condensate reservoirs, composed by carbon dioxide,
nitrogen, light, intermediate and heavy hydrocarbons. It exhibits a large
window of retrograde condensation behavior in the temperature range from
$12^{\circ}C$ to $207^{\circ}C$, as shown in the phase diagram in Figure 5.
For our analyses, a temperature of $60^{\circ}C$ was chosen as the reservoir
temperature.

Table 1: Gas-condensate mixture composition. Component | Molar Fraction
---|---
$CO_{2}$ | 0.05
$N_{2}$ | 0.02
$C_{1}$ | 0.65
$C_{2}$ | 0.13
$C_{3}$ | 0.07
$C_{6}$ | 0.05
$C_{10}$ | 0.025
$C_{16}$ | 0.005
Figure 5: Fluid mixture phase envelope

After condensate accumulation, the gases chosen as candidates for injection
aiming condensate enhanced recovery were produced gas obtained during the flow
of the mixture presented in Table 1, methane ($C_{1}$), ethane ($C_{2}$),
carbon dioxide ($CO_{2}$) and nitrogen ($N_{2}$). For the sake of appraising
preliminarily the ability of the four last gases to re-vaporize condensate,
Figure 6 illustrates the liquid dropout as a function of the pressure obtained
by mixing $C_{1}$, $C_{2}$, $CO_{2}$ or $N_{2}$ with the composition presented
in Table 1, at different molar fractions and a temperature of $60^{\circ}C$.
It can be noticed from Figure 6(a) that the addition of $25\%$ in moles of the
tested gases to the reservoir composition leads to a substantial reduction in
maximum liquid dropout, varying from $31.8\%$, in the case of $C_{2}$, to
$49.6\%$, in the case of $C_{1}$. It is also noticeable that the effects of
the tested gases on the mixture’s dew point pressure are very distinct.
$C_{2}$ and $CO_{2}$ reduce the saturation pressure, thus allowing the mixture
to remain as a single gas phase at lower pressures. The opposite is verified
as the reservoir fluid is mixed with $C_{1}$ and $N_{2}$. While the negative
effect of $C_{1}$ on the dew point pressure is mild, $N_{2}$ increases
significantly its value, and can lead to very early condensation in porous
media. The same trend of effects is observed as the content of the tested
gases is increased to $50\%$, Figure 6(b), and $75\%$, Figure 6(c). The
maximum liquid dropout is progressively reduced, with $C_{1}$ and $CO_{2}$
being the most, and $C_{2}$ and $N_{2}$ the least effective gases to lower the
maximum volume of liquid formed. As for the effects on the dew point pressure,
at higher molar fractions $C_{2}$ and $CO_{2}$ gradually reduce it, while $N2$
increases it. $C_{1}$, on the other hand, slightly augments and then reduces
the mixture’s dew point pressure, at the molar fractions of $50\%$ and $75\%$,
respectively. These results suggest that mixtures rich in $C_{1}$, $C_{2}$,
$CO_{2}$ and $N_{2}$ have a positive prospect to work as condensate re-
vaporization agents and, consequently, improve gas flow in porous media.

(a) 75% of the gas-condensate mixture and 25% of the injected gases

(b) 50% of the gas-condensate mixture and 50% of the injected gases

(c) 25% of the gas-condensate mixture and 75% of the injected gases

Figure 6: Liquid dropout at a temperature of $60^{\circ}C$

Additionally to the injection of $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ or
produced gas individually, mixtures of $C_{1}$, $C_{2}$, $CO_{2}$ and $N_{2}$
with the produced gas, at molar fractions of $50\%$, were also injected in the
networks to assess their ability to recover condensate after liquid banking.
These tests were devised considering two different scenarios. First, to take
advantage of the produced gas availability, we wanted to test whether directly
injecting mixtures of the produced gas with $C_{1}$, $C_{2}$, $CO_{2}$ or
$N_{2}$ could satisfactorily enhance condensate recovery. Secondly, even
during the injection of pure gases, the flow of mixtures of injected and
reservoir fluids will take place, especially at regions far from the injection
wells and/or at small injection volumes. These test, therefore, serve also to
appraise the effects of in situ mixing with the reservoir fluid in the
performance of $C_{1}$, $C_{2}$, $CO_{2}$ and $N_{2}$ to recover condensate.

###  3.3 Flow conditions

Two steps of injection were performed in the networks for every evaluated
case. First, the reservoir fluid was injected so that condensate accumulated
in the porous medium, mimicking the process of condensate banking in the near
wellbore regions. Then, the injected composition in the networks was altered,
so that the liquid buildup was followed by gas injection and condensate
recovery.

The process of condensate accumulation was reproduced by injecting 25 pore
volumes (PV) of the mixture presented in Table 1 in the networks at different
pressures. As the timing for gas injection is one of the most relevant
parameters for the method’s success [15], we wanted to test its performance at
various reservoir depletion stages. Six pressure values were used, 22, 21.5,
21, 19.5, 17.75 and 14.75 MPa, which corresponded to liquid dropout
saturations of $2.3\%$, $9.7\%$, $13.3\%$, $17.1\%$, $17.3\%$ and $14.7\%$,
respectively. With that, a broad range of depletion scenarios could be
covered, from the early stages of condensate formation to pressures below the
maximum liquid dropout.

After the first 25 PV of injection, $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$,
produced gas or their mixtures were injected in the the networks. This step
also lasted until 25 PV were injected. Under certain conditions, the recovery
of heavy components by gas injection can be a very slow process, requiring
tens to hundreds of pore volumes injected [8]. For this reason, in some cases,
we could not reach a steady state flow.

The prescribed boundary conditions for both steps were molar flow rate at the
network inlet and pressure at the outlet. For all tested injection scenarios,
a molar flow rate was imposed so that a gas flowing velocity of $35m/day$ was
achieved. It has been reported in the literature [8] that injection rate does
not influence significantly the performance of gas injection as a gas-
condensate EOR method and we chose, therefore, not to explore the effects of
this parameter.

###  3.4 Condensate saturation reduction

Figures 7(a)-(f) illustrate the evolution of the condensate saturation in the
networks with time, during the injection of the gas-condensate mixture
followed by $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ or produced gas, at different
pressures. The first 25 PV of injection represent the condensate buildup
during the injection of the reservoir fluid, and the last 25 PV represent the
condensate recovery during gas injection. At the pressure of $22MPa$ (Fig.
7(a)), just below the dew point pressure, $C_{2}$ and $CO_{2}$ produce almost
identical effects on the accumulated condensate. Both reduced the liquid
content faster and to a lower saturation than the other tested gases. $C_{1}$
also reduced significantly the condensate saturation, from $13.93\%$ to
$1.08\%$, while $N_{2}$ and the produced gas displayed the lowest capacity to
recover condensate, leaving $3.33\%$ and $2.74\%$ of liquid in the network,
respectively. Similar results were obtained with gas injection at the
pressures of $21.5$ (Fig. 7(b)) and $21MPa$ (Fig. 7(c)). Even though the
amount of accumulated condensate in the networks increased significantly,
corresponding to $27.65\%$ at $21.5MPa$, and $32.06\%$ at $21MPa$, the
injected gases cleared the liquid damage efficiently. This is a good
indicative that, at high pressures, all tested gases could be used to support
reservoir pressure and mobilize accumulated liquid. As the injection pressure
is lowered below $20MPa$, however, the ability of the gases to recover
condensate is clearly reduced, as shown in Figures 7(d) to (f). The rate of
liquid re-vaporization is progressively slowed down during the injection of
all tested gases, indicating that more injected volume is required for the
same volume of recovered condensate, as the reservoir becomes more depleted.
At the lowest tested pressure of $14.75MPa$, all injection scenarios perform
similarly, leaving most of the condensate in the porous medium after 25 PV.
These findings support the hypothesis that the timing for injection of gases
in gas-condensate reservoirs is crucial for effective condensate recovery and
gas flow improvement.

(a) $P=22MPa$

(b) $P=21.5MPa$

(c) $P=21MPa$

(d) $P=19.5MPa$

(e) $P=17.75MPa$

(f) $P=14.75MPa$

Figure 7: Condensate saturation evolution with time in the networks. Injection
of the reservoir mixture followed by $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ or
produced gas.

To further illustrate the effect of pressure in the ability of the tested
gases to recover condensate, Figure 8 displays the $CO_{2}$ molar
concentration in the networks following the 25 PV of $CO_{2}$ injection. From
$22MPa$ down to $19.5MPa$, $CO_{2}$ injection’s micro-sweep efficiency is
gradually reduced, but most of the pore space is cleared from its liquid
content. At $17.75MPa$, the reduction in condensate recovery becomes more
prominent and, at $14.75MPa$, most of the pore space is bypassed by $CO_{2}$
injection.

(a) $P=22MPa$

(b) $P=21.5MPa$

(c) $P=21MPa$

(d) $P=19.5MPa$

(e) $P=17.75MPa$

(f) $P=14.75MPa$

Figure 8: $CO_{2}$ concentration in the networks after 25 PV of $CO_{2}$
injection at different pressures

The condensate saturation evolution during the injection of the gas mixtures
containing $50\%$ in moles of the produced gas and $50\%$ of $C_{1}$, $C_{2}$,
$CO_{2}$ or $N_{2}$ is presented in Figure 9. Additionally, a comparison
between the liquid saturation reductions achieved with the injection of pure
$C_{1}$, $C_{2}$, $CO_{2}$ or $N_{2}$ and their mixtures with the produced gas
is presented in Figure 10.

(a) $P=22MPa$

(b) $P=21.5MPa$

(c) $P=21MPa$

(d) $P=19.5MPa$

(e) $P=17.75MPa$

(f) $P=14.75MPa$

Figure 9: Condensate saturation evolution with time in the networks. Injection
of gases at $50\%$ molar fraction.

When injected mixed with the produced gas, $C_{1}$ and $N_{2}$ were rather
ineffective to reduce the accumulated condensate in the porous medium. At the
highest evaluated pressure of $22MPa$, both gases even enlarged the liquid
content in the networks, as seen in Figure 9(a). While the injection of the
$C_{1}$ mixture led to an increase of $25.69\%$ in the liquid saturation in
the network, the $N_{2}$ mixture enlarged it by $76.26\%$. This result is
related to the negative effect of $C_{1}$ and $N_{2}$ on the reservoir
mixture’s liquid dropout at high pressures, as shown in Figure 6(b). Yet, even
at lower pressures, where the effects of $C_{1}$ and $N_{2}$ on condensate
dropout become positive, very high liquid saturations remained in the networks
after gas injection. This suggests that, to be successful as a condensate
recovery method, the concentration of $C_{1}$ and $N_{2}$ flowing in porous
media should be kept high, as the effects of mixing with the reservoir fluids
can significantly reduce their ability to recover liquid. Contrarily, the
results obtained with the injection of the $C_{2}$ and $CO_{2}$ mixtures
showed very little sensitivity to the dilution with $50\%$ of the produced
gas, especially at high pressures.

(a) $C_{1}$ Injection

(b) $C_{2}$ Injection

(c) $CO_{2}$ Injection

(d) $N_{2}$ Injection

Figure 10: Condensate saturation reductions after the injection of different
gases, both pure or mixed with $50\%$ in moles with the produced gas.

###  3.5 Recovery of heavy components

(a) $C_{6}$ Recovery

(b) $C_{10}$ Recovery

(c) $C_{16}$ Recovery

Figure 11: Recovery per mixture component. Injection of gases at $99\%$ molar
fraction.

Ideally, a gas mixture injected in a gas-condensate reservoir should be able
to re-vaporize not only light, but also medium and heavy components
accumulated in the porous medium during condensate banking. For this reason,
the recoveries of the three heaviest components of the gas-condensate fluid
used in this study, namely hexane, $C_{6}$, decane, $C_{10}$ and hexadecane,
$C_{16}$, were quantified for all tested gas injection scenarios. The
recoveries obtained with the injection of $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$
or produced gas are presented in Figure 11, while the recoveries obtained by
injecting mixtures of these gases with the produced gas are shown in Figure
12.

The results in Figure 11 suggest that, when flowing at sufficiently high
concentrations, $C_{1}$, $C_{2}$ and $CO_{2}$ are able to recover the three
heaviest components evenly. For these cases, the recoveries of each analyzed
component can be directly related to the reduction in liquid saturation shown
in Figure 10. During the injection of $N_{2}$, however, the recovery of
hexadecane was considerably lower than the recoveries of hexane and decane.
Therefore, in this case, during the re-vaporization of condensate, the
remaining liquid in the porous medium becomes particularly rich in $C_{16}$,
impeding its recovery. This indicates that, even when flowing at high
concentrations, nitrogen may be unable to retrieve the heaviest components
accumulated in gas-condensate reservoirs. Finally, among the tested gases, the
lowest recoveries of $C_{6}$, $C_{10}$ and $C_{16}$ were obtained with the
injection of produced gas. In this case, even at high pressures, when the
liquid saturations could be significantly reduced after gas flooding, valuable
heavier components remained trapped in the porous medium.

(a) $C_{6}$ Recovery

(b) $C_{10}$ Recovery

(c) $C_{16}$ Recovery

Figure 12: Recovery per mixture component. Injection of gases at $50\%$ molar
fraction.

Mixing $C_{1}$, $C_{2}$, $CO_{2}$ or $N_{2}$ with the produced gas at equal
molar fractions before injection impacted significantly the recovery of heavy
components, as shown in Figure 12. The injection of the mixture of $C_{1}$ and
produced gas resulted in superior recoveries of hexane then those obtained
with the injection of produced gas only, at all tested pressures. The same
effect was observed for the recovery of decane, at pressures above $22MPa$,
and hexadecane, at pressures above $21MPa$. When compared to the injection of
pure $C_{1}$, however, the mixture of $C_{1}$ and produced gas led to
significant lower recoveries of heavy components. Mixing $N_{2}$ with the
produced gas led to even less favorable results. In these cases, the
recoveries of heavy components were not only significantly lower than those
obtained with the injection of pure $N_{2}$, but also, in most cases, than
those achieved after the injection of produced gas only. Additionally, at high
pressures, injecting the mixture of $N_{2}$ and produced gas led to an
increased amount of trapped $C_{10}$ and $C{16}$ in the porous medium, even
for the cases where the liquid saturation was reduced ($P=21.5MPa$ and
$P=21MPa$). Therefore, in these cases, instead of re-vaporizing heavier
components, the presence of $N_{2}$ in the injected mixture prompted a
migration of $C_{10}$ and $C{16}$ from the gas to the liquid phase.

The injection of the gas mixtures rich in $C_{2}$ and $CO_{2}$ displayed
overall the best results among the tested cases. With the exception of the
injection at the lowest tested pressure, adding $CO_{2}$ and $C_{2}$ to the
produced gas increase significantly the recovery of the heavier components
accumulated in the network. This positive outcome reinforces the hypothesis
arisen in section 3.4 that $C_{2}$ and $CO_{2}$ are more suitable candidates
for condensate enhanced recovery than $C_{1}$ and $N_{2}$.

###  3.6 Gas Relative Permeability

As a way of directly quantifying the gas flow improvement after gas injection,
in this section, we compared the gas relative permeabilities before and after
the treatment. Figure 13 presents the gas relative permeabilities calculated
after the flow of 25 PV of the gas-condensate fluid (labeled as GC Mixture)
through the network and also after the 25 PV of gas flooding, using the
different tested compositions.

As can be seen in Figure 13(a), after the injection of pure $C_{1}$, $C_{2}$,
$CO_{2}$, $N_{2}$ or produced gas, the gas relative permeability through the
network increased significantly, at all testes pressures. At $P\geq 19.5MPa$,
injecting $C_{2}$ and $CO_{2}$ almost recovered completely the medium’s damage
due to condensate banking. As for the injection of $C_{1}$, $N_{2}$ and
produced gas, very high values of gas relative permeability, superior to
$k_{rg}=0.8$, were achieved at $P\geq 21MPa$. Even for the two lowest tested
pressures, in which cases the liquid saturation could not be reduced to values
below $20\%$, the gas relative permeabilities underwent at least a tenfold
increase after gas flooding.

(a) Gas relative permeabilities before and after the injection of 25 PV of
$C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ or produced gas.

(b) Gas relative permeabilities before and after the injection of 25 PV of
$C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$ mixed with the produced gas (50% in
moles).

Figure 13: Gas relative permeabilities before and after gas injection.

In Figure 13(b) we can see the gas relative permeabilities before and after
the injection of $C_{1}$, $C_{2}$, $CO_{2}$ or $N_{2}$ mixed with he produced
gas. At high injection pressures, mixtures of $C_{2}$ or $CO_{2}$ with the
produced gas could still revert successfully the damage caused by condensate
accumulation, leading to $k_{rg}\approx 1$. Yet, at the same high pressures,
the injection of mixtures of $C_{1}$ and $N_{2}$ with the produced gas was not
as beneficial to the gas flow. At $P=22MPa$, the $C_{1}$ and $N_{2}$ mixtures
reduced the gas relative permeability from $k_{rg}=0.757$ to $k_{rg}=0.435$
and $k_{rg}=0.0417$, respectively. These negative results were expected, as
the injection of these mixtures increased the liquid saturation in the porous
medium at the highest tested pressure (seen in Fig. 10b and d), due to their
tendency to prompt an early condensation when mixed with the reservoir fluid
(seen in Fig. 6). At $P=21.5MPa$ and $P=21MPa$, however, the mixture of
$N_{2}$ and produced gas still impacted negatively the gas flow, even though
the condensate saturation in the porous medium was reduced after gas flooding.
We attribute this gas flow hindrance to the effect of $N_{2}$ on the
interfacial tension of the reservoir fluid’s gas and liquid phases, explained
next.

Injecting gas on the reservoir after condensate accumulation affects not only
the phases bulk properties (e.g. saturation, as seen in section 3.4, and
composition, section 3.5), but also interfacial properties. Considering the
potential impacts on relative permeability curves [33, 34, 10], these effects
should be taken into account when choosing a gas mixture composition for
injection in gas-condensate reservoirs. Figure 14 presents the interfacial
tension for the mixtures containing $50\%$ $C_{1}$, $C_{2}$, $CO_{2}$ and
$N_{2}$ and $50\%$ the reservoir fluid, in moles, calculated with the
correlation proposed by Weinaug and Katz [30]. The results obtained with the
correlation indicate that, for all evaluated pressures, mixing $C_{1}$ and
$N2$ with the reservoir fluid leads to an increase in the interfacial tension
between gas and condensate. $C_{2}$, on the other hand, reduces the
interfacial tension, while $CO_{2}$ has little effect on it.

Figure 14: Calculated interfacial tension for mixtures of the reservoir fluid
with $C_{1}$, $C_{1}$, $CO_{2}$ or $N_{2}$

Therefore, the significant increase in interfacial tension observed with the
addition of nitrogen to the reservoir fluid mixture can strengthen the
mechanism behind condensate blockage, overriding the positive effects on flow
achieved with the reduction of liquid saturation. These results are
particularly relevant considering that most field scale numerical studies
neglect the possible effects of the injected gas composition on the micro-
scale gas-condensate displacement. With this oversight, estimated recoveries
obtained with $C_{1}$ and, especially, $N_{2}$ injection in gas-condensate
reservoirs [16, 17, 18] might have been overestimated. For this reason,
further studies concerning the effects of interfacial tension alteration
during gas injection for condensate enhanced recovery are highly recommended.

##  4 Discussion

The presented results emphasize the relevance of the injected gas composition
and depletion stage in the reservoir for the success of gas injection as a
condensate enhanced method. The findings suggest that $CO_{2}$ and $C_{2}$ are
promising candidates for the studied method. Both gases displayed the overall
highest liquid reduction potential, recovered efficiently heavy components and
induced improvement in the gas-condensate relative permeability. Therefore,
the use of $CO_{2}$ and $C_{2}$, either pure or mixed with another components
due to cost and availability restrictions, should be further investigated for
applications in partial pressure maintenance or huff-n-puff injection schemes.
Pure $C_{1}$ injection also presented rather positive results and should be
considered for injection aiming full or partial pressure maintenance.
Relatively slow re-vaporization rate of condensate and potential increase in
both interfacial tension and dew point pressure could, however, be a concern
for the applications of $C_{1}$ flooding. Produced gas flooding led to similar
liquid saturation reductions and gas relative permeability improvements as
observed for the $C_{1}$ injection, but could not efficiently recover trapped
heavy components from the porous medium, being considered, therefore, a
slightly inferior option for condensate recovery. As for $N_{2}$ injection,
given the feeble condensate re-vaporization tendency when mixed with the
produced gas, low ability to recover heavy components and significant increase
in interfacial tension, it should be carefully evaluated before being
considered as a candidate for gas-condensate enhanced recovery. The relatively
low injection costs [16] might be canceled out by early condensation, poor
two-phase flow performance and entrapment of heavy components in the porous
medium. The level of depletion in the reservoir previous to the gas injection
also affected notably the results. For the evaluated gas-condensate fluid,
performing gas injection below the pressure associated with the maximum liquid
dropout in the reservoir could not enhance significantly condensate recovery,
for any tested gas composition. For this reason, the implementation of the EOR
method in reservoirs at advanced stages of production might result in limited
production improvement.

The generalizability of the presented results is, however, bounded by some
modeling limitations and particularities of the studied cases. First, the
observed recoveries stemmed from mixing injected and reservoir fluids may
depend upon the composition chosen to represent gas-condensate fluids. The
extension of the analyses to different gas-condensate compositions could lead
to variation in the results and, therefore, should be carried out.
Additionally, porous media morphology can impact the performance of gas
injection for condensate recovery, especially at low pressures, where the
injected gas miscibility is limited. More heterogeneous media than the
evaluated sandstone, for instance, may provide preferential paths for gas flow
that deter the method’s success, while very regular media could lead to
increased recoveries. Finally, further studies should be performed to
investigate the effects of gas injection rate on condensate recovery. As for
the model limitations, the used compositional pore-network model does not take
into account molecular diffusion. It has been pointed out in the literature
[17, 35] that the effects of molecular diffusion during gas injection in gas-
condensate reservoirs are not significant. Yet, for the analyses at low
pressures, in which condensate recovery was low, implementing molecular
diffusion in the model could produce more favorable results. Also, the
implemented compositional formulation assumes thermodynamic equilibrium in
each pore at every simulation time step. For this reason, non-equilibrium
effects in the gas-condensate flow are not taken into consideration. Further
investigation should be carried out on gas injection parameters, e.g. high
flow rates, that may call for the incorporation of non-equilibrium effects in
the model.

##  5 Conclusions

The present study provided a pore-scale evaluation of gas injection as a
condensate enhanced recovery method. By using a compositional pore-network
model, the effects of injected composition and pressure could be investigated
by quantifying liquid saturation reduction, heavy component recovery and the
impact on gas relative permeabilities. With these criteria, the performances
of injecting $C_{1}$, $C_{2}$, $CO_{2}$, $N_{2}$, produced gas or their
mixtures to recover the condensate accumulated in a sandstone based pore-
network at different depletion stages were compared.

The results indicated that gas injection can produce substantial flow
improvement in damaged gas-condensate reservoirs, given that the injected gas
composition and pressure for the treatment are adequately determined.
Injecting $C_{2}$ and $CO_{2}$ produced the most favorable results among the
tested gases, while $C_{1}$ and produced gas exhibited moderate positive
effects on flow and $N_{2}$ displayed overall a less beneficial performance
for condensate recovery.

As a preliminary investigation of gas injection in gas-condensate reservoirs
at the pore-scale, this work exposed how changing the phases bulk and
interfacial properties can impact significantly the two-phase flow during
condensate recovery and, consequently, the method’s efficiency. Interfacial
tension effects, often overlooked in the field literature, showed particularly
relevance in gas flow performance after gas injection. By neglecting the
potential negative effects of increase in interfacial tension during the
injection of $C_{1}$ and $N_{2}$, condensate recovery predictions obtained
with reservoir-scale models may be overestimated. As future work, a systematic
pore-scale investigation of the parameters affecting condensate recovery with
gas injection, as flow rate, gas-condensate fluid composition and porous media
heterogeneity should be carried out. With that, data for up-scaling pore-scale
effects pertaining to particular gas injection scenarios could be generated
for reservoir-scale modeling. This could lead to more realistic recovery
estimations and benefit gas-condensate fields development planning.

## References

  * Afidick et al. [1994] D. Afidick, N. Kaczorowski, S. Bette, et al.,  Production performance of a retrograde gas reservoir: a case study of the arun field,  in: SPE Asia pacific oil and gas conference, Society of Petroleum Engineers, 1994.
  * Vo et al. [2015] H. X. Vo, R. N. Horne, et al.,  Experimental study of composition variation during flow of gas-condensate,  in: SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, 2015\.
  * Sayed et al. [2016] M. A. Sayed, G. A. Al-Muntasheri, et al.,  Mitigation of the effects of condensate banking: a critical review,  SPE Production & Operations 31 (2016) 85–102.
  * Ganie et al. [2019] K. Ganie, A. K. Idris, D. F. Mohshim, W. R. W. Sulaiman, I. M. Saaid, A. A. Malik,  A review on the wettability alteration mechanism in condensate banking removal,  Journal of Petroleum Science and Engineering 183 (2019) 106431.
  * Hassan et al. [2019] A. Hassan, M. Mahmoud, A. Al-Majed, M. B. Alawi, S. Elkatatny, M. BaTaweel, A. Al-Nakhli,  Gas condensate treatment: A critical review of materials, methods, field applications, and new solutions,  Journal of Petroleum Science and Engineering 177 (2019) 602–613.
  * Luo et al. [2001] K. Luo, S. Li, X. Zheng, G. Chen, Z. Dai, N. Liu, et al.,  Experimental investigation into revaporization of retrograde condensate by lean gas injection,  in: SPE Asia Pacific Oil and Gas Conference and Exhibition, Society of Petroleum Engineers, 2001\.
  * Ahmed et al. [1998] T. Ahmed, J. Evans, R. Kwan, T. Vivian, et al.,  Wellbore liquid blockage in gas-condensate reservoirs,  in: SPE eastern regional meeting, Society of Petroleum Engineers, 1998.
  * Al-Anazi et al. [2004] H. A. Al-Anazi, M. M. Sharma, G. A. Pope, et al.,  Revaporization of condensate with methane flood,  in: SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, 2004\.
  * Gachuz-Muro et al. [2011] H. Gachuz-Muro, B. E. Gonzalez Valtierra, E. E. Luna, B. Aguilar Lopez, et al.,  Laboratory tests with co2, n2 and lean natural gas in a naturally fractured gas-condensate reservoir under hp/ht conditions,  in: SPE Enhanced Oil Recovery Conference, Society of Petroleum Engineers, 2011.
  * Al-Abri et al. [2012] A. Al-Abri, H. Sidiq, R. Amin,  Mobility ratio, relative permeability and sweep efficiency of supercritical co2 and methane injection to enhance natural gas and condensate recovery: Coreflooding experimentation,  Journal of Natural Gas Science and Engineering 9 (2012) 166–171.
  * Sheng [2015] J. J. Sheng,  Increase liquid oil production by huff-n-puff of produced gas in shale gas condensate reservoirs,  Journal of Unconventional Oil and Gas Resources 11 (2015) 19–26.
  * Meng et al. [2016] X. Meng, J. J. Sheng, et al.,  Experimental study on revaporization mechanism of huff-n-puff gas injection to enhance condensate recovery in shale gas condensate reservoirs,  in: SPE improved oil recovery conference, Society of Petroleum Engineers, 2016.
  * Sharma et al. [2018] S. Sharma, J. J. Sheng, Z. Shen,  A comparative experimental study of huff-n-puff gas injection and surfactant treatment in shale gas-condensate cores,  Energy & fuels 32 (2018) 9121–9131.
  * Meng et al. [2018] X. Meng, Z. Meng, J. Ma, T. Wang, et al.,  Performance evaluation of co 2 huff-n-puff gas injection in shale gas condensate reservoirs,  Energies 12 (2018) 1–18.
  * Marokane et al. [2002] D. Marokane, A. Logmo-Ngog, R. Sarkar, et al.,  Applicability of timely gas injection in gas condensate fields to improve well productivity,  in: SPE/DOE Improved Oil Recovery Symposium, Society of Petroleum Engineers, 2002.
  * Linderman et al. [2008] J. T. Linderman, F. S. Al-Jenaibi, S. G. Ghori, K. Putney, J. Lawrence, M. Gallat, K. Hohensee, et al.,  Feasibility study of substituting nitrogen for hydrocarbon in a gas recycle condensate reservoir,  in: Abu Dhabi International Petroleum Exhibition and Conference, Society of Petroleum Engineers, 2008\.
  * Taheri et al. [2013] A. Taheri, L. Hoier, O. Torsaeter, et al.,  Miscible and immiscible gas injection for enhancing of condensate recovery in fractured gas condensate reservoirs,  in: EAGE Annual Conference & Exhibition incorporating SPE Europec, Society of Petroleum Engineers, 2013.
  * Fath et al. [2016] A. H. Fath, A.-R. Pouranfard, R. Parandvar, S. Pourhadi,  An investigation of different gas injection scenarios as enhanced condensate recovery method in a naturally fractured gas-condensate reservoir,  Petroleum Science and Technology 34 (2016) 295–301.
  * Sharma and Sheng [2017] S. Sharma, J. J. Sheng,  A comparative study of huff-n-puff gas and solvent injection in a shale gas condensate core,  Journal of Natural Gas Science and Engineering 38 (2017) 549–565.
  * Zhang et al. [2020] A. Zhang, Z. Fan, L. Zhao,  An investigation on phase behaviors and displacement mechanisms of gas injection in gas condensate reservoir,  Fuel 268 (2020) 117373\.
  * Reis and Carvalho [2021] P. Reis, M. Carvalho,  Pore-scale compositional modeling of gas-condensate flow: Effects of interfacial tension and flow velocity on relative permeability,  Journal of Petroleum Science and Engineering (2021) 108454.
  * Reis and Carvalho [2020] P. K. Reis, M. S. Carvalho,  Pore-scale analysis of condensate blockage mitigation by wettability alteration,  Energies 13 (2020) 4673.
  * Santos and Carvalho [2020] M. Santos, M. Carvalho,  Pore network model for retrograde gas flow in porous media,  Journal of Petroleum Science and Engineering 185 (2020) 106635.
  * CoŞkuner [1997] G. CoŞkuner,  Microvisual study of multiphase gas condensate flow in porous media,  Transport in Porous Media 28 (1997) 1–18.
  * Al Harrasia et al. [2009] M. Al Harrasia, C. Grattonib, Q. Fishera, S. Al-Hinaic,  Condensate displacement mechanisms in low permeability rocks,  Proceedings, Society of Core Analysts (2009) 15.
  * Jamiolahmady et al. [2003] M. Jamiolahmady, A. Danesh, D. Tehrani, D. Duncan,  Positive effect of flow velocity on gas–condensate relative permeability: network modelling and comparison with experimental results,  Transport in Porous Media 52 (2003) 159–183.
  * Al-Gharbi and Blunt [2005] M. S. Al-Gharbi, M. J. Blunt,  Dynamic network modeling of two-phase drainage in porous media,  Physical Review E 71 (2005) 016308.
  * Sochi [2013] T. Sochi,  Newtonian flow in converging-diverging capillaries,  International Journal of Modeling, Simulation, and Scientific Computing 4 (2013) 1350011\.
  * Beresnev et al. [2009] I. A. Beresnev, W. Li, R. D. Vigil,  Condition for break-up of non-wetting fluids in sinusoidally constricted capillary channels,  Transport in porous media 80 (2009) 581–604.
  * Weinaug and Katz [1943] C. F. Weinaug, D. L. Katz,  Surface tensions of methane-propane mixtures,  Industrial & Engineering Chemistry 35 (1943) 239–246. doi:10.1021/ie50398a028.
  * Collins et al. [1992] D. Collins, L. Nghiem, Y. Li, J. Grabonstotter, et al.,  An efficient approach to adaptive-implicit compositional simulation with an equation of state,  SPE reservoir engineering 7 (1992) 259–264.
  * Dong [2008] H. Dong, Micro-CT imaging and pore network extraction, Ph.D. thesis, Department of Earth Science and Engineering, Imperial College London, 2008.
  * Al-Abri and Amin [2010] A. Al-Abri, R. Amin,  Phase behaviour, fluid properties and recovery efficiency of immiscible and miscible condensate displacements by scco 2 injection: experimental investigation,  Transport in porous media 85 (2010) 743–756.
  * Al-Abri and Amin [2011] A. Al-Abri, R. Amin,  Effect of hydrocarbon and non-hydrocarbon gas injection on the interfacial tension of a gas condensate system,  Chemical engineering & technology 34 (2011) 127–133.
  * Shtepani et al. [2006] E. Shtepani, et al.,  Co2 sequestration in depleted gas/condensate reservoirs,  in: SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, 2006\.